# SaffmanMeiLiftForce

 Register Blogs Members List Search Today's Posts Mark Forums Read

 October 13, 2013, 12:16 SaffmanMeiLiftForce #1 Member   Evangelos Join Date: Sep 2011 Posts: 85 Rep Power: 10 Hi Can anyone provide an example of SaffmanMeiLiftForce dictionary ? thanks a lot !

 October 13, 2013, 15:51 #2 Super Moderator   Bruno Santos Join Date: Mar 2009 Location: Lisbon, Portugal Posts: 10,956 Blog Entries: 43 Rep Power: 121 Hi Evangelos, I had no idea there was such a thing as "SaffmanMeiLiftForce" in OpenFOAM... But if you can provide an example case, something simple based on one of OpenFOAM's tutorials, in which this dictionary can be used, then it should be pretty simple to reverse engineer it! Actually, OpenFOAM helps with this in two ways: The "banana" way: http://openfoamwiki.net/index.php/Op...de/Use_bananas The source code way, namely by looking at the source code and or the code documentation: http://www.openfoam.org/docs/cpp/ Best regards, Bruno __________________ OpenFOAM: FAQ | Getting started Forum: How to get help, to post code/output and forum guide Read this before sending me PM

 October 13, 2013, 17:00 #3 Member   Evangelos Join Date: Sep 2011 Posts: 85 Rep Power: 10 Hello Bruno ! I try to calculate the orbit of a particle when injected into a domain without fluid motion.I use icoUncoupledKinematicFoam but i am a little bit confused! All the solvers take into account dragforce and gravity without calculate the lift from the fluid ! In openfoam 2.2.0 particle lift is available So in the kinematicTransportproperties/particleforces under the drag force and gravity i add particle force and the program show me the valid particle forces There are two " lift methods " TomiyamaLift and SaffmanMeiLiftForce http://www.cfd-online.com/Forums/ope...rce-model.html

 October 14, 2013, 05:06 #4 Member   Evangelos Join Date: Sep 2011 Posts: 85 Rep Power: 10 I need something like that http://www.cfd-online.com/Forums/ope...ictionary.html but for saffmanmeiliftforce!!!

 October 14, 2013, 17:47 #5 Super Moderator   Bruno Santos Join Date: Mar 2009 Location: Lisbon, Portugal Posts: 10,956 Blog Entries: 43 Rep Power: 121 Hi Evangelos, Do you mean "icoUncoupledKinematicParcelFoam"? Do you think that the tutorial "lagrangian/icoUncoupledKinematicParcelFoam/hopper" can be used to test this dictionary you are asking about? If so, I can give it a try next weekend. Best regards, Bruno __________________ OpenFOAM: FAQ | Getting started Forum: How to get help, to post code/output and forum guide Read this before sending me PM

 October 5, 2014, 00:34 #6 Senior Member     Join Date: Jan 2010 Posts: 347 Blog Entries: 2 Rep Power: 13 Hi all, Do you find any dictionary for saffmanMeiLiftForce?

October 5, 2014, 05:20
#7
Super Moderator

Bruno Santos
Join Date: Mar 2009
Location: Lisbon, Portugal
Posts: 10,956
Blog Entries: 43
Rep Power: 121
Greetings Maysam,

Quote:
 Originally Posted by maysmech Do you find any dictionary for saffmanMeiLiftForce?
At least I haven't, because Evangelos didn't answer back

Therefore my questions and offer still stand:
Quote:
 Originally Posted by wyldckat Do you mean "icoUncoupledKinematicParcelFoam"? Do you think that the tutorial "lagrangian/icoUncoupledKinematicParcelFoam/hopper" can be used to test this dictionary you are asking about? If so, I can give it a try next weekend.
Best regards,
Bruno
__________________

 October 5, 2014, 05:33 #8 Senior Member     Join Date: Jan 2010 Posts: 347 Blog Entries: 2 Rep Power: 13 For example for DPMFoam

October 5, 2014, 07:58
#9
Super Moderator

Bruno Santos
Join Date: Mar 2009
Location: Lisbon, Portugal
Posts: 10,956
Blog Entries: 43
Rep Power: 121
OK, so as I mentioned in one of the posts above, having an example case as a basis and using the "banana" trick, here's what I've done (used OpenFOAM 2.3.x for this example):
1. Made a copy of the tutorial "lagrangian/DPMFoam/Goldschmidt" and ran the utility blockMesh.
2. Then I had to look at which file needed to modified. All pointed to the file "constant/kinematicCloudProperties".
3. Now I had to figure out where the Lift comes in. After going around looking at the other tutorials and a few threads, the conclusion was that this block of dictionary is where we can stack up the various forces at work:
Code:
particleForces
{
ErgunWenYuDrag
{
alphac alpha.air;
}
gravity;
}
4. And here's where the "banana" comes in:
Code:
particleForces
{
ErgunWenYuDrag
{
alphac alpha.air;
}
gravity;
banana;
}
5. As soon as I run DPMFoam, it complains about "banana" not being valid:
Code:
--> FOAM FATAL ERROR:
Unknown particle force type banana, constructor not in hash table

Valid particle force types are:

13
(
ErgunWenYuDrag
PlessisMasliyahDrag
SRF
SaffmanMeiLiftForce
TomiyamaLift
WenYuDrag
gravity
nonInertialFrame
nonSphereDrag
paramagnetic
sphereDrag
virtualMass
)
There's our pretty lift model "SaffmanMeiLiftForce".
6. OK, rename "banana" to "SaffmanMeiLiftForce" and added brackets with a "banana" inside :
Code:
particleForces
{
ErgunWenYuDrag
{
alphac alpha.air;
}
gravity;
SaffmanMeiLiftForce
{
banana;
}
}
mmm... the solver crashed. So much for the "banana".
7. OK, if we look at the source code for this model:
Code:
template<class CloudType>
Foam::SaffmanMeiLiftForce<CloudType>::SaffmanMeiLiftForce
(
CloudType& owner,
const fvMesh& mesh,
const dictionary& dict,
const word& forceType
)
:
LiftForce<CloudType>(owner, mesh, dict, forceType)
{}
there is no coefficient necessary.
8. Therefore, this should be enough:
Code:
particleForces
{
ErgunWenYuDrag
{
alphac alpha.air;
}
gravity;
SaffmanMeiLiftForce;
}
9. Er, nope... it complains it should be a dictionary itself. Then back to this:
Code:
particleForces
{
ErgunWenYuDrag
{
alphac alpha.air;
}
gravity;
SaffmanMeiLiftForce
{
}
}
10. So we go back to the crash. The error given is this:
Code:
--> FOAM FATAL ERROR:

request for volVectorField U from objectRegistry region0 failed
available objects of type volVectorField are
1(U.air)
OK, it's looking for the "U" field and it only found "U.air". So, either the model is just wrong, or something else is missing.
11. Looking at the class from which this model derives, the constructor has a vital clue to what's missing:
Code:
template<class CloudType>
Foam::LiftForce<CloudType>::LiftForce
(
CloudType& owner,
const fvMesh& mesh,
const dictionary& dict,
const word& forceType
)
:
ParticleForce<CloudType>(owner, mesh, dict, forceType, true),
UName_(this->coeffs().template lookupOrDefault<word>("U", "U")),
curlUcInterpPtr_(NULL)
{}
It's looking for a new name for "U", something like:
Code:
U U.air;
12. So let's try this:
Code:
particleForces
{
ErgunWenYuDrag
{
alphac alpha.air;
}
gravity;
SaffmanMeiLiftForce
{
U U.air;
}
}
Side note: curiously enough, Maysam got a similar answer a few minutes ago in post #4 at http://www.cfd-online.com/Forums/ope...tml#post512965
13. Running DPMFoam once again gives this error message:
Code:
--> FOAM FATAL IO ERROR:
keyword curlUcDt is undefined in dictionary "/home/ofuser/OpenFOAM/ofuser-2.3.x/run/tutorials/lagrangian/DPMFoam/Goldschmidt/constant/kinematicCloudProperties.solution.interpolationSchemes"

file: /home/ofuser/OpenFOAM/ofuser-2.3.x/run/tutorials/lagrangian/DPMFoam/Goldschmidt/constant/kinematicCloudProperties.solution.interpolationSchemes from line 27 to line 29.
OK, we need to look at this this block of dictionary:
Code:
solution
{
active          true;
coupled         true;
transient       yes;
cellValueSourceCorrection off;

interpolationSchemes
{
rho.air         cell;
U.air           cellPoint;
mu.air          cell;
}

integrationSchemes
{
U               Euler;
}

sourceTerms
{
schemes
{
U semiImplicit 1;
}
}
}
14. Let's try:
Code:
interpolationSchemes
{
rho.air         cell;
U.air           cellPoint;
mu.air          cell;
curlUcDt        banana;
}
The error message stated by DPMFoam:
Code:
--> FOAM FATAL ERROR:
Unknown interpolation type banana for field curlUcDt

Valid interpolation types :

6
(
cell
cellPatchConstrained
cellPoint
cellPointFace
cellPointWallModified
pointMVC
)
15. I have no idea. Let's try the most conservative/complex one "cellPointFace"... nope, it crashes. Let's try all others... nope, all of them crash.
16. The error message of the crash is something like this:
Code:
Solving 3-D cloud kinematicCloud
#0  Foam::error::printStack(Foam::Ostream&) in "/home/ofuser/OpenFOAM/OpenFOAM-2.3.x/platforms/linux64GccDPOpt/lib/libOpenFOAM.so"
#1  Foam::sigFpe::sigHandler(int) in "/home/ofuser/OpenFOAM/OpenFOAM-2.3.x/platforms/linux64GccDPOpt/lib/libOpenFOAM.so"
#2   in "/lib/x86_64-linux-gnu/libc.so.6"
#3  Foam::SaffmanMeiLiftForce<Foam::KinematicCloud<Foam::Cloud<Foam::CollidingParcel<Foam::KinematicParcel<Foam::particle> > > > >::Cl(Foam::CollidingParcel<Foam::KinematicParcel<Foam::particle> > const&, Foam::Vector<double> const&, double, double) const in "/home/ofuser/OpenFOAM/OpenFOAM-2.3.x/platforms/linux64GccDPOpt/lib/liblagrangianIntermediate.so"
#4  Foam::LiftForce<Foam::KinematicCloud<Foam::Cloud<Foam::CollidingParcel<Foam::KinematicParcel<Foam::particle> > > > >::calcCoupled(Foam::CollidingParcel<Foam::KinematicParcel<Foam::particle> > const&, double, double, double, double) const in "/home/ofuser/OpenFOAM/OpenFOAM-2.3.x/platforms/linux64GccDPOpt/lib/liblagrangianIntermediate.so"
#5
in "/home/ofuser/OpenFOAM/OpenFOAM-2.3.x/platforms/linux64GccDPOpt/bin/DPMFoam"
#6
So there was a "sigFpe", hence division by zero or infinite or something like that, in the method "::Cl" (it's in line #3).
17. OK, looking at the source code for this method:
Quote:
 https://github.com/OpenFOAM/OpenFOAM...MeiLiftForce.C Code: template Foam::scalar Foam::SaffmanMeiLiftForce::SaffmanMeiLiftForce::Cl ( const typename CloudType::parcelType& p, const vector& curlUc, const scalar Re, const scalar muc ) const { scalar Rew = p.rhoc()*mag(curlUc)*sqr(p.d())/(muc + ROOTVSMALL); scalar beta = 0.5*(Rew/(Re + ROOTVSMALL)); scalar alpha = 0.3314*sqrt(beta); scalar f = (1.0 - alpha)*exp(-0.1*Re) + alpha; scalar Cld = 0.0; if (Re < 40) { Cld = 6.46*f; } else { Cld = 6.46*0.0524*sqrt(beta*Re); } return 3.0/(mathematical::twoPi*sqrt(Rew))*Cld; }
there are a few possible locations where it might have broken.
18. OK, let's go old-school debugging:
Code:
template<class CloudType>
Foam::scalar Foam::SaffmanMeiLiftForce<CloudType>::SaffmanMeiLiftForce::Cl
(
const typename CloudType::parcelType& p,
const vector& curlUc,
const scalar Re,
const scalar muc
) const
{
scalar Rew = p.rhoc()*mag(curlUc)*sqr(p.d())/(muc + ROOTVSMALL);
Info << "Rew: " << Rew << endl;

scalar beta = 0.5*(Rew/(Re + ROOTVSMALL));
Info << "Re: " << Re << endl;
Info << "beta: " << beta << endl;

scalar alpha = 0.3314*sqrt(beta);
Info << "alpha: " << alpha << endl;

scalar f = (1.0 - alpha)*exp(-0.1*Re) + alpha;
Info << "f: " << f << endl;

scalar Cld = 0.0;
if (Re < 40)
{
Cld = 6.46*f;
}
else
{
Cld = 6.46*0.0524*sqrt(beta*Re);
}

Info << "Cld: " << Cld << endl;

return 3.0/(mathematical::twoPi*sqrt(Rew))*Cld;
}
we print out everything so that we can have a look at each value. It's good to know what's happening under-the-hood, instead of having to guess: http://openfoamwiki.net/index.php/HowTo_debugging
19. And we have to build the modified library "src/lagrangian/intermediate":
Code:
wmake libso \$FOAM_SRC/lagrangian/intermediate
20. Let's have another go with DPMFoam... the output we now get:
Code:
Solving 3-D cloud kinematicCloud
Rew: 0
Re: 71.5702
beta: 0
alpha: 0
f: 0.00077937
Cld: 0
that's waaaay too many zeros... this leads to this expression to go haywire:
Code:
return 3.0/(mathematical::twoPi*sqrt(Rew))*Cld;
Because it becomes:
Code:
3.0/(2*Pi*sqrt(0.0))*0.0
where it divides by zero.
21. OK, we can try and insert a division protection like it's done for the others:
Code:
return 3.0/(mathematical::twoPi*sqrt(Rew+SMALL))*Cld;
22. Build the library again and gave DPMFoam another run... didn't crash, but gave a lot of junk output. So I Ctrl+C to cancel the run and went back to the code and removed the "Info" lines I had added. And built the library again.
23. Another run with DPMFoam and... and it's finally working! As to whether it's giving good or bad values, I have no idea.
In the meantime I'll report this bug, although I'm not sure if the fix I presented is correct or not .

-------------------
edit: Bug reported here: http://www.openfoam.org/mantisbt/view.php?id=1408

Additional note: Letting the case run first without any particles should initialize the "U.air" field, should no longer require you to modify the source code, since the flow field would all not be zero.
-------------------

Best regards,
Bruno

Last edited by wyldckat; October 5, 2014 at 08:16. Reason: see "edit:" and "Additional note:"

 October 5, 2014, 21:25 #10 Senior Member     Join Date: Jan 2010 Posts: 347 Blog Entries: 2 Rep Power: 13 Thank you very much Bruno. I used a lot from your step by step debugging. Especially Banana and Old-School debugging were fantastic. This problem is for other forces too such as pressure gradient force. However thanks for submit of the bug. wyldckat likes this.

 October 18, 2014, 14:38 #11 Super Moderator   Bruno Santos Join Date: Mar 2009 Location: Lisbon, Portugal Posts: 10,956 Blog Entries: 43 Rep Power: 121 FYI: Bug has been fixed: http://www.openfoam.org/mantisbt/view.php?id=1408 The official fix can be seen here: https://github.com/OpenFOAM/OpenFOAM...1d371bd9be3000 Mojtaba.a, lithos and CUnsworth like this.

November 21, 2014, 08:54
#12
New Member

Michael
Join Date: Dec 2011
Location: Geneva
Posts: 28
Rep Power: 10
Quote:
 Originally Posted by wyldckat FYI: Bug has been fixed: http://www.openfoam.org/mantisbt/view.php?id=1408 The official fix can be seen here: https://github.com/OpenFOAM/OpenFOAM...1d371bd9be3000
Dear Bruno, you made my day. Thank you very much for this little detail, explained very understandable. Cheers, Michael

P.S.: I used this with 2.2.x.

Last edited by lithos; November 21, 2014 at 10:11.

 June 27, 2018, 03:05 Brownian motion force not running when using uncoupledKinematicParcelFoam #13 Member   AJAY BHANDARI Join Date: Jul 2015 Location: INDIA Posts: 57 Rep Power: 6 Hi all, I think my post best suits here. Basically i am solving a steady state continuity and darcy equation to get steady state pressure and velocity fields, which go into a transient convection diffusion equation to solve time dependent concentration of a species in a computational domain. Now i want to couple my developed eulerain solver with lagrangian one. For that i am using uncoupledKinematicParcelFoam and coupling it with my eulerian solver. The solver is compiled succesfully. But when i am running the case it is giving following error. Code: Constructing kinematicCloud kinematicCloud --> FOAM Warning : From function entry::getKeyword(keyType&, Istream&) in file db/dictionary/entry/entryIO.C at line 80 Reading /home/ajay/OpenFOAM/ajay-4.1/run/mytutorials/srmparticlenew/constant/kinematicCloudProperties found on line 64 the punctuation token '{' expected either } or EOF Constructing particle forces Selecting particle force sphereDrag Selecting particle force gravity Selecting particle force BrownianMotion --> FOAM FATAL ERROR: Unknown particle force type BrownianMotion, constructor not in hash table Valid particle force types are: 13 ( ErgunWenYuDrag PlessisMasliyahDrag SRF SaffmanMeiLiftForce TomiyamaLift WenYuDrag gravity nonInertialFrame nonSphereDrag paramagnetic pressureGradient sphereDrag virtualMass ) From function static Foam::autoPtr > Foam::ParticleForce::New(CloudType&, const Foam::fvMesh&, const Foam::dictionary&, const Foam::word&) [with CloudType = Foam::KinematicCloud > > >] in file lagrangian/intermediate/lnInclude/ParticleForceNew.C at line 47. FOAM exiting I know the error is coming because of addition of Brownian motion force in the case. Can anybody please help why the brownian motion force is not running for the case or where i am doing the mistake. I have compiled the Lagrangian library properly. But in that there is no Brownian motion force in that. Brownian motion force is situated in different folder named turbulence. I compiled it as well but not running. Any help will be appreciated. Thanks Ajay

 August 22, 2018, 10:30 #14 Super Moderator   Bruno Santos Join Date: Mar 2009 Location: Lisbon, Portugal Posts: 10,956 Blog Entries: 43 Rep Power: 121 Greetings Ajay, Sorry for the late reply, but I only managed to look into this now. Here is a quick notes answer: That force is deployed with this macro: https://github.com/OpenFOAM/OpenFOAM...ulenceForces.H That macro is only used in a few cloud types: https://github.com/OpenFOAM/OpenFOAM...rbulenceForces Essentially if your clouds are not of type spray, reactive nor thermal, then you cannot use that force, given that it depends on thermodynamic characteristics of the fluid(s). From your description, it's not clear if there is even temperature being used in your fluid transport, therefore this force may very likely not be usable in the first place. The warning given above the error message is of some concern, because a misplaced bracket pair "{ }" is enough to break the way the dictionary is read... Best regards, Bruno __________________ OpenFOAM: FAQ | Getting started Forum: How to get help, to post code/output and forum guide Read this before sending me PM