CFD Online Logo CFD Online URL
www.cfd-online.com
[Sponsors]
Home > Forums > OpenFOAM Running, Solving & CFD

PimpleFoam: Solver Crashes for simple laminar flow

Register Blogs Members List Search Today's Posts Mark Forums Read

Reply
 
LinkBack Thread Tools Display Modes
Old   May 1, 2014, 20:53
Default PimpleFoam: Solver Crashes for simple laminar flow
  #1
Member
 
India
Join Date: Oct 2012
Posts: 84
Rep Power: 4
mayank.dce2k7 is on a distinguished road
Hi Foamers,

I am trying to calculate transient flow field in a single bifurcating tube(using OF-2.1.1 and swak4Foam) with inlet boundary condition varying with space and time:

I will start illustrating my case by giving 'checkMesh' report:

Create time

Create polyMesh for time = 0

Time = 0

Mesh stats
points: 140311
faces: 1127971
internal faces: 1083988
cells: 521531
boundary patches: 4
point zones: 0
face zones: 1
cell zones: 1

Overall number of cells of each type:
hexahedra: 0
prisms: 125835
wedges: 0
pyramids: 0
tet wedges: 0
tetrahedra: 395696
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 910 612 ok (non-closed singly connected)
OUTLET1 562 390 ok (non-closed singly connected)
OUTLET2 566 394 ok (non-closed singly connected)
WALL 41945 21065 ok (non-closed singly connected)

Checking geometry...
Overall domain bounding box (-0.108847 -0.161672 -0.0190498) (0.108851 0.127 0.0190498)
Mesh (non-empty, non-wedge) directions (1 1 1)
Mesh (non-empty) directions (1 1 1)
Boundary openness (2.61084e-17 -1.48172e-16 3.15641e-17) OK.
Max cell openness = 2.91142e-16 OK.
Max aspect ratio = 7.00156 OK.
Minumum face area = 1.85341e-07. Maximum face area = 3.64501e-06. Face area magnitudes OK.
Min volume = 5.00262e-11. Max volume = 2.11586e-09. Total volume = 0.000347759. Cell volumes OK.
Mesh non-orthogonality Max: 63.6578 average: 14.781
Non-orthogonality check OK.
Face pyramids OK.
Max skewness = 2.94579 OK.
Coupled point location match (average 0) OK.

Mesh OK.

End

Boundary Conditions:

0/U

Code:

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

internalField  uniform (0 0 0);

boundaryField
{
    INLET
    {
                  
  
     type            groovyBC;

     variables "vavg=(-0.000163-0.00336*cos(time()*0.3254)-0.3416*sin(time()*0.3254));r2=(pow(pos().x,2)+pow(pos().z,2));R2=sum(area())/pi;para=-((R2-r2)/R2)*normal();";
     valueExpression "2*vavg*para";
     value        uniform (0 0 0);
 
    }

    OUTLET1
    {
       type  zeroGradient;
    }

    OUTLET2
    {
        type            zeroGradient;
     
}

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

   
}
0/p

Code:
dimensions      [0 2 -2 0 0 0 0];

internalField   uniform 0;

boundaryField
{
    INLET
    {
        type            zeroGradient;
    }

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

    WALL
    {
        type            zeroGradient;
    }

    
}
controlDict

Code:
libs ( "libOpenFOAM.so"
      "libsimpleSwakFunctionObjects.so"
      "libswakFunctionObjects.so"
      "libgroovyBC.so" ) ;
application     pimpleFoam;

startFrom       startTime;

startTime       0;

stopAt          endTime;

endTime         20;

deltaT          0.02;

writeControl    adjustableRunTime;

writeInterval   0.4;

purgeWrite      0;

writeFormat     ascii;

writePrecision  6;

writeCompression uncompressed;

timeFormat      general;

timePrecision   6;

runTimeModifiable true;

adjustTimeStep  no;

maxCo           100;

functions {
    initSwak {
        // needed to allow DataEntry to work
        type initSwakFunctionObject;
        region region0;
    }
};

functions
{
    volFlowInlet
    {
        type                swakExpression;
        valueType           patch;
        patchName           INLET;
        
        verbose             true;
        expression          "phi";
        accumulations       ( sum );
    }
    
    volFlowOutlet1
    {
        type                swakExpression;
        valueType           patch;
        patchName           OUTLET1;
        
        verbose             true;
        expression          "phi";
        accumulations       ( sum );
    }

    volFlowOutlet2
    {
        type                swakExpression;
        valueType           patch;
        patchName           OUTLET2;
        
        verbose             true;
        expression          "phi";
        accumulations       ( sum );
    }
}
fvSchemes

Code:
ddtSchemes
{
    default         CrankNicholson 1;
}

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

divSchemes
{
    default         none;
    div(phi,U)      Gauss linearUpwind grad(U);
    div(phi,k)      Gauss limitedLinear 1;
    div(phi,epsilon) Gauss limitedLinear 1;
    div(phi,R)      Gauss limitedLinear 1;
    div(R)          Gauss linear;
    div(phi,nuTilda) Gauss limitedLinear 1;
    div((nuEff*dev(T(grad(U))))) Gauss linear;
}

laplacianSchemes
{
    default         none;
    laplacian(nuEff,U) Gauss linear corrected;
    laplacian((1|A(U)),p) Gauss linear corrected;
    laplacian((1|A(U)),p_rgh) Gauss linear corrected;
    laplacian(rAU,p_rgh) Gauss linear corrected;
    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;
    interpolate(U)  linear;
}

snGradSchemes
{
    default         corrected;
}

fluxRequired
{
    default         no;
    p               ;
    p_rgh           ;
}
fvSolution

Code:
solvers
{
    p
    {
        solver          GAMG;
        tolerance       1e-06;
        relTol          0.01;
        smoother        GaussSeidel;
        cacheAgglomeration true;
        nCellsInCoarsestLevel 10;
        agglomerator    faceAreaPair;
        mergeLevels     1;
    }

    pFinal
    {
        solver          GAMG;
        tolerance       1e-06;
        relTol          0.01;
        smoother        GaussSeidel;
        cacheAgglomeration true;
        nCellsInCoarsestLevel 10;
        agglomerator    faceAreaPair;
        mergeLevels     1;
    }

    "(U|kl|kt|omega)"
    {
        solver          PBiCG;
        preconditioner  DILU;
        tolerance       1e-06;
        relTol          0;
    }

    "(U|kl|kt|omega)Final"
    {
        $U;
        tolerance       1e-06;
        relTol          0;
    }
}

PIMPLE
{
    nOuterCorrectors    50;
    nCorrectors         1;
    nNonOrthogonalCorrectors 2;
    pRefCell            0;
    pRefValue           0;
    residualControl
{
"(p|U)"
{
tolerance 1e-4;
relTol 0;
}
}
}

relaxationFactors
{
    fields
    {
    }
    equations
    {
        "U.*"           0.5;
        "kl.*"          1;
        "kt.*"          1;
        "omega.*"       1;
    }
}
The simulation blows blows at time=0.9s or little higher. The pulse lasts 20s and the reynolds number peaks(Re=2200) at 5s and 15s. If I choose adjustable time step then the it blows by showing time step in range of 10^-13 or so. If I choose a fixed time step of 0.01 or 0.02 or 0.005 it blows showing high courant number or unusually high flow rates. I have employed meshes of various kinds.....even pure tetrahedral mesh doesn't seem to work.

Here is the link to download the case directory along with the mesh:

https://www.dropbox.com/sh/u7sqlfvqtc6uynq/1MMKNryk7Z

Following are the error message snippets:

Code:
Expression volFlowInlet : sum=0.000155198
Expression volFlowOutlet1 : sum=-1.01936e+12
Expression volFlowOutlet2 : sum=1.01938e+12
Courant Number mean: 4.53531e+15 max: 9.33679e+17
Time = 1.23
Also the time step continuity errors goes exponentially high:

Code:
PIMPLE: iteration 50
DILUPBiCG: Solving for Ux, Initial residual = 0.0124649, Final residual = 3.77817e-07, No Iterations 5
DILUPBiCG: Solving for Uy, Initial residual = 0.00973289, Final residual = 3.76167e-07, No Iterations 5
DILUPBiCG: Solving for Uz, Initial residual = 0.0198487, Final residual = 4.315e-08, No Iterations 6
GAMG: Solving for p, Initial residual = 0.0537055, Final residual = 0.00027507, No Iterations 3
GAMG: Solving for p, Initial residual = 0.00744695, Final residual = 6.84879e-05, No Iterations 3
GAMG: Solving for p, Initial residual = 0.00143221, Final residual = 8.22524e-06, No Iterations 3
time step continuity errors : sum local = 3.47333e+09, global = -7.64618e+08, cumulative = -2.18106e+10
PIMPLE: not converged within 50 iterations
ExecutionTime = 55306.7 s ClockTime = 88055 s
Velocity fields have exploded and there is a flow reversal too. I tried avoiding flow reversal using outletinlet boudary condition since my outlets should behave as inlets in the beginning:

Code:
Code:
type  outletInlet;
outletValue uniform 0;
value uniform 0;
but it didn't help either.
Please suggest me how to make it work .
mayank.dce2k7 is offline   Reply With Quote

Reply

Thread Tools
Display Modes

Posting Rules
You may not post new threads
You may not post replies
You may not post attachments
You may not edit your posts

BB code is On
Smilies are On
[IMG] code is On
HTML code is Off
Trackbacks are On
Pingbacks are On
Refbacks are On


Similar Threads
Thread Thread Starter Forum Replies Last Post
Maximum words supported by the flow solver exceeded rae CFX 3 June 6, 2013 06:49
Compiled library vs. inInclude Files, DSMC solver crashes after run GPesch OpenFOAM Programming & Development 8 April 18, 2013 07:17
Laminar flow mech FLUENT 4 February 1, 2007 15:18
Laminar flow or Turbulent flow mech FLUENT 0 January 27, 2007 19:51
modelling laminar and turbulent flow in the same Ol FLUENT 2 November 25, 2005 03:52


All times are GMT -4. The time now is 01:16.