CFD Online Discussion Forums

CFD Online Discussion Forums (http://www.cfd-online.com/Forums/)
-   OpenFOAM Running, Solving & CFD (http://www.cfd-online.com/Forums/openfoam-solving/)
-   -   Diverging Simulation using rhoSimplecFoam, kOmegaSST (http://www.cfd-online.com/Forums/openfoam-solving/111455-diverging-simulation-using-rhosimplecfoam-komegasst.html)

AndreasP. January 7, 2013 16:35

Diverging Simulation using rhoSimplecFoam, kOmegaSST
 
Hello everybody,

I have been trying for quite some time to run a compressible simulation of a NACA 0012 profile, using rhoSimplecFoam and kOmegaSST for turbulence modeling.
My mesh is unstructured, checkMesh indicates no problems.
Running a simulation at Mach 0.7 generates following output:

Starting time loop

Time = 0.1

GAMG: Solving for Ux, Initial residual = 1, Final residual = 0.0615838, No Iterations 1
GAMG: Solving for Uy, Initial residual = 1, Final residual = 0.0598451, No Iterations 1
GAMG: Solving for p, Initial residual = 1, Final residual = 27.4647, No Iterations 1000
time step continuity errors : sum local = 0.793546, global = 0.47889, cumulative = 0.47889
rho max/min : 0.1 0.1
GAMG: Solving for h, Initial residual = 0.627704, Final residual = 0.906171, No Iterations 1000
GAMG: Solving for omega, Initial residual = 0.00118687, Final residual = 5.36443e-06, No Iterations 1
smoothSolver: Solving for k, Initial residual = 1, Final residual = 0.0494026, No Iterations 2
ExecutionTime = 7.73 s ClockTime = 7 s

Time = 0.2

GAMG: Solving for Ux, Initial residual = 0.335287, Final residual = 0.00814846, No Iterations 2
GAMG: Solving for Uy, Initial residual = 0.0892067, Final residual = 0.00223789, No Iterations 2
#0 Foam::error::printStack(Foam::Ostream&) in "/opt/openfoam201/platforms/linux64GccDPOpt/lib/libOpenFOAM.so"
#1 Foam::sigFpe::sigHandler(int) in "/opt/openfoam201/platforms/linux64GccDPOpt/lib/libOpenFOAM.so"
#2 in "/lib/libc.so.6"
(...
...)
in "/opt/openfoam201/platforms/linux64GccDPOpt/bin/rhoSimplecFoam"
#10 __libc_start_main in "/lib/libc.so.6"
#11
in "/opt/openfoam201/platforms/linux64GccDPOpt/bin/rhoSimplecFoam"
Floating point exception

I presume there's something off when calculating for p and h. I have no clue what the problem is though, an incompressible simulation, using simpleFoam and Spallart-Allmaras works fine.

Switching transonic off in fvSolution seems to stabilize the simulate at least so much to allow a minimal amount of postprocessing; every single variable diverges massively. Changing the relaxation factors has no significant effect.

I am quite new to OpenFoam, and have generated my files according to various threats in this forum as well as to the tutorial cases.
Here are my files fvSolution, fvSchemes, p and U.

Code:

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

solvers
{
   
    p0
    {
        solver          GAMG;
        tolerance      1e-10;
        relTol          0.1;
        smoother        GaussSeidel;
        nPreSweeps      0;
        nPostSweeps    2;
        nFinestSweeps  2;
        cacheAgglomeration true;
        nCellsInCoarsestLevel 20;
        agglomerator    faceAreaPair;
        mergeLevels    1;
    }
    p1
    {
        solver          GAMG;
        tolerance      1e-10;
        relTol          0.1;
        smoother        GaussSeidel;
        nPreSweeps      0;
        nPostSweeps    2;
        nFinestSweeps  2;
        cacheAgglomeration true;
        nCellsInCoarsestLevel 20;
        agglomerator    faceAreaPair;
        mergeLevels    1;
    }
    p
    {
        solver          GAMG;
        tolerance      1e-10;
        relTol          0.1;
        smoother        GaussSeidel;
        nPreSweeps      0;
        nPostSweeps    2;
        nFinestSweeps  2;
        cacheAgglomeration true;
        nCellsInCoarsestLevel 20;
        agglomerator    faceAreaPair;
        mergeLevels    1;
    }
    omega0
    {
        solver          PBiCG;
        preconditioner  DILU;
        tolerance      1e-10;
        relTol          0.1;
    }

    omega1
    {
        solver          smoothSolver;
        smoother        GaussSeidel;
        nSweeps        1;
        tolerance      1e-10;
        relTol          0.1;
    }

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

    U0
    {
        solver          PBiCG;
        preconditioner  DILU;
        tolerance      1e-10;
        relTol          0.1;
    }

    U1
    {
        solver          smoothSolver;
        smoother        GaussSeidel;
        nSweeps        1;
        tolerance      1e-10;
        relTol          0.1;
    }

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

    h0
    {
        solver          PBiCG;
        preconditioner  DILU;
        tolerance      1e-10;
        relTol          0.1;
    }

    h1
    {
        solver          smoothSolver;
        smoother        GaussSeidel;
        nSweeps        1;
        tolerance      1e-10;
        relTol          0.1;
    }

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

    k0
    {
        solver          PBiCG;
        preconditioner  DILU;
        tolerance      1e-10;
        relTol          0.1;
    }

    k1
    {
        solver          PBiCG;
        preconditioner  DILU;
        tolerance      1e-10;   
        relTol          0.1;
    }

    k
    {
        solver          smoothSolver;
        smoother        GaussSeidel;
        tolerance        1e-10;
        relTol          0.1;
        nSweeps          1;
    }

   

 
}

SIMPLE
{
    nNonOrthogonalCorrectors 0;
    rhoMin          rhoMin [1 -3 0 0 0] 0.1;
    rhoMax          rhoMax [1 -3 0 0 0] 1.0;
    transonic      yes;                                     
    pRefPoint        (100 0 0);
    pRefValue      11000;
}

relaxationFactors
{
    p              0.5;
    U              0.7;
    k              0.5;
    omega          0.5;           
    rho            1;                     
    h            0.95;
}


relaxationFactors0
{
    p            0.3;
    rho          0.4;               
    U            0.7;
    h            0.7;
    k            0.7;
    omega        0.7;
}

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

Code:

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

ddtSchemes
{
    default            steadyState;
}

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

divSchemes
{
    default            none;

    div(phi,U)          Gauss upwind; //limitedLinearV 1; //upwind;
    div((muEff*dev2(T(grad(U)))))      Gauss linear;
    div(phi,h)          Gauss upwind; //limitedLinear 1; //upwind;

    div(phi,k)          Gauss upwind;
    div(phi,omega)    Gauss upwind;
    div(phid,p)        Gauss upwind;
    div(U,p)            Gauss linear;
}

laplacianSchemes
{
    default        Gauss linear limited 0.5;
}

interpolationSchemes
{
    default                none;
    interpolate(rho)        linear;
    div(U,p)                upwind phi;
    interpolate((psi*U))    linear;
    interpolate(U)          linear;
    UD                      upwind phid;
    interpolate(p)          linear;
    interpolate(((rho|(A(U)-H(1)))-(rho|A(U)))) linear;
    interpolate((rho*U))    linear;
}

snGradSchemes
{
    default                corrected;
}

fluxRequired
{
    default                no;
    p;
    pCorr;
}

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

Code:

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

internalField  uniform 11000;

boundaryField
{
    NONE1      //Nearfield1
    {
        type            empty;
       
    }

    NONE2    //Nearfield2
    {
        type            empty;
    }
 
    NONE3    //Pressure Side
    {
        type            zeroGradient;
    }

    NONE4    //Nearfield3
    {
        type            empty;
    }
   
    NONE5    //Nearfield4
    {
        type            empty;
    }

    NONE6    //Suction Side
    {
        type            zeroGradient;
    }

    NONE7    //Wake resolution1
    {
        type            empty;
    }

    NONE8    //Wake resolution2
    {
        type            empty;
    }
   
    NONE9    //Trailing Edge
    {
        type            zeroGradient;
    }

    NONE10    //Farfield1
    {
        type            empty;
   
    }
    NONE11    //Farfield2
    {
        type            empty;
     
    }
   
    NONE12    //Outer rim, lower half
    {
        type            freestreamPressure;   
        //type            fixedValue;
        value          uniform 11000;
    }
   
    NONE13    //Outer rim, upper half
    {
        type            freestreamPressure;
        //type            fixedValue;
        value          uniform 11000;
    }
}

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

Code:

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

internalField  uniform (239.266 20.933 0);


boundaryField
{
    NONE1      //Nearfield1
    {
        type            empty;
       
    }

    NONE2    //Nearfield2
    {
        type            empty;
    }
 
    NONE3    //Pressure Side
    {
        type            fixedValue;
        value          uniform (0 0 0);
    }

    NONE4    //Nearfield3
    {
        type            empty;
    }
   
    NONE5    //Nearfield4
    {
        type            empty;
    }

    NONE6    //Suction Side
    {
      type            fixedValue;
      value          uniform (0 0 0);
    }

    NONE7    //Wake resolution1
    {
        type            empty;
    }

    NONE8    //Wake resolution2
    {
        type            empty;
    }
   
    NONE9    //Trailing Edge
    {
        type            fixedValue;
        value          uniform (0 0 0);
    }

    NONE10    //Farfield1
    {
        type            empty;
   
    }
    NONE11    //Farfield2
    {
        type            empty;
     
    }
   
    NONE12    //Outer rim, lower half
    {
        type            freestream;
        freestreamValue uniform (239.266 20.933 0); 
    }
   
    NONE13    //Outer rim, upper half
    {
      type            freestream;
      freestreamValue uniform (239.266 20.933 0);
    }
}

Thanks in advance for your help.

AndreasP. January 10, 2013 02:16

*Bump*

Any help would be greatly appreciated :)

hannes January 27, 2013 06:19

Quote:

Originally Posted by AndreasP. (Post 400610)
Hello everybody,
rho max/min : 0.1 0.1

Looks like you always hit the lower limit of rho (given in fvSolution/SIMPLE-section).
I think you should adapt "rhoMin" to your quite low abient pressure of 11kPa.

AndreasP. January 27, 2013 06:37

Hey Hannes,
thanks for the reply!
Unfortunately, decreasing rhoMin seems to have no effect...

fredo490 March 30, 2013 11:12

It's a bit late to answer but I often see thread with problems over rhoMin and rhoMax.

Those two values are some "limiters" to avoid unphysical values of rho. This limiter has to be temporary (often at the beginning of the computation to avoid any divergence).

You need to know first the freestream value of rho. If you are lazy you can use Wolfram: http://www.wolframalpha.com/input/?i...=&equal=Submit

Then you have to choose your rhoMin and rhoMax arround the freestream value. And here you need to be smart and use your experience because the min and max often depend of the geometry. Until Mach 0.5 the density usually change of less than +/- 0.2.


All times are GMT -4. The time now is 06:09.