zhyang-dev/OpenFOAM-12
1
0
Fork 0
Code Issues Pull Requests Actions Packages Projects Releases Wiki Activity

Renamed phiAlpha -P alphaPhi for consistency with Euler-Euler solvers

Browse Source
This commit is contained in:
Henry Weller
2015-09-11 17:52:43 +01:00
parent cc5f67a0ff
commit 1b835f64f0
11 changed files with 105 additions and 105 deletions
Download Patch File Download Diff File Expand all files Collapse all files

6
applications/solvers/multiphase/compressibleInterFoam/alphaEqns.H
View File

@ -45,7 +45,7 @@
} }
surfaceScalarField phiAlpha1 surfaceScalarField alphaPhi1
( (
fvc::flux fvc::flux
( (
@ -66,7 +66,7 @@
geometricOneField(), geometricOneField(),
alpha1, alpha1,
phi, phi,
phiAlpha1, alphaPhi1,
Sp, Sp,
Su, Su,
1, 1,
@ -75,7 +75,7 @@
surfaceScalarField rho1f(fvc::interpolate(rho1)); surfaceScalarField rho1f(fvc::interpolate(rho1));
surfaceScalarField rho2f(fvc::interpolate(rho2)); surfaceScalarField rho2f(fvc::interpolate(rho2));
rhoPhi = phiAlpha1*(rho1f - rho2f) + phi*rho2f; rhoPhi = alphaPhi1*(rho1f - rho2f) + phi*rho2f;
alpha2 = scalar(1) - alpha1; alpha2 = scalar(1) - alpha1;
} }

20
applications/solvers/multiphase/compressibleMultiphaseInterFoam/multiphaseMixtureThermo/multiphaseMixtureThermo.C
View File

@ -942,14 +942,14 @@ void Foam::multiphaseMixtureThermo::solveAlphas
surfaceScalarField phic(mag(phi_/mesh_.magSf())); surfaceScalarField phic(mag(phi_/mesh_.magSf()));
phic = min(cAlpha*phic, max(phic)); phic = min(cAlpha*phic, max(phic));
PtrList<surfaceScalarField> phiAlphaCorrs(phases_.size()); PtrList<surfaceScalarField> alphaPhiCorrs(phases_.size());
int phasei = 0; int phasei = 0;
forAllIter(PtrDictionary<phaseModel>, phases_, phase) forAllIter(PtrDictionary<phaseModel>, phases_, phase)
{ {
phaseModel& alpha = phase(); phaseModel& alpha = phase();
phiAlphaCorrs.set alphaPhiCorrs.set
( (
phasei, phasei,
new surfaceScalarField new surfaceScalarField
@ -964,7 +964,7 @@ void Foam::multiphaseMixtureThermo::solveAlphas
) )
); );
surfaceScalarField& phiAlphaCorr = phiAlphaCorrs[phasei]; surfaceScalarField& alphaPhiCorr = alphaPhiCorrs[phasei];
forAllIter(PtrDictionary<phaseModel>, phases_, phase2) forAllIter(PtrDictionary<phaseModel>, phases_, phase2)
{ {
@ -974,7 +974,7 @@ void Foam::multiphaseMixtureThermo::solveAlphas
surfaceScalarField phir(phic*nHatf(alpha, alpha2)); surfaceScalarField phir(phic*nHatf(alpha, alpha2));
phiAlphaCorr += fvc::flux alphaPhiCorr += fvc::flux
( (
-fvc::flux(-phir, alpha2, alpharScheme), -fvc::flux(-phir, alpha2, alpharScheme),
alpha, alpha,
@ -988,7 +988,7 @@ void Foam::multiphaseMixtureThermo::solveAlphas
geometricOneField(), geometricOneField(),
alpha, alpha,
phi_, phi_,
phiAlphaCorr, alphaPhiCorr,
zeroField(), zeroField(),
zeroField(), zeroField(),
1, 1,
@ -999,7 +999,7 @@ void Foam::multiphaseMixtureThermo::solveAlphas
phasei++; phasei++;
} }
MULES::limitSum(phiAlphaCorrs); MULES::limitSum(alphaPhiCorrs);
rhoPhi_ = dimensionedScalar("0", dimensionSet(1, 0, -1, 0, 0), 0); rhoPhi_ = dimensionedScalar("0", dimensionSet(1, 0, -1, 0, 0), 0);
@ -1025,8 +1025,8 @@ void Foam::multiphaseMixtureThermo::solveAlphas
{ {
phaseModel& alpha = phase(); phaseModel& alpha = phase();
surfaceScalarField& phiAlpha = phiAlphaCorrs[phasei]; surfaceScalarField& alphaPhi = alphaPhiCorrs[phasei];
phiAlpha += upwind<scalar>(mesh_, phi_).flux(alpha); alphaPhi += upwind<scalar>(mesh_, phi_).flux(alpha);
volScalarField::DimensionedInternalField Sp volScalarField::DimensionedInternalField Sp
( (
@ -1096,12 +1096,12 @@ void Foam::multiphaseMixtureThermo::solveAlphas
( (
geometricOneField(), geometricOneField(),
alpha, alpha,
phiAlpha, alphaPhi,
Sp, Sp,
Su Su
); );
rhoPhi_ += fvc::interpolate(alpha.thermo().rho())*phiAlpha; rhoPhi_ += fvc::interpolate(alpha.thermo().rho())*alphaPhi;
Info<< alpha.name() << " volume fraction, min, max = " Info<< alpha.name() << " volume fraction, min, max = "
<< alpha.weightedAverage(mesh_.V()).value() << alpha.weightedAverage(mesh_.V()).value()

32
applications/solvers/multiphase/driftFluxFoam/alphaEqn.H
View File

@ -23,32 +23,32 @@
<< " Max(" << alpha1.name() << ") = " << max(alpha1).value() << " Max(" << alpha1.name() << ") = " << max(alpha1).value()
<< endl; << endl;
tmp<surfaceScalarField> tphiAlphaUD(alpha1Eqn.flux()); tmp<surfaceScalarField> talphaPhiUD(alpha1Eqn.flux());
phiAlpha = tphiAlphaUD(); alphaPhi = talphaPhiUD();
if (alphaApplyPrevCorr && tphiAlphaCorr0.valid()) if (alphaApplyPrevCorr && talphaPhiCorr0.valid())
{ {
Info<< "Applying the previous iteration correction flux" << endl; Info<< "Applying the previous iteration correction flux" << endl;
MULES::correct MULES::correct
( (
alpha1, alpha1,
phiAlpha, alphaPhi,
tphiAlphaCorr0(), talphaPhiCorr0(),
mixture.alphaMax(), mixture.alphaMax(),
0 0
); );
phiAlpha += tphiAlphaCorr0(); alphaPhi += talphaPhiCorr0();
} }
// Cache the upwind-flux // Cache the upwind-flux
tphiAlphaCorr0 = tphiAlphaUD; talphaPhiCorr0 = talphaPhiUD;
} }
for (int aCorr=0; aCorr<nAlphaCorr; aCorr++) for (int aCorr=0; aCorr<nAlphaCorr; aCorr++)
{ {
tmp<surfaceScalarField> tphiAlphaUn tmp<surfaceScalarField> talphaPhiUn
( (
fvc::flux fvc::flux
( (
@ -66,14 +66,14 @@
if (MULESCorr) if (MULESCorr)
{ {
tmp<surfaceScalarField> tphiAlphaCorr(tphiAlphaUn() - phiAlpha); tmp<surfaceScalarField> talphaPhiCorr(talphaPhiUn() - alphaPhi);
volScalarField alpha10("alpha10", alpha1); volScalarField alpha10("alpha10", alpha1);
MULES::correct MULES::correct
( (
alpha1, alpha1,
tphiAlphaUn(), talphaPhiUn(),
tphiAlphaCorr(), talphaPhiCorr(),
mixture.alphaMax(), mixture.alphaMax(),
0 0
); );
@ -81,23 +81,23 @@
// Under-relax the correction for all but the 1st corrector // Under-relax the correction for all but the 1st corrector
if (aCorr == 0) if (aCorr == 0)
{ {
phiAlpha += tphiAlphaCorr(); alphaPhi += talphaPhiCorr();
} }
else else
{ {
alpha1 = 0.5*alpha1 + 0.5*alpha10; alpha1 = 0.5*alpha1 + 0.5*alpha10;
phiAlpha += 0.5*tphiAlphaCorr(); alphaPhi += 0.5*talphaPhiCorr();
} }
} }
else else
{ {
phiAlpha = tphiAlphaUn; alphaPhi = talphaPhiUn;
MULES::explicitSolve MULES::explicitSolve
( (
alpha1, alpha1,
phi, phi,
phiAlpha, alphaPhi,
mixture.alphaMax(), mixture.alphaMax(),
0 0
); );
@ -106,7 +106,7 @@
if (alphaApplyPrevCorr && MULESCorr) if (alphaApplyPrevCorr && MULESCorr)
{ {
tphiAlphaCorr0 = phiAlpha - tphiAlphaCorr0; talphaPhiCorr0 = alphaPhi - talphaPhiCorr0;
} }
alpha2 = 1.0 - alpha1; alpha2 = 1.0 - alpha1;

16
applications/solvers/multiphase/driftFluxFoam/alphaEqnSubCycle.H
View File

@ -1,9 +1,9 @@
{ {
surfaceScalarField phiAlpha surfaceScalarField alphaPhi
( (
IOobject IOobject
( (
"phiAlpha", "alphaPhi",
runTime.timeName(), runTime.timeName(),
mesh mesh
), ),
@ -19,11 +19,11 @@
if (nAlphaSubCycles > 1) if (nAlphaSubCycles > 1)
{ {
dimensionedScalar totalDeltaT = runTime.deltaT(); dimensionedScalar totalDeltaT = runTime.deltaT();
surfaceScalarField phiAlphaSum surfaceScalarField alphaPhiSum
( (
IOobject IOobject
( (
"phiAlphaSum", "alphaPhiSum",
runTime.timeName(), runTime.timeName(),
mesh mesh
), ),
@ -38,10 +38,10 @@
) )
{ {
#include "alphaEqn.H" #include "alphaEqn.H"
phiAlphaSum += (runTime.deltaT()/totalDeltaT)*phiAlpha; alphaPhiSum += (runTime.deltaT()/totalDeltaT)*alphaPhi;
} }
phiAlpha = phiAlphaSum; alphaPhi = alphaPhiSum;
} }
else else
{ {
@ -59,7 +59,7 @@
alpha1Eqn.solve(mesh.solver("alpha1Diffusion")); alpha1Eqn.solve(mesh.solver("alpha1Diffusion"));
phiAlpha += alpha1Eqn.flux(); alphaPhi += alpha1Eqn.flux();
alpha2 = 1.0 - alpha1; alpha2 = 1.0 - alpha1;
Info<< "Phase-1 volume fraction = " Info<< "Phase-1 volume fraction = "
@ -69,6 +69,6 @@
<< endl; << endl;
} }
rhoPhi = phiAlpha*(rho1 - rho2) + phi*rho2; rhoPhi = alphaPhi*(rho1 - rho2) + phi*rho2;
rho = mixture.rho(); rho = mixture.rho();
} }

2
applications/solvers/multiphase/driftFluxFoam/createFields.H
View File

@ -135,4 +135,4 @@ mesh.setFluxRequired(p_rgh.name());
mesh.setFluxRequired(alpha1.name()); mesh.setFluxRequired(alpha1.name());
// MULES Correction // MULES Correction
tmp<surfaceScalarField> tphiAlphaCorr0; tmp<surfaceScalarField> talphaPhiCorr0;

38
applications/solvers/multiphase/interFoam/alphaEqn.H
View File

@ -98,19 +98,19 @@
<< " Max(" << alpha1.name() << ") = " << max(alpha1).value() << " Max(" << alpha1.name() << ") = " << max(alpha1).value()
<< endl; << endl;
tmp<surfaceScalarField> tphiAlphaUD(alpha1Eqn.flux()); tmp<surfaceScalarField> talphaPhiUD(alpha1Eqn.flux());
phiAlpha = tphiAlphaUD(); alphaPhi = talphaPhiUD();
if (alphaApplyPrevCorr && tphiAlphaCorr0.valid()) if (alphaApplyPrevCorr && talphaPhiCorr0.valid())
{ {
Info<< "Applying the previous iteration compression flux" << endl; Info<< "Applying the previous iteration compression flux" << endl;
MULES::correct(alpha1, phiAlpha, tphiAlphaCorr0(), 1, 0); MULES::correct(alpha1, alphaPhi, talphaPhiCorr0(), 1, 0);
phiAlpha += tphiAlphaCorr0(); alphaPhi += talphaPhiCorr0();
} }
// Cache the upwind-flux // Cache the upwind-flux
tphiAlphaCorr0 = tphiAlphaUD; talphaPhiCorr0 = talphaPhiUD;
alpha2 = 1.0 - alpha1; alpha2 = 1.0 - alpha1;
@ -122,7 +122,7 @@
{ {
surfaceScalarField phir(phic*mixture.nHatf()); surfaceScalarField phir(phic*mixture.nHatf());
tmp<surfaceScalarField> tphiAlphaUn tmp<surfaceScalarField> talphaPhiUn
( (
fvc::flux fvc::flux
( (
@ -141,33 +141,33 @@
// Calculate the Crank-Nicolson off-centred alpha flux // Calculate the Crank-Nicolson off-centred alpha flux
if (ocCoeff > 0) if (ocCoeff > 0)
{ {
tphiAlphaUn = talphaPhiUn =
cnCoeff*tphiAlphaUn + (1.0 - cnCoeff)*phiAlpha.oldTime(); cnCoeff*talphaPhiUn + (1.0 - cnCoeff)*alphaPhi.oldTime();
} }
if (MULESCorr) if (MULESCorr)
{ {
tmp<surfaceScalarField> tphiAlphaCorr(tphiAlphaUn() - phiAlpha); tmp<surfaceScalarField> talphaPhiCorr(talphaPhiUn() - alphaPhi);
volScalarField alpha10("alpha10", alpha1); volScalarField alpha10("alpha10", alpha1);
MULES::correct(alpha1, tphiAlphaUn(), tphiAlphaCorr(), 1, 0); MULES::correct(alpha1, talphaPhiUn(), talphaPhiCorr(), 1, 0);
// Under-relax the correction for all but the 1st corrector // Under-relax the correction for all but the 1st corrector
if (aCorr == 0) if (aCorr == 0)
{ {
phiAlpha += tphiAlphaCorr(); alphaPhi += talphaPhiCorr();
} }
else else
{ {
alpha1 = 0.5*alpha1 + 0.5*alpha10; alpha1 = 0.5*alpha1 + 0.5*alpha10;
phiAlpha += 0.5*tphiAlphaCorr(); alphaPhi += 0.5*talphaPhiCorr();
} }
} }
else else
{ {
phiAlpha = tphiAlphaUn; alphaPhi = talphaPhiUn;
MULES::explicitSolve(alpha1, phiCN, phiAlpha, 1, 0); MULES::explicitSolve(alpha1, phiCN, alphaPhi, 1, 0);
} }
alpha2 = 1.0 - alpha1; alpha2 = 1.0 - alpha1;
@ -177,7 +177,7 @@
if (alphaApplyPrevCorr && MULESCorr) if (alphaApplyPrevCorr && MULESCorr)
{ {
tphiAlphaCorr0 = phiAlpha - tphiAlphaCorr0; talphaPhiCorr0 = alphaPhi - talphaPhiCorr0;
} }
if if
@ -186,18 +186,18 @@
== fv::EulerDdtScheme<vector>::typeName == fv::EulerDdtScheme<vector>::typeName
) )
{ {
rhoPhi = phiAlpha*(rho1 - rho2) + phiCN*rho2; rhoPhi = alphaPhi*(rho1 - rho2) + phiCN*rho2;
} }
else else
{ {
if (ocCoeff > 0) if (ocCoeff > 0)
{ {
// Calculate the end-of-time-step alpha flux // Calculate the end-of-time-step alpha flux
phiAlpha = (phiAlpha - (1.0 - cnCoeff)*phiAlpha.oldTime())/cnCoeff; alphaPhi = (alphaPhi - (1.0 - cnCoeff)*alphaPhi.oldTime())/cnCoeff;
} }
// Calculate the end-of-time-step mass flux // Calculate the end-of-time-step mass flux
rhoPhi = phiAlpha*(rho1 - rho2) + phi*rho2; rhoPhi = alphaPhi*(rho1 - rho2) + phi*rho2;
} }
Info<< "Phase-1 volume fraction = " Info<< "Phase-1 volume fraction = "

6
applications/solvers/multiphase/interFoam/createFields.H
View File

@ -121,11 +121,11 @@ mesh.setFluxRequired(p_rgh.name());
mesh.setFluxRequired(alpha1.name()); mesh.setFluxRequired(alpha1.name());
// MULES flux from previous time-step // MULES flux from previous time-step
surfaceScalarField phiAlpha surfaceScalarField alphaPhi
( (
IOobject IOobject
( (
"phiAlpha", "alphaPhi",
runTime.timeName(), runTime.timeName(),
mesh, mesh,
IOobject::READ_IF_PRESENT, IOobject::READ_IF_PRESENT,
@ -135,4 +135,4 @@ surfaceScalarField phiAlpha
); );
// MULES Correction // MULES Correction
tmp<surfaceScalarField> tphiAlphaCorr0; tmp<surfaceScalarField> talphaPhiCorr0;

42
applications/solvers/multiphase/interFoam/interMixingFoam/alphaEqns.H
View File

@ -40,7 +40,7 @@
// Create the complete convection flux for alpha1 // Create the complete convection flux for alpha1
surfaceScalarField phiAlpha1 surfaceScalarField alphaPhi1
( (
fvc::flux fvc::flux
( (
@ -63,13 +63,13 @@
); );
// Create the bounded (upwind) flux for alpha1 // Create the bounded (upwind) flux for alpha1
surfaceScalarField phiAlpha1BD surfaceScalarField alphaPhi1BD
( (
upwind<scalar>(mesh, phi).flux(alpha1) upwind<scalar>(mesh, phi).flux(alpha1)
); );
// Calculate the flux correction for alpha1 // Calculate the flux correction for alpha1
phiAlpha1 -= phiAlpha1BD; alphaPhi1 -= alphaPhi1BD;
// Calculate the limiter for alpha1 // Calculate the limiter for alpha1
if (LTS) if (LTS)
@ -83,8 +83,8 @@
rDeltaT, rDeltaT,
geometricOneField(), geometricOneField(),
alpha1, alpha1,
phiAlpha1BD, alphaPhi1BD,
phiAlpha1, alphaPhi1,
zeroField(), zeroField(),
zeroField(), zeroField(),
1, 1,
@ -99,8 +99,8 @@
1.0/runTime.deltaT().value(), 1.0/runTime.deltaT().value(),
geometricOneField(), geometricOneField(),
alpha1, alpha1,
phiAlpha1BD, alphaPhi1BD,
phiAlpha1, alphaPhi1,
zeroField(), zeroField(),
zeroField(), zeroField(),
1, 1,
@ -109,7 +109,7 @@
} }
// Create the complete flux for alpha2 // Create the complete flux for alpha2
surfaceScalarField phiAlpha2 surfaceScalarField alphaPhi2
( (
fvc::flux fvc::flux
( (
@ -126,13 +126,13 @@
); );
// Create the bounded (upwind) flux for alpha2 // Create the bounded (upwind) flux for alpha2
surfaceScalarField phiAlpha2BD surfaceScalarField alphaPhi2BD
( (
upwind<scalar>(mesh, phi).flux(alpha2) upwind<scalar>(mesh, phi).flux(alpha2)
); );
// Calculate the flux correction for alpha2 // Calculate the flux correction for alpha2
phiAlpha2 -= phiAlpha2BD; alphaPhi2 -= alphaPhi2BD;
// Further limit the limiter for alpha2 // Further limit the limiter for alpha2
if (LTS) if (LTS)
@ -146,8 +146,8 @@
rDeltaT, rDeltaT,
geometricOneField(), geometricOneField(),
alpha2, alpha2,
phiAlpha2BD, alphaPhi2BD,
phiAlpha2, alphaPhi2,
zeroField(), zeroField(),
zeroField(), zeroField(),
1, 1,
@ -162,8 +162,8 @@
1.0/runTime.deltaT().value(), 1.0/runTime.deltaT().value(),
geometricOneField(), geometricOneField(),
alpha2, alpha2,
phiAlpha2BD, alphaPhi2BD,
phiAlpha2, alphaPhi2,
zeroField(), zeroField(),
zeroField(), zeroField(),
1, 1,
@ -172,32 +172,32 @@
} }
// Construct the limited fluxes // Construct the limited fluxes
phiAlpha1 = phiAlpha1BD + lambda*phiAlpha1; alphaPhi1 = alphaPhi1BD + lambda*alphaPhi1;
phiAlpha2 = phiAlpha2BD + lambda*phiAlpha2; alphaPhi2 = alphaPhi2BD + lambda*alphaPhi2;
// Solve for alpha1 // Solve for alpha1
solve(fvm::ddt(alpha1) + fvc::div(phiAlpha1)); solve(fvm::ddt(alpha1) + fvc::div(alphaPhi1));
// Create the diffusion coefficients for alpha2<->alpha3 // Create the diffusion coefficients for alpha2<->alpha3
volScalarField Dc23(D23*max(alpha3, scalar(0))*pos(alpha2)); volScalarField Dc23(D23*max(alpha3, scalar(0))*pos(alpha2));
volScalarField Dc32(D23*max(alpha2, scalar(0))*pos(alpha3)); volScalarField Dc32(D23*max(alpha2, scalar(0))*pos(alpha3));
// Add the diffusive flux for alpha3->alpha2 // Add the diffusive flux for alpha3->alpha2
phiAlpha2 -= fvc::interpolate(Dc32)*mesh.magSf()*fvc::snGrad(alpha1); alphaPhi2 -= fvc::interpolate(Dc32)*mesh.magSf()*fvc::snGrad(alpha1);
// Solve for alpha2 // Solve for alpha2
fvScalarMatrix alpha2Eqn fvScalarMatrix alpha2Eqn
( (
fvm::ddt(alpha2) fvm::ddt(alpha2)
+ fvc::div(phiAlpha2) + fvc::div(alphaPhi2)
- fvm::laplacian(Dc23 + Dc32, alpha2) - fvm::laplacian(Dc23 + Dc32, alpha2)
); );
alpha2Eqn.solve(); alpha2Eqn.solve();
// Construct the complete mass flux // Construct the complete mass flux
rhoPhi = rhoPhi =
phiAlpha1*(rho1 - rho3) alphaPhi1*(rho1 - rho3)
+ (phiAlpha2 + alpha2Eqn.flux())*(rho2 - rho3) + (alphaPhi2 + alpha2Eqn.flux())*(rho2 - rho3)
+ phi*rho3; + phi*rho3;
alpha3 = 1.0 - alpha1 - alpha2; alpha3 = 1.0 - alpha1 - alpha2;

22
applications/solvers/multiphase/interPhaseChangeFoam/alphaEqn.H
View File

@ -10,7 +10,7 @@
const volScalarField& vDotvAlphal = vDotAlphal[1](); const volScalarField& vDotvAlphal = vDotAlphal[1]();
const volScalarField vDotvmcAlphal(vDotvAlphal - vDotcAlphal); const volScalarField vDotvmcAlphal(vDotvAlphal - vDotcAlphal);
tmp<surfaceScalarField> tphiAlpha; tmp<surfaceScalarField> talphaPhi;
if (MULESCorr) if (MULESCorr)
{ {
@ -37,14 +37,14 @@
<< " Max(" << alpha1.name() << ") = " << max(alpha1).value() << " Max(" << alpha1.name() << ") = " << max(alpha1).value()
<< endl; << endl;
tphiAlpha = alpha1Eqn.flux(); talphaPhi = alpha1Eqn.flux();
} }
volScalarField alpha10("alpha10", alpha1); volScalarField alpha10("alpha10", alpha1);
for (int aCorr=0; aCorr<nAlphaCorr; aCorr++) for (int aCorr=0; aCorr<nAlphaCorr; aCorr++)
{ {
tmp<surfaceScalarField> tphiAlphaCorr tmp<surfaceScalarField> talphaPhiCorr
( (
fvc::flux fvc::flux
( (
@ -62,7 +62,7 @@
if (MULESCorr) if (MULESCorr)
{ {
tphiAlphaCorr() -= tphiAlpha(); talphaPhiCorr() -= talphaPhi();
volScalarField alpha100("alpha100", alpha10); volScalarField alpha100("alpha100", alpha10);
alpha10 = alpha1; alpha10 = alpha1;
@ -71,8 +71,8 @@
( (
geometricOneField(), geometricOneField(),
alpha1, alpha1,
tphiAlpha(), talphaPhi(),
tphiAlphaCorr(), talphaPhiCorr(),
vDotvmcAlphal, vDotvmcAlphal,
( (
divU*(alpha10 - alpha100) divU*(alpha10 - alpha100)
@ -85,12 +85,12 @@
// Under-relax the correction for all but the 1st corrector // Under-relax the correction for all but the 1st corrector
if (aCorr == 0) if (aCorr == 0)
{ {
tphiAlpha() += tphiAlphaCorr(); talphaPhi() += talphaPhiCorr();
} }
else else
{ {
alpha1 = 0.5*alpha1 + 0.5*alpha10; alpha1 = 0.5*alpha1 + 0.5*alpha10;
tphiAlpha() += 0.5*tphiAlphaCorr(); talphaPhi() += 0.5*talphaPhiCorr();
} }
} }
else else
@ -100,20 +100,20 @@
geometricOneField(), geometricOneField(),
alpha1, alpha1,
phi, phi,
tphiAlphaCorr(), talphaPhiCorr(),
vDotvmcAlphal, vDotvmcAlphal,
(divU*alpha1 + vDotcAlphal)(), (divU*alpha1 + vDotcAlphal)(),
1, 1,
0 0
); );
tphiAlpha = tphiAlphaCorr; talphaPhi = talphaPhiCorr;
} }
alpha2 = 1.0 - alpha1; alpha2 = 1.0 - alpha1;
} }
rhoPhi = tphiAlpha()*(rho1 - rho2) + phi*rho2; rhoPhi = talphaPhi()*(rho1 - rho2) + phi*rho2;
Info<< "Liquid phase volume fraction = " Info<< "Liquid phase volume fraction = "
<< alpha1.weightedAverage(mesh.V()).value() << alpha1.weightedAverage(mesh.V()).value()

20
applications/solvers/multiphase/multiphaseInterFoam/multiphaseMixture/multiphaseMixture.C
View File

@ -573,14 +573,14 @@ void Foam::multiphaseMixture::solveAlphas
surfaceScalarField phic(mag(phi_/mesh_.magSf())); surfaceScalarField phic(mag(phi_/mesh_.magSf()));
phic = min(cAlpha*phic, max(phic)); phic = min(cAlpha*phic, max(phic));
PtrList<surfaceScalarField> phiAlphaCorrs(phases_.size()); PtrList<surfaceScalarField> alphaPhiCorrs(phases_.size());
int phasei = 0; int phasei = 0;
forAllIter(PtrDictionary<phase>, phases_, iter) forAllIter(PtrDictionary<phase>, phases_, iter)
{ {
phase& alpha = iter(); phase& alpha = iter();
phiAlphaCorrs.set alphaPhiCorrs.set
( (
phasei, phasei,
new surfaceScalarField new surfaceScalarField
@ -595,7 +595,7 @@ void Foam::multiphaseMixture::solveAlphas
) )
); );
surfaceScalarField& phiAlphaCorr = phiAlphaCorrs[phasei]; surfaceScalarField& alphaPhiCorr = alphaPhiCorrs[phasei];
forAllIter(PtrDictionary<phase>, phases_, iter2) forAllIter(PtrDictionary<phase>, phases_, iter2)
{ {
@ -605,7 +605,7 @@ void Foam::multiphaseMixture::solveAlphas
surfaceScalarField phir(phic*nHatf(alpha, alpha2)); surfaceScalarField phir(phic*nHatf(alpha, alpha2));
phiAlphaCorr += fvc::flux alphaPhiCorr += fvc::flux
( (
-fvc::flux(-phir, alpha2, alpharScheme), -fvc::flux(-phir, alpha2, alpharScheme),
alpha, alpha,
@ -619,7 +619,7 @@ void Foam::multiphaseMixture::solveAlphas
geometricOneField(), geometricOneField(),
alpha, alpha,
phi_, phi_,
phiAlphaCorr, alphaPhiCorr,
zeroField(), zeroField(),
zeroField(), zeroField(),
1, 1,
@ -630,7 +630,7 @@ void Foam::multiphaseMixture::solveAlphas
phasei++; phasei++;
} }
MULES::limitSum(phiAlphaCorrs); MULES::limitSum(alphaPhiCorrs);
rhoPhi_ = dimensionedScalar("0", dimensionSet(1, 0, -1, 0, 0), 0); rhoPhi_ = dimensionedScalar("0", dimensionSet(1, 0, -1, 0, 0), 0);
@ -652,19 +652,19 @@ void Foam::multiphaseMixture::solveAlphas
{ {
phase& alpha = iter(); phase& alpha = iter();
surfaceScalarField& phiAlpha = phiAlphaCorrs[phasei]; surfaceScalarField& alphaPhi = alphaPhiCorrs[phasei];
phiAlpha += upwind<scalar>(mesh_, phi_).flux(alpha); alphaPhi += upwind<scalar>(mesh_, phi_).flux(alpha);
MULES::explicitSolve MULES::explicitSolve
( (
geometricOneField(), geometricOneField(),
alpha, alpha,
phiAlpha, alphaPhi,
zeroField(), zeroField(),
zeroField() zeroField()
); );
rhoPhi_ += phiAlpha*alpha.rho(); rhoPhi_ += alphaPhi*alpha.rho();
Info<< alpha.name() << " volume fraction, min, max = " Info<< alpha.name() << " volume fraction, min, max = "
<< alpha.weightedAverage(mesh_.V()).value() << alpha.weightedAverage(mesh_.V()).value()

6
applications/solvers/multiphase/twoLiquidMixingFoam/alphaEqn.H
View File

@ -1,7 +1,7 @@
{ {
word alphaScheme("div(phi,alpha)"); word alphaScheme("div(phi,alpha)");
surfaceScalarField phiAlpha surfaceScalarField alphaPhi
( (
phi.name() + alpha1.name(), phi.name() + alpha1.name(),
fvc::flux fvc::flux
@ -12,9 +12,9 @@
) )
); );
MULES::explicitSolve(alpha1, phi, phiAlpha, 1, 0); MULES::explicitSolve(alpha1, phi, alphaPhi, 1, 0);
rhoPhi = phiAlpha*(rho1 - rho2) + phi*rho2; rhoPhi = alphaPhi*(rho1 - rho2) + phi*rho2;
Info<< "Phase-1 volume fraction = " Info<< "Phase-1 volume fraction = "
<< alpha1.weightedAverage(mesh.Vsc()).value() << alpha1.weightedAverage(mesh.Vsc()).value()