[WaterHammer1985]

Hello,

I am running a simple open channel flow simulation using interFoam. At the outlet, the max phase 1 volume fraction is greater than 1. After some time of running (+1 min flow time) the solution goes unstable and begins to diverge. Any ideas on what causes that and if there is a way to keep the alpha.water bounded between 0 and 1?

Boundary Conditions and fv files are below:

U

boundaryField
{
    inlet
    {
        type            variableHeightFlowRateInletVelocity;
        flowRate        constant 100;
        alpha           alpha.water;
        value           uniform (0 0 0);
    }
    outlet
    {
        type            pressureInletOutletVelocity;
        value           uniform (0 0 0);
    }
    concrete
    {
        type            fixedValue;
        value           uniform (0 0 0);
    }
    bed
    {
        type            fixedValue;
        value           uniform (0 0 0);
    }
    wall
    {
        type            slip;
    }
    top
    {
        type            pressureInletOutletVelocity;
        value           uniform (0 0 0);
    }
}

p_rgh:

boundaryField
{
    inlet       {type            zeroGradient;}

    outlet
    {
        type            fixedValue;
        value           uniform 2673221;
    }

    wall        {type            zeroGradient;}
    bed         {type            zeroGradient;}
    concrete    {type            zeroGradient;}

    top
    {
        type            totalPressure;
        gamma           1;
        p0              uniform 0;
        value           uniform 0;
    }

    defaultFaces
    {
        type            empty;
    }
}

alpha.water:

boundaryField
{
    inlet
    {
        type            variableHeightFlowRate;
        lowerBound      0;
        upperBound      1;
        value           uniform 0;
    }

    outlet
    {
        type            inletOutlet;
        inletValue      uniform 0;
        value           uniform 1;
    }

    wall
    {
        type            zeroGradient;
    }

    bed
    {
        type            zeroGradient;
    }

    concrete
    {
        type            zeroGradient;
    }


    top
    {
        type            inletOutlet;
        inletValue      uniform 0;
        value           uniform 0;
    }

    defaultFaces
    {
        type            empty;
    }
}

fvSolution:

solvers
{
    "alpha.water.*"
    {
        nAlphaCorr      2;
        nAlphaSubCycles 2;
        alphaOuterCorrectors yes;
        cAlpha          1;

        MULESCorr       yes;
        nLimiterIter    3;

        solver          smoothSolver;
        smoother        symGaussSeidel;
        tolerance       1e-8;
        relTol          0;
    }

    pcorr GAMG
    {
        tolerance        1e-7;
        relTol           0.0;

        smoother         DIC;//GaussSeidel;
        nPreSweeps       0;
        nPostSweeps      2;
        nFinestSweeps    2;

        cacheAgglomeration true;
        nCellsInCoarsestLevel 10;
        agglomerator     faceAreaPair;
        mergeLevels      1;
        maxIter 10;
    };

    p_rgh GAMG
    {
        tolerance        1e-7;
        relTol           0.0;

        smoother         DIC;//GaussSeidel;
        nPreSweeps       0;
        nPostSweeps      2;
        nFinestSweeps    2;

        cacheAgglomeration true;
        nCellsInCoarsestLevel 10;
        agglomerator     faceAreaPair;
        mergeLevels      1;
    };
    p_rghFinal GAMG
    {
        tolerance        1e-8;
        relTol           0.0;

        smoother         DIC;//GaussSeidel;
        nPreSweeps       0;
        nPostSweeps      2;
        nFinestSweeps    2;

        cacheAgglomeration true;
        nCellsInCoarsestLevel 10;
        agglomerator     faceAreaPair;
        mergeLevels      1;
    };



    "(U|k|omega).*"
    {
        solver          smoothSolver;
        smoother        symGaussSeidel;
        nSweeps         1;
        tolerance       1e-6;
        relTol          0.1;
    };
}


PIMPLE
{
    momentumPredictor no;
    nOuterCorrectors 1;
    nCorrectors     1;
    nNonOrthogonalCorrectors 2;
}

relaxationFactors
{
    fields
    {
    }
    equations
    {
        ".*" 0.1;
    }
}

fvSchemes:

ddtSchemes
{
    default         Euler;
}

gradSchemes
{
  default cellMDLimited Gauss linear 1;
  grad(k) faceLimited leastSquares 1;
  grad(omega) faceLimited leastSquares 1;
}

divSchemes
{
    div(rhoPhi,U)  Gauss upwind;
    div(phi,alpha)  Gauss vanLeer;
    div(phirb,alpha) Gauss interfaceCompression 1;
    div(phi,p_rgh)  Gauss upwind;
    div(phi,k)      Gauss upwind;//vanLeer;
    div(((rho*nuEff)*dev2(T(grad(U))))) Gauss linear;
    div(phi,omega) Gauss upwind;
}

laplacianSchemes
{
    default         Gauss linear uncorrected;
}

interpolationSchemes
{
    default         linear;
}

snGradSchemes
{
    default         uncorrected;
}

fluxRequired
{
    default         no;
    p_rgh;
    pcorr;
    gamma;
}

wallDist
{
    method meshWave;
}

[WaterHammer1985]

Output Log:

Courant Number mean: 0.017532 max: 0.618729
Interface Courant Number mean: 0.00174773 max: 0.371229
deltaT = 0.00198997
Time = 600

PIMPLE: iteration 1
MULES: Solving for alpha.water
Phase-1 volume fraction = 0.315577  Min(alpha.water) = -2.47033e-323  Max(alpha.water) = 1.99015
MULES: Solving for alpha.water
Phase-1 volume fraction = 0.315545  Min(alpha.water) = -4.12902e-13  Max(alpha.water) = 1.98918
DICPCG:  Solving for p_rgh, Initial residual = 0.0021341, Final residual = 5.86447e-05, No Iterations 2
DICPCG:  Solving for p_rgh, Initial residual = 5.86446e-05, Final residual = 2.66273e-06, No Iterations 84
DICPCG:  Solving for p_rgh, Initial residual = 2.59791e-06, Final residual = 1.22832e-07, No Iterations 45
time step continuity errors : sum local = 2.96715e-08, global = -1.10132e-08, cumulative = 0.00124093
DICPCG:  Solving for p_rgh, Initial residual = 5.60846e-05, Final residual = 2.76946e-06, No Iterations 83
DICPCG:  Solving for p_rgh, Initial residual = 2.87651e-06, Final residual = 1.39745e-07, No Iterations 46
DICPCG:  Solving for p_rgh, Initial residual = 1.39763e-07, Final residual = 9.93079e-08, No Iterations 9
time step continuity errors : sum local = 2.30964e-08, global = 1.20769e-08, cumulative = 0.00124095
DICPCG:  Solving for p_rgh, Initial residual = 8.74883e-06, Final residual = 4.11144e-07, No Iterations 12
DICPCG:  Solving for p_rgh, Initial residual = 4.11166e-07, Final residual = 9.70768e-08, No Iterations 73
DICPCG:  Solving for p_rgh, Initial residual = 9.71144e-08, Final residual = 9.71144e-08, No Iterations 0
time step continuity errors : sum local = 2.25762e-08, global = -1.12547e-09, cumulative = 0.00124094
ExecutionTime = 12325.9 s  ClockTime = 12533 s

[WaterHammer1985]

Results output is attached showing location in domain where the issue is occurring.

[Bodo1993]

Hi,

I am wondering if you figured out a solution to your problem.

Thanks.

[WaterHammer1985]

Originally Posted by Bodo1993

Hi,

I am wondering if you figured out a solution to your problem.

Thanks.

Use a "totalPressure" boundary condition for the outlet.

[Bodo1993]

Thanks for the reply