CFD Online Discussion Forums

CFD Online Discussion Forums (https://www.cfd-online.com/Forums/)
-   OpenFOAM Running, Solving & CFD (https://www.cfd-online.com/Forums/openfoam-solving/)
-   -   Simulation adaption for turbulent flow (https://www.cfd-online.com/Forums/openfoam-solving/201644-simulation-adaption-turbulent-flow.html)

Triggin May 7, 2018 07:33

Simulation adaption for turbulent flow
 
1 Attachment(s)
Hello to all.

I'm currently working on the following case:
I have a domain (done with sHM - checkMesh states OK), where i would like to simulate the flow behaviour with simpleFOAM.

The simulation converges for laminar flow (switching constant/turbulenceProperties to laminar). So far so good.

As soon as switch to RAS/kOmegaSST, it doesn't. Means: the simulation works, but the residuals (e.g. velocities) don't go below 1e-2 (playing with relaxations-factors doesn't lead to better results).

From my point of view, the reasons are limited:
  • Either k, omega or nut- BC's
  • or solver settings
I have the feeling, that when I'm simulating with turbulences, some unsteady-vortexes appears (please have a look at the probes in the lower subfigures > attachment).


My question: Has anyone experienced anything like this before? Do additional (special) settings have to be made in the solver, when the turbulence behaviour changes?
Should I run it with pimpleFoam and use the current solution as initial-file? I hope the simulation converges then...(?)



Here are my BC's:
U:
Code:

internalField  uniform (0 0 0);

boundaryField
{
    mainWalls
    {  type            fixedValue;
        value          uniform (0 0 0);
    }
    symmetry
    {  type            symmetryPlane;
    }
    inOut
    {  type            inletOutlet;
        inletValue      uniform (0 0 0);
        value          $internalField;
    }
    walls
    {  type            fixedValue;
        value          uniform (0 0 0);
    }
    rotWall
    {  type            rotatingWallVelocity;
        origin          (0 0 0.110835);
        axis            (1 0 0);
        omega          $rotRoll;
    }
    lowerInlet
    {  type            fixedValue;
        value          uniform (0 0 0);
    }
    upperInlet
    {  type            surfaceNormalFixedValue;
        refValue        uniform -5;
    }
}


p:
Code:

internalField  uniform $ambPress;

boundaryField
{
    mainWalls
    {  type            zeroGradient;
    }
    symmetry
    {  type            symmetryPlane;
    }
    inOut
    {  type            fixedValue;
        value          $internalField;
    }
    walls
    {  type            zeroGradient;
    }
    rotWall
    {  type            zeroGradient;
    }
    lowerInlet
    {  type            zeroGradient;
    }
    upperInlet
    {  type            zeroGradient;
    }
}


k:
Code:

internalField        uniform $myK;

boundaryField
{
    mainWalls
    {  type            kqRWallFunction;
        value          $internalField;
    }
    symmetry
    {  type            symmetryPlane;
    }
    inOut
    {  type            inletOutlet;
        inletValue      $internalField;
        value          $internalField;
    }
    walls
    {  type            kqRWallFunction;
        value          $internalField;
    }
    rotWall
    {  type            kqRWallFunction;
        value          $internalField;
    }
    lowerInlet
    {  type            fixedValue;
        value          $internalField;
    }
    upperInlet
    {  type            fixedValue;
        value          $internalField;
    }
}




omega:
Code:

internalField  uniform $myOmega;

boundaryField
{
    mainWalls
    {  type            omegaWallFunction;
        value          $internalField;
    }
    symmetry
    {  type            symmetryPlane;
    }
    inOut
    {  type            inletOutlet;
        inletValue      $internalField;
    value          $internalField;
    }
    walls
    {  type            omegaWallFunction;
        value          $internalField;
    }
    rotWall
    {  type            omegaWallFunction;
        value          $internalField;
    }
    lowerInlet
    {  type            fixedValue;
        value          $internalField;
    }
    upperInlet
    {  type            fixedValue;
        value          $internalField;
    }
}




nut:

Code:

internalField  uniform 0;

boundaryField
{
    mainWalls
    {  type            nutkWallFunction;
        value          uniform 0;
    }
    symmetry
    {  type            symmetryPlane;
    }
    inOut
    {  // type            inletOutlet;
    // inletValue      uniform 0;
    // value          $internalField;
        type            calculated;
        value          $internalField;
    }
    walls
    {  type            nutkWallFunction;
        value          uniform 0;
    }
    rotWall
    {  type            nutkWallFunction;
        value          uniform 0;
    }
    lowerInlet
    {  type            calculated;
        value          uniform 0;
    }
    upperInlet
    {  type            calculated;
        value          uniform 0;
    }
}


And my fvSolution:
Code:

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

    "(U|C|k|epsilon|omega|f|v2)"
    {  solver          smoothSolver;
        smoother        symGaussSeidel;
        tolerance      1e-05;
        relTol          0.1;
    }

    Phi
    { $p; }
}

SIMPLE
{  nCorrectors              1;
    nNonOrthogonalCorrectors 1;
    consistent            yes;
    residualControl
    {  p              1e-4;
        U              1e-5;
    C              1e-4;
        "(k|epsilon|omega|f|v2)" 1e-3;
    }
}

potentialFlow
{ nNonOrthogonalCorrectors 10; }

relaxationFactors
{
    equations
    {
        U              0.95;
        ".*"            0.95;
    }
}


piu58 May 7, 2018 08:07

> the flow behaviour with simpleFOAM.

If you use a steady state solver you have to ensure that a steady state solution exist. Otherwise you don't get a stable solution.

I would try pimpleFoam and look how the flow evolves. May be, there are problematic areas in your geometry, like inconsistent b.c. or problems with the mesh. With pimpleFoam it is easier to detect wehter the solution moves in a nonphysical direction and where this starts.

It may be, however, that you case is transient per se, and a steady state solver cannot be used.

Triggin May 7, 2018 13:13

2 Attachment(s)
Dear Uwe

Thank you for your reply. You are right. It is possible that I could have some at least some kind of taylor-couette vortexes (see: https://youtu.be/BNiocOsgxW8)

I've started a transient calculation. But somehow the initial residuals won't go below 1e-3 (see the pictures) - there is a jump in the residuals, where i changes the relaxations-factors to see their influence (from 0.8 to 0.6).

Ok, I have to admit, that I'm not very familiar with transient calculations. But none of the initial results are going below the mentioned 1e-3. Where does this come from?

Here is the controlDict:
Code:

application    pimpleFoam;
startFrom      latestTime;
startTime      0;
stopAt          endTime;
endTime        1;
deltaT          0.001;
writeControl    adjustableRunTime;
writeInterval  0.005;
purgeWrite      0;
writeFormat    ascii;
writePrecision  6;
writeCompression off;
timeFormat      general;
timePrecision  6;
runTimeModifiable yes;
adjustTimeStep  yes;
maxCo          50;

and the fvSolution file:
Code:

{  p
    {  solver          GAMG;
        tolerance      1e-7;
        relTol          0.01;
        smoother        GaussSeidel;
        nPreSweeps      0;
        nPostSweeps    2;
        cacheAgglomeration on;
        agglomerator    faceAreaPair;
        nCellsInCoarsestLevel 100;
        mergeLevels    1;
    }
    pFinal
    {  $p;
        relTol          0;
    }

    "(U|k|epsilon|omega)"
    {  solver          smoothSolver;
        smoother        symGaussSeidel;
        tolerance      1e-8;
        relTol          0.1;
    }

    "(p|U|k|epsilon|omega)Final"
    {  $U;
        tolerance      1e-05;
        relTol          0;
    }

    Phi
    { $p; }
}

PIMPLE
{  nNonOrthogonalCorrectors 1; // 10;
    nOuterCorrectors 40;
    nCorrectors    2;

    residualControl
    { U
      { tolerance  1e-4;
        relTol      0;
      }
      p
      { tolerance  1e-4;
        relTol      0;
      }
    }
}

potentialFlow
{ nNonOrthogonalCorrectors 10; }

relaxationFactors
{  fields
    { p                  0.6;
      pFinal                1;
    }
    equations
    { "U|k|epsilon|omega"        0.6;
      "(U|k|epsilon|omega)Final" 1;
    }
}

Every help is very welcome.
Best, Triggin

Triggin May 17, 2018 03:38

1 Attachment(s)
Dear all

I finally managed it to converge the calculations...well at least to a residuals below 5e-4.
The main error was the use of a second order scheme, which lead to oscillations. At the moment I'm using upwind.

Nonetheless I'm not very satisfied with the results at the moment. On one hand because of the high residuals with oscillations. At the other hand because of the high numerical diffusion of the upwind scheme.

Therefore I tried to run a second simulation, with the results from the the upwind-calculation as initial conditions.
But: No matter what kind of scheme I'm going to use, the simulation always starts to diverge. What could be the reason for this?

Any help, would be appreciated.

My mesh should be fine:
  • Max skewness: 3.99
  • Max non-ortho: 59.99 / av. 10.3
  • Max aspect ratio: 15.4
  • No collapsing boundary layers
fvScheme:
Code:

ddtSchemes
{ default        steadyState; }
gradSchemes
{ default        Gauss linear; }
divSchemes
{
    default        none;
    div(phi,U)      bounded Gauss upwind;
    div(phi,C)      bounded Gauss linearUpwind grad(C);
    div(phi,k)      bounded Gauss limitedLinear 1;
    div(phi,omega)  bounded Gauss limitedLinear 1;
    div((nuEff*dev2(T(grad(U))))) Gauss linear;
}
laplacianSchemes
{ default        Gauss linear corrected; }
interpolationSchemes
{ default        linear; }
snGradSchemes
{ default        corrected; }
wallDist
{ method meshWave; }



All times are GMT -4. The time now is 17:49.