# TwoPhaseEulerFOAM application

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

 November 16, 2006, 07:25 Description: The application #1 Senior Member   Rasmus Hemph Join Date: Mar 2009 Location: Sweden Posts: 108 Rep Power: 13 Description: The application gives erroneous solution for the Ua-field when the particle-particle force is included (g0 > 0). If settling particles are simulated, the velocity of the particles does not approach zero as they reach the bottom of the domain. Solution. The ppMagf-term in alphaEqn.H: should be divided by alphaf on row 26. Old: ppMagf = g0*rUaAf*min(exp(preAlphaExp*(alphaf - alphaMax)), expMax); New: ppMagf = 1.0/(rhoa*max(alphaf,SMALL))*g0*rUaAf*min(exp(preAlpha Exp*(alphaf - alphaMax)), expMax); Also row 32 in alphaEqn.H needs to be changed. Old: alphaEqn -= fvm::laplacian(ppMagf, alpha); New: alphaEqn -= fvm::laplacian(alphaf*ppMagf, alpha); Solver/Application: twoPhaseEulerFoam Source file: alphaEqn.H Testcase: Platform: All Version: All Notes: The dimension of g0 is in Pascal in most references. In twoPhaseEulerFoam the dimensions of g0 is Pa/(kg/m3). To aid comparison between models, the ppMagf term should be divided by rho, as ppMagf = 1.0/(rhoa*max(alphaf,SMALL))*g0*rUaAf*min(exp(preAlpha Exp*(alphaf - alphaMax)), expMax); with a corresonding change to the dimensions of g0 in constant/ppProperties

 November 16, 2006, 07:33 Not a bug, they were left out #2 Super Moderator     Niklas Nordin Join Date: Mar 2009 Location: Stockholm, Sweden Posts: 694 Rep Power: 25 Not a bug, they were left out deliberately. ppMagf is only important around and below alphaMax and the g0 constant can be modified to account for the division/multiplication of alpha. Niklas

 November 16, 2006, 12:09 Sorry, but why is this term le #3 Senior Member   Rasmus Hemph Join Date: Mar 2009 Location: Sweden Posts: 108 Rep Power: 13 Sorry, but why is this term left out? According to the equations they should be there. With the proposed fix included I get expected results. Without it, physics are wrong and results suspicious. Try the supplied test case for an example. I also disagree about ppMagf. It is certainly important over a wide range of volume fractions in packing of for instance a soft material. //Rasmus

 November 16, 2006, 13:18 To expand the discussion. Alth #4 Senior Member   Rasmus Hemph Join Date: Mar 2009 Location: Sweden Posts: 108 Rep Power: 13 To expand the discussion. Although you could exclude the division of alpha in the ppMagf-calculation (for hard materials), the "bug", I believe, is in the laplacian-term. There ppMagf needs to be multiplied by alpha to give correct results.

 November 22, 2006, 15:45 Does the term in pEqn.H stay u #5 Senior Member   Alberto Passalacqua Join Date: Mar 2009 Location: Ames, Iowa, United States Posts: 1,912 Rep Power: 32 Does the term in pEqn.H stay unchanged? Regards, Alberto __________________ Alberto Passalacqua GeekoCFD - A free distribution based on openSUSE 64 bit with CFD tools, including OpenFOAM. Available as in both physical and virtual formats (current status: http://albertopassalacqua.com/?p=1541) OpenQBMM - An open-source implementation of quadrature-based moment methods. To obtain more accurate answers, please specify the version of OpenFOAM you are using.

 January 19, 2007, 07:01 Sorry, missed your message Alb #6 Senior Member   Rasmus Hemph Join Date: Mar 2009 Location: Sweden Posts: 108 Rep Power: 13 Sorry, missed your message Alberto. I believe that pEqn.H remains unchanged. The change in the code applies to the ppMagf-term, which is calculated in alphaEqn.H.

 January 19, 2007, 07:43 OK. Thanks for your confirmati #7 Senior Member   Alberto Passalacqua Join Date: Mar 2009 Location: Ames, Iowa, United States Posts: 1,912 Rep Power: 32 OK. Thanks for your confirmation. :-) __________________ Alberto Passalacqua GeekoCFD - A free distribution based on openSUSE 64 bit with CFD tools, including OpenFOAM. Available as in both physical and virtual formats (current status: http://albertopassalacqua.com/?p=1541) OpenQBMM - An open-source implementation of quadrature-based moment methods. To obtain more accurate answers, please specify the version of OpenFOAM you are using.

 May 26, 2008, 03:12 It seems to me the equations t #8 Member   Juho Peltola Join Date: Mar 2009 Location: Finland Posts: 89 Rep Power: 13 It seems to me the equations that Rasmus reports as erroneous are still used in the twoPhaseEulerFoam in the 1.4.1 release. Has there been chenges to the solver after the bug report? Is the bug report still valid?

 May 26, 2008, 05:59 Are we all agreed that the pro #9 Senior Member   Join Date: Mar 2009 Posts: 854 Rep Power: 18 Are we all agreed that the proposed change is correct? If so I will put it in 1.5. Do we have a reference or derivation which covers this? Henry

 May 26, 2008, 23:09 Expressions for the solids nor #10 Senior Member   Alberto Passalacqua Join Date: Mar 2009 Location: Ames, Iowa, United States Posts: 1,912 Rep Power: 32 Expressions for the solids normal stress modulus are empirical. You can find a summary and a list of references here http://www.scielo.br/scielo.php?script=sci_arttext&pid=S0100-73862001000200006&l ng=&nrm=iso&tlng=#tab01 together with a short explanation of how they're used inside the particle phase momentum equation. I have almost all the original references, if necessary. Alberto __________________ Alberto Passalacqua GeekoCFD - A free distribution based on openSUSE 64 bit with CFD tools, including OpenFOAM. Available as in both physical and virtual formats (current status: http://albertopassalacqua.com/?p=1541) OpenQBMM - An open-source implementation of quadrature-based moment methods. To obtain more accurate answers, please specify the version of OpenFOAM you are using.

 May 27, 2008, 03:04 From these do you agree with R #11 Senior Member   Join Date: Mar 2009 Posts: 854 Rep Power: 18 From these do you agree with Rasmus' change? Henry

 May 27, 2008, 07:26 Which change is under discussi #12 Senior Member   Rasmus Hemph Join Date: Mar 2009 Location: Sweden Posts: 108 Rep Power: 13 Which change is under discussion here? I have reported five different issues. I do not have the 1.41 source heres, so I can not tell which ones have been corrected. First, there is a missing alpha in the laplacian-expression of the alphaEqn. The inclusion of alpha can be shown from the equations and was found to be important to the physics. Second, ppMagf should be divided by alpha*rho. Niklas left this out and left it to the user to modify the expression for G. However, for soft materials, there is a particle-particle force over a wide range of volume fractions, which falsifies the assumption of alpha at packing being a constant. (Also, as alpha has no dimensions, the user would need to read through the source code to understand that she needed to multiply her constant by alpha..) Thirdly, there is (was?) a grave bug in the interpolation of particle-particle forces to the faces. It has been reported here http://www.cfd-online.com/OpenFOAM_D...tml?1173086743 Fourthly, there is a problem with certain boundary conditions for particles. It was reported here http://www.cfd-online.com/OpenFOAM_D...tml?1170065229 Rasmus

 May 27, 2008, 07:50 Hi Rasmus, I assumed that i #13 Senior Member   Join Date: Mar 2009 Posts: 854 Rep Power: 18 Hi Rasmus, I assumed that it was the change at the top of this list being discussed. I am happy to include any changes you think appropriate but as I am not currently developing or running this solver I can't test or verify it. If you could send me the latest code I would be happy to merge it into the 1.5 release if you would like me to. Thanks Henry

 May 27, 2008, 07:51 As far as I can see at least i #14 Member   Juho Peltola Join Date: Mar 2009 Location: Finland Posts: 89 Rep Power: 13 As far as I can see at least issues 1-3 exist in the 1.4.1.

 May 27, 2008, 09:00 Hi Henry and Juho, I think tha #15 Senior Member   Rasmus Hemph Join Date: Mar 2009 Location: Sweden Posts: 108 Rep Power: 13 Hi Henry and Juho, I think that it would be good if the changes came into 1.5. They have been found to increase stability quite a lot. The code I have been running is largely modified for a possible article and is hard to submit for inclusion as is. However, here is a suggestion of a fix for the bug in interpolation (third bug in my post) and particle forces calculation (first and second). The changes needs to go into alphaEqn.H. Warning, I do not have access to a linux computer at the moment, so there is possibly some syntatic errors. Note that the dimensions of g0 need to change in the tutorial test case as well. I would be very happy if someone was willing to test it on the tutorials. Old code: if (g0.value() > 0.0) { surfaceScalarField alphaf = fvc::interpolate(alpha); // Correct the particle-particle force magnitude ppMagf = g0*rUaAf*min(exp(preAlphaExp*(alphaf - alphaMax)), expMax); alphaEqn -= fvm::laplacian(ppMagf, alpha); } New code: if (g0.value() > 0.0) { volScalarField G = (1.0/(rhoa*(alpha+scalar(0.0001))))*g0*min(exp(preAlpha Exp*(alpha - alphaMax)), expMax)); ppMagf = rUaAf*fvc::interpolate(G); surfaceScalarField alphaf = fvc::interpolate(alpha); alphaEqn -= fvm::laplacian(((alphaf+scalar(0.0001))*ppMagf), alpha); } Rasmus

 May 27, 2008, 09:17 Thank you. It compiled fine. R #16 Member   Juho Peltola Join Date: Mar 2009 Location: Finland Posts: 89 Rep Power: 13 Thank you. It compiled fine. Requires an ugly entry to fvScemes. I can run tests on the tutorials. At least it seems to run fine on another test case I have.

 May 27, 2008, 09:55 Good. Thanks for trying it out #17 Senior Member   Rasmus Hemph Join Date: Mar 2009 Location: Sweden Posts: 108 Rep Power: 13 Good. Thanks for trying it out. To get rid of the ugly entry, I believe you can give the term a name, which is then called in fvSchemes. In Ueqn.H for example, it is fvm::div(phia, Ua, "div(phia,Ua)").

 May 27, 2008, 10:26 Is there a reason not to use #18 Member   Juho Peltola Join Date: Mar 2009 Location: Finland Posts: 89 Rep Power: 13 Is there a reason not to use max(alphaf,SMALL) max(alpha,SMALL) instead of (alphaf+scalar(0.0001)) (alpha+scalar(0.0001)) ? The latter seems a lot less stable.

 May 27, 2008, 11:09 False alarm. The instability w #19 Member   Juho Peltola Join Date: Mar 2009 Location: Finland Posts: 89 Rep Power: 13 False alarm. The instability wasn't related to the max(alphaf,SMALL) vs. (alphaf+scalar(0.0001)). But is there a reason to choose on over the other? Does max() compare and replace values of each cell in of any field? A difference between the approaches is of course that the max() filters out even negative values of magnitude larger than 0.0001 and for debugging reasons it migth be nice that such values cause an error.

 May 27, 2008, 14:32 The reason for going with 0.00 #20 Senior Member   Rasmus Hemph Join Date: Mar 2009 Location: Sweden Posts: 108 Rep Power: 13 The reason for going with 0.0001 here was an attempt to make the division a little better numerically, as SMALL is 1e-15. The same term is multiplied in the laplacian equation later. I do not know if it made any practical difference, but I agree that is a bit ugly with a free scalar like that. If I recall correctly, you can not do max(alphaf, SMALL) unless both terms are of the same dimension (because of the replacement), so you need to create a new dimensioned scalar, but principally it should be ok.

 Posting Rules You may not post new threads You may not post replies You may not post attachments You may not edit your posts BB code is On Smilies are On [IMG] code is On HTML code is OffTrackbacks are Off Pingbacks are On Refbacks are On Forum Rules

 Similar Threads Thread Thread Starter Forum Replies Last Post sara OpenFOAM Running, Solving & CFD 2 November 6, 2008 19:26 alberto OpenFOAM Bugs 2 May 20, 2008 21:25 alondono OpenFOAM Bugs 1 February 19, 2008 20:01 hemph OpenFOAM Bugs 0 November 16, 2006 07:27 newbee OpenFOAM 0 March 27, 2006 08:41

All times are GMT -4. The time now is 01:45.