CFD Online Discussion Forums

CFD Online Discussion Forums (
-   OpenFOAM Programming & Development (
-   -   Which procedure to solve Ueqn without momentum predictor. (

sharonyue July 2, 2013 21:04

Which procedure to solve Ueqn without momentum predictor.
This is icoFoam's code,:


fvVectorMatrix UEqn
          + fvm::div(phi, U)
          - fvm::laplacian(nu, U)

        solve(UEqn == -fvc::grad(p));//this term is only executed when the momentumPredictor flag is set to true

        // --- PISO loop

        for (int corr=0; corr<nCorr; corr++)
            volScalarField rAU(1.0/UEqn.A());

            volVectorField HbyA("HbyA", U);
            HbyA = rAU*UEqn.H();
            surfaceScalarField phiHbyA
                (fvc::interpolate(HbyA) & mesh.Sf())
              + fvc::ddtPhiCorr(rAU, U, phi)

            adjustPhi(phiHbyA, U, p);

            for (int nonOrth=0; nonOrth<=nNonOrthCorr; nonOrth++)
                fvScalarMatrix pEqn
                    fvm::laplacian(rAU, p) == fvc::div(phiHbyA)

My confusion is if I set the momentum predictor is off, which line is to solve U?

santiagomarquezd July 2, 2013 23:37

Hi, in the PISO method, U is predicted and then corrected in the pressure-velocity coupling loop, so that, if you don't predict the velocity, the velocity used to start the coupling loop is a derivative of the velocity in the last time step as is done in interFoam and other solvers.

In the code you've cited the matrix for the momentum equation is assembled so that you have the A and H methods available for the use in the PISO loop, the line:

HbyA = rAU*UEqn.H();

calculates a first approximation for U using the H operator applied to the velocity of the last-timestep or which was predicted when you solved the UEqn (if you did it) and the A operator. If you need it, please check the definition of these operators in chapter 2.7 of my Ph.D. thesis

The link to my PhD thesis is provided in my user page at

and Jasak's Ph.D. and other's thesis.


All times are GMT -4. The time now is 18:24.