[pedroamaral]

Hi Foamers,

i'm quite new on OpenFoam, and i'm having some troubles to run a model with simpleFoam.

my simulation runs until 5-15 iterations in average (depending on my relaxation factors) and crashes after that. I'm getting this error message:

Error message:

 
Time = 6

smoothSolver:  Solving for Ux, Initial residual = 0.999149, Final residual = 0.0996028, No Iterations 16
smoothSolver:  Solving for Uy, Initial residual = 0.999832, Final residual = 0.0874885, No Iterations 17
smoothSolver:  Solving for Uz, Initial residual = 0.999298, Final residual = 0.0873607, No Iterations 17
GAMG:  Solving for p, Initial residual = 1, Final residual = nan, No Iterations 1000
time step continuity errors : sum local = nan, global = -nan, cumulative = -nan
smoothSolver:  Solving for epsilon, Initial residual = nan, Final residual = nan, No Iterations 1000
smoothSolver:  Solving for k, Initial residual = nan, Final residual = nan, No Iterations 1000
ExecutionTime = 685.74 s  ClockTime = 744 s

[1] [2] [3] [0]
[3]
[3] --> FOAM FATAL IO ERROR:
[3] wrong token type - expected Scalar, found on line 0 the word 'nan'


[0]
[0] --> FOAM FATAL IO ERROR:
[0] wrong token type - expected Scalar, found on line 0 the word 'nan'
[0]
[0] file: /home/colaborador/Maxzone_v7/processor0/system/data.solverPerformance.p at line [3]
[3] file: /home/colaborador/Maxzone_v7/processor3/system/data.solverPerformance.p at line 0.
[3]
[3]     From function 0.
[0]
[0]     From function Foam::Istream& Foam::operator>>(Foam::Istream&, Foam::doubleScalar&)
[0]     in file lnInclude/Scalar.C at line 93.
[0]
FOAM parallel run exiting
[0]
Foam::Istream& Foam::operator>>(Foam::Istream&, Foam::doubleScalar&)
[3]     in file lnInclude/Scalar.C at line 93.
[3]
FOAM parallel run exiting
[3]

[2]
[2] --> FOAM FATAL IO ERROR:
[2] wrong token type - expected Scalar, found on line 0 the word 'nan'
[2]
[2] file: /home/colaborador/Maxzone_v7/processor2/system/data.solverPerformance.p at line 0.
[2]
[2]     From function Foam::Istream& Foam::operator>>(Foam::Istream&, Foam::doubleScalar&)
[2]     in file lnInclude/Scalar.C at line 93.
[2]
FOAM parallel run exiting
[2]
[1]
[1] --> FOAM FATAL IO ERROR:
[1] wrong token type - expected Scalar, found on line 0 the word 'nan'
[1]
[1] file: /home/colaborador/Maxzone_v7/processor1/system/data.solverPerformance.p at line 0.
[1]
[1]     From function Foam::Istream& Foam::operator>>(Foam::Istream&, Foam::doubleScalar&)
[1]     in file lnInclude/Scalar.C at line 93.
[1]
FOAM parallel run exiting
[1]
--------------------------------------------------------------------------
MPI_ABORT was invoked on rank 1 in communicator MPI COMMUNICATOR 3 SPLIT FROM 0
with errorcode 1.

NOTE: invoking MPI_ABORT causes Open MPI to kill all MPI processes.
You may or may not see output from other processes, depending on
exactly when Open MPI kills them.
--------------------------------------------------------------------------
[KOT81:03349] 3 more processes have sent help message help-mpi-api.txt / mpi-abort
[KOT81:03349] Set MCA parameter "orte_base_help_aggregate" to 0 to see all help / error messages

And before this message, the values of k and epsilon are try to be bounded but they explodes to number like #e^150 and became a NAN. I've tried to change relaxation factors, fvSchemes, change solver and still nothing works. I've read quite thread in this forum about nan and the foamers usually recommend to check BC and IC, i've checked like 10 times and still nothing. Anyone have any suggestion of what could i do? I'm attaching thumbs with my geometry and my fvSolution and fvSchemes.
My boundary conditions are:


Velocity:

dimensions      [0 1 -1 0 0 0 0];

internalField   uniform (0 0 0);

boundaryField
{
    inlet
    {
        type            fixedValue;
        value           uniform (0.31 0 0);
    }

    outlet
    {
        type            zeroGradient;
    }

    upperWall
    {
        type            fixedValue;
        value           uniform (0 0 0);
    }
   
    lowerWall
    {
        type            fixedValue;
        value           uniform (2.92 0 0);
    }

    frontAndBack
    {
        type            zeroGradient;
    }
}

Pressure:

dimensions      [0 2 -2 0 0 0 0];

internalField   uniform 0;

boundaryField
{
    inlet
    {
        type            zeroGradient;
    }

    outlet
    {
        type            fixedValue;
        value           uniform 0;
    }

    upperWall
    {
        type            zeroGradient;
    }

    lowerWall
    {
        type            zeroGradient;
    }

    frontAndBack
    {
        type            fixedValue;
        value           uniform 0;
    }
}

Iteration right before the crash:

Time = 5

smoothSolver:  Solving for Ux, Initial residual = 0.995088, Final residual = 0.0958079, No Iterations 13
smoothSolver:  Solving for Uy, Initial residual = 0.984021, Final residual = 0.0909928, No Iterations 13
smoothSolver:  Solving for Uz, Initial residual = 0.999731, Final residual = 0.0981715, No Iterations 13
GAMG:  Solving for p, Initial residual = 1, Final residual = 0.0801382, No Iterations 3
time step continuity errors : sum local = 2.13492e+39, global = -1.70627e+38, cumulative = -1.70627e+38
smoothSolver:  Solving for epsilon, Initial residual = 1, Final residual = 0.0326797, No Iterations 3
bounding epsilon, min: -6.68363e+72 max: 2.29307e+73 average: 4.19343e+66
smoothSolver:  Solving for k, Initial residual = 1, Final residual = 0.00813179, No Iterations 1
bounding k, min: -1.94819e+57 max: 3.8381e+59 average: 5.65189e+53
ExecutionTime = 105.32 s  ClockTime = 119 s

fvSolution:

solvers
{
    p
    {
        solver          GAMG;
        tolerance       1e-06;
        relTol          0.1;
        smoother        GaussSeidel;
    }

    "(U|k|epsilon|omega|f|v2)"
    {
        solver          smoothSolver;
        smoother        symGaussSeidel;
        tolerance       1e-05;
        relTol          0.1;
    }
}

SIMPLE
{
    nNonOrthogonalCorrectors 0;
    consistent      yes;

    residualControl
    {
        p               1e-2;
        U               1e-3;
        "(k|epsilon|omega|f|v2)" 1e-2;
    }
}

relaxationFactors
{
    equations
    {
        U               0.95; // 0.9 is more stable but 0.95 more convergent
        ".*"            0.95; // 0.9 is more stable but 0.95 more convergent
    }
}

fvSchemes:

ddtSchemes
{
    default         steadyState;
}

gradSchemes
{
    default         Gauss linear;
}

divSchemes
{
    default         none;
    div(phi,U)      bounded Gauss linearUpwind grad(U);
    div(phi,k)      bounded Gauss limitedLinear 1;
    div(phi,epsilon) bounded Gauss limitedLinear 1;
    div(phi,omega)  bounded Gauss limitedLinear 1;
    div(phi,v2)     bounded Gauss limitedLinear 1;
    div((nuEff*dev2(T(grad(U))))) Gauss linear;
    div(nonlinearStress) Gauss linear;
}

laplacianSchemes
{
    default         Gauss linear corrected;
}

interpolationSchemes
{
    default         linear;
}

snGradSchemes
{
    default         corrected;
}

wallDist
{
    method meshWave;
 }

[pete20r2]

Hi Pedro,
A few things:
1. Can you share all of your boundary conditions, especially for k and epsilon.
2. What y+ values are you getting?
3. try changing the div schemes for k and eps back to upwind for now.
4. I think your issue is actually with 'p' field, try these URF (under relaxation factor) settings

relaxationFactors
{
    fields
    {
        p               0.3;
    }
    equations
    {
        "U.*"               0.7;
        k               0.7;
        epsilon     0.7;
    }
}

[pedroamaral]

Hi Peter,

Thanks for the replying.

I've made a few changes in the simulations after i submit this thread and add your suggestions, and it works! Thanks a lot!
I've removed the big rectangles from above (image attached), remove the initial velocity from the inlet (let the velocity from lowerwall input the initial velocity) and made your suggestions.









For the questions you made:
1) Sure, it's attached right bellow:

k

dimensions      [0 2 -2 0 0 0 0];

internalField   uniform 0.000360375;

boundaryField
{
    inlet
    {
        type            zeroGradient;
    }
    outlet
    {
        type            zeroGradient;
    }
    upperWall
    {
        type            kqRWallFunction;
        value           uniform 0.000360375;
    }
    lowerWall
    {
        type            kqRWallFunction;
        value           uniform 0.000360375;
    }
    frontAndBack
    {
        type            zeroGradient;
    }
}

Epsilon

dimensions      [0 2 -3 0 0 0 0];

internalField   uniform 0.0000138269;

boundaryField
{
    inlet
    {
        type            zeroGradient;
    }
    outlet
    {
        type            zeroGradient;
    }
    upperWall
    {
        type            epsilonWallFunction;
        value           uniform 0.0000138269;
    }
    lowerWall
    {
        type            epsilonWallFunction;
        value           uniform 0.0000138269;
    }
    frontAndBack
    {
        type            zeroGradient;
    }
}

U

dimensions      [0 1 -1 0 0 0 0];

internalField   uniform (0 0 0);

boundaryField
{
    inlet
    {
        type            zeroGradient;
    }

    outlet
    {
        type            zeroGradient;
    }

    upperWall
    {
        type            fixedValue;
        value           uniform (0 0 0);
    }
    
    lowerWall
    {
        type            fixedValue;
        value           uniform (2.92 0 0);
    }

    frontAndBack
    {
        type            zeroGradient;
    }
}

Pressure

dimensions      [0 2 -2 0 0 0 0];

internalField   uniform 0;

boundaryField
{
    inlet
    {
        type            zeroGradient;
    }

    outlet
    {
        type            fixedValue;
        value           uniform 0;
    }

    upperWall
    {
        type            zeroGradient;
    }

    lowerWall
    {
        type            zeroGradient;
    }

    frontAndBack
    {
        type            zeroGradient;
    }
}

I'm using this equations to obtain the values of k and epsilon.



--------------------------

2) Y+? I dont get it, what do you mean? (Sorry for the question, i'm really new to OF and CFD in general)

--------------------------

3) Ok! But can you explain why is better run in this schemes? I'm trying to understand the logic behind.

--------------------------

4) Ok! Changed.

--------------------------

I'm still getting "time step continuity errors" but i don't understand why/what is this. I set my time step based on my length and initial velocity (l=10m; u=0,31m/s; Delta t = 0,001s). I think if Delta t is really small, this shouldn't been happening, do you know why is this is happening?

[Rasmusiwersen]

For anyone that might need a second possible solution:

I had the same issue, epsilon exploding which ultimately lead to expoloding velocities.

I've tried updating the mesh several times, although in the end the quality of the mesh wasn't really affected by increasing refinements of surfaces and regions. Truth be told, i didn't try the above suggested solutions, so do try them first.

But what i did was to change the turbulence model from k-epsilon to a k-omega SST which uses both k-epsilon and k-omega essentially. Look more into if changing turbulence model might work for you before you do it though.

Happy foaming.