[Elol]

Hello Foamers,

Right now , I am trying to develope my own incombressible multiphase slover starting from icoFoam solver. So, I bring the first step which I will add the alphaEqn form driftFLuxFoam to the solver (I use driftFluxFoam Solver as a refrence). However there is some questions regarding the solution method of the alphaEqn.

the alphaEqn.H file solve this eqn ddt(alpha,t)+div(alpha, U) = 0 as follow:

{
word alphaScheme("div(phi,alpha)");
word alpharScheme("div(phirb,alpha)");

if (MULESCorr)
{
fvScalarMatrix alpha1Eqn
(
fv::EulerDdtScheme<scalar>(mesh).fvmDdt(alpha1)
+ fv::gaussConvectionScheme<scalar>
(
mesh,
phi,
upwind<scalar>(mesh, phi)
).fvmDiv(phi, alpha1)
);

solve(alpha1Eqn);

Info<< "Phase-1 volume fraction = "
<< alpha1.weightedAverage(mesh.Vsc()).value()
<< " Min(" << alpha1.name() << ") = " << min(alpha1).value()
<< " Max(" << alpha1.name() << ") = " << max(alpha1).value()
<< endl;

tmp<surfaceScalarField> talphaPhiUD(alpha1Eqn.flux());
alphaPhi = talphaPhiUD();

if (alphaApplyPrevCorr && talphaPhiCorr0.valid())
{
Info<< "Applying the previous iteration correction flux" << endl;

MULES::correct
(
geometricOneField(),
alpha1,
alphaPhi,
talphaPhiCorr0.ref(),
UniformField<scalar>(mixture.alphaMax()),
zeroField()
);

alphaPhi += talphaPhiCorr0();
}

// Cache the upwind-flux
talphaPhiCorr0 = talphaPhiUD;
}

for (int aCorr=0; aCorr<nAlphaCorr; aCorr++)
{
tmp<surfaceScalarField> talphaPhiUn
(
fvc::flux
(
phi,
alpha1,
alphaScheme
)
+ fvc::flux
(
phir,
alpha1,
alpharScheme
)
);

if (MULESCorr)
{
tmp<surfaceScalarField> talphaPhiCorr(talphaPhiUn() - alphaPhi);
volScalarField alpha10("alpha10", alpha1);

MULES::correct
(
geometricOneField(),
alpha1,
talphaPhiUn(),
talphaPhiCorr.ref(),
UniformField<scalar>(mixture.alphaMax()),
zeroField()
);

// Under-relax the correction for all but the 1st corrector
if (aCorr == 0)
{
alphaPhi += talphaPhiCorr();
}
else
{
alpha1 = 0.5*alpha1 + 0.5*alpha10;
alphaPhi += 0.5*talphaPhiCorr();
}
}
else
{
alphaPhi = talphaPhiUn;

MULES::explicitSolve
(
geometricOneField(),
alpha1,
phi,
alphaPhi,
UniformField<scalar>(mixture.alphaMax()),
zeroField()
);
}
}

if (alphaApplyPrevCorr && MULESCorr)
{
talphaPhiCorr0 = alphaPhi - talphaPhiCorr0;
}

alpha2 = 1.0 - alpha1;

Info<< "Phase-1 volume fraction = "
<< alpha1.weightedAverage(mesh.Vsc()).value()
<< " Min(" << alpha1.name() << ") = " << min(alpha1).value()
<< " Max(" << alpha1.name() << ") = " << max(alpha1).value()
<< endl;
}




What I understand till this moment is that it starts to solve the equation using low order scheme to generate talphaPhiUD() in this section

fvScalarMatrix alpha1Eqn
(
fv::EulerDdtScheme<scalar>(mesh).fvmDdt(alpha1)
+ fv::gaussConvectionScheme<scalar>
(
mesh,
phi,
upwind<scalar>(mesh, phi)
).fvmDiv(phi, alpha1)
);

solve(alpha1Eqn);

Info<< "Phase-1 volume fraction = "
<< alpha1.weightedAverage(mesh.Vsc()).value()
<< " Min(" << alpha1.name() << ") = " << min(alpha1).value()
<< " Max(" << alpha1.name() << ") = " << max(alpha1).value()
<< endl;

tmp<surfaceScalarField> talphaPhiUD(alpha1Eqn.flux());
alphaPhi = talphaPhiUD();

then and *if condition* is applied to perform MULES

if (alphaApplyPrevCorr && talphaPhiCorr0.valid())
{
Info<< "Applying the previous iteration correction flux" << endl;

MULES::correct
(
geometricOneField(),
alpha1,
alphaPhi,
talphaPhiCorr0.ref(),
UniformField<scalar>(mixture.alphaMax()),
zeroField()
);

alphaPhi += talphaPhiCorr0();
}

alphaApplyPrevCorr is an inlet from the user and here is my question :

1. why they apply MULES correction to the low order fluxes ?
2. what this condition (talphaPhiCorr0.valid()) mean ?

Thank you guys in advance.

[SHUBHAM9595]

Hi Elol,


I'm also facing bit difficulty in understanding specific meaning of these terms clearly. Can you share if you found an answer or some resource where it is mentioned clearly.

[Elol]

Hi MNM,

Sorry for the late answer, I don't know if you find the answer or not. After some research, I really advise you to take a look on the PhD thesis of Santiago Marquez Damian. It is entitled as " An Extended Mixture Model for the Simultaneous Treatment of Short and Long Scale Interfaces". He explained in chapter 2 the whole idea of weighting factors and how it is implemnted in openFoam.

Enjoy reading it. It is really valuable work.

Best regards
Elol