[mayank.dce2k7]

Hi Foamers,

I am facing SimpleFoam convergence problem with a really basic simulation i.e steady state flow simulation in a single bifurcation tube by specifying fully developed velocity profile at inlet(using swak4Foam toolbox).

After digging through various threads related to convergence problem, I tried changing divergence schemes, played with relaxation factors and tried grid refinement but nothing helped. Situation became so bad the case is not converging using first order upwind scheme let one second order schemes. Below are some details of my case set up:

1.) Mesh is generated on ANSYS ICEM using tetrahedral and prism elements. My model contains about 0.6 million elements.

2.) checkMesh result:

Create time

Create polyMesh for time = 0

Time = 0

Mesh stats
points: 195781
faces: 1354300
internal faces: 1310124
cells: 603415
boundary patches: 4
point zones: 0
face zones: 1
cell zones: 1

Overall number of cells of each type:
hexahedra: 0
prisms: 250764
wedges: 0
pyramids: 0
tet wedges: 0
tetrahedra: 352651
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 41794 20988 ok (non-closed singly connected)
INLET 1047 787 ok (non-closed singly connected)
OUTLET1 666 523 ok (non-closed singly connected)
OUTLET2 669 528 ok (non-closed singly connected)

Checking geometry...
Overall domain bounding box (-0.108855 -0.161689 -0.01905) (0.108855 0.127 0.01905)
Mesh (non-empty, non-wedge) directions (1 1 1)
Mesh (non-empty) directions (1 1 1)
Boundary openness (6.33376e-17 2.25628e-17 -1.46543e-16) OK.
Max cell openness = 3.10378e-16 OK.
Max aspect ratio = 18.7479 OK.
Minumum face area = 1.85691e-07. Maximum face area = 4.16397e-06. Face area magnitudes OK.
Min volume = 1.38121e-11. Max volume = 2.46223e-09. Total volume = 0.000347757. Cell volumes OK.
Mesh non-orthogonality Max: 76.7538 average: 13.2348
*Number of severely non-orthogonal faces: 1.
Non-orthogonality check OK.
<<Writing 1 non-orthogonal faces to set nonOrthoFaces
Face pyramids OK.
Max skewness = 2.96438 OK.
Coupled point location match (average 0) OK.

Mesh OK.

End

3.) 0/U file

FoamFile
{
    version     2.0;
    format      ascii;
    class       volVectorField;
    object      U;
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

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

internalField  uniform (0 -0.172 0);

boundaryField
{
    INLET
    {

       type            groovyBC;
       variables "r2=(pow(pos().x,2)+pow(pos().z,2));R2=pow(0.019,2);para=-((R2-r2)/R2)*normal();";
       valueExpression "2*0.172*para"; //where 0.351 is the average velocity
       value        uniform (0 0 0);      
     
    }

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

    OUTLET1
    {

        type            flowRateInletVelocity;
        flowRateExpression "0";
        flowRate swak { 
            variables "phi1{INLET}=-0.52*0.172*sum(area());"; 
            expression "phi1"; // volume flow going out of domain from outlet 1 which 0.52 times that at inlet
            valueType patch;
            patchName OUTLET1;};

     value           uniform ( 0 0 0 );
    }

    OUTLET2
    {

      type  zeroGradient;
       
    }
}

2.) 0/p

FoamFile
{
    version     2.0;
    format      ascii;
    class       volScalarField;
    object      p;
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

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

internalField   uniform 0;

boundaryField
{
    INLET
    {
        type            zeroGradient;
    }

    OUTLET1
    {
      type            zeroGradient;
    }
   
    OUTLET2
    {
      type            fixedValue;
      value           uniform 0;
    }  

    WALL
    {
        type            zeroGradient;
    }

    
}

3.)

fvSchemes

FoamFile
{
    version     2.0;
    format      ascii;
    class       dictionary;
    location    "system";
    object      fvSchemes;
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

ddtSchemes
{
    default         steadyState;
}

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

divSchemes
{
    default         none;
    div(phi,U)      Gauss upwind;
    div(phi,nuTilda) Gauss upwind grad(nuTilda);
    div((nuEff*dev(T(grad(U))))) Gauss linear;
}

laplacianSchemes
{
    default         none;
    laplacian(nuEff,U) Gauss linear limited 0.333;
    laplacian((1|A(U)),p) Gauss linear limited 0.333;
    laplacian(DnuTildaEff,nuTilda) Gauss linear limited 0.333;
    laplacian(1,p)  Gauss linear limited 0.333;
}

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

snGradSchemes
{
    default         limited 0.333;
}

fluxRequired
{
    default         no;
    p               ;
}


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

4.)

/*--------------------------------*- C++ -*----------------------------------*\
| =========                 |                                                 |
| \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox           |
|  \\    /   O peration     | Version:  2.1.1                                 |
|   \\  /    A nd           | Web:      www.OpenFOAM.org                      |
|    \\/     M anipulation  |                                                 |
\*---------------------------------------------------------------------------*/
FoamFile
{
    version     2.0;
    format      ascii;
    class       dictionary;
    location    "system";
    object      fvSolution;
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

solvers
{
    p
    {
        solver          GAMG;
        tolerance       1e-02;
        relTol          0.1;
        smoother        GaussSeidel;
        nPreSweeps      0;
        nPostSweeps     2;
        cacheAgglomeration true;
        nCellsInCoarsestLevel 10;
        agglomerator    faceAreaPair;
        mergeLevels     1;
    }

    U
    {
        solver          smoothSolver;
        smoother        GaussSeidel;
        nSweeps         2;
        tolerance       1e-02;
        relTol          0.1;
    }

    nuTilda
    {
        solver          smoothSolver;
        smoother        GaussSeidel;
        nSweeps         2;
        tolerance       1e-02;
        relTol          0.1;
    }
}

SIMPLE
{
    nNonOrthogonalCorrectors 2;
    pRefCell        0;
    pRefValue       0;

    residualControl
    {
        p               1e-4;
        U               1e-4;
        nuTilda         1e-4;
    }
}

relaxationFactors
{
    fields
    {
        p               0.9;
    }
    equations
    {
        U               0.1;
        nuTilda         0.2;
    }
}


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

Any help will be appreciated!

[mayank.dce2k7]

ok, one improvement on my part. I got it converged for first order upwind scheme but still no luck with second schemes(SFDCV, linearUpwind etc). I took first order result and used it for my simulation with second order divergence scheme but with no luck.

[RodriguezFatz]

Hi, did you try to use some limiters on the gradients? Sounds like this could be a problem here, since first order works and second order doesn't.