applications/solvers/compressible/rhoSimpleFoam/rhoPorousMRFSimpleFoam/pEqn.H

   1 if (pressureImplicitPorosity)
   2 {
   3     U = trTU()&UEqn().H();
   4 }
   5 else
   6 {
   7     U = trAU()*UEqn().H();
   8 }
   9
  10 UEqn.clear();
  11
  12 bool closedVolume = false;
  13
  14 if (simple.transonic())
  15 {
  16     surfaceScalarField phid
  17     (
  18         "phid",
  19         fvc::interpolate(psi)*(fvc::interpolate(U) & mesh.Sf())
  20     );
  21     mrfZones.relativeFlux(fvc::interpolate(psi), phid);
  22
  23     for (int nonOrth=0; nonOrth<=simple.nNonOrthCorr(); nonOrth++)
  24     {
  25         tmp<fvScalarMatrix> tpEqn;
  26
  27         if (pressureImplicitPorosity)
  28         {
  29             tpEqn = (fvc::div(phid, p) - fvm::laplacian(rho*trTU(), p));
  30         }
  31         else
  32         {
  33             tpEqn = (fvc::div(phid, p) - fvm::laplacian(rho*trAU(), p));
  34         }
  35
  36         tpEqn().setReference(pRefCell, pRefValue);
  37
  38         tpEqn().solve();
  39
  40         if (nonOrth == simple.nNonOrthCorr())
  41         {
  42             phi == tpEqn().flux();
  43         }
  44     }
  45 }
  46 else
  47 {
  48     phi = fvc::interpolate(rho*U) & mesh.Sf();
  49     mrfZones.relativeFlux(fvc::interpolate(rho), phi);
  50
  51     closedVolume = adjustPhi(phi, U, p);
  52
  53     for (int nonOrth=0; nonOrth<=simple.nNonOrthCorr(); nonOrth++)
  54     {
  55         tmp<fvScalarMatrix> tpEqn;
  56
  57         if (pressureImplicitPorosity)
  58         {
  59             tpEqn = (fvm::laplacian(rho*trTU(), p) == fvc::div(phi));
  60         }
  61         else
  62         {
  63             tpEqn = (fvm::laplacian(rho*trAU(), p) == fvc::div(phi));
  64         }
  65
  66         tpEqn().setReference(pRefCell, pRefValue);
  67
  68         tpEqn().solve();
  69
  70         if (nonOrth == simple.nNonOrthCorr())
  71         {
  72             phi -= tpEqn().flux();
  73         }
  74     }
  75 }
  76
  77 #include "incompressible/continuityErrs.H"
  78
  79 // Explicitly relax pressure for momentum corrector
  80 p.relax();
  81
  82 if (pressureImplicitPorosity)
  83 {
  84     U -= trTU()&fvc::grad(p);
  85 }
  86 else
  87 {
  88     U -= trAU()*fvc::grad(p);
  89 }
  90
  91 U.correctBoundaryConditions();
  92
  93 // For closed-volume cases adjust the pressure and density levels
  94 // to obey overall mass continuity
  95 if (closedVolume)
  96 {
  97     p += (initialMass - fvc::domainIntegrate(psi*p))
  98         /fvc::domainIntegrate(psi);
  99 }
 100
 101 rho = thermo.rho();
 102 rho = max(rho, rhoMin);
 103 rho = min(rho, rhoMax);
 104 rho.relax();
 105 Info<< "rho max/min : " << max(rho).value() << " " << min(rho).value() << endl;