applications/solvers/multiphase/MRFInterFoam/pEqn.H

   1 {
   2     volScalarField rUA = 1.0/UEqn.A();
   3     surfaceScalarField rUAf = fvc::interpolate(rUA);
   4
   5     U = rUA*UEqn.H();
   6
   7     surfaceScalarField phiU
   8     (
   9         "phiU",
  10         (fvc::interpolate(U) & mesh.Sf())
  11     //+ fvc::ddtPhiCorr(rUA, rho, U, phi)
  12     );
  13     mrfZones.relativeFlux(phiU);
  14
  15     phi = phiU +
  16         (
  17             fvc::interpolate(interface.sigmaK())*
  18             fvc::snGrad(alpha1)*mesh.magSf()
  19           + fvc::interpolate(rho)*(g & mesh.Sf())
  20         )*rUAf;
  21     adjustPhi(phi, U, p);
  22
  23     for(int nonOrth=0; nonOrth<=nNonOrthCorr; nonOrth++)
  24     {
  25         fvScalarMatrix pdEqn
  26         (
  27             fvm::laplacian(rUAf, pd) == fvc::div(phi)
  28         );
  29
  30         pdEqn.setReference(pdRefCell, pdRefValue);
  31
  32         if (corr == nCorr-1 && nonOrth == nNonOrthCorr)
  33         {
  34             pdEqn.solve(mesh.solutionDict().solver(pd.name() + "Final"));
  35         }
  36         else
  37         {
  38             pdEqn.solve(mesh.solutionDict().solver(pd.name()));
  39         }
  40
  41         if (nonOrth == nNonOrthCorr)
  42         {
  43             phi -= pdEqn.flux();
  44         }
  45     }
  46
  47     U += rUA*fvc::reconstruct((phi - phiU)/rUAf);
  48     U.correctBoundaryConditions();
  49 }