CFD Online Discussion Forums - How to Include particle implicit and EXPLICIT source terms in VOF momentum equation

CFD Online Discussion Forums (https://www.cfd-online.com/Forums/)

- OpenFOAM Programming & Development (https://www.cfd-online.com/Forums/openfoam-programming-development/)

- - How to Include particle implicit and EXPLICIT source terms in VOF momentum equation (https://www.cfd-online.com/Forums/openfoam-programming-development/158263-how-include-particle-implicit-explicit-source-terms-vof-momentum-equation.html)

How to Include particle implicit and EXPLICIT source terms in VOF momentum equation

hello
i intend to couple VOF (interFoam) & DPMFOAM but have some problems in including particle implicit and explicit terms in VOF momentum equation. i have ask my questions in bold and underined style within the codes. i really appreciate annnnny help

my momentum question is (4.53):
Attachment 41544

so i wrote the equation as:

Code:

1. fvVectorMatrix UEqn

2.  (

3.      fvm::ddt(rhoAlphac, U) + fvm::div(alphacf*rhoPhi, U)  //the 1st &  2nd terms

5.    - fvm::Sp((fvc::ddt(rhoAlphac) + fvc::div(alphacf*rhoPhi)), U) //for numerical stability

7.      - turbulence->divDevRhoReff(rhoAlphac,U)

8.     + cloudSU //add the effect of particle on momentum equation

9. which i dont know whether i insert it corectly or not!

10. );



//and in momentumPredictor:



11. surfaceScalarField phiForces

12. (

13.    (fvc::interpolate(rAU*cloudVolSUSu) & mesh.Sf()) 

14. );





15. solve

16.     (

17.         UEqn

18.      ==

19.         fvc::reconstruct

20.         (

21.                 phiForces/rAUf  ///the 4th term on the right hand of eq. (4.53) as defined above

22.                 +

23.               (

24.                 fvc::interpolate(interface.sigmaK())*fvc::snGrad(alpha1)   //the 3rd term

25.               - ghf*fvc::snGrad(rhoAlphac)   //the 2nd term

26.               - alphacf*fvc::snGrad(p_rgh)   //the 1st term

27.            ) * mesh.magSf()

28.                 // dont know how to insert  5th term on right hand 

29.             )

30.        );

furthermore, i appreciate any other point (forexpalme one may suggest to split the 1st term on left hand side of equation into implicit and explicit term as: Attachment 41545
with regards

annyy reply!
oh pleassssse! :(

OK, here is any reply: I am stuck, too.
However, may I show how far I got. I introduced the kinematicCollidingCloud in the beginning of the
interFoam.C and added the necessary fields in createFields.H:

tmp<surfaceScalarField> tphiAlphaCorr0;

word kinematicCloudName("kinematicCloud");
args.optionReadIfPresent("cloudName", kinematicCloudName);

Info<< "Constructing kinematicCloud " << kinematicCloudName << endl;
basicKinematicCollidingCloud kinematicCloud //basicKinematicCollidingCloud
(
kinematicCloudName,
rhoInf,
U,
mu,
g
);

IOobject Hheader
(
"H",
runTime.timeName(),
mesh,
IOobject::MUST_READ,
IOobject::AUTO_WRITE
);

autoPtr<volVectorField> HPtr;

if (Hheader.headerOk())
{
Info<< "\nReading field H\n" << endl;

HPtr.reset(new volVectorField (Hheader, mesh));
}

IOobject HdotGradHheader
(
"HdotGradH",
runTime.timeName(),
mesh,
IOobject::MUST_READ,
IOobject::AUTO_WRITE
);

autoPtr<volVectorField> HdotGradHPtr;

if (HdotGradHheader.headerOk())
{
Info<< "Reading field HdotGradH" << endl;

HdotGradHPtr.reset(new volVectorField(HdotGradHheader, mesh));
}

#include "createNonInertialFrameFields.H"

For my first approach I neglected the alpha phase reduction due to particles present. I added in interFoam.C simply after the alpha transport equations:

mu = threePhaseProperties.nu()*rhoInfValue; // needed for cloud, adapt to phases
kinematicCloud.evolve();
fvVectorMatrix cloudSU(kinematicCloud.SU(U));//rate of momentum exchange per volume
volVectorField cloudVolSUSu
(
IOobject
(
"cloudVolSUSu",
runTime.timeName(),
mesh
),
mesh,
dimensionedVector
(
"0",
cloudSU.dimensions()/dimVolume,
vector::zero
),
zeroGradientFvPatchVectorField::typeName
);

cloudVolSUSu.internalField() = -cloudSU.source()/mesh.V();
cloudVolSUSu.correctBoundaryConditions();
cloudSU.source() = vector::zero;

Info<<min(kinematicCloud.UTrans())<<endl;
Info<<min(kinematicCloud.UCoeff())<<endl;
Info<<max(kinematicCloud.UTrans())<<endl;
Info<<max(kinematicCloud.UCoeff())<<endl;

and in the UEqun the corresponding momentum as ' == cloudSU'.
The solver runs with particles colliding at the walls but crashes when particles collide with each other. However, that might be due to my little knowledge about lagrangian cloud stability...