CFD Online Discussion Forums

CFD Online Discussion Forums (https://www.cfd-online.com/Forums/)
-   OpenFOAM Running, Solving & CFD (https://www.cfd-online.com/Forums/openfoam-solving/)
-   -   simpleFoam no convergence (https://www.cfd-online.com/Forums/openfoam-solving/132404-simplefoam-no-convergence.html)

ilpaso March 31, 2014 07:48
simpleFoam no convergence
 
Hi,
I'm trying to set up a case with simpleFoam solver without turbolence.

The geometry is a pipe with a diameter of 3mm and a length of 117mm.


U:
The volume flow rate at inflow is 5 ml/s, so the velocity field is
Code:
inflow
{
  type fixedValue;
  value uniform (0.707355 0 0)
}
outflow
{
  type zeroGradient;
}
wall
{
  type fixedValue;
  value uniform (0 0 0);
}
nu:
The density rho is 1065kg/m^3 (blood)
The dynamic viscosity is 3.4e-3 Pa*s (blood)
So my cinematic viscosity is:
Code:
nu    nu [0  2 -1 0 0 0 0 ] 3.2e-6
p:
pressure boundary condictions are
Code:
inflow
{
  type zeroGradient;
}
outflow
{
  type fixedValue;
  value uniform 0;
}
wall
{
  type zeroGradient;
}
turbolenceProperties:
I'd like to use a laminar model, so I set:
Code:
simulationType laminar;
turbolenceModel laminar;
tubolence off;

The solution doesn't converge but if I change the inflow velocity to
Code:
inflow
{
  type fixedValue;
  value uniform (0.707355e-3 0 0)
}
and the nu to:
Code:
nu    nu [0  2 -1 0 0 0 0 ] 3.2e-3
the solution converges after 43 iterations


What is wrong in my case?

alexeym March 31, 2014 07:57
Hi,

you've reduced inlet velocity and increased viscosity, so Re went from 663.15 to 0.00066315, that's why it converged quickly.

What's the number of time steps in your original simulation (with higher Re), what's in your SIMPLE dictionary?

ilpaso March 31, 2014 08:10
thank you for the reply

after 1000 timesteps the residuals are horizontal but there is no convergence

this is the SIMPLE dictionary
Code:
SIMPLE
{
    nNonOrthogonalCorrectors 0;
    residualControl
    {
        p              1e-2;
        U              1e-3;
        "(k|epsilon)"  1e-3;
    }

alexeym March 31, 2014 08:19
Well...

1. what are the values of residuals for pressure and velocity (flat ones)?
2. I suppose there's no error in the direction of the velocity (i.e. it really goes along the pipe not in normal direction)
3. is it cylindrical pipe? can you show checkMesh output?

ilpaso March 31, 2014 08:27
1) residuals:
Code:
smoothSolver:  Solving for Ux, Initial residual = 0.000598572, Final residual = 4.88836e-05, No Iterations 3
smoothSolver:  Solving for Uy, Initial residual = 0.00176527, Final residual = 0.000138654, No Iterations 3
smoothSolver:  Solving for Uz, Initial residual = 0.00174763, Final residual = 0.000137999, No Iterations 3
GAMG:  Solving for p, Initial residual = 0.00861638, Final residual = 0.000187439, No Iterations 2
time step continuity errors : sum local = 2.47772e-08, global = 1.5928e-10, cumulative = -8.88454e-05
2) yes, the direction is ok

3) checkMesh output:
Code:
Mesh stats
    points:          2387
    faces:            9722
    internal faces:  8578
    cells:            3660
    faces per cell:  5
    boundary patches: 3
    point zones:      0
    face zones:      0
    cell zones:      0

Overall number of cells of each type:
    hexahedra:    0
    prisms:        3660
    wedges:        0
    pyramids:      0
    tet wedges:    0
    tetrahedra:    0
    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                 
    wall                900      930      ok (non-closed singly connected) 
    inflow              122      77      ok (non-closed singly connected) 
    outflow            122      77      ok (non-closed singly connected) 

Checking geometry...
    Overall domain bounding box (0 -1.4918 -1.5) (50 1.4918 1.5)
    Mesh (non-empty, non-wedge) directions (1 1 1)
    Mesh (non-empty) directions (1 1 1)
    Boundary openness (6.36006e-18 3.39489e-18 -1.83353e-18) OK.
    Max cell openness = 1.57376e-16 OK.
    Max aspect ratio = 16.779 OK.
    Minimum face area = 0.0297778. Maximum face area = 1.11068.  Face area magnitudes OK.
    Min volume = 0.0496099. Max volume = 0.217021.  Total volume = 350.847.  Cell volumes OK.
    Mesh non-orthogonality Max: 24.1209 average: 6.15059
    Non-orthogonality check OK.
    Face pyramids OK.
    Max skewness = 0.333255 OK.
    Coupled point location match (average 0) OK.

Mesh OK.

alexeym March 31, 2014 08:34
What's in your fvSolution? What if you switch GAMG/smoothSolver to PCG/PBiCG?

ilpaso March 31, 2014 08:48
What means change GAMG/smoothSolver to PCG/PBiCG is obscure for me (I've to study the solver before run the simulation!! I know).
Please can you help me to set up the dictionary?

The second question is: I've another geometry with a variable diameter. The lower value is 1.2mm and the reynolds at that section is about 1600. I think I've to use a turbolence model. How change the dictionaries in order to activate this feature?

Thank you very much

this is my fvSolution.
Code:
solvers
{
    p
    {
        solver          GAMG;
        smoother        GaussSeidel;
        cacheAgglomeration true;
        nCellsInCoarsestLevel 10;
        agglomerator    faceAreaPair;
        mergeLevels    1;

        tolerance      1e-06;
        relTol          0.05;
    }

    pFinal
    {
        solver          GAMG;
        smoother        GaussSeidel;
        cacheAgglomeration true;
        nCellsInCoarsestLevel 10;
        agglomerator    faceAreaPair;
        mergeLevels    1;

        tolerance      1e-06;
        relTol          0;
    }

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

    "(U|k|epsilon)Final"
    {
        solver          PBiCG;
        preconditioner  DILU;

        tolerance      1e-05;
        relTol          0;
    }
}

PIMPLE
{
    nOuterCorrectors 4;
    nCorrectors    1;
    nNonOrthogonalCorrectors 0;
    pRefCell        0;
    pRefValue      0;
}

SIMPLE
{
    nNonOrthogonalCorrectors 0;
    residualControl
    {
        p              1e-2;
        U              1e-3;
        "(k|epsilon)"  1e-3;
    }
}

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

cache
{
    grad(U);
}

alexeym March 31, 2014 09:00
Well,

you can just take fvSchemes/fvSolution dictionaries from $FOAM_TUTORIALS/incompressible/simpleFoam/pitzDaily.

In the same tutorial you'll find all necessary modifications to run case with turbulence.

by switching from GAMG/smoothSolver to PCG/PBiCGI meant changing:

Code:
solvers
{
    p
    {
        solver          GAMG;
        smoother        GaussSeidel;
        cacheAgglomeration true;
        nCellsInCoarsestLevel 10;
        agglomerator    faceAreaPair;
        mergeLevels    1;

        tolerance      1e-06;
        relTol          0.05;
    }

    "(U|k|epsilon)"
    {
        solver          smoothSolver;
        smoother        GaussSeidel;
        tolerance      1e-05;
        relTol          0.1;
    }
...
}
to

Code:
solvers
{
    p
    {
        solver          PCG;
        preconditioner  DIC;
        tolerance      1e-06;
        relTol          0.01;
    }

    U
    {
        solver          PBiCG;
        preconditioner  DILU;
        tolerance      1e-05;
        relTol          0.1;
    }
...
}

ilpaso March 31, 2014 09:19
It doesn't converge.
But the residuals are low! or not? :(
Do you think the problem is in the mesh?

these are the residuals:

Code:
DILUPBiCG:  Solving for Ux, Initial residual = 0.000595124, Final residual = 1.49812e-07, No Iterations 1
DILUPBiCG:  Solving for Uy, Initial residual = 0.0018561, Final residual = 3.85916e-07, No Iterations 1
DILUPBiCG:  Solving for Uz, Initial residual = 0.00181595, Final residual = 4.51476e-07, No Iterations 1
DICPCG:  Solving for p, Initial residual = 0.00817729, Final residual = 7.80714e-05, No Iterations 27
time step continuity errors : sum local = 1.06637e-08, global = -2.70428e-10, cumulative = 9.34916e-06

alexeym March 31, 2014 09:33
Are you sure your mesh is 3 mm in diameter? checkMesh thinks it is 3 m in diameter

Code:
Overall domain bounding box (0 -1.4918 -1.5) (50 1.4918 1.5)
(AFAIK these number are in meters)

in this case, you've got Re = UD/Nu = 663145.3125 and it's obviously non-laminar case and therefore you've got no convergence as you're running without turbulence.

ilpaso March 31, 2014 09:47
now the mesh is scaled
Code:
Overall domain bounding box (0 -0.0014918 -0.0015) (0.05 0.0014918 0.0015)
but there is no convergence. I'm desperate!

Code:
DILUPBiCG:  Solving for Ux, Initial residual = 0.000565348, Final residual = 2.30843e-05, No Iterations 1
DILUPBiCG:  Solving for Uy, Initial residual = 0.00543364, Final residual = 0.000228817, No Iterations 1
DILUPBiCG:  Solving for Uz, Initial residual = 0.0052694, Final residual = 0.000218604, No Iterations 1
DICPCG:  Solving for p, Initial residual = 0.00364692, Final residual = 3.63526e-05, No Iterations 55
time step continuity errors : sum local = 0.000746084, global = 5.79239e-06, cumulative = -0.000574256

ilpaso March 31, 2014 09:58
1 Attachment(s)
attached here there is the plot of the residuals

alexeym March 31, 2014 10:08
1 Attachment(s)
I don't know what I'm doing wrong but attached case converged in 143 iterations. I've taken velocity and viscosity from your post. Two points that are different from your case:

1. Fully hexagonal mesh.
2. van Leer scheme for velocity discretisation.

(if you'd like to run attached case, you'll need Gmsh to create mesh)

ilpaso March 31, 2014 10:16
1 Attachment(s)
thank you Alexey!
I'll try with gmsh.

Attached here you can find my polymesh directory.

alexeym March 31, 2014 10:22
With attached polyMesh directory case converged in 42 iterations.

ilpaso March 31, 2014 10:22
1 Attachment(s)
and here the case

alexeym March 31, 2014 10:30
There's not convergence with "limitedLinear 1.0" scheme but if I change it to upwind it converges in 42 iteration (64 iterations with GammaV, 68 iteration with vanLeerV, 65 with linear).

ilpaso March 31, 2014 10:37
with your fvSolution and fvScheme it converges in 67 iterations.
Where can I change the scheme from "limitedLinear 1.0" scheme to upwind in order to reduce the number o iterations?
A lower number of iterations is very important because I need to run the case recursively because I'll couple it with a 0D model.

alexeym March 31, 2014 10:41
Your current fvSchemes:

Code:
divSchemes
{
    default        none;
    div(phi,U)      bounded Gauss limitedLinearV 1;
    div(phi,k)      bounded Gauss limitedLinear 1;
    div(phi,epsilon) bounded Gauss limitedLinear 1;
    div(phi,R)      bounded Gauss limitedLinear 1;
    div(R)          Gauss linear;
    div(phi,nuTilda) bounded Gauss limitedLinear 1;
    div((nuEff*dev(T(grad(U))))) Gauss linear;
}
change it to

Code:
divSchemes
{
    default        none;
    div(phi,U)      bounded Gauss upwind;
    div(phi,k)      bounded Gauss upwind;
    div(phi,epsilon) bounded Gauss upwind;
    div(phi,R)      bounded Gauss upwind;
    div(R)          Gauss linear;
    div(phi,nuTilda) bounded Gauss upwind;
    div((nuEff*dev(T(grad(U))))) Gauss linear;
}

ilpaso March 31, 2014 10:52
Ok. I'll try this evening with the different divSchemes.
Now It converges.

Thank you very much Alexey for your help!! I've lost a lot of hours with this problem!


All times are GMT -4. The time now is 22:36.