[JasonG]

I have finally began exploring OpenFOAM, and trying to evaluate if the software will meet my needs in terms of functionality. I plan to attend one of the future training courses, but in the mean time I am trying to see how far I can get on my own.

The current model is a section of internal passages in a structure. I attempted to place a prismatic boundary layer on the surfaces that have flow across them, which admittedly I am still learning and not sure about the appropriate layer dimensions that this model calls for. I was attempting to look at a laminar incompressible case using the icoFoam solver.

The outlet velocity is known, and I wish to determine what the pressure drop is through the passage. My main goal is to view the steady state condition pressure drop, and I can accept some loss of resolution to improve computing time.

I found that I can not maintain a courant number less than 1 unless my timestep is at least 1e-7. I let this iterate over the course of 2 days, and it appeared to be converging. Then, the pressure residual started oscillating around 8e-6 and 1e-5. Eventually, it finally started to diverge and the courant number went very high.

Here are my initial conditions:

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

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

internalField uniform 0;

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

outlet
{

type zeroGradient;
}

boundary
{
type zeroGradient;
}

}

// ************************************************** *********************** //"


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 -.998);
}

inlet
{
type zeroGradient;
}

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


}

// ************************************************** *********************** //"

Attached is a picture of the mesh and pressure/velocity plots prior to divergence.

I am not sure if I have specified a boundary condition incorrectly, an improper mesh, or if I am using the wrong solver. I greatly appreciate any help.

Thanks.

[JasonG]

Also, is there a simple way to specify boundary conditions to a face based on a local coordinate system. The inlet face is not perfectly aligned with the global Cartesian, but I suppose I could transform the velocity vector if I knew the angle it was mis-aligned.

[MartinB]

Hi Jason,
you can try simpleFoam. It's much faster and more stable...
In the following thread are some hints how to speed up simpleFoam:
http://www.cfd-online.com/Forums/ope...nvergence.html

Martin

[JasonG]

Originally Posted by MartinB

Hi Jason,
you can try simpleFoam. It's much faster and more stable...
In the following thread are some hints how to speed up simpleFoam:
http://www.cfd-online.com/Forums/ope...nvergence.html

Martin

Thank you for the suggestion. I was thinking about simpleFoam, but I know my current flow regime is certainly laminar. If I set the turbulence to "off", will I still need to specify the other initial conditions such as: "epsilon, nut, nuTilda, and R" ?

[MartinB]

No, you don't need to specify other values/Parameters/fields than U and p.

Just select

RASModel        laminar;
turbulence      off;

in "constant/RASProperties".

[JasonG]

Thanks again, I am going to give that a try on simplified passage.

[MartinB]

Try this "system/fvSolution":

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

solvers
{
    U                                   // linear equation system solver for U
    {
        solver          smoothSolver;   // solver type
        smoother        GaussSeidel;    // smoother type
        tolerance       1e-06;          // solver finishes if either absolute
        relTol          0.01;           // tolerance is reached or the relative
                                        // tolerance here
        nSweeps         1;              // setting for smoothSolver
        maxIter         100;            // limitation of iterations number
    }
    p                                  // linear equation system solver for p
    {
        solver          GAMG;           // very efficient multigrid solver
        tolerance       1e-07;          // solver finishes if either absolute
        relTol          0.001;          // tolerance is reached or the relative
                                        // tolerance here
        minIter         3;              // a minimum number of iterations
        maxIter         100;            // limitation of iterions number
        smoother        DIC;            // setting for GAMG
        nPreSweeps      1;              // 1 for p, set to 0 for all other!
        nPostSweeps     2;              // 2 is fine
        nFinestSweeps   2;              // 2 is fine
        scaleCorrection true;           // true is fine
        directSolveCoarsestLevel false; // false is fine
        cacheAgglomeration on;          // on is fine; set to off, if dynamic
                                        // mesh refinement is used!
        nCellsInCoarsestLevel 500;      // 500 is fine,
                                        // otherwise sqrt(number of cells)
        agglomerator    faceAreaPair;   // faceAreaPair is fine
        mergeLevels     1;              // 1 is fine
    }
}

SIMPLE
{
    nNonOrthogonalCorrectors 1;         // 0 is fine for perfect mesh; increase
                                        // to 1 or 2 for non orthogonal mesh

    convergence         1e-5;           // convergence criteria steady state
}

relaxationFactors
{
    p              0.3;                // 0.3 is stable, decrease for bad mesh
    U               0.7;                // 0.7 is stable, decrease for bad mesh
}

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

And this for "system/fvSchemes":

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

ddtSchemes
{
    default         steadyState;    // no time dependence here
}

gradSchemes
{
    default         Gauss linear;               // Gauss linear or leastSquares
    grad(p)        faceLimited leastSquares 0 1;   // these entries overwrite
    grad(U)         Gauss linear;               // the "default" setting
}

divSchemes
{
    default         none;
    div(phi,U)      Gauss upwind;               // stable: Gauss upwind
    //div(phi,U)      Gauss linearUpwind Gauss; // faster linearUpwind Gauss
//    div(phi,k)      Gauss upwind;             // for turbulent flow only
//    div(phi,epsilon) Gauss upwind;            // for turbulent flow only
//    div(phi,R)      Gauss upwind;             // for turbulent flow only
//    div(R)          Gauss linear;             // for turbulent flow only
//    div(phi,nuTilda) Gauss upwind;            // for turbulent flow only
    div((nuEff*dev(grad(U).T()))) Gauss linear; // this scheme is not used for
                                                // steady and laminar flow, but
                                                // it is read by the solver, so
                                                // it must be defined.
}

laplacianSchemes
{
    default         none;
    laplacian(nuEff,U) Gauss linear corrected;          // for good meshes
    //laplacian(nuEff,U) Gauss linear limited 0.7;      // for bad meshes

    laplacian((1|A(U)),p) Gauss linear corrected;      // for good meshes
    //laplacian((1|A(U)),p) Gauss linear limited 0.7;  // for bad meshes

//  for turbulent flow only:
//    laplacian(DkEff,k) Gauss linear corrected;
//    laplacian(DepsilonEff,epsilon) Gauss linear corrected;
//    laplacian(DREff,R) Gauss linear corrected;
//    laplacian(DnuTildaEff,nuTilda) Gauss linear corrected;
}

interpolationSchemes
{
    default         linear;     // linear is fine
    U               linear;
}

snGradSchemes
{
    default         corrected;  // limited 0.7; for bad meshes;
                                // must fit to laplacianSchemes
}

fluxRequired
{
    default         no;
    p              ;
    phi             ;
}

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

Here comes the "system/controlDict":

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

application     simpleFoam;   //  which solver (for documentation)

startFrom       startTime;              //  firstTime, startTime, latestTime

startTime       0;                      //  set > 0 to continue

stopAt          endTime;                //  writeNow, endTime, nextWrite, ...

endTime         600;                    //  Latest timestep allowed

deltaT          1;                      //  simple counter for steadyState

writeControl    adjustableRunTime;      //  or uncommon: cpuTime, clockTime

writeInterval   50;                     //  time step write interval

purgeWrite      0;                      //  1 recycles time steps storage
                                        //  0 keeps all time steps on disk

writeFormat     ascii;                  //  ascii: readable, binary: smaller

writePrecision  6;                      //  precision for ascii format

writeCompression uncompressed;          //  compressed for gzipped files

timeFormat      general;                //  fixed, scientific or general

timePrecision   6;                      //  precision for time handling

runTimeModifiable yes;                  //  yes: OF reads dictionaries each
                                        //  time step

graphFormat     raw;                    //  raw, gnuplot, xmgr, jplot

//libs()        for user libraries, p.e. user boundary conditions
//functions()   for special functions

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

If your simulation does not start converging properly after some 75 iterations, change the fvScheme entries according to the comments for lower quality meshes, and increase the nNonOrthogonalCorrectors in fvSolution.

I suppose simpleFoam will solve your simulation in a couple of minutes ;-)

Martin

[JasonG]

Wow, thanks alot! Yes, I know I am only running on a core 2 duo with 5gb of ram, but felt 2 days was a bit much!

[JasonG]

Right now it is at Time= 180, but the residuals are slowly dropping. I think for now I'll let it solve through in the background for a bit and see if it reaches convergence, attached are some initial plots.

[JasonG]

I left the solution to run over night at work, but decided to restart it on my home pc. So far, it is very very slowly dropping the residual; much like the icoFoam did just before it went unstable.

I feel that some of the sharp corners with overly dense mesh concentration may be causing some of the flow instability. I am at least learning, that I have much to learn when switching from a standard structural approach to fluids .



Here is an excerpt of the log file:


Time = 1596

smoothSolver: Solving for Ux, Initial residual = 5.04376e-05, Final residual = 9.87483e-07, No Iterations 6
smoothSolver: Solving for Uy, Initial residual = 5.51604e-05, Final residual = 5.5295e-07, No Iterations 7
smoothSolver: Solving for Uz, Initial residual = 2.065e-05, Final residual = 7.26189e-07, No Iterations 5
GAMG: Solving for p, Initial residual = 0.00138633, Final residual = 8.14647e-07, No Iterations 5
GAMG: Solving for p, Initial residual = 0.0016717, Final residual = 8.46325e-07, No Iterations 5
GAMG: Solving for p, Initial residual = 0.00174926, Final residual = 6.99528e-07, No Iterations 5
time step continuity errors : sum local = 0.000138876, global = -2.97148e-07, cumulative = -6.34467e-05
ExecutionTime = 1416.44 s ClockTime = 1423 s





Time = 1603

smoothSolver: Solving for Ux, Initial residual = 4.989e-05, Final residual = 9.76406e-07, No Iterations 6
smoothSolver: Solving for Uy, Initial residual = 5.40616e-05, Final residual = 5.39811e-07, No Iterations 7
smoothSolver: Solving for Uz, Initial residual = 2.02852e-05, Final residual = 7.14789e-07, No Iterations 5
GAMG: Solving for p, Initial residual = 0.00136729, Final residual = 5.67102e-07, No Iterations 5
GAMG: Solving for p, Initial residual = 0.00167996, Final residual = 1.14624e-06, No Iterations 4
GAMG: Solving for p, Initial residual = 0.00162179, Final residual = 7.83107e-07, No Iterations 5
time step continuity errors : sum local = 0.000155448, global = -1.93752e-07, cumulative = -6.60641e-05
ExecutionTime = 1455.49 s ClockTime = 1463 s

[alberto]

Hi,

may I ask you to run a checkMesh on your case? What is the number of cells you are using and what is the size of the system?

Best,

[MartinB]

I think, Alberto is right:
run a "checkMesh -allTopology -allGeometry"
I suppose that there are cell determinant problems...
Martin

[JasonG]

Here is the output generated from the above command:

Build  : 1.7.0-279cc8e8233b
Exec   : checkMesh -allTopology -allGeometry
Date   : Jul 23 2010
Time   : 08:21:33
Host   : jason-desktop
PID    : 7142
Case   : /home/jason/OpenFOAM/jason-1.7.0/run/testing/pipe_test_2
nProcs : 1
SigFpe : Enabling floating point exception trapping (FOAM_SIGFPE).

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

Create polyMesh for time = 0

Time = 0

Mesh stats
    points:           60110
    internal points:  53412
    edges:            294764
    internal edges:   275060
    internal edges using one boundary point:   7452
    internal edges using two boundary points:  0
    faces:            425684
    internal faces:   412676
    cells:            191029
    boundary patches: 3
    point zones:      0
    face zones:       0
    cell zones:       0

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

Checking topology...
    Boundary definition OK.
    Point usage OK.
    Upper triangular ordering OK.
    Topological cell zip-up check OK.
    Face vertices OK.
    Face-face connectivity OK.
    Number of regions: 1 (OK).

Checking patch topology for multiply connected surfaces ...
    Patch               Faces    Points   Surface topology                   Bounding box
    boundary            12374    6219     ok (non-closed singly connected)   (-0.03277 0.00680315 -0.0723646) (-0.0184584 0.0211148 -0.0500063)
    inlet               284      248      ok (non-closed singly connected)   (-0.03277 0.00680315 -0.0563375) (-0.0294114 0.0101618 -0.0516125)
    outlet              350      295      ok (non-closed singly connected)   (-0.0267988 0.0144653 -0.0723646) (-0.0214612 0.0198247 -0.0723646)

Checking geometry...
    Overall domain bounding box (-0.03277 0.00680315 -0.0723646) (-0.0184584 0.0211148 -0.0500063)
    Mesh (non-empty, non-wedge) directions (1 1 1)
    Mesh (non-empty) directions (1 1 1)
    Boundary openness (-3.64353e-19 -2.55572e-19 -5.58326e-19) OK.
    Max cell openness = 2.04311e-16 OK.
    Max aspect ratio = 14.7891 OK.
    Minumum face area = 5.09213e-12. Maximum face area = 3.86097e-07.  Face area magnitudes OK.
    Min volume = 1.08791e-17. Max volume = 7.96135e-11.  Total volume = 8.08761e-07.  Cell volumes OK.
    Mesh non-orthogonality Max: 60.2802 average: 19.3085
    Non-orthogonality check OK.
    Face pyramids OK.
    Max skewness = 3.16024 OK.
    Min/max edge length = 1.10157e-06 0.001065 OK.
    All angles in faces OK.
    Face flatness (1 = flat, 0 = butterfly) : average = 0.999739  min = 0.895423
    All face flatness OK.
    Cell determinant (wellposedness) : minimum: 0.0383356 average: 1.54342
    Cell determinant check OK.

Mesh OK.

End

Is the non-orthogonality rather high?

[JasonG]

I was also wondering if it might be appropriate to loosen up on my pressure and maybe velocity tolerance. For these initial stages I am attempting to show comparative pressure drops between designs, and 5-10% within the converged result may be adequate.

[alberto]

Hi,

the mesh passes all the tests without errors or warnings. About tolerances, try to plot the residuals, and check if they became flat or the keep dropping, and stop the simulation after they stayed flat for a while. That should be the converged solution.

Best,

[JasonG]

I thought the solution couldn't be assumed converged until the tolerance set is met, is a steady mean residual over time more important than hitting the tolerance?

[alberto]

The tolerance is an arbitrary value, and it is not necessarily something you can reach, especially if very low.

Usually residuals are an indication of convergence that has to me considered together with other elements, both numerical and physical.
For example, your residuals can be high, but never change for a high number of iterations, and that means your solution is not changing at all.

[JasonG]

Thanks again, this has been very helpful. I am relieved to know this case will not require 3 days to run .

[alberto]

It should really take less thank 1 hour! :-)

[JasonG]

Yes, looking at the plot I think I could have killed it after 200-500 iterations which would be 20-40mins .