[popovich]

Hi all,
i made a solver based on interPhaseChangFoam and compressibleInterfoam.
i have some question about the alpha equation. Can you give me some suggestion?

alpha Equation in compressibleInterFoam is as followed
{
word alphaScheme("div(phi,alpha)");
word alpharScheme("div(phirb,alpha)");

surfaceScalarField phir(phic*interface.nHatf());

for (int gCorr=0; gCorr<nAlphaCorr; gCorr++)
{
volScalarField:imensionedInternalField Sp
(
IOobject
(
"Sp",
runTime.timeName(),
mesh
),
mesh,
dimensionedScalar("Sp", dgdt.dimensions(), 0.0)
);

volScalarField:imensionedInternalField Su
(
IOobject
(
"Su",
runTime.timeName(),
mesh
),
// Divergence term is handled explicitly to be
// consistent with the explicit transport solution
divU*min(alpha1, scalar(1))
);

forAll(dgdt, celli)
{
if (dgdt[celli] > 0.0 && alpha1[celli] > 0.0)
{
Sp[celli] -= dgdt[celli]*alpha1[celli];
Su[celli] += dgdt[celli]*alpha1[celli];
}
else if (dgdt[celli] < 0.0 && alpha1[celli] < 1.0)
{
Sp[celli] += dgdt[celli]*(1.0 - alpha1[celli]);
}
}


surfaceScalarField phiAlpha1
(
fvc::flux
(
phi,
alpha1,
alphaScheme
)
+ fvc::flux
(
-fvc::flux(-phir, alpha2, alpharScheme),
alpha1,
alpharScheme
)
);

MULES::explicitSolve
(
geometricOneField(),
alpha1,
phi,
phiAlpha1,
Sp,
Su,
1,
0
);

surfaceScalarField rho1f(fvc::interpolate(rho1));
surfaceScalarField rho2f(fvc::interpolate(rho2));
rhoPhi = phiAlpha1*(rho1f - rho2f) + phi*rho2f;

alpha2 = scalar(1) - alpha1;
}

Info<< "Liquid phase volume fraction = "
<< alpha1.weightedAverage(mesh.V()).value()
<< " Min(" << alpha1.name() <<") = " << min(alpha1).value()
<< " Min(" << alpha2.name() <<") = " << min(alpha2).value()
<< endl;
}



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

surfaceScalarField phir("phir", phic*interface.nHatf());

Pair<tmp<volScalarField> > vDotAlphal =
twoPhaseProperties->vDotAlphal();
const volScalarField& vDotcAlphal = vDotAlphal[0]();
const volScalarField& vDotvAlphal = vDotAlphal[1]();
const volScalarField vDotvmcAlphal(vDotvAlphal - vDotcAlphal);
for (int gCorr=0; gCorr<nAlphaCorr; gCorr++)
{
volScalarField:imensionedInternalField Sp
(
IOobject
(
"Sp",
runTime.timeName(),
mesh
),
vDotvmcAlphal
);

volScalarField:imensionedInternalField Su
(
IOobject
(
"Su",
runTime.timeName(),
mesh
),
// Divergence term is handled explicitly to be
// consistent with the explicit transport solution
vDotcAlphal //divU*min(alpha1, scalar(1))
);

forAll(dgdt, celli)
{
if (dgdt[celli] >= 0.0 )
{
Su[celli] += (dgdt[celli]*alpha1[celli]*alpha2[celli]);
}
else
{
Sp[celli] += (dgdt[celli]*alpha2[celli]);
}
if (divU[celli] >= 0.0 )
{
Su[celli] += alpha1[celli]*divU[celli];
}
else
{
Sp[celli] += divU[celli];
}
}


surfaceScalarField phiAlpha1
(
fvc::flux
(
phi,
alpha1,
alphaScheme
)
+ fvc::flux
(
-fvc::flux(-phir, alpha2, alpharScheme),
alpha1,
alpharScheme
)
);

MULES::explicitSolve
(
geometricOneField(),
alpha1,
phi,
phiAlpha1,
Sp,
Su,
1,
0
);

surfaceScalarField rho1f(fvc::interpolate(rho1));
surfaceScalarField rho2f(fvc::interpolate(rho2));
rhoPhi = phiAlpha1*(rho1f - rho2f) + phi*rho2f;

alpha1 = min(max(alpha1, scalar(0)), scalar(1));
alpha2 = scalar(1) - alpha1;
}

Info<< "Liquid phase volume fraction = "
<< alpha1.weightedAverage(mesh.V()).value()
<< " Min(" << alpha1.name() <<") = " << min(alpha1).value()
<< " Min(" << alpha2.name() <<") = " << min(alpha2).value()
<< endl;
}

[elmo555]

Hey,

first of all, it's easier to read if you put your code in the CODE-tags, then the formatting will be nice.

And what exactly is your question?