CFD Online Discussion Forums

CFD Online Discussion Forums (http://www.cfd-online.com/Forums/)
-   OpenFOAM (http://www.cfd-online.com/Forums/openfoam/)
-   -   Bug in the PilchErdman Breakup Model (http://www.cfd-online.com/Forums/openfoam/91795-bug-pilcherdman-breakup-model.html)

xck1986 August 23, 2011 04:57

Bug in the PilchErdman Breakup Model
 
Hi everyone,

I have notice, there is a bug in the PilchErdman Breakup Model in OpenFOAM.

In the source code of the Model, the following two equ.,

scalar Vd = Urmag*rho12*(B1_*taubBar * B2_*taubBar*taubBar);
......
scalar Ds = 2.0*Wec*sigma*Vd1/(Vd1*rhoc*sqr(Urmag));

which are not the same as in the paper: Use of breakup time
data and velocity history date to predict the maximun size of
stable fragments for acceleration
induced breakup of a liquid
drop
.

I think it should be rewritten like following:

scalar Vd = Urmag*rho12*(B1_*taubBar + B2_*taubBar*taubBar);
......
scalar Ds = 2.0*Wec*sigma/(Vd1*rhoc*sqr(Urmag));

Does anyone also notice this error.

Thanks

Chenkai

xck1986 August 30, 2011 04:59

Nobody has ever used this model?

The problem is that, when I change the code of this model as in the paper. My Simulation results are even worse than before. I don't if here is really a bug or not.
Hope someone can help me?

Thanks a lot!

bigeddy January 16, 2013 08:08

A year late but this post already exits so im updating it here.

Im looking at this model currently for my own purposes. You are correct, the model should be what you suggest as a fix i.e.

scalar Vd = Urmag*rho12*(B1_*taubBar + B2_*taubBar*taubBar);
......
scalar Ds = 2.0*Wec*sigma/(Vd1*rhoc*sqr(Urmag));

Also,

// update the droplet diameter according to the rate eq. (implicitly)
d = (d + frac*Ds)/(1.0 + frac);

makes no sense whatsoever. It should be

d = d*(1.0-frac)+Ds*frac

when im less busy ill update the git repository if someone else don't get to it first.

niklas January 16, 2013 14:28

Quote:

Originally Posted by bigeddy (Post 402173)
A year late but this post already exits so im updating it here.

Im looking at this model currently for my own purposes. You are correct, the model should be what you suggest as a fix i.e.

scalar Vd = Urmag*rho12*(B1_*taubBar + B2_*taubBar*taubBar);
......
scalar Ds = 2.0*Wec*sigma/(Vd1*rhoc*sqr(Urmag));

Also,

// update the droplet diameter according to the rate eq. (implicitly)
d = (d + frac*Ds)/(1.0 + frac);

makes no sense whatsoever. It should be

d = d*(1.0-frac)+Ds*frac

when im less busy ill update the git repository if someone else don't get to it first.

BUG!!! Inconcievable :)
I cant belieave I've missed this one.
I will talk to the guys and correct it.
About the 'makes no sense whatsoever' part. You could do it like that. Its not wrong, but neither is the implemented update. (as shown below)

f = dt/tau

we're trying to solve this eq:
dD/dt = -(D-Ds)/tau

1) explicit, approximate D with D0 in the RHS
(D1-D0)/dt = -(D0-Ds)/tau
D1 – D0 = -(D0 – Ds)*f
D1 = D0 - D0*f + Ds*f
D1 = (1-f)*D0 + Ds*f


2) implicit, approximate D with D1 in the RHS
(D1-D0)/dt = -(D1 – Ds)/tau
(D1-D0)= -D1*f + Ds*f
(1+f)*D1 - D0 = Ds*f
(1+f)*D1 = D0 + Ds*f
D1 = (D0 + f*Ds)/(1+f)

niklas January 16, 2013 15:29

now...you may wonder 'why is method 2 to prefer'?

think of large timesteps (f >> 1)

first method yields

D1 -> f*(Ds-D0) (makes no physical sense)

second method yields

D1 -> f*Ds/f = Ds (makes perfect sense)

bigeddy January 16, 2013 18:07

Ahhh, I see, thanks for reminding me of the two methods!


All times are GMT -4. The time now is 17:37.