[MrFourier]

Hi everyone,

I'm new to OpenFOAM, and in preparation for a more complex simulation I'm currently working on evaluating some turbulence models on a simple mesh.

I have a 2D-flow between two smooth no-slip walls. Boundary conditions on inlet and outlet are periodic, except for a pressure jump. I'm using simpleFoam and tried out both k-omega SST and Spalart-Allmaras turbulence models, with all default coefficients, and compared the results. This plot is for a mesh with a wall y+ of 0.5 (I get the same results using different mesh resolutions).

comparison_turb_small.png

I'm not sure what I'm supposed to get, but I figure the k-omega SST results are quite good, whereas the Spalart-Allmaras results look unusable. The Spalart results can't be the way they're supposed to be, can they? I must be doing something wrong!

Except for the nuTilda file, the settings are exactly the same, so here's my nuTilda file:

dimensions      [ 0 2 -1 0 0 0 0 ];
internalField   uniform 1.e-5;//0.05;
boundaryField
{
    inlet
    {
        type            cyclicAMI;
    }
    outlet
    {
        type            cyclicAMI;
    }
    wall
    {
        type            fixedValue;
        value           uniform 0;
    }
    "sym(A|B)"
    {
        type            symmetryPlane;
    }
}

Does anyone have an idea how to fix this?
Are there any other configuration files that I should post here?

Thank you!

EDIT: Btw, this is OF 2.1

[vwibaut]

Hi MrFourier,

Do you find something? I have a similar problem

Thanks for your help

[kiddmax]

Hi Valentin,

Did you find something?

Best
Ye

[kovamaniac]

Hi.

As I am having trouble in making the SST run, could you please post the rest of your files for the setup? My solution does not converge and I cannot understand why.

[RodriguezFatz]

Anastasios what is the problem? Post your relevant settings.

[kovamaniac]

Originally Posted by RodriguezFatz

Anastasios what is the problem? Post your relevant settings.

Hi Philipp.

My settings are the following:
Initial conditions for U

dimensions [0 1 -1 0 0 0 0];
internalField uniform ( 68.4 0. 0. );
　
boundaryField
{
farfield2
{
type empty;
}
farfield1
{
type empty;
}
inflow
{
type fixedValue;
value uniform (68.4 0. 0.);
}
outflow
{
type zeroGradient;
}
symmetry
{
type symmetryPlane;
}
wall
{
type fixedValue;
value uniform (0. 0. 0.);
}
top
{
type slip;
}
　
}

Initial conditions for P

dimensions [0 2 -2 0 0 0 0];
internalField uniform 0.;
　
boundaryField
{
farfield2
{
type empty;
}
farfield1
{
type empty;
}
inflow
{
type /*freestreamPressure;*/zeroGradient;
}
outflow
{
type zeroGradient;
}
symmetry
{
type symmetryPlane;
}
wall
{
type zeroGradient;
}
top
{
type slip;
}
　
}

Initial conditions for k

dimensions [0 2 -2 0 0 0 0];
internalField uniform 0.001052676;
boundaryField
{
farfield2
{
type empty;
}
farfield1
{
type empty;
}
inflow
{
type fixedValue;
value uniform 0.001052676;
}
outflow
{
type zeroGradient;
}
symmetry
{
type symmetryPlane;
}
wall
{
type fixedValue;
value uniform 0;
}
top
{
type slip;
}
}

Initial conditions for omega

 
dimensions [0 0 -1 0 0 0 0];
internalField uniform 8550;
boundaryField
{
farfield2
{
type empty;
}
farfield1
{
type empty;
}
inflow
{
type fixedValue;
value uniform 8550;
}
outflow
{
type zeroGradient;
}
symmetry
{
type symmetryPlane;
}
wall
{
//type fixedValue;
type omegaWallFunction;
//value uniform 1.581e+08;
value uniform 1;
}
top
{
type slip;
}
}

Initial conditions for nut

dimensions [0 2 -1 0 0 0 0];
internalField uniform 0;
boundaryField
{
farfield2
{
type empty;
}
farfield1
{
type empty;
}
inflow
{
type calculated;
value uniform 0;
}
outflow
{
type calculated;
value uniform 0;
}
symmetry
{
type symmetryPlane;
}
wall
{
type nutUSpaldingWallFunction;
value uniform 0;
}
top
{
type slip;
}
}

The system files are:

solvers
{
p GAMG
{
preconditioner DIC;
mergeLevels 1;
smoother GaussSeidel;
agglomerator faceAreaPair;
nCellsInCoarsestLevel 100;
cacheAgglomeration on; 
tolerance 1e-07;
relTol 0.01;
};
U GAMG
{
preconditioner DILU;
mergeLevels 1;
smoother GaussSeidel;
agglomerator faceAreaPair;
nCellsInCoarsestLevel 100;
cacheAgglomeration on; 
tolerance 1e-07;
relTol 0.1;
};
k GAMG
{
preconditioner DILU;
mergeLevels 1; 
smoother GaussSeidel;
agglomerator faceAreaPair;
nCellsInCoarsestLevel 100;
cacheAgglomeration on; 
tolerance 1e-07;
relTol 0.1;
};
omega GAMG
{
preconditioner DILU;
mergeLevels 1; 
smoother GaussSeidel;
agglomerator faceAreaPair;
nCellsInCoarsestLevel 100;
cacheAgglomeration on; 
tolerance 1e-07;
relTol 0.1;
};
}
SIMPLE
{
nNonOrthogonalCorrectors 0;
pRefCell 0;
pRefValue 0;
}
relaxationFactors
{
p 0.3; //0.3 original!!
U 0.7;
k 0.7;
omega 0.7;
}

ddtSchemes
{
default steadyState;
}
gradSchemes
{
default Gauss linear;
grad(U) Gauss linear;
grad(p) Gauss linear;
}
divSchemes
{
default none;
div(phi,U) bounded Gauss linearUpwindV grad(U);
div(phi,k) bounded Gauss upwind;
div(phi,omega) bounded Gauss upwind;
div((nuEff*dev(T(grad(U))))) Gauss linear;
}
laplacianSchemes
{
default Gauss linear corrected;
}
interpolationSchemes
{
default linear;
}
snGradSchemes
{
default corrected;
}
fluxRequired
{
default no;
p;
}

My meshes come from NASA website and the y+ value is below 1 everywhere, that's why I do not use wall functions. I am not sure about the top boundary condition for which I put "slip".

With the above settings the solution converges but the values are strange; I have values in the flow field greater than the input velocity.

[RodriguezFatz]

Originally Posted by kovamaniac

With the above settings the solution converges but the values are strange; I have values in the flow field greater than the input velocity.

Ok, so you have a flat plate on the bottom, inlet on the left and outlet on the right?
At your inlet you put in some constant velocity profile. Now, at the bottom wall a boundary layer develops with increasing size along the flow direction. In this boundary layer the velocity is low. But the mass flow through each vertical line (integrated) must be constant. Thus, if you have a lower velocity at the bottom you must get a higher value at some other places.

[kovamaniac]

Originally Posted by RodriguezFatz

Ok, so you have a flat plate on the bottom, inlet on the left and outlet on the right?
At your inlet you put in some constant velocity profile. Now, at the bottom wall a boundary layer develops with increasing size along the flow direction. In this boundary layer the velocity is low. But the mass flow through each vertical line (integrated) must be constant. Thus, if you have a lower velocity at the bottom you must get a higher value at some other places.

You are correct. Mass conservation. Thank you very much.

[RodriguezFatz]

Does it work if you set inlet and outlet to zerogradient and the top wall to some fixed velocity value?

[kovamaniac]

Originally Posted by RodriguezFatz

Does it work if you set inlet and outlet to zerogradient and the top wall to some fixed velocity value?

This reminds me the way I was doing it with fluent where the top was also an inflow bc. Let me try this and I will let you know.

[kovamaniac]

Originally Posted by RodriguezFatz

Does it work if you set inlet and outlet to zerogradient and the top wall to some fixed velocity value?

Exactly the same situation.

My main problem seems to be the pressure convergence. In both the very fine mesh (~2100000 cells) and the coarse presents a very oscillatory behavior. In the fine mesh after some iterations it seems to somehow blow up.
Please find attached the residuals for both cases. I tried to reduced the order of the Laplacian but the same.

[kovamaniac]

Originally Posted by RodriguezFatz

Does it work if you set inlet and outlet to zerogradient and the top wall to some fixed velocity value?

Also I am not sure about the accurace as the results differ significantly form NASA's. I do not so much care about the values but from the fact that my results do not follow well the law of the wall.

[RodriguezFatz]

I would try to run it with gradient limiters. Go for

default cellMDLimited Gauss linear 1.0;

and see what happens.
In my experience, the convergence of the k-omega-sst model strongly depends on the gradient limiters.

[kovamaniac]

Originally Posted by RodriguezFatz

I would try to run it with gradient limiters. Go for

default cellMDLimited Gauss linear 1.0;

and see what happens.
In my experience, the convergence of the k-omega-sst model strongly depends on the gradient limiters.

The same behavior and convergence in higher values.

In any case, I solved the issue with the inaccuracy.

FYI, the wallShearStress DOES NOT WORK FOR LOW-RE GRIDS.

Calculating u_tau with the wallGradU utility resulted in correct results.

[RodriguezFatz]

1) Did you deactivate ("//") the gradients for U and p and set all the gradients as I suggested?
2) Are you sure, you use the correct syntax for relaxation factors? I use them like
http://www.openfoam.org/docs/user/fvSolution.php

[kovamaniac]

Originally Posted by RodriguezFatz

1) Did you deactivate ("//") the gradients for U and p and set all the gradients as I suggested?
2) Are you sure, you use the correct syntax for relaxation factors? I use them like
http://www.openfoam.org/docs/user/fvSolution.php

I use the setup for the fvSolution file suggested for the simpleFoam at the link you gave; however this is for a laminar case -this is what it says, despite having k and epsilon defined.

If I use the gradscheme you propose for default, the result is convergence in higher values again with this small oscillation. I cannot understand why. At the moment I compromise with this as the results at least for the flat plate flow are satisfactory. To be honest I can understand the difference of the solver schemes, but I do not have the intuition. It is too mathematical for me and I have not been taught such advanced linear algebra topics.

I will use some more complicated models and let you know.

[RodriguezFatz]

Can you upload the complete setup somewhere?

[kovamaniac]

Originally Posted by RodriguezFatz

Can you upload the complete setup somewhere?

Finally,

I followed your advice and changed the gradient schemes. That behaved correctly.

However, the most important issue for the accuracy results was to understand that the wallShearStress utility outputs in m2/s2, which is Pa/rho and not Pa. That was a big mistake which took me a lot to understand -to be honest I hadn't paid attention.

Thank you very much for your help anyway.