[cheng1988sjtu]

I'm right now trying to understand how twoPhaseEulerFoam works, and get what they are erally solving. However, There seems something odd in the code.

Following is the source code for Ua and Ub in twoPhaseEulerFoam:

fvVectorMatrix UaEqn(Ua, Ua.dimensions()*dimVol/dimTime);
00002 fvVectorMatrix UbEqn(Ub, Ub.dimensions()*dimVol/dimTime);
00003
00004 {
00005 {
00006 volTensorField gradUaT = fvc::grad(Ua)().T();
00007 volTensorField Rca
00008 (
00009 "Rca",
00010 ((2.0/3.0)*I)*(sqr(Ct)*k + nuEffa*tr(gradUaT)) - nuEffa*gradUaT
00011 );
00012
00013 if (kineticTheory.on())
00014 {
00015 Rca -= ((kineticTheory.lambda()/rhoa)*tr(gradUaT))*tensor(I);
00016 }
00017
00018 surfaceScalarField phiRa =
00019 -fvc::interpolate(nuEffa)*mesh.magSf()*fvc::snGrad( alpha)
00020 /fvc::interpolate(alpha + scalar(0.001));
00021
00022 UaEqn =
00023 (
00024 (scalar(1) + Cvm*rhob*beta/rhoa)*
00025 (
00026 fvm::ddt(Ua)
00027 + fvm::div(phia, Ua, "div(phia,Ua)")
00028 - fvm::Sp(fvc::div(phia), Ua)
00029 )
00030
00031 - fvm::laplacian(nuEffa, Ua)
00032 + fvc::div(Rca)
00033
00034 + fvm::div(phiRa, Ua, "div(phia,Ua)")
00035 - fvm::Sp(fvc::div(phiRa), Ua)
00036 + (fvc::grad(alpha)/(fvc::average(alpha) + scalar(0.001)) & Rca)
00037 ==
00038 // g // Buoyancy term transfered to p-equation
00039 - fvm::Sp(beta/rhoa*K, Ua)
00040 //+ beta/rhoa*K*Ub // Explicit drag transfered to p-equation
00041 - beta/rhoa*(liftCoeff - Cvm*rhob*DDtUb)
00042 );
00043
00044 UaEqn.relax();
00045 }
00046
00047 {
00048 volTensorField gradUbT = fvc::grad(Ub)().T();
00049 volTensorField Rcb
00050 (
00051 "Rcb",
00052 ((2.0/3.0)*I)*(k + nuEffb*tr(gradUbT)) - nuEffb*gradUbT
00053 );
00054
00055 surfaceScalarField phiRb =
00056 -fvc::interpolate(nuEffb)*mesh.magSf()*fvc::snGrad( beta)
00057 /fvc::interpolate(beta + scalar(0.001));
00058
00059 UbEqn =
00060 (
00061 (scalar(1) + Cvm*rhob*alpha/rhob)*
00062 (
00063 fvm::ddt(Ub)
00064 + fvm::div(phib, Ub, "div(phib,Ub)")
00065 - fvm::Sp(fvc::div(phib), Ub)
00066 )
00067
00068 - fvm::laplacian(nuEffb, Ub)
00069 + fvc::div(Rcb)
00070
00071 + fvm::div(phiRb, Ub, "div(phib,Ub)")
00072 - fvm::Sp(fvc::div(phiRb), Ub)
00073
00074 + (fvc::grad(beta)/(fvc::average(beta) + scalar(0.001)) & Rcb)
00075 ==
00076 // g // Buoyancy term transfered to p-equation
00077 - fvm::Sp(alpha/rhob*K, Ub)
00078 //+ alpha/rhob*K*Ua // Explicit drag transfered to p-equation
00079 + alpha/rhob*(liftCoeff + Cvm*rhob*DDtUa)
00080 );
00081
00082 UbEqn.relax();
00083 }
00084 }

in the code it seems very different for solving Ua and Ub (Noted as red in the quoted codes, compare line 10 and 52, line 24 and 61, line 41 and 79), However, We can define which phase is the dispersed phase in the "constant" files. Which means that we Ua and Ub should be solved in the similar way.

Or does openfoam has assigned use Ua to calculate dispersed phase and Ub for continuous phase? for example I define b as dispersed phase, then openfoam will solve Ua for b?

Also, note that only in Ueqn for a we can find the kinetic theory option. Is that true?

What's the momentum equation to guide this code? I'm relatively new to twoPhaseEulerFoam. Any help will be highly appreciated!

[ulisses]

Hi Charlie,

I am using kinetic theory in twoPhaseEulerFoam but my solution diverges. Maybe it could be because of that you showed up. Have you achieved a conclusion about that?

Thanks

[cheng1988sjtu]

Hi,

There are two options to simulate the particle phase turbulence (In the code).

If we turn on the kinetic theory, then use kinetic theory to calculate the effective viscosity of particle phase.

If we turn the kinetic theory off, then we use turbulent response coefficient Ct to model the turbulent viscosity of particle phase.

As far as I'm concerned, the alphaEqn is more of significance leading to numerical divergence, since the alphaEqn in twoPhaseEulerFoam is nonlinear and also need to be bounded on both ends. try to refer to PhD thesis of Henrik Rusche, 2002.

About the phase selection, I just confirm to the convention, that is, use default setting, a as disperse phase, b as continuous phase.

Any question are welcomed.

Charlie

Originally Posted by ulisses

Hi Charlie,

I am using kinetic theory in twoPhaseEulerFoam but my solution diverges. Maybe it could be because of that you showed up. Have you achieved a conclusion about that?

Thanks