applications/solvers/multiphase/interPhaseChangeFoam/pEqn.H

   1 {
   2     volScalarField rAU(1.0/UEqn.A());
   3     surfaceScalarField rAUf(fvc::interpolate(rAU));
   4
   5     U = rAU*UEqn.H();
   6
   7     surfaceScalarField phiU
   8     (
   9         "phiU",
  10         (fvc::interpolate(U) & mesh.Sf())
  11       + fvc::ddtPhiCorr(rAU, rho, U, phi)
  12     );
  13
  14     adjustPhi(phiU, U, p_rgh);
  15
  16     phi = phiU +
  17         (
  18             fvc::interpolate(interface.sigmaK())*fvc::snGrad(alpha1)
  19           - ghf*fvc::snGrad(rho)
  20         )*rAUf*mesh.magSf();
  21
  22     Pair<tmp<volScalarField> > vDotP = twoPhaseProperties->vDotP();
  23     const volScalarField& vDotcP = vDotP[0]();
  24     const volScalarField& vDotvP = vDotP[1]();
  25
  26     for (int nonOrth=0; nonOrth<=pimple.nNonOrthCorr(); nonOrth++)
  27     {
  28         fvScalarMatrix p_rghEqn
  29         (
  30             fvc::div(phi) - fvm::laplacian(rAUf, p_rgh)
  31           - (vDotvP - vDotcP)*(pSat - rho*gh) + fvm::Sp(vDotvP - vDotcP, p_rgh)
  32         );
  33
  34         p_rghEqn.setReference(pRefCell, pRefValue);
  35
  36         p_rghEqn.solve
  37         (
  38             mesh.solver(p_rgh.select(pimple.finalInnerIter(corr, nonOrth)))
  39         );
  40
  41         if (nonOrth == pimple.nNonOrthCorr())
  42         {
  43             phi += p_rghEqn.flux();
  44         }
  45     }
  46
  47     U += rAU*fvc::reconstruct((phi - phiU)/rAUf);
  48     U.correctBoundaryConditions();
  49
  50     #include "continuityErrs.H"
  51
  52     p == p_rgh + rho*gh;
  53
  54     if (p_rgh.needReference())
  55     {
  56         p += dimensionedScalar
  57         (
  58             "p",
  59             p.dimensions(),
  60             pRefValue - getRefCellValue(p, pRefCell)
  61         );
  62         p_rgh = p - rho*gh;
  63     }
  64 }