CFD Online Discussion Forums

CFD Online Discussion Forums (
-   OpenFOAM Running, Solving & CFD (
-   -   interFoam algorithm (

Pedro24 February 18, 2013 06:46

interFoam algorithm
Hi :),

For each time step, I found that the order of the computations is a little bit strange.

1) It computes the new surface force term using the alpha field at the time "n"

2) Compute the alpha field "n+1" using the velocity field U at the time "n"

        #include "alphaEqnSubCycle.H"
3) Compute the pressure-velocity field at the time "n+1" using the surface term force at the time "n"

        // --- Pressure-velocity PIMPLE corrector loop
        while (pimple.loop())
            #include "UEqn.H"

            // --- Pressure corrector loop
            while (pimple.correct())
                #include "pEqn.H"

            if (pimple.turbCorr())

Why don't we inverse the step 1 and 2. It will allows to use the surface term force at the time "n+1" in the navier-stokes equations ?

It seems more logical to update the surface force term (step 1) just after computing the new alpha field. In my mind, the algorithm should be :

1) update properties time "n"
2) compute implicitly velocity-pressure field "n+1"
3) compute alpha field "n+1" using the velocity-pressure previously computed

Am I wrong ?

Thank You,


Pedro24 February 19, 2013 07:20

I confirm what I said before, this is an error in the "interFoam" solver (and the derivated solvers).

In the H. Rusche thesis, page 162, we can see the following solution procedure :

steps 1-4 : refer to the moving frame

step 5 : transport "alpha1"

step 6 : update properties (smooth gamma => curvature computation => surface force term)

step 7 and 8 : PISO-loop to compute pressure-velocity fields.

Maybe this error should be reported to the OpenFOAM Foundation to correct it for the following versions.


akidess February 19, 2013 08:51

How are step 5-8 different than what is already implemented?

Pedro24 February 19, 2013 08:58

In the current version of "interFoam", the step 6 is computed before the step 5.

First it computes the surface force term :


and then it computes the new alpha1 field :


#include "alphaEqnSubCycle.H"
So, the surface force term used in the PISO loop is computed using the old-time alpha field, which is different from the usual solution procedure.

akidess February 19, 2013 10:50

The interface curvature for the surface force term is calculated by calling interface.correct() after the alphaEqnSubCycle.

Pedro24 February 19, 2013 11:03

Ok, I had inverted the "correct" function of twoPhaseMixture for the viscosity with the "correct" function of interfaceProperties for the curvature.

Thank you for your reply

All times are GMT -4. The time now is 00:52.