simpleFoam not converging?

Urchin · July 27, 2022, 10:50

Hi,

I've been trying to get a custom case running using simpleFoam - the case is in theory quite simple, one fixed pressure inlet and one fixed velocity outlet (that I will later swap in for a fixed volumetric flow rate outlet).

The simpleFoam solver seems to end early without convering, output below.

Code:

// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
Create time

Create mesh for time = 0


SIMPLE: convergence criteria
    field p     tolerance 0.01
    field U     tolerance 0.001
    field "(k|epsilon|omega|f|v2)"     tolerance 0.001

Reading field p

Reading field U

Reading/calculating face flux field phi

Selecting incompressible transport model Newtonian
Selecting turbulence model type laminar
Selecting laminar stress model Stokes
No MRF models present

No finite volume options present

Starting time loop

Time = 1

smoothSolver:  Solving for Ux, Initial residual = 1, Final residual = 0.0294685, No Iterations 2
smoothSolver:  Solving for Uy, Initial residual = 1, Final residual = 0.0668456, No Iterations 1
smoothSolver:  Solving for Uz, Initial residual = 1, Final residual = 0.0747838, No Iterations 2
GAMG:  Solving for p, Initial residual = 1, Final residual = 0.0184127, No Iterations 2
GAMG:  Solving for p, Initial residual = 0.452291, Final residual = 0.0186042, No Iterations 1
GAMG:  Solving for p, Initial residual = 0.0913773, Final residual = 0.00562931, No Iterations 1
time step continuity errors : sum local = 8.32867e-07, global = -2.72535e-14, cumulative = -2.72535e-14
ExecutionTime = 0.058 s  ClockTime = 0 s

Time = 2

smoothSolver:  Solving for Ux, Initial residual = 0.695821, Final residual = 0.0251763, No Iterations 2
smoothSolver:  Solving for Uy, Initial residual = 0.632267, Final residual = 0.0188064, No Iterations 2
smoothSolver:  Solving for Uz, Initial residual = 0.360745, Final residual = 0.0117506, No Iterations 2
GAMG:  Solving for p, Initial residual = 0.567775, Final residual = 0.0469497, No Iterations 1
GAMG:  Solving for p, Initial residual = 0.140186, Final residual = 0.0123325, No Iterations 1
GAMG:  Solving for p, Initial residual = 0.034633, Final residual = 0.00218083, No Iterations 2
time step continuity errors : sum local = 2.41008e-07, global = 5.57529e-12, cumulative = 5.54804e-12
ExecutionTime = 0.069 s  ClockTime = 0 s

Time = 3

smoothSolver:  Solving for Ux, Initial residual = 0.567759, Final residual = 0.0244518, No Iterations 2
smoothSolver:  Solving for Uy, Initial residual = 0.491631, Final residual = 0.027497, No Iterations 2
smoothSolver:  Solving for Uz, Initial residual = 0.281927, Final residual = 0.00766034, No Iterations 2
GAMG:  Solving for p, Initial residual = 0.136629, Final residual = 0.0119827, No Iterations 1
GAMG:  Solving for p, Initial residual = 0.0586256, Final residual = 0.004022, No Iterations 1
GAMG:  Solving for p, Initial residual = 0.0136292, Final residual = 0.000861313, No Iterations 2
time step continuity errors : sum local = 1.42179e-07, global = 1.59204e-11, cumulative = 2.14684e-11
ExecutionTime = 0.081 s  ClockTime = 0 s

Time = 4

smoothSolver:  Solving for Ux, Initial residual = 0.424208, Final residual = 0.00987025, No Iterations 2
smoothSolver:  Solving for Uy, Initial residual = 0.397932, Final residual = 0.0125944, No Iterations 2
smoothSolver:  Solving for Uz, Initial residual = 0.120651, Final residual = 0.00666165, No Iterations 2
GAMG:  Solving for p, Initial residual = 0.0942463, Final residual = 0.00665727, No Iterations 1
GAMG:  Solving for p, Initial residual = 0.0289666, Final residual = 0.00104403, No Iterations 2
GAMG:  Solving for p, Initial residual = 0.0095302, Final residual = 0.000736445, No Iterations 1
time step continuity errors : sum local = 1.34236e-07, global = 2.93937e-11, cumulative = 5.08622e-11
ExecutionTime = 0.093 s  ClockTime = 0 s


*** continues a bit ***


Time = 113

smoothSolver:  Solving for Ux, Initial residual = 0.643398, Final residual = 0.0500197, No Iterations 12
smoothSolver:  Solving for Uy, Initial residual = 0.65007, Final residual = 0.0483867, No Iterations 5
smoothSolver:  Solving for Uz, Initial residual = 0.71876, Final residual = 0.0614276, No Iterations 4
GAMG:  Solving for p, Initial residual = 0.732957, Final residual = 0.0433888, No Iterations 2
GAMG:  Solving for p, Initial residual = 0.608711, Final residual = 0.0317925, No Iterations 1
GAMG:  Solving for p, Initial residual = 0.499677, Final residual = 0.0171388, No Iterations 1
time step continuity errors : sum local = 1.66446, global = -0.0479148, cumulative = -0.0450688
ExecutionTime = 1.369 s  ClockTime = 1 s

Time = 114

smoothSolver:  Solving for Ux, Initial residual = 0.656116, Final residual = 0.0442092, No Iterations 28
smoothSolver:  Solving for Uy, Initial residual = 0.644418, Final residual = 0.0452629, No Iterations 34
smoothSolver:  Solving for Uz, Initial residual = 0.631059, Final residual = 0.0511848, No Iterations 31
GAMG:  Solving for p, Initial residual = 0.9998, Final residual = 0.0843136, No Iterations 3
GAMG:  Solving for p, Initial residual = 0.994701, Final residual = 0.000595515, No Iterations 1
GAMG:  Solving for p, Initial residual = 0.899573, Final residual = 7.32625e-05, No Iterations 1
time step continuity errors : sum local = 3.97051e+06, global = -10147.3, cumulative = -10147.4
ExecutionTime = 1.387 s  ClockTime = 1 s

Time = 115

smoothSolver:  Solving for Ux, Initial residual = 0.740435, Final residual = 0.0571045, No Iterations 47
smoothSolver:  Solving for Uy, Initial residual = 0.740483, Final residual = 0.0687116, No Iterations 47
smoothSolver:  Solving for Uz, Initial residual = 0.741827, Final residual = 0.0724437, No Iterations 52
GAMG:  Solving for p, Initial residual = 1, Final residual = 0.0998332, No Iterations 12
GAMG:  Solving for p, Initial residual = 1, Final residual = 4.85116e-10, No Iterations 1
GAMG:  Solving for p, Initial residual = 0.439633, Final residual = 9.28016e-11, No Iterations 1
time step continuity errors : sum local = 5.78855e+18, global = 7.63848e+17, cumulative = 7.63848e+17
ExecutionTime = 1.413 s  ClockTime = 1 s

Time = 116

smoothSolver:  Solving for Ux, Initial residual = 0.614144, Final residual = 2.64685e+06, No Iterations 1000
smoothSolver:  Solving for Uy, Initial residual = 0.612062, Final residual = 3.58646e+06, No Iterations 1000
smoothSolver:  Solving for Uz, Initial residual = 0.76164, Final residual = 8.99608e+06, No Iterations 1000

I am looking to start with laminar flow so I have made the following modifications to turbulenceProperties:

- Change simulation type to Laminar
- Change tubulence to off

Code:

/*--------------------------------*- C++ -*----------------------------------*\
| =========                 |                                                 |
| \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox           |
|  \\    /   O peration     | Version:  v2206                                 |
|   \\  /    A nd           | Website:  www.openfoam.com                      |
|    \\/     M anipulation  |                                                 |
\*---------------------------------------------------------------------------*/
FoamFile
{
    version     2.0;
    format      ascii;
    class       dictionary;
    object      turbulenceProperties;
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

simulationType      laminar;

RAS
{
    // Tested with kEpsilon, realizableKE, kOmega, kOmegaSST,
    // ShihQuadraticKE, LienCubicKE.
    RASModel        kEpsilon;

    turbulence      off;

    printCoeffs     on;
}


// ************************************************************************* //

And the output of checkMesh is as follows (unfortunately I am unable to publicly post the mesh)

Code:


// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
Create time

Create mesh for time = 0

Time = 0

Mesh stats 
    points:           2052
    faces:            13345
    internal faces:   9483
    cells:            5707
    faces per cell:   4
    boundary patches: 3
    point zones:      0
    face zones:       0
    cell zones:       1

Overall number of cells of each type:
    hexahedra:     0
    prisms:        0
    wedges:        0
    pyramids:      0
    tet wedges:    0
    tetrahedra:    5707
    polyhedra:     0

Checking topology...
    Boundary definition OK.
    Cell to face addressing OK.
    Point usage OK.
    Upper triangular ordering OK.
    Face vertices OK.
    Number of regions: 1 (OK).

Checking patch topology for multiply connected surfaces...
    Patch               Faces    Points   Surface topology                  
    inlet               12       11       ok (non-closed singly connected)  
    outlet              30       25       ok (non-closed singly connected)  
    defaultFaces        3820     1921     ok (non-closed singly connected)  

Checking faceZone topology for multiply connected surfaces...
    No faceZones found.

Checking basic cellZone addressing...
    CellZone            Cells        Points       Volume       BoundingBox
    fluid               5707         2052         631.946      (-8.84167 -6.3597 -5.27608) (4.55833 56.7699 3.56623)

Checking geometry...
    Overall domain bounding box (-8.84167 -6.3597 -5.27608) (4.55833 56.7699 3.56623)
    Mesh has 3 geometric (non-empty/wedge) directions (1 1 1)
    Mesh has 3 solution (non-empty) directions (1 1 1)
    Boundary openness (-2.57052e-18 1.1789e-17 9.14779e-18) OK.
    Max cell openness = 1.86801e-15 OK.
    Max aspect ratio = 346.237 OK.
    Minimum face area = 5.12208e-05. Maximum face area = 1.52766.  Face area magnitudes OK.
    Min volume = 3.85096e-06. Max volume = 0.580758.  Total volume = 631.946.  Cell volumes OK.
    Mesh non-orthogonality Max: 87.5162 average: 28.9304
   *Number of severely non-orthogonal (> 70 degrees) faces: 151.
    Non-orthogonality check OK.
  <<Writing 151 non-orthogonal faces to set nonOrthoFaces
    Face pyramids OK.
    Max skewness = 2.47004 OK.
    Coupled point location match (average 0) OK.

Mesh OK.

End

I have tried things like:
- Increasing/decreasing nOrthogonal
- Reducing the flow rate

And none seem to have it converge. I wonder if someone has pointers about what I should be looking into to start learning how to solve this?

Thank you!

Here are my initial conditions

p

Code:

/*--------------------------------*- C++ -*----------------------------------*\
| =========                 |                                                 |
| \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox           |
|  \\    /   O peration     | Version:  v2206                                 |
|   \\  /    A nd           | Website:  www.openfoam.com                      |
|    \\/     M anipulation  |                                                 |
\*---------------------------------------------------------------------------*/
FoamFile
{
    version     2.0;
    format      ascii;
    class       volScalarField;
    object      p;
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

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

internalField   uniform 0;

boundaryField
{
    outlet
    {
        type            zeroGradient;
    }

    inlet
    {
        type            fixedValue;
        value           uniform 0;
    }

    defaultFaces
    {
        type            zeroGradient;
    }

}


// ************************************************************************* //

U

Code:

/*--------------------------------*- C++ -*----------------------------------*\
| =========                 |                                                 |
| \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox           |
|  \\    /   O peration     | Version:  v2206                                 |
|   \\  /    A nd           | Website:  www.openfoam.com                      |
|    \\/     M anipulation  |                                                 |
\*---------------------------------------------------------------------------*/
FoamFile
{
    version     2.0;
    format      ascii;
    class       volVectorField;
    object      U;
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

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

internalField   uniform (0 0 0);

boundaryField
{


    outlet
    {
        type            fixedValue;
        value            uniform (0 0 0.005);
    }

    inlet
    {
        type            zeroGradient;
    }

    defaultFaces
    {
        type            noSlip;
    }

}


// ************************************************************************* //

k

Code:

/*--------------------------------*- C++ -*----------------------------------*\
| =========                 |                                                 |
| \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox           |
|  \\    /   O peration     | Version:  v2206                                 |
|   \\  /    A nd           | Website:  www.openfoam.com                      |
|    \\/     M anipulation  |                                                 |
\*---------------------------------------------------------------------------*/
FoamFile
{
    version     2.0;
    format      ascii;
    class       volScalarField;
    object      k;
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

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

internalField   uniform 0.375;

boundaryField
{
    outlet
    {
        type            fixedValue;
        value           uniform 0.375;
    }

    inlet
    {
        type            zeroGradient;
    }

    defaultFaces
    {
        type            kqRWallFunction;
        value           uniform 0.375;
    }

}


// ************************************************************************* //

nut

Code:

/*--------------------------------*- C++ -*----------------------------------*\
| =========                 |                                                 |
| \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox           |
|  \\    /   O peration     | Version:  v2206                                 |
|   \\  /    A nd           | Website:  www.openfoam.com                      |
|    \\/     M anipulation  |                                                 |
\*---------------------------------------------------------------------------*/
FoamFile
{
    version     2.0;
    format      ascii;
    class       volScalarField;
    object      nut;
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

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

internalField   uniform 0;

boundaryField
{
    inlet
    {
        type            calculated;
        value           uniform 0;
    }

    outlet
    {
        type            calculated;
        value           uniform 0;
    }

    defaultFaces
    {
        type            nutkWallFunction;
        value           uniform 0;
    }

}


// ************************************************************************* //

nuTilda

Code:

/*--------------------------------*- C++ -*----------------------------------*\
| =========                 |                                                 |
| \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox           |
|  \\    /   O peration     | Version:  v2206                                 |
|   \\  /    A nd           | Website:  www.openfoam.com                      |
|    \\/     M anipulation  |                                                 |
\*---------------------------------------------------------------------------*/
FoamFile
{
    version     2.0;
    format      ascii;
    class       volScalarField;
    object      nuTilda;
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

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

internalField   uniform 0;

boundaryField
{
    outlet
    {
        type            fixedValue;
        value           uniform 0;
    }

    inlet
    {
        type            zeroGradient;
    }

    defaultFaces
    {
        type            zeroGradient;
    }

}


// ************************************************************************* //

omega

Code:

/*--------------------------------*- C++ -*----------------------------------*\
| =========                 |                                                 |
| \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox           |
|  \\    /   O peration     | Version:  v2206                                 |
|   \\  /    A nd           | Website:  www.openfoam.com                      |
|    \\/     M anipulation  |                                                 |
\*---------------------------------------------------------------------------*/
FoamFile
{
    version     2.0;
    format      ascii;
    class       volScalarField;
    object      omega;
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

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

internalField   uniform 440.15;

boundaryField
{
    outlet
    {
        type            fixedValue;
        value           $internalField;
    }

    inlet
    {
        type            zeroGradient;
    }

    defaultFaces
    {
        type            omegaWallFunction;
        value           $internalField;
    }

}


// ************************************************************************* //

epsilon

Code:

/*--------------------------------*- C++ -*----------------------------------*\
| =========                 |                                                 |
| \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox           |
|  \\    /   O peration     | Version:  v2206                                 |
|   \\  /    A nd           | Website:  www.openfoam.com                      |
|    \\/     M anipulation  |                                                 |
\*---------------------------------------------------------------------------*/
FoamFile
{
    version     2.0;
    format      ascii;
    class       volScalarField;
    object      epsilon;
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

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

internalField   uniform 14.855;

boundaryField
{
    outlet
    {
        type            fixedValue;
        value           uniform 14.855;
    }

    inlet
    {
        type            zeroGradient;
    }

    defaultFaces
    {
        type            epsilonWallFunction;
        value           uniform 14.855;
    }

}


// ************************************************************************* //

Urchin · July 27, 2022, 23:08

Here are two helpful videos:
https://www.youtube.com/watch?v=UM8I6CsVl4Y
https://www.youtube.com/watch?v=zDFk0GCZG68

After watching those I tried the following:

Removed all except p and U files from 0 (not necessary since turbulence model is laminar)

And commented out the RAS portion of turbulenceProperties

Code:

simulationType      laminar;

/*
RAS
{
    // Tested with kEpsilon, realizableKE, kOmega, kOmegaSST,
    // ShihQuadraticKE, LienCubicKE.
    RASModel        kEpsilon;

    turbulence      off;

    printCoeffs     on;
}
*/

This didn't converge either.

I then changed the viscosity to 1e-3 and it converged in 36 steps. Since the modelled fluid viscosity is closer to 1e-6, I wonder if this means I need to refine the mesh further. I'll try that and see how it goes.

Pappelau · July 28, 2022, 01:16

Hello Urchin your simple Simulation fails because of the very bad mesh : Orthogonality > 70 °

Try to remesh it you will not get a good convergence regardless of the solver settings.
Another point is why you use a theth mesh?

Urchin · July 28, 2022, 06:29

Hi Pappelau,

Thank you for the tip and you were indeed correct, the solution converges much better when I have a more refined mesh.

I used Gmsh to import my CAD file and generate the mesh and it has defaulted to tets I believe.

What would you recommend instead of a tetrahedral mesh? I picked it for no good reason to be honest and am happy to try and switch.

Thank you for your help.

geth03 · July 28, 2022, 16:59

Non-orthog. Quality larger than 85 is absolutely a no go. You should definitely remesh your geometry.

Generally, your mesh related error will be lower when your mesh cells are aligned with the flow. That is hardly the case for terahedral cells, but that should not mean that you can not use them, you need to adjust your fvSchemes entries accordingly.

You can use snappyHexMesh to create a hexahedral dominant mesh, that will lower the cell count overall and speedup the simulation.

Pappelau · July 29, 2022, 07:19

Try this workflow for further simulations :

CAD -> STEP (If u draw the CAD yourself u can directly use Salome)

SALOME -> STL
STL + blockMesh +SnappyHexMesh = Good Hexaedral Mesh

This Workflow can be a rocky road, but there are very good Tutorials how to Mesh things from Wolfdynamics or Tobias Holzmann.

Greetings Pappelau

Urchin · July 29, 2022, 07:34

Thank you geth03 and Pappelau, looks like I have some more learning to do. I appreciate your help pointing me in the right direction.

Similar Threads
Thread	Thread Starter	Forum	Replies	Last Post
SimpleFoam Simulation for windtunnel is not converging	Wanting	OpenFOAM	2	July 29, 2022 11:08
simpleFoam is not converging at high flow rate	minzhang	OpenFOAM Running, Solving & CFD	0	July 10, 2020 14:58
kEspilon / simpleFoam - p non converging	KYPCK444	OpenFOAM Running, Solving & CFD	5	February 27, 2018 09:25
SimpleFoam Case not converging	Doflo	OpenFOAM Running, Solving & CFD	2	December 18, 2014 12:09
simpleFoam parallel solver & Fluent polyhedral mesh	Zlatko	OpenFOAM Running, Solving & CFD	3	September 26, 2014 06:53

July 27, 2022, 23:08		#2
Urchin New Member Join Date: Jul 2022 Posts: 4 Rep Power: 3	Here are two helpful videos: https://www.youtube.com/watch?v=UM8I6CsVl4Y https://www.youtube.com/watch?v=zDFk0GCZG68 After watching those I tried the following: Removed all except p and U files from 0 (not necessary since turbulence model is laminar) And commented out the RAS portion of turbulenceProperties Code: simulationType laminar; /* RAS { // Tested with kEpsilon, realizableKE, kOmega, kOmegaSST, // ShihQuadraticKE, LienCubicKE. RASModel kEpsilon; turbulence off; printCoeffs on; } */ This didn't converge either. I then changed the viscosity to 1e-3 and it converged in 36 steps. Since the modelled fluid viscosity is closer to 1e-6, I wonder if this means I need to refine the mesh further. I'll try that and see how it goes.

July 28, 2022, 01:16		#3
Pappelau Member Lukas Join Date: Sep 2021 Posts: 36 Rep Power: 4	Hello Urchin your simple Simulation fails because of the very bad mesh : Orthogonality > 70 ° Try to remesh it you will not get a good convergence regardless of the solver settings. Another point is why you use a theth mesh?

July 28, 2022, 06:29		#4
Urchin New Member Join Date: Jul 2022 Posts: 4 Rep Power: 3	Hi Pappelau, Thank you for the tip and you were indeed correct, the solution converges much better when I have a more refined mesh. I used Gmsh to import my CAD file and generate the mesh and it has defaulted to tets I believe. What would you recommend instead of a tetrahedral mesh? I picked it for no good reason to be honest and am happy to try and switch. Thank you for your help.

July 28, 2022, 16:59		#5
geth03 Senior Member Join Date: Dec 2019 Location: Cologne, Germany Posts: 355 Rep Power: 8	Non-orthog. Quality larger than 85 is absolutely a no go. You should definitely remesh your geometry. Generally, your mesh related error will be lower when your mesh cells are aligned with the flow. That is hardly the case for terahedral cells, but that should not mean that you can not use them, you need to adjust your fvSchemes entries accordingly. You can use snappyHexMesh to create a hexahedral dominant mesh, that will lower the cell count overall and speedup the simulation.

July 29, 2022, 07:19		#6
Pappelau Member Lukas Join Date: Sep 2021 Posts: 36 Rep Power: 4	Try this workflow for further simulations : CAD -> STEP (If u draw the CAD yourself u can directly use Salome) SALOME -> STL STL + blockMesh +SnappyHexMesh = Good Hexaedral Mesh This Workflow can be a rocky road, but there are very good Tutorials how to Mesh things from Wolfdynamics or Tobias Holzmann. Greetings Pappelau

July 29, 2022, 07:34		#7
Urchin New Member Join Date: Jul 2022 Posts: 4 Rep Power: 3	Thank you geth03 and Pappelau, looks like I have some more learning to do. I appreciate your help pointing me in the right direction.