Commit 25040ad0 authored by Henry Weller's avatar Henry Weller
Browse files

reactingTwoPhaseEulerFoam: Update cell-based UEqn

to use the transport matrix construction in MovingPhaseModel
parent 557c8adb
...@@ -31,5 +31,6 @@ PtrList<fvVectorMatrix> UEqns(phases.size()); ...@@ -31,5 +31,6 @@ PtrList<fvVectorMatrix> UEqns(phases.size());
UEqns[phasei].relax(); UEqns[phasei].relax();
fvOptions.constrain(UEqns[phasei]); fvOptions.constrain(UEqns[phasei]);
fvOptions.correct(U);
} }
} }
...@@ -3,61 +3,36 @@ Info<< "Constructing momentum equations" << endl; ...@@ -3,61 +3,36 @@ Info<< "Constructing momentum equations" << endl;
fvVectorMatrix U1Eqn(U1, rho1.dimensions()*U1.dimensions()*dimVol/dimTime); fvVectorMatrix U1Eqn(U1, rho1.dimensions()*U1.dimensions()*dimVol/dimTime);
fvVectorMatrix U2Eqn(U2, rho2.dimensions()*U2.dimensions()*dimVol/dimTime); fvVectorMatrix U2Eqn(U2, rho2.dimensions()*U2.dimensions()*dimVol/dimTime);
volScalarField Kd(fluid.Kd());
{ {
volScalarField Vm(fluid.Vm()); autoPtr<phaseSystem::momentumTransferTable>
momentumTransferPtr(fluid.momentumTransfer());
const volScalarField dmdt12(posPart(fluid.dmdt())); phaseSystem::momentumTransferTable&
const volScalarField dmdt21(negPart(fluid.dmdt())); momentumTransfer(momentumTransferPtr());
{ {
U1Eqn = U1Eqn =
( (
fvm::ddt(alpha1, rho1, U1) + fvm::div(alphaRhoPhi1, U1) phase1.UEqn()
- fvm::Sp(phase1.continuityError() + dmdt21, U1) + dmdt21*U2
+ MRF.DDt(alpha1*rho1, U1)
+ phase1.turbulence().divDevRhoReff(U1)
== ==
- Vm *momentumTransfer[phase1.name()]
*(
fvm::ddt(U1)
+ fvm::div(phi1, U1)
- fvm::Sp(fvc::div(phi1), U1)
- phase2.DUDt()
)
- MRF.DDt(Vm, U1 - U2)
+ fvOptions(alpha1, rho1, U1) + fvOptions(alpha1, rho1, U1)
); );
U1Eqn.relax(); U1Eqn.relax();
U1Eqn += fvm::Sp(Kd, U1);
fvOptions.constrain(U1Eqn); fvOptions.constrain(U1Eqn);
U1.correctBoundaryConditions();
fvOptions.correct(U1); fvOptions.correct(U1);
} }
{ {
U2Eqn = U2Eqn =
( (
fvm::ddt(alpha2, rho2, U2) + fvm::div(alphaRhoPhi2, U2) phase2.UEqn()
- fvm::Sp(phase2.continuityError() - dmdt12, U2) - dmdt12*U1
+ MRF.DDt(alpha2*rho2, U2)
+ phase2.turbulence().divDevRhoReff(U2)
== ==
- Vm *momentumTransfer[phase2.name()]
*(
fvm::ddt(U2)
+ fvm::div(phi2, U2)
- fvm::Sp(fvc::div(phi2), U2)
- phase1.DUDt()
)
- MRF.DDt(Vm, U2 - U1)
+ fvOptions(alpha2, rho2, U2) + fvOptions(alpha2, rho2, U2)
); );
U2Eqn.relax(); U2Eqn.relax();
U2Eqn += fvm::Sp(Kd, U2);
fvOptions.constrain(U2Eqn); fvOptions.constrain(U2Eqn);
U2.correctBoundaryConditions();
fvOptions.correct(U2); fvOptions.correct(U2);
} }
} }
...@@ -31,6 +31,8 @@ surfaceScalarField alpharAUf2 ...@@ -31,6 +31,8 @@ surfaceScalarField alpharAUf2
fvc::interpolate(max(alpha2, phase2.residualAlpha())*rAU2) fvc::interpolate(max(alpha2, phase2.residualAlpha())*rAU2)
); );
volScalarField Kd(fluid.Kd());
// Turbulent diffusion, particle-pressure, lift and wall-lubrication fluxes // Turbulent diffusion, particle-pressure, lift and wall-lubrication fluxes
tmp<surfaceScalarField> phiF1; tmp<surfaceScalarField> phiF1;
tmp<surfaceScalarField> phiF2; tmp<surfaceScalarField> phiF2;
......
Supports Markdown
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment