wrong(?) equation for Capillary Force in gran cohesion easo/capillary/viscous

giuraso's picture
Submitted by giuraso on Wed, 02/17/2016 - 17:18

Hi,
I have checked the capillary force model in cohesion easo/capillary/viscous.
I think the equation you give as reference is

B = (-0.34*LN(volBondScaled)-0.96)*contactAngleEff*contactAngleEff - 0.019*LN(volBondScaled) + 0.48
C = 0.0042LN(volBondScaled)+0.078;

but in the code you use LOG instead of LN, and converted LN into LOG.
So you divided the coefficient 0.34, 0.0019 and 0.0042 by 2.303.

B = (-0.148*log(volBondScaled)-0.96)*contactAngleEff*contactAngleEff - 0.0082*log(volBondScaled) + 0.48
C = 0.0018*log(volBondScaled)+0.078;

But you need to multiply these numbers by 2.303 in order to convert LOG to LN instead.
The correct expression for B and C should be

B = (-0.783*log(volBondScaled)-0.96)*contactAngleEff*contactAngleEff - 0.0438*log(volBondScaled) + 0.48
C = 0.00967*log(volBondScaled)+0.078;

Can you please check if what I say is correct?

Regards,
Giuseppe

Forums:

ckloss's picture

ckloss | Mon, 02/22/2016 - 11:45

Hi Giuseppe,

I contacted Mrs. Easo - let's see what she says!

best wishes
Christoph

Daniel Queteschiner | Mon, 02/22/2016 - 15:14

Where did you get the reference equation from?
I think the reference for the equation in question is
Soulie, Cherblanc, Youssoufi, Saix, Intl. J Numerical and Analytical Methods in Geomechanics, p213, 30 (2006)
The equation on page 218 is the same as written in the LIGGGHTS manual at
http://www.cfdem.com/media/DEM/docu/gran_cohesion_easo_capillary_viscous....
And that's exactly what's in the code. Note that in C++ the log function returns the natural logarithm.

giuraso's picture

giuraso | Mon, 02/22/2016 - 17:01

I missed the fact that "in C++ the log function returns the natural logarithm". Thank you. That's why I tried to find the reason of discrepancy on the conversion.
BY THE WAY:
On Soulie's Thesis ( https://tel.archives-ouvertes.fr/tel-00010079/document ) page 23 you find the following expressions for B and C:

B = (-0.34*LN(V/R^3)-0.96)*contactAngle^2 - 0.019*LN(V/R^3) + 0.48
C = 0.0042LN(V/R^3)+0.078;
Together with the reference Mikami 1998 (Numerical simulation of cohesive powder behavior in a fluidized bed).

Later on at page 40 of Soulie's Thesis you find the expression for particles of different size
B = (-0.148*LN(volBondStar)-0.96)*contactAngleEff*contactAngleEff - 0.0082*LN(volBondStar) + 0.48
C = 0.0018*LN(volBondStar)+0.078;

and in this case the expression of B and C are correct in LIGGGHTS.

I think I need to ask to the author of the Thesis for clarification.

Thank you Christoph and Daniel for your answers

giuraso's picture

giuraso | Mon, 02/29/2016 - 15:43

Hi all,

I had an email conversation with Liza Easo and Christoph Kloss.
Liza and I agree that Soulié's model in the case of two particles of the same size must match Mikami's model which expression for parameters B and C are

B = (-0.34*LN(volBondScaled)-0.96)*contactAngleEff*contactAngleEff - 0.019*LN(volBondScaled) + 0.48
C = 0.0042LN(volBondScaled)+0.078;

For this reason we believe that the coefficient of the model easo capillary viscous, implemented in LIGGGHTS must be changed.

I think the error is somehow a typo in the reference paper and it is not related to a bug or a mistake of DCS team. My apologies.

Kind regards,
Giuseppe

Additional info:
Mikami 1998 (Numerical simulation of cohesive powder behavior in a fluidized bed), also in Soulie's Thesis ( https://tel.archives-ouvertes.fr/tel-00010079/document ) page 23 you find the following expressions for B and C:
B = (-0.34*LN(V/R^3)-0.96)*contactAngle^2 - 0.019*LN(V/R^3) + 0.48
C = 0.0042LN(V/R^3)+0.078;

Later on at page 40 of Soulie's Thesis you find the expression for particles of different size
B = (-0.148*LN(volBondStar)-0.96)*contactAngleEff*contactAngleEff - 0.0082*LN(volBondStar) + 0.48
C = 0.0018*LN(volBondStar)+0.078;

alex.z | Wed, 09/12/2018 - 19:07

I've been looking into this model and there are apparently some other mistakes:
- Average contact angle is defined as:

contactAngleEff = 0.5 * contactAngle[itype] * contactAngle[jtype] instead of contactAngleEff = 0.5 * (contactAngle[itype] + contactAngle[jtype])

- Derjaguin average radius is defined as:
rEff = radi*radj / (radi+radj), should be rEff = 2 * radi*radj / (radi+radj)

Regards,
Alex

Jeff_DONG | Sat, 02/27/2021 - 02:16

Hi Alex,

I just found the original one 0.5 * contactAngle[itype] * contactAngle[jtype] makes the Fcap results less sensitive when the contact angle ranges from 0-30 degrees, which means Fcap does not have much difference in this range. However the modified one 0.5 * (contactAngle[itype] + contactAngle[jtype]) resolve this issue.

In Bhageshvar and Christoph 2014 Powder tech, they use contact angle 0 which basically deliminate the influence of the equation above.

So, have you figure out the reason?

Thanks,
Jeff

syedUConn | Fri, 06/21/2019 - 23:05

For my project, I need to make a model of fluid bed drying. I have hot air coming in from the bottom and heating wet particles and eventually drying the particles. For this, I could run the simulation up to the heat transfer from the fluid to the particle, which is pretty straightforward and a lot of tutorials helped me. But, now I am really stuck as I need to incorporate moisture in particles and simulate drying. My first idea was to use cfdemSolverPisoSTM and consider the moisture as species but I do not know what is speciesCCapacity and how can I relate this parameter to the drying rate. So I could not go anywhere.
Then I thought of using the easo capillary model. I incorporated the model as specified by the documentation. But I receive the following error (as I need to know the liquid content and how it varies with time),

dump dmp all custom 20000 ../DEM/post/dump*.liggghts_run id type x y z ix iy iz vx vy vz fx fy fz omegax omegay omegaz radius f_Temp[0] f_heatFlux[0] f_surfaceLiquidContent[0]
ERROR: Could not find dump custom fix ID (../dump_custom.cpp:1384)
(../dump_custom.cpp:1384)

Also before the error, I receive the following warning,

WARNING:
The contact model you specified is not located in any whitelist.
Because of this the model will be run in an unoptimized version (increasing runtime by up to 20%).
In order to optimize this model you have the following options:
(i) Run the genAutoExamplesWhitelist.sh script in your LIGGGHTS(R) source folder to automatically parse the input script
(ii) Add the model combination by hand to your style_contact_model_user.whitelist that can be found in your LIGGGHTS(R) source folder
If you perform one of the steps above LIGGGHTS(R) needs to be recompiled to generate the optimized code.

(../contact_models.h:364)

I will really appreciate if someone can give some input for this.

Thanks in advance