pEqn.H

Go to the documentation of this file.

{
    volScalarField rAU(1.0/UEqn.A());
    volVectorField HbyA("HbyA", U);
    HbyA = rAU*UEqn.H();
 
 
    // Define coefficients and pseudo-velocities for RCM interpolation
    // M[U] = AU - H = -grad(p)
    // U = H/A - 1/A grad(p)
    // H/A = U + 1/A grad(p)
    surfaceScalarField rhorAUf
    (
        "rhorAUf",
        fvc::interpolate(rho)/fvc::interpolate(UEqn.A())
    );
 
    surfaceVectorField rhoHbyAf
    (
        "rhoHbyAf",
        fvc::interpolate(rho)*fvc::interpolate(U)
      + rhorAUf*fvc::interpolate(fvc::grad(p))
    );
 
    #include "resetBoundaries.H"
 
    if (pimple.nCorrPISO() <= 1)
    {
        tUEqn.clear();
    }
 
    if (pimple.transonic())
    {
         FatalError
             << "\nTransonic option not available for " << args.executable()
             << exit(FatalError);
    }
    else
    {
        // Rhie & Chow interpolation (part 1)
        surfaceScalarField phiHbyA
        (
            "phiHbyA",
            (
                (rhoHbyAf & mesh.Sf())
              + rhorAUf*fvc::interpolate(rho)*fvc::ddtCorr(U, phiByRho)
              + fvc::interpolate(rho)
              * fvc::alphaCorr(U, phiByRho, pimple.finalInnerIter())
            )
        );
 
        MRF.makeRelative(fvc::interpolate(rho), phiHbyA);
 
        // Non-orthogonal pressure corrector loop
        while (pimple.correctNonOrthogonal())
        {
            // Pressure corrector
            fvScalarMatrix pEqn
            (
                fvm::ddt(psi, p)
              + fvc::div(phiHbyA)
              - fvm::laplacian(rhorAUf, p)
              ==
                fvOptions(psi, p, rho.name())
            );
 
            pEqn.solve(mesh.solver(p.select(pimple.finalInnerIter())));
 
            // Rhie & Chow interpolation (part 2)
            if (pimple.finalNonOrthogonalIter())
            {
                phi = phiHbyA + pEqn.flux();
            }
        }
    }
 
    phiByRho = phi/fvc::interpolate(rho);
 
    #include "rhoEqn.H"
    #include "compressibleContinuityErrs.H"
 
    // Explicitly relax pressure for momentum corrector
    p.relax();
 
    U = HbyA - rAU*fvc::grad(p);
    U.correctBoundaryConditions();
    fvOptions.correct(U);
}
 
rho = thermo.rho();