applications/solvers/combustion/fireFoam/pEqn.H

   1 bool closedVolume = false;
   2
   3 rho = thermo.rho();
   4
   5 volScalarField rUA = 1.0/UEqn.A();
   6 surfaceScalarField rhorUAf("(rho*(1|A(U)))", fvc::interpolate(rho*rUA));
   7 U = rUA*UEqn.H();
   8
   9 surfaceScalarField phiU
  10 (
  11     fvc::interpolate(rho)
  12    *(
  13         (fvc::interpolate(U) & mesh.Sf())
  14       + fvc::ddtPhiCorr(rUA, rho, U, phi)
  15     )
  16 );
  17
  18 phi = phiU + rhorUAf*fvc::interpolate(rho)*(g & mesh.Sf());
  19
  20 for (int nonOrth=0; nonOrth<=nNonOrthCorr; nonOrth++)
  21 {
  22        surfaceScalarField rhorUAf = fvc::interpolate(rho*rUA);
  23
  24         fvScalarMatrix pEqn
  25         (
  26             fvm::ddt(psi,p)
  27           + fvc::div(phi)
  28           - fvm::laplacian(rhorUAf, p)
  29         );
  30
  31         closedVolume = p.needReference();
  32
  33         if (corr == nCorr-1 && nonOrth == nNonOrthCorr)
  34         {
  35             pEqn.solve(mesh.solver(p.name() + "Final"));
  36         }
  37         else
  38         {
  39             pEqn.solve(mesh.solver(p.name()));
  40         }
  41
  42         if (nonOrth == nNonOrthCorr)
  43         {
  44             phi += pEqn.flux();
  45         }
  46 }
  47
  48 DpDt = fvc::DDt(surfaceScalarField("phiU", phi/fvc::interpolate(rho)), p);
  49
  50 #include "rhoEqn.H"
  51 #include "compressibleContinuityErrs.H"
  52
  53 U += rUA*fvc::reconstruct((phi - phiU)/rhorUAf);
  54 U.correctBoundaryConditions();
  55
  56 // For closed-volume cases adjust the pressure and density levels
  57 // to obey overall mass continuity
  58 if (closedVolume)
  59 {
  60     p +=
  61         (initialMass - fvc::domainIntegrate(thermo.psi()*p))
  62        /fvc::domainIntegrate(thermo.psi());
  63     rho = thermo.rho();
  64 }