[1988]

hello everyone
I am somehow new in openfoam please help me to solve the problem.
I have a simulated a pipe which water is flowing in and I am using icoFoam solver the check mesh is here

Create time

Create polyMesh for time = 0

Time = 0

Mesh stats
    points:           96136
    faces:            267065
    internal faces:   246979
    cells:            85674
    faces per cell:   6
    boundary patches: 3
    point zones:      0
    face zones:       1
    cell zones:       1

Overall number of cells of each type:
    hexahedra:     85674
    prisms:        0
    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                  
    body                19650    19700    ok (non-closed singly connected)  
    inlet               218      244      ok (non-closed singly connected)  
    outlet              218      244      ok (non-closed singly connected)  

Checking geometry...
    Overall domain bounding box (-0.0142294 0 -0.0151572) (0.0129588 0.0472016 0.0022)
    Mesh (non-empty, non-wedge) directions (1 1 1)
    Mesh (non-empty) directions (1 1 1)
    Boundary openness (8.88725e-17 9.11711e-17 1.4754e-17) OK.
    Max cell openness = 2.45267e-16 OK.
    Max aspect ratio = 5.50297 OK.
    Minimum face area = 5.30557e-09. Maximum face area = 1.31635e-07.  Face area magnitudes OK.
    Min volume = 1.66462e-12. Max volume = 3.77399e-11.  Total volume = 1.09182e-06.  Cell volumes OK.
    Mesh non-orthogonality Max: 20.6887 average: 5.3841
    Non-orthogonality check OK.
    Face pyramids OK.
    Max skewness = 0.435615 OK.
    Coupled point location match (average 0) OK.

Mesh OK.

End

I can not understand why the solution is diverged and why the courant number suddenly explodes to unbelievable amounts.
here I put fvsolution and fv schemes

ddtSchemes
{
    default         Euler;
}

gradSchemes
{
    default         Gauss linear;
    grad(p)         Gauss linear;
}

divSchemes
{
    default         none;
    div(phi,U)      Gauss linear;
}

laplacianSchemes
{
    default         none;
    laplacian(nu,U) Gauss linear orthogonal;
    laplacian((1|A(U)),p) Gauss linear orthogonal;
}

interpolationSchemes
{
    default         linear;
    interpolate(HbyA) linear;
}

snGradSchemes
{
    default         orthogonal;
}

fluxRequired
{
    default         no;
    p               ;
}

solvers
{
    p
    {
        solver          PCG;
        preconditioner  DIC;
        tolerance       1e-06;
        relTol          0;
    }

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

PISO
{
    nCorrectors     2;
    nNonOrthogonalCorrectors 2;
    pRefCell        0;
    pRefValue       0;
}

but this is courant numbers after just one time step(time steps are 0.001)

Time = 0.004

Courant Number mean: 2.82932e+17 max: 1.79062e+21

Is there any idea?

[Pj.]

I'm not an expert either with icoFoam.

You could look at the number of iterations for each field (u,p) to see which one diverges.
Sometimes the problem is also caused by bad cells (even few can create big problems). Try run a checkMesh to see if the cells are all ok.
If the mesh is ok you can try lowering the time-step to have a more stable problem.

The question "why my simulation diverges?" is not an easy one to answer to, there could be tons of reasons...

[ArathoN]

How did you compute the Courant number? remember to choose the smallest grid space when you'll calculate it. And try to change the solver in fvschemes or fvsolutions.

[Pj.]

Originally Posted by ArathoN

How did you compute the Courant number? remember to choose the smallest grid space when you'll calculate it. And try to change the solver in fvschemes or fvsolutions.

You does not compute it: OpenFOAM automatically return the mean and the max value at each iteration.

I think that such a rapid "explosion" is a mesh related problem.

However you can also try to plot the value of Courant using the "Co" utility (just tipe "Co" in the prompt) and visualise it with paraView. This can help to understand where the Courant is exploded and to find the corresponding cells. This however works when you have "fisical" Co. Such an huge Co suggest the that U field is already (excuse me the technical term ) fucked up.

[ArathoN]

I said computed because you need to set the time step in controldict, and to do so you need the grid spacing (better the smallest one) and a desired Co (i usually set it to 0.5 for stability). You can also try to limit the Max Co in controldict and see if it works.

I had a similar problem and i solved it by lowering the time step and changing the schemes.

[olivierG]

hello,
You should add the 0/U and 0/p file if you want some help here.
However, common mistake are:
- bad scalling for your mesh see bounding box, here ~4 cm box is ok ?
- bad BC and initial value.

regards,
Olivier

[1988]

thanks all of you for your attention
the mesh was made in gambit and scaled to mm in openfoam.I have put check mesh result and my time step as you can see is some thing about 0.001 and I have to decrease it to 0.00001 but it is really time consuming and it is not logically true in this case.
these are my p and u boundary conditions.

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

internalField   uniform (0 0 0);

boundaryField
{
    inlet
    {
      type            groovyBC;
      valueExpression "vector(0,0.15-0.15*cos(2.5*pi*time()),0)";
      gradientExpression "vector(0,0,0)";
      fractionExpression "1";
      value           uniform (0 0 0);
    }

    outlet
    {
        type            zeroGradient;
    }

    body
    {
        type            fixedValue;
        value           uniform (0 0 0);
    }
}

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

internalField   uniform 0;

boundaryField
{
    inlet
    {
        type            zeroGradient;
    }

    outlet
    {
       type            groovyBC;
      valueExpression "12000-2667*cos(2.5*pi*time())";
      gradientExpression "0";
      fractionExpression "1";
      value           uniform 0;
    }

Is there any idea?

[Pj.]

as Olivier pointed out: is it correct for the mesh to bee THAT SMALL?

Your overall domain is about 4cm large! The cells side are less than 1mm, is that correct?

Your smallest cell is 1.66462e-12 m^3. With a 0.3 maximum inlet velocity 0.001s time step is not very small.

[1988]

yes that is true but let me explain it maybe you have better idea.
I am working on a vein which is in mm scale and it needs 100000 mesh but I ran in in meter scale with no problem and the courant number even with 1 second time steps(not 0.001 or small ones) was less than 1 by IcoFoam solver however I have got the wrong answer for wall shear stress .By scaling the geometry to mm my wall shear stress come closer to fact nevertheless the courant number explode and I prefer to change my solver to PimpleFoam because I think that it lets me to define a courant number and it adjusts the time steps but I am not sure this is the right decision and I need to help me to solve the problem.
thank you

[ArathoN]

Originally Posted by 1988

yes that is true but let me explain it maybe you have better idea.
I am working on a vein which is in mm scale and it needs 100000 mesh but I ran in in meter scale with no problem and the courant number even with 1 second time steps(not 0.001 or small ones) was less than 1 by IcoFoam solver however I have got the wrong answer for wall shear stress .By scaling the geometry to mm my wall shear stress come closer to fact nevertheless the courant number explode and I prefer to change my solver to PimpleFoam because I think that it lets me to define a courant number and it adjusts the time steps but I am not sure this is the right decision and I need to help me to solve the problem.
thank you

Are you working with blood flow? because you need to consider a non-newtonian flow with the appropriate coefficients. Unfortunately I have no knowledge on such flows but in transport properties you can set the constants.

To calculate the time step you need to choose a Co number (less than 0.5 for LES or less than 1 generally, I usually tend to choose a given Co less than 0.5 for stability even if sometimes it isn't needed), then you need to see the smallest grid spacing possibly in all the 3 direction.

Now you can calculate the time step needed and choose the smallest, obviously if the mesh is really small the time step will be really small.

As a partially OT matter what is the smallest mesh/geometry than openfoam can study because i'm studying a channel with a height of 0.025m and a relative timestep of 2e-6 (Max Co 0.5). The solver has some problems especially if i use cyclic BCs with LES+pimplefoam (or even RANS+pimpleFoam).

[1988]

I don't use now newtonian methods because of some special biological things and you can consider it as laminar, unsteady,newtonian flow but the dimensions are the same as things which you saw in check mesh part.