About underrelaxation of pressure in OpenFOAM solvers

I am curious about when p should be relaxed after solved in OpenFOAM's fluid solvers.

pEqn.H exists in each solvers' source code directory contains solving and relaxing pressure, and correcting velocity. So let's take a look at pEqn.H of simpleFoam because it is simplest.

pEqn.H of simpleFoam (OF6)

 1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
{
    volScalarField rAU(1.0/UEqn.A());
    volVectorField HbyA(constrainHbyA(rAU*UEqn.H(), U, p));
    surfaceScalarField phiHbyA("phiHbyA", fvc::flux(HbyA));
    MRF.makeRelative(phiHbyA);
    adjustPhi(phiHbyA, U, p);

    tmp<volScalarField> rAtU(rAU);

    if (simple.consistent())
    {
        rAtU = 1.0/(1.0/rAU - UEqn.H1());
        phiHbyA +=
            fvc::interpolate(rAtU() - rAU)*fvc::snGrad(p)*mesh.magSf();
        HbyA -= (rAU - rAtU())*fvc::grad(p);
    }

    tUEqn.clear();

    // Update the pressure BCs to ensure flux consistency
    constrainPressure(p, U, phiHbyA, rAtU(), MRF);

    // Non-orthogonal pressure corrector loop
    while (simple.correctNonOrthogonal())
    {
        fvScalarMatrix pEqn
        (
            fvm::laplacian(rAtU(), p) == fvc::div(phiHbyA)
        );

        pEqn.setReference(pRefCell, pRefValue);

        pEqn.solve();

        if (simple.finalNonOrthogonalIter())
        {
            phi = phiHbyA - pEqn.flux();
        }
    }

    #include "continuityErrs.H"

    // Explicitly relax pressure for momentum corrector
    p.relax();

    // Momentum corrector
    U = HbyA - rAtU()*fvc::grad(p);
    U.correctBoundaryConditions();
    fvOptions.correct(U);
}

Because I already explained the algorithm in a past article, I don't mention here about Line. 1-18 except that if simple.consistent() is true, it is SIMPLEC, otherwise it is SIMPLE.

Pressure is solved in Line. 26-33, and phi is corrected using unrelaxed pressure in Line. 37.

And then pressure is relaxed in Line. 44 before velocity is corrected using relaxed pressure in Line. 47.

What I'm curious about is, why velocity is corrected based on the relaxed pressure. In order to meet the continuity, phi must be corrected by unrelaxed pressure. So Line. 37 is correct.

And also, phi is theoretically the product of rho and U. So it seems preferable to me to correct U with unrelaxed pressure as well.

It might be for stability, but I haven't seen any blow up issue as far as I tested with Ahmed body and racing car simulations.

Comments