
[Sponsors] 
April 9, 2018, 04:24 
Spatial source term

#1 
Member
Join Date: Dec 2012
Posts: 33
Rep Power: 13 
Hi all,
I've been trying to add a spatial source term to the pimpleFoam application but struggling with the compilation. Hoping someone here can help me out. I need to add a source that ramps from an initial value to a final value. Here's what I did so far: Defining information I will read from transportProperties: Code:
// Read in ramp data from transport properties Info<< "Reading ramp properties\n" << endl; IOdictionary transportProperties ( IOobject ( "transportProperties", runTime.constant(), mesh, IOobject::MUST_READ_IF_MODIFIED, IOobject::NO_WRITE ) ); dimensionedVector ramp_beg ( transportProperties.lookup("ramp_beg") ); dimensionedVector ramp_end ( transportProperties.lookup("ramp_end") ); dimensionedVector ramp_loc ( transportProperties.lookup("ramp_loc") ); // End of ramp reading Code:
volVectorField argument = (mesh.C()  ramp_loc); volVectorField force ( ramp_beg + 0.5*(1.0+Foam::tanh(argument.value()))*(ramp_end  ramp_beg) ); // Solve the Momentum equation tmp<fvVectorMatrix> UEqn ( fvm::ddt(U) + fvm::div(phi, U) + turbulence>divDevReff(U) + force // This is the ramped up pressure gradient: +dp/dx == fvOptions(U) );
Any suggestions how to fix it? Thanks! 

April 9, 2018, 06:13 

#2 
Senior Member
Cyprien
Join Date: Feb 2010
Location: Stanford University
Posts: 299
Rep Power: 18 
argument.value() only works if argument is a dimensionedVector.
Can you compile if you remove .value() ? 

April 9, 2018, 07:15 

#3  
Member
Join Date: Dec 2012
Posts: 33
Rep Power: 13 
Quote:
Code:
UEqn.H:3:40: error: conversion from ‘Foam::tmp<Foam::GeometricField<Foam::Vector<double>, Foam::fvPatchField, Foam::volMesh> >’ to nonscalar type ‘Foam::dimensionedVector {aka Foam::dimensioned<Foam::Vector<double> >}’ requested dimensionedVector argument = (mesh.C()  ramp_loc); ^ UEqn.H:6:48: error: no matching function for call to ‘tanh(Foam::Vector<double>&)’ ramp_beg + 0.5*(1.0+Foam::tanh(argument.value()))*(ramp_end  ramp_beg) ^ Code:
UEqn.H:6:40: error: no matching function for call to ‘tanh(Foam::volVectorField&)’ ramp_beg + 0.5*(1.0+Foam::tanh(argument))*(ramp_end  ramp_beg) ^ 

April 9, 2018, 08:46 

#4 
Senior Member
Cyprien
Join Date: Feb 2010
Location: Stanford University
Posts: 299
Rep Power: 18 
of course, because "(mesh.C()  ramp_loc) " returns a volVectorField and not a dimensionedVector.
The other problem in your code is the function Foam::tanh( ) that is not defined for a volVectorField (what is an hyperbolic tangent of a vector ?). Try Code:
volVectorField argument = (mesh.C()  ramp_loc); volVectorField force ( "force", ramp_beg + 0.5*(1.0+Foam::tanh(mag(argument)))*(ramp_end  ramp_beg) ); 

April 9, 2018, 09:58 

#5  
Member
Join Date: Dec 2012
Posts: 33
Rep Power: 13 
Quote:
Code:
> FOAM FATAL ERROR: Argument of trancendental function not dimensionless From function trans(const dimensionSet&) in file dimensionSet/dimensionSet.C at line 430. 

April 9, 2018, 10:02 

#6 
Senior Member
Cyprien
Join Date: Feb 2010
Location: Stanford University
Posts: 299
Rep Power: 18 
it probably means that your formula is not correct. "argument" should be dimensionless to be apply to tanh().


April 10, 2018, 21:27 

#7  
Member
Join Date: Dec 2012
Posts: 33
Rep Power: 13 
Quote:
I was thinking if it is possible to absorb this pressure gradient into the fvc::grad(p) term to calculate for an effective p. The current approach seems like the obtained pressure does not reflect the actual pressure in the system. Would it also affect what the U values would be? 

April 11, 2018, 05:42 

#8 
Senior Member
Cyprien
Join Date: Feb 2010
Location: Stanford University
Posts: 299
Rep Power: 18 
If you remove fvc::grad(p) from the momentum equation, you can not define an equation for the pressure that satisfies the conservation of the mass (see the description of the PISO algorithm for the pressurevelocity coupling). Basically, in what you did, your velocity is not divergence free.
If you want a ramp of pressure, I will suggest two strategies:  first, you apply a pressure difference at the boundaries of your domain and you apply a ramp to one end of the domain (this is what people do usually, and also what make sense if you want to compare with some lab experiments)  Or, you decompose the pressure gradient into a perturbation and a mean value. The mean value is your source term and the pressure gradient correspond to gradient of the perturbation. In this case, you keep fvc::grad(p) (corresponding to the gradient of the perturbation) in the momentum equation, you add your body force (the mean pressure gradient that applies to every cells of the domain) and you solve the pressure equation (you don't have to modify anything). This strategy is usually used for simulations with periodic boundary conditions. 

April 12, 2018, 00:47 

#9  
Member
Join Date: Dec 2012
Posts: 33
Rep Power: 13 
Quote:
Quote:
For the second point, is it the same as I have been doing i.e. adding force to UEqn and solving with fvc::grad(p) term using the PISO algorithm? Or you are suggesting to add the force on the right side of the solve() such that solve (UEqn == fvc::grad(p) + force) and modify the pEqn to be consistent, like the treatment of gravity in buoyantPimpleFoam? 

April 12, 2018, 04:43 

#10 
Senior Member
Cyprien
Join Date: Feb 2010
Location: Stanford University
Posts: 299
Rep Power: 18 
Your right for the first point. You consider the ramp only at the inlet boundary.
For the second point, it is important to combine this approach with cyclic boundary condition (so, may be it is not relevant in your case, meaning that you should probably deal with the first option). Here, again, there are two ways for adding a source term in the momentum equation:  you add your source term when you build the matrix: Code:
fvVectorMatrix UEqn ( fvm::ddt(U) + ... + force ); solve (UEqn == fvc::grad(p));  you write Code:
fvVectorMatrix UEqn ( fvm::ddt(U) + ... ); solve (UEqn == fvc::grad(p)+force); The first option is probably the simplest. Cheers, 

April 23, 2021, 15:06 

#11 
New Member
irwin
Join Date: Apr 2014
Location: New Delhi, India
Posts: 19
Rep Power: 12 
I am starting with 2d pipe flow simulation . I have added source(pressure gradient) with cyclic inlet and outlet in icoFoam solver . I checked it with adding same pressure gradient without cyclic boundary condition at the inlet . Both of them gives completely different results . In fact solver with added source term (pressure gradient ) gives completely different results compared to analytical solution. Can you please suggest something in case I am doing something wrong .
thanks 

Tags 
pimplefoam, source term, spatial 
Thread Tools  Search this Thread 
Display Modes  


Similar Threads  
Thread  Thread Starter  Forum  Replies  Last Post 
Source Term due to evaporation in energy transport equation  styleworker  OpenFOAM Programming & Development  3  September 7, 2022 04:09 
[OpenFOAM.org] Patches to compile OpenFOAM 2.2 on Mac OS X  gschaider  OpenFOAM Installation  136  October 10, 2017 18:25 
polynomial BC  srv537  OpenFOAM PreProcessing  4  December 3, 2016 10:07 
"parabolicVelocity" in OpenFoam 2.1.0 ?  sawyer86  OpenFOAM Running, Solving & CFD  21  February 7, 2012 12:44 
UDFs for Scalar Eqn  Fluid/Solid HT  Greg Perkins  FLUENT  0  October 11, 2000 04:43 