k-omega SST flow over a cylinder

GP_FDU · May 30, 2013, 08:47

I simulated flow over a cylinder with k-omega SST model. But the result was very wrong.
Here is some configuration of my case. D=3.81cm, U=0.19m/s, fluid is set to be water.
U and p B.C. is set as usual.
B.C. of k is fixedValue=1e-12 at the cylinder.
B.C. of nut is fixedValue=nutUSpaldingWallFunciton at the cylinder or fixedValue=1e-12.
B.C. of omega is omegaWallFunction.
y+is around 1 according to y+=0.172*delta_y/D*Re^0.9.

With this configuration, there is vortex shedding and the frequency seems to be Ok, but the pressure at the cylinder is totally wrong, which leads to a wrong drag and lift force.

Anybody know what is wrong with my case, please give me some advice.

OpenF · May 30, 2013, 09:20

Quote:

Originally Posted by GP_FDU

I simulated flow over a cylinder with k-omega SST model. But the result was very wrong.
Here is some configuration of my case. D=3.81cm, U=0.19m/s, fluid is set to be water.
U and p B.C. is set as usual.
B.C. of k is fixedValue=1e-12 at the cylinder.
B.C. of nut is fixedValue=nutUSpaldingWallFunciton at the cylinder or fixedValue=1e-12.
B.C. of omega is omegaWallFunction.
y+is around 1 according to y+=0.172*delta_y/D*Re^0.9.

With this configuration, there is vortex shedding and the frequency seems to be Ok, but the pressure at the cylinder is totally wrong, which leads to a wrong drag and lift force.

Anybody know what is wrong with my case, please give me some advice.

It's difficult to say, maybe you can put fvScheme file, fvSolution file and your boundaries conditions online.

GP_FDU · May 30, 2013, 11:01

Hi Anselme,

the main boundaries are as below:
omega boundary:

Code:

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

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

internalField   uniform 0;

boundaryField
{
    right
    {
        type            zeroGradient;
    }
    up
    {
        type            symmetryPlane;
    }
    left
    {
        type            fixedValue;
        value           uniform 0.2;
    }
    cylinder
    {
        type            omegaWallFunction;
        Cmu             0.09;
        kappa           0.41;
        E               9.8;
        beta1           0.075;
        value           uniform 1e7;
    }
    defaultFaces
    {
        type            empty;
    }
}


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

k boundary

Code:

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

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

internalField   uniform 0;

boundaryField
{
    right
    {
        type            zeroGradient;
    }
    up
    {
        type            symmetryPlane;
    }
    left
    {
        type            fixedValue;
        value           uniform 0;
    }
    cylinder
    {
        type            fixedValue;
        value           uniform 1e-12;
    }
    defaultFaces
    {
        type            empty;
    }
}


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

nut boundary:

Code:

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

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

internalField   uniform 0;

boundaryField
{
    right
    {
        type            zeroGradient;
    }
    up
    {
        type            symmetryPlane;
    }
    left
    {
        type            fixedValue;
        value           uniform 0;
    }
    cylinder
    {
        type            nutUSpaldingWallFunction;
        Cmu             0.09;
        kappa           0.41;
        E               9.8;
        value           uniform 1e-12;
    }
    defaultFaces
    {
        type            empty;
    }
}


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

GP_FDU · May 30, 2013, 11:05

The fvscheme and fvsolution file are :

Code:

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

ddtSchemes
{
    default Euler;
}

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

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

laplacianSchemes
{
    default         Gauss linear limited corrected 0.5;
}

interpolationSchemes
{
    default         linear;
}

snGradSchemes
{
    default         corrected;
}

fluxRequired
{
    default         no;
    pcorr           ;
    p;
}

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

Code:

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

solvers
{
    pcorr
    {
        solver           GAMG;
        tolerance        0.02;
        relTol           0;
        smoother         GaussSeidel;
        nPreSweeps       0;
        nPostSweeps      2;
        cacheAgglomeration on;
        agglomerator     faceAreaPair;
        nCellsInCoarsestLevel 10;
        mergeLevels      1;
    }

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

    pFinal
    {
        $p;
        tolerance        1e-9;
        relTol           0;
    }

    "(U|k|omega)"
    {
        solver          PBiCG;
        preconditioner  DILU;
        tolerance       1e-07;
        relTol          0.1;
    }

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

    cellDisplacement
    {
        solver          GAMG;
        tolerance       1e-5;
        relTol          0;
        smoother        GaussSeidel;
        cacheAgglomeration true;
        nCellsInCoarsestLevel 10;
        agglomerator    faceAreaPair;
        mergeLevels     1;
    }
}

PIMPLE
{
    correctPhi          yes;
    nOuterCorrectors    2;
    nCorrectors         1;
    nNonOrthogonalCorrectors 0;
}

relaxationFactors
{
    fields
    {
        p               0.3;
    }
    equations
    {
        "(U|k|omega)"   0.7;
        "(U|k|omega)Final" 1.0;
    }
}

cache
{
    grad(U);
}

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

GP_FDU · May 31, 2013, 00:10

Does anybody have some ideas?

fredo490 · June 2, 2013, 13:40

1) you can try to set 1 or 2 non orthogonal correction (nNonOrthogonalCorrectors ) depending of your mesh.

2) Is it normal that your pcorr solver has only a tolerance of 0.02 ?!

3) Your pressure relaxation number looks a bit low for an unsteady simulation. Ideally it has to be close to 1 to really catch the unsteadiness. 0.8 would be much better I think.

4) Where did you get your k and omega value for the inlet ? Those two values are really important to predict the boundary layer separation.

5) use the "yPlusRAS -latestTime" to really compute the yplus of your case

...

GP_FDU · June 4, 2013, 08:46

hi Heckmann,

I modified the files according to your suggestion. But It seemed to have little influence on the result.
I tried two cases of which the yplus is 1.1 and 0.8 computed by the utility.

Something interesting is that the pressure near the cylinder seems to be OK if you ignore the sudden change of them. It suddenly change to a very high or very low value. it looked like this:

fredo490 · June 4, 2013, 08:48

how is your mesh ? Can you give a screenshot ?
Also, you can run a "checkMesh" utility and give the report

GP_FDU · June 4, 2013, 10:40

I checked mesh, it showed the mesh is all OK.
here is the report:

Code:

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

Create polyMesh for time = 0

Time = 0

Mesh stats
    points:           115260
    internal points:  0
    faces:            228930
    internal faces:   113670
    cells:            57100
    faces per cell:   6
    boundary patches: 5
    point zones:      0
    face zones:       0
    cell zones:       0

Overall number of cells of each type:
    hexahedra:     57100
    prisms:        0
    wedges:        0
    pyramids:      0
    tet wedges:    0
    tetrahedra:    0
    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                  
    right               130      262      ok (non-closed singly connected)  
    up                  540      1084     ok (non-closed singly connected)  
    left                130      262      ok (non-closed singly connected)  
    cylinder            260      520      ok (non-closed singly connected)  
    defaultFaces        114200   115260   ok (non-closed singly connected)  

Checking geometry...
    Overall domain bounding box (-0.2286 -0.2667 -0.01905) (0.5715 0.2667 0.01905)
    Mesh (non-empty, non-wedge) directions (1 1 0)
    Mesh (non-empty) directions (1 1 0)
    All edges aligned with or perpendicular to non-empty directions.
    Boundary openness (5.842462579e-18 -3.1197547e-17 1.227081693e-15) OK.
    Max cell openness = 2.093623056e-16 OK.
    Max aspect ratio = 10.99434892 OK.
    Minimum face area = 1.690836995e-08. Maximum face area = 0.0004056543995.  Face area magnitudes OK.
    Min volume = 6.442088951e-10. Max volume = 1.940056445e-06.  Total volume = 0.01621663178.  Cell volumes OK.
    Mesh non-orthogonality Max: 44.11247936 average: 9.104769412
    Non-orthogonality check OK.
    Face pyramids OK.
    Max skewness = 0.5934534386 OK.
    Coupled point location match (average 0) OK.

Mesh OK.

I know the mesh is a little coarse. But the I tried some fine mesh, the result showed no big difference.

fredo490 · June 4, 2013, 11:04

From the picture on post 7, it looks like you have a divergence...

1) Personally I would not set k = 0 at the inlet since it doesn't have any physical meaning. You can use the tool http://www.cfd-online.com/Tools/turbulence.php to estimate it. Roughly, you can set it as equal to 1e-6.

2) I'm not familiar with pimpleFoam (I work more on rhoPimpleFoam) but I think that you can set the inlet boundary of nut as "calculated" or "zeroGradient" since nut is derived from k and omega (If I remember well my classics).

3) maybe you have a problem with the time step. What is your maximum Courant Number (given at each iteration) ? It should be lower than 0.2 to avoid divergence or some sporadic crazy behavior. You can also try to use a second order scheme for ddt or maybe a cranknicholson.

4) try to change the turbulent div scheme to linearUpwind. It's often more stable.

GP_FDU · June 5, 2013, 05:10

Hi Heckmann,

Thanks very much for your patience and kindness. Following your suggestions, I finally succeed to get a reasonable result.

I changed the maxCo to 0.2 the pressure relaxation number to 0.8, and the time scheme to CrankNicolson 0.5.

I am not good at CFD. For many parameters in the case files, sometimes I don't know how to set them -- I don't know what value is appropriate. I am trying to learn. Could you give me some advise?

fredo490 · June 5, 2013, 05:15

CFD is most of the time a matter of experience. I've never found any book or class able to really prepare a student to deal with "advanced" cfd. The only thing you have to do is to keep working and trying

In French we say: "c'est en forgeant qu'on devient forgeron", that you could roughly translate to "practice makes perfect".

OpenF · June 5, 2013, 05:54

Quote:

Originally Posted by GP_FDU

Hi Heckmann,

Thanks very much for your patience and kindness. Following your suggestions, I finally succeed to get a reasonable result.

I changed the maxCo to 0.2 the pressure relaxation number to 0.8, and the time scheme to CrankNicolson 0.5.

I am not good at CFD. For many parameters in the case files, sometimes I don't know how to set them -- I don't know what value is appropriate. I am trying to learn. Could you give me some advise?

I completely agree with Frederick, particularly with OpenFOAM, where you set all the parameters of the computation. For example, with Fluent, you still don't understand what you are doing ..

I use OpenFOAM since two month and I agree with you, it is not always easy ..

I can only advise you to learn about CFD online and trust people who have written hundreds of posts..

I am particularly interested in your blockMeshDict your cylinder. I would be very grateful if you can give me your file.

Thank you,
Anselme

GP_FDU · June 5, 2013, 06:04

Hi, Bonnet

My blockmeshDict is copied from the tutorial case in potentialFoam. You can download my file from the attachment, or go to the tutorial to find the original file.

OpenF · June 5, 2013, 07:07

Thank you, I didn't see that ! But It's only a semi-cylinder. So I duplicated it but It seems to have a problem with my mesh regarding the connections in the symmetry plane. If you look on the picture, there is abnormal features in the symmetry plane

Code:

/*--------------------------------*- C++ -*----------------------------------*\
| =========                 |                                                 |
| \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox           |
|  \\    /   O peration     | Version:  2.2.0                                 |
|   \\  /    A nd           | Web:      www.OpenFOAM.org                      |
|    \\/     M anipulation  |                                                 |
\*---------------------------------------------------------------------------*/
FoamFile
{
      version         2.0; 
      format          ascii; 
      root            ""; 
      case            ""; 
      instance        ""; 
      local           ""; 
   
      class           dictionary; 
      object          blockMeshDict; 
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
convertToMeters 0.0381;
vertices 
(
    
 (0.5 0 -0.5) //0
        (2 0 -0.5)
        (15 0 -0.5)
        (15 1.414214 -0.5)
        (1.414214 1.414214 -0.5)
        (0.353553 0.353553 -0.5)
        (15 7 -0.5)
        (1.414214 7 -0.5)
        (0 7 -0.5)
        (0 2 -0.5)
        (0 0.5 -0.5) //10
        (-0.5 0 -0.5)
        (-2 0 -0.5)
        (-6 0 -0.5)
        (-6 1.414214 -0.5)
        (-1.414214 1.414214 -0.5)
        (-0.353553 0.353553 -0.5)
        (-6 7 -0.5)
        (-1.414214 7 -0.5)
 (0.5 0 0.5)
        (2 0 0.5) //20
        (15 0 0.5)
        (15 1.414214 0.5)
        (1.414214 1.414214 0.5)
        (0.353553 0.353553 0.5)
        (15 7 0.5)
        (1.414214 7 0.5)
        (0 7 0.5)
        (0 2 0.5)
        (0 0.5 0.5)
        (-0.5 0 0.5) //30
        (-2 0 0.5)
        (-6 0 0.5)
        (-6 1.414214 0.5)
        (-1.414214 1.414214 0.5)
        (-0.353553 0.353553 0.5)
        (-6 7 0.5)
        (-1.414214 7 0.5) //37

 (15 -1.414214 -0.5)  //38
        (1.414214 -1.414214 -0.5)
        (0.353553 -0.353553 -0.5)
        (15 -7 -0.5)
        (1.414214 -7 -0.5)
        (0 -7 -0.5)
        (0 -2 -0.5)
        (0 -0.5 -0.5) //45
 (-6 -1.414214 -0.5) //46
        (-1.414214 -1.414214 -0.5)
        (-0.353553 -0.353553 -0.5)
        (-6 -7 -0.5)
        (-1.414214 -7 -0.5) //50
 (15 -1.414214 0.5) //51
        (1.414214 -1.414214 0.5)
        (0.353553 -0.353553 0.5)
        (15 -7 0.5)
        (1.414214 -7 0.5)
        (0 -7 0.5)
        (0 -2 0.5)
        (0 -0.5 0.5) //58
 (-6 -1.414214 0.5) //59
        (-1.414214 -1.414214 0.5)
        (-0.353553 -0.353553 0.5)
        (-6 -7 0.5)
        (-1.414214 -7 0.5) //63
);

blocks
(
    hex (5 4 9 10 24 23 28 29) (100 25 1) simpleGrading (60 1 1)  //0
    hex (0 1 4 5 19 20 23 24) (100 25 1) simpleGrading (60 1 1)   //1
    hex (1 2 3 4 20 21 22 23) (150 25 1) simpleGrading (2 1 1)   //2
    hex (4 3 6 7 23 22 25 26) (150 40 1) simpleGrading (2 5 1)   //3
    hex (9 4 7 8 28 23 26 27) (25 40 1) simpleGrading (1 5 1)   //4
    hex (15 16 10 9 34 35 29 28) (100 25 1) simpleGrading (0.0166667 1 1) //5
    hex (12 11 16 15 31 30 35 34) (100 25 1) simpleGrading (0.0166667 1 1) //6
    hex (13 12 15 14 32 31 34 33) (40 25 1) simpleGrading (1 1 1)  //7
    hex (14 15 18 17 33 34 37 36) (40 40 1) simpleGrading (1 5 1)  //8
    hex (15 9 8 18 34 28 27 37) (25 40 1) simpleGrading (1 5 1)   //9
      hex (40 39 44 45 53 52 57 58) (100 25 1) simpleGrading (60 1 1) //10
      hex (0 1 39 40 19 20 52 53) (100 25 1) simpleGrading (60 1 1)  //11
      hex (1 2 38 39 20 21 51 52) (150 25 1) simpleGrading (2 1 1)   //12
      hex (39 38 41 42 52 51 54 55) (150 40 1) simpleGrading (2 5 1) //13
      hex (44 39 42 43 57 52 55 56) (25 40 1) simpleGrading (1 5 1)  //14
      hex (47 48 45 44 60 61 58 57) (100 25 1) simpleGrading (0.0166667 1 1) //15 
      hex (12 11 48 47 31 30 61 60) (100 25 1) simpleGrading (0.0166667 1 1) //16
      hex (13 12 47 46 32 31 60 59) (40 25 1) simpleGrading (1 1 1) //17
      hex (46 47 50 49 59 60 63 62) (40 40 1) simpleGrading (1 5 1) //18
      hex (47 44 43 50 60 57 56 63) (25 40 1) simpleGrading (1 5 1) //19
 
 
);
edges
(
    arc 0 5 (0.469846 0.17101 -0.5)
    arc 5 10 (0.17101 0.469846 -0.5)
    arc 1 4 (1.879386 0.68404 -0.5)
    arc 4 9 (0.68404 1.879386 -0.5)
    arc 19 24 (0.469846 0.17101 0.5)
    arc 24 29 (0.17101 0.469846 0.5)
    arc 20 23 (1.879386 0.68404 0.5)
    arc 23 28 (0.68404 1.879386 0.5)
    arc 11 16 (-0.469846 0.17101 -0.5)
    arc 16 10 (-0.17101 0.469846 -0.5)
    arc 12 15 (-1.879386 0.68404 -0.5)
    arc 15 9 (-0.68404 1.879386 -0.5)
    arc 30 35 (-0.469846 0.17101 0.5)
    arc 35 29 (-0.17101 0.469846 0.5)
    arc 31 34 (-1.879386 0.68404 0.5)
    arc 34 28 (-0.68404 1.879386 0.5)
    arc 0 40 (0.469846 -0.17101 -0.5)
    arc 40 45 (0.17101 -0.469846 -0.5)
    arc 1 39 (1.879386 -0.68404 -0.5)
    arc 39 44 (0.68404 -1.879386 -0.5)
    arc 19 53 (0.469846 -0.17101 0.5)
    arc 53 58 (0.17101 -0.469846 0.5)
    arc 20 52 (1.879386 -0.68404 0.5)
    arc 52 57 (0.68404 -1.879386 0.5)
    arc 11 48 (-0.469846 -0.17101 -0.5)
    arc 48 45 (-0.17101 -0.469846 -0.5)
    arc 12 47 (-1.879386 -0.68404 -0.5)
    arc 47 44 (-0.68404 -1.879386 -0.5)
    arc 30 61 (-0.469846 -0.17101 0.5)
    arc 61 58 (-0.17101 -0.469846 0.5)
    arc 31 60 (-1.879386 -0.68404 0.5)
    arc 60 57 (-0.68404 -1.879386 0.5)
);
boundary
(
    inlet
    {
        type patch;
        faces
        (
          (14 13 32 33) 
          (17 14 33 36) 
          (46 13 32 59) 
          (49 46 59 62) 
        );
    }
    outlet
    {
        type patch;
        faces
        (
          (2 3 22 21) 
          (3 6 25 22) 
          (2 38 51 21) 
          (38 41 54 51) 
        );
    }

    cylinder
    {
        type wall;
        faces
        (
 (30 11 16 35)
 (35 16 10 29)
 (10 29 24 5)
 (24 5 0 19)
 (0 19 53 40)
 (53 40 45 58)
 (45 58 61 48)
 (48 61 30 11)
        );
    }
  frontAndBack
  {
    type empty;
    faces
       (
 (14 15 18 17)
 (15 18 8 9)
 (8 9 4 7)
 (7 6 3 4)
 (14 15 12 13)
 (12 15 16 11)
 (15 16 10 9)
 (9 10 5 4)
 (4 5 0 1)
 (4 1 2 3)
 (36 37 34 33)
 (37 34 28 27)
 (27 28 23 26)
 (26 23 22 25)
 (33 34 31 32)
 (34 31 30 35)
 (34 35 29 28)
 (28 29 24 23)
 (23 24 19 20)
 (23 20 21 22)

 (46 47 50 49)
 (47 50 43 44)
 (43 44 39 42)
 (42 41 38 39)
 (46 47 12 13)
 (12 47 48 11)
 (47 48 45 44)
 (44 45 40 39)
 (39 40 0 1)
 (39 1 2 38)
 (62 63 60 59)
 (63 60 57 56)
 (56 57 52 55)
 (55 52 51 54)
 (59 60 31 32)
 (60 31 30 61)
 (60 61 58 57)
 (57 58 53 52)
 (52 53 19 20)
 (52 20 21 51)
       );
  }
);
   
// ************************************************************************* //

GP_FDU · June 5, 2013, 10:04

Use mirrorMesh utility to create a mirror mesh. Set control parameters in mirrorMeshDict in system directory.

OpenF · June 5, 2013, 11:08

thank you!
I did what you said, I put the file "mirrorMeshDict" in the system directory, but before and after the computation, I can't generate the mesh, I still have a half-cylinder .. Maybe I need to put something in the controlDict or somewhere?

GP_FDU · June 5, 2013, 11:35

which version do you use? if of1.6ext, please check the 0 directory.

MDB · June 8, 2013, 13:29

Dear all,

I have been following this thread because I am starting to use OpenFOAM and trying to simulate the case of a circular cylinder, starting from a simulation in steady and 2D conditions.

Not extent of some trouble I have finnally built a case which runs properly, but I get values of drag around 0.012 only! My fvSchemes and fvSolution scripts are below these lines. Any thoughts will be more than welcome

Code:

/*--------------------------------*- 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;
    object      fvSchemes;
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
ddtSchemes
{
    default         steadyState;
}
gradSchemes
{
    default         Gauss linear;
    grad(p)     Gauss linear;
    grad(U)     Gauss linear;
}
divSchemes
{
    default         none;
    div(phi,U)      Gauss linearUpwind grad(U);
    div(phi,k)      Gauss linearUpwind grad(k);
    div(phi,omega)  Gauss linearUpwind grad(omega);
    div((nuEff*dev(T(grad(U))))) Gauss linear;
}
laplacianSchemes
{
    default         Gauss linear corrected;
    laplacian(nuEff,U) Gauss linear corrected;
    laplacian((1|A(U)),p) Gauss linear corrected;
    laplacian(DkEff,k) Gauss linear corrected;
    laplacian(DomegaEff,omega) Gauss linear corrected;
}
interpolationSchemes
{
    default         linear;
}
snGradSchemes
{
    default         corrected;
}
fluxRequired
{
    default         no;
    p;
}
// ************************************************************************* //

Code:

/*--------------------------------*- 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;
    object      fvSolution;
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
solvers
{
    p
    {
        solver           GAMG;
        tolerance        1e-16;
        relTol           0.1;
        smoother         GaussSeidel;
        nPreSweeps       0;
        nPostSweeps      2;
        cacheAgglomeration on;
        agglomerator     faceAreaPair;
        nCellsInCoarsestLevel 10;
        mergeLevels      1;
    }
    pFinal
    {
        solver          GAMG;
        tolerance       1e-16;
        relTol          0;
        smoother        GaussSeidel;
        cacheAgglomeration true;
        nCellsInCoarsestLevel 10;
        agglomerator    faceAreaPair;
        mergeLevels     1;
    }

    U
    {
        solver           PBiCG;
        preconditioner   DILU;
        tolerance        1e-16;
        relTol           0.1;
    }
    k
    {
        solver           PBiCG;
        preconditioner  DILU;
        tolerance        1e-16;
        relTol           0.1;
 minIter   1;
    }
   omega
    {
        solver           PBiCG;
        preconditioner         DILU;
        tolerance        1e-16;
        relTol           0.1;
 minIter   1;
    }
}
SIMPLE
{
    nNonOrthogonalCorrectors 0;
    
    residualControl
    {
        p                1e-16;
        U                1e-16;
        k                1e-16;
        omega            1e-16;
    }
}
potentialFlow
{
    nNonOrthogonalCorrectors 10;
}
relaxationFactors
{
    fields
    {
        p               0.5;
    }
    equations
    {
        U               0.5;
        k               0.5;
        omega           0.5;
    }
}
cache
{
    grad(U);
}
// ************************************************************************* //

fredo490 · June 8, 2013, 13:54

Did you set the good value for the force function (in controldict) ? I mean, does your reference length and all other arguments are ok ?

May 30, 2013, 08:47	k-omega SST flow over a cylinder	#1
GP_FDU New Member Gao Pan Join Date: Mar 2013 Posts: 13 Rep Power: 13	I simulated flow over a cylinder with k-omega SST model. But the result was very wrong. Here is some configuration of my case. D=3.81cm, U=0.19m/s, fluid is set to be water. U and p B.C. is set as usual. B.C. of k is fixedValue=1e-12 at the cylinder. B.C. of nut is fixedValue=nutUSpaldingWallFunciton at the cylinder or fixedValue=1e-12. B.C. of omega is omegaWallFunction. y+is around 1 according to y+=0.172delta_y/DRe^0.9. With this configuration, there is vortex shedding and the frequency seems to be Ok, but the pressure at the cylinder is totally wrong, which leads to a wrong drag and lift force. Anybody know what is wrong with my case, please give me some advice.

June 4, 2013, 11:04		#10
fredo490 Senior Member HECKMANN Frédéric Join Date: Jul 2010 Posts: 249 Rep Power: 17	From the picture on post 7, it looks like you have a divergence... 1) Personally I would not set k = 0 at the inlet since it doesn't have any physical meaning. You can use the tool http://www.cfd-online.com/Tools/turbulence.php to estimate it. Roughly, you can set it as equal to 1e-6. 2) I'm not familiar with pimpleFoam (I work more on rhoPimpleFoam) but I think that you can set the inlet boundary of nut as "calculated" or "zeroGradient" since nut is derived from k and omega (If I remember well my classics). 3) maybe you have a problem with the time step. What is your maximum Courant Number (given at each iteration) ? It should be lower than 0.2 to avoid divergence or some sporadic crazy behavior. You can also try to use a second order scheme for ddt or maybe a cranknicholson. 4) try to change the turbulent div scheme to linearUpwind. It's often more stable. Rafael_Coelho likes this.

June 5, 2013, 05:15		#12
fredo490 Senior Member HECKMANN Frédéric Join Date: Jul 2010 Posts: 249 Rep Power: 17	CFD is most of the time a matter of experience. I've never found any book or class able to really prepare a student to deal with "advanced" cfd. The only thing you have to do is to keep working and trying In French we say: "c'est en forgeant qu'on devient forgeron", that you could roughly translate to "practice makes perfect". louisgag, JR22, OpenF and 2 others like this.

Similar Threads
Thread	Thread Starter	Forum	Replies	Last Post
k-omega SST simulation of turbulent flow around a circular cylinder	DanM	OpenFOAM Running, Solving & CFD	17	October 13, 2016 14:29
benchmark: flow over a circular cylinder	goodegg	Main CFD Forum	12	January 22, 2013 12:47
Flow induced vibration of a mobile cylinder	Hooman	Main CFD Forum	0	December 31, 2010 09:48
Simulation Flow Around cylinder 3D	Jwolf	CFX	19	November 25, 2009 15:21
Flow over a cylinder	Anna	Main CFD Forum	9	March 24, 2006 15:32

May 31, 2013, 00:10		#5
GP_FDU New Member Gao Pan Join Date: Mar 2013 Posts: 13 Rep Power: 13	Does anybody have some ideas?

June 2, 2013, 13:40		#6
fredo490 Senior Member HECKMANN Frédéric Join Date: Jul 2010 Posts: 249 Rep Power: 17	1) you can try to set 1 or 2 non orthogonal correction (nNonOrthogonalCorrectors ) depending of your mesh. 2) Is it normal that your pcorr solver has only a tolerance of 0.02 ?! 3) Your pressure relaxation number looks a bit low for an unsteady simulation. Ideally it has to be close to 1 to really catch the unsteadiness. 0.8 would be much better I think. 4) Where did you get your k and omega value for the inlet ? Those two values are really important to predict the boundary layer separation. 5) use the "yPlusRAS -latestTime" to really compute the yplus of your case ...

June 4, 2013, 08:48		#8
fredo490 Senior Member HECKMANN Frédéric Join Date: Jul 2010 Posts: 249 Rep Power: 17	how is your mesh ? Can you give a screenshot ? Also, you can run a "checkMesh" utility and give the report

June 5, 2013, 05:10		#11
GP_FDU New Member Gao Pan Join Date: Mar 2013 Posts: 13 Rep Power: 13	Hi Heckmann, Thanks very much for your patience and kindness. Following your suggestions, I finally succeed to get a reasonable result. I changed the maxCo to 0.2 the pressure relaxation number to 0.8, and the time scheme to CrankNicolson 0.5. I am not good at CFD. For many parameters in the case files, sometimes I don't know how to set them -- I don't know what value is appropriate. I am trying to learn. Could you give me some advise?

June 5, 2013, 11:08		#17
OpenF New Member Anselme Bonnet Join Date: Apr 2013 Location: France Posts: 14 Rep Power: 13	thank you! I did what you said, I put the file "mirrorMeshDict" in the system directory, but before and after the computation, I can't generate the mesh, I still have a half-cylinder .. Maybe I need to put something in the controlDict or somewhere?

June 5, 2013, 11:35		#18
GP_FDU New Member Gao Pan Join Date: Mar 2013 Posts: 13 Rep Power: 13	which version do you use? if of1.6ext, please check the 0 directory.

June 8, 2013, 13:54		#20
fredo490 Senior Member HECKMANN Frédéric Join Date: Jul 2010 Posts: 249 Rep Power: 17	Did you set the good value for the force function (in controldict) ? I mean, does your reference length and all other arguments are ok ?