[MIchael_P]

Hi,
I am starting to use Openfoam and want to validate some physical model test results.

Water flows with 0.28m³/s through an upstream flume around an obstacle out through a pipe. At the end of the pipe there should be a drop towards a hydro turbine but is not modeled here. The goal is to derive the velocity distribution at the cross section.

https://www.dropbox.com/s/uqq6168osg...sults.jpg?dl=0




The shaded results are from a CFD Software COMET "kepsilon model" from a dissertation while the white lines are the results of a physical model test presented in another PhD thesis. The results are shown as velocities divided by the mean velocity of the section.


My results show the peak velocities on the right side of the cross section, being completely wrong.

I have tried to lengthen the approach channel with no effect.
I have changed the initial conditions without effect.
I have changed the background mesh to 2*2*2. There was also only a small change. Peaks moved to the lower right of the cross sections.
I do not what I am doing wrong.



The case files can be found her under this link.

https://www.dropbox.com/s/254zfs9u0o4cyee/CASE.zip?dl=0

I am using ESI OpenFoamv1712 in a Virtual Box

I prepared the model in Bricscad. I created the surface mesh of the model as well the background mesh (5*5*5cm blocks) in Salome Meca.
For meshing I used SnappyHexMesh.
CheckMesh tells me that the mesh is fine.
Create mesh for time = 0

Time = 0

Mesh stats
points: 547336
faces: 1327985
internal faces: 1164729
cells: 392144
faces per cell: 6.356629
boundary patches: 4
point zones: 0
face zones: 0
cell zones: 0

Overall number of cells of each type:
hexahedra: 318486
prisms: 18311
wedges: 0
pyramids: 0
tet wedges: 8
tetrahedra: 0
polyhedra: 55339
Breakdown of polyhedra by number of faces:
faces number of cells
4 1252
5 1183
6 6458
7 52
8 47
9 39071
11 12
12 6979
15 283
18 2

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 1596 1964 ok (non-closed singly connected)
Outlet 588 729 ok (non-closed singly connected)
Wall 156669 160189 ok (non-closed singly connected)
Air_Surface 4403 5215 ok (non-closed singly connected)

Checking faceZone topology for multiply connected surfaces...
No faceZones found.

Checking basic cellZone addressing...
No cellZones found.

Checking geometry...
Overall domain bounding box (36.35 39.6 10.4) (43.35 41.25 11.4)
Mesh has 3 geometric (non-empty/wedge) directions (1 1 1)
Mesh has 3 solution (non-empty) directions (1 1 1)
Boundary openness (-9.481594e-17 -3.323153e-14 7.796728e-14) OK.
Max cell openness = 3.511362e-16 OK.
Max aspect ratio = 6.143295 OK.
Minimum face area = 1.04191e-05. Maximum face area = 0.002756876. Face area
magnitudes OK.
Min volume = 2.065489e-08. Max volume = 0.0001394037. Total volume = 7.5558
05. Cell volumes OK.
Mesh non-orthogonality Max: 53.08211 average: 11.39594
Non-orthogonality check OK.
Face pyramids OK.
Max skewness = 2.758769 OK.
Coupled point location match (average 0) OK.

Mesh OK.

ControlDict

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

application simpleFoam;

startFrom startTime;

startTime 0;

stopAt endTime;

endTime 2000;

deltaT 1;

writeControl timeStep;

writeInterval 10;

purgeWrite 0;

writeFormat ascii;

writePrecision 7;

writeCompression off;

timeFormat general;

timePrecision 6;

runTimeModifiable true;

functions
{
// #includeFunc streamlines
}

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

fvSchemes

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

solvers
{
p
{
solver GAMG;
tolerance 1e-06;
relTol 1e-06;
smoother GaussSeidel;
}

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

SIMPLE
{
nNonOrthogonalCorrectors 0;
consistent yes;
residualControl
{
p 1e-2;
U 1e-3;
"(k|epsilon|omega|f|v2)" 1e-3;
}
}

relaxationFactors
{
fields
{
p 0.3; // from Nagy tutorial
}
equations
{
U 0.7; // 0.9 is more stable but 0.95 more convergent
k 0.7; // 0.9 is more stable but 0.95 more convergent
epsilon 0.7; // 0.9 is more stable but 0.95 more convergent
//not used here R 0.7; // 0.9 is more stable but 0.95 more convergent
// not used here nuTilda 0.7; // 0.9 is more stable but 0.95 more convergent
}
}


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

TransportProperties

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

transportModel Newtonian;

nu [0 2 -1 0 0 0 0] 1e-06;


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

turbulenceProperties

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

simulationType RAS;

RAS
{
// Tested with kEpsilon, realizableKE, kOmega, kOmegaSST, v2f,
// ShihQuadraticKE, LienCubicKE.
RASModel kEpsilon;

turbulence on;

printCoeffs on;
}


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

U at boundary

/*--------------------------------*- C++ -*----------------------------------*\
| ========= | |
| \\ / F ield | OpenFOAM: The Open Source CFD Toolbox |
| \\ / O peration | Version: 4.x |
| \\ / A nd | Web: www.OpenFOAM.org |
| \\/ M anipulation | |
\*---------------------------------------------------------------------------*/
FoamFile
{
version 2.0;
format ascii;
class volVectorField;
object U;
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
#include "include/initialConditions"
dimensions [0 1 -1 0 0 0 0];

internalField uniform (0 0 0);

boundaryField
{
Inlet
{
type fixedValue;
value uniform ($inletFlowRate 0 0);
}

Outlet
{
type zeroGradient;
}

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

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

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

p at the boundaries

/*--------------------------------*- C++ -*----------------------------------*\
| ========= | |
| \\ / F ield | OpenFOAM: The Open Source CFD Toolbox |
| \\ / O peration | Version: 4.x |
| \\ / A nd | Web: www.OpenFOAM.org |
| \\/ M anipulation | |
\*---------------------------------------------------------------------------*/
FoamFile
{
version 2.0;
format ascii;
class volScalarField;
object p;
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
#include "include/initialConditions"
dimensions [0 2 -2 0 0 0 0];

internalField uniform $pressure;

boundaryField
{
Inlet
{
type zeroGradient;
}

Outlet
{
type fixedValue;
value uniform $pressure;
}

Wall
{
type zeroGradient;
}
Air_Surface
{
type zeroGradient;
}
}

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

epsilon at the boundaries

/*--------------------------------*- C++ -*----------------------------------*\
| ========= | |
| \\ / F ield | OpenFOAM: The Open Source CFD Toolbox |
| \\ / O peration | Version: 4.x |
| \\ / A nd | Web: www.OpenFOAM.org |
| \\/ M anipulation | |
\*---------------------------------------------------------------------------*/
FoamFile
{
version 2.0;
format ascii;
class volScalarField;
location "0";
object epsilon;
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
#include "include/initialConditions"
dimensions [0 2 -3 0 0 0 0];

internalField uniform $turbulentEpsilon;

boundaryField
{
Inlet
{
type fixedValue;
value uniform $turbulentEpsilon;
}
Outlet
{
type zeroGradient;
}
Wall
{
type epsilonWallFunction;
value uniform $turbulentEpsilon;
}
Air_Surface
{
type epsilonWallFunction;
value uniform $turbulentEpsilon;
}

}


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

k at the boundaries

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

internalField uniform $turbulentKE;

boundaryField
{
Inlet
{
type fixedValue;
value uniform $turbulentKE;
}
Outlet
{
type zeroGradient;
}
Wall
{
type kqRWallFunction;
value uniform $turbulentKE;
}
Air_Surface
{
type kqRWallFunction;
value uniform $turbulentKE
;
}
}


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

nut at boundaries

/*--------------------------------*- C++ -*----------------------------------*\
| ========= | |
| \\ / F ield | OpenFOAM: The Open Source CFD Toolbox |
| \\ / O peration | Version: 4.x |
| \\ / A nd | Web: www.OpenFOAM.org |
| \\/ M anipulation | |
\*---------------------------------------------------------------------------*/
FoamFile
{
version 2.0;
format ascii;
class volScalarField;
location "0";
object nut;
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
#include "include/initialConditions"
dimensions [0 2 -1 0 0 0 0];

internalField uniform $initial_nut;

boundaryField
{
Inlet
{
type calculated;
value uniform $initial_nut;
}
Outlet
{
type calculated;
value uniform 0;
}
Wall
{
type nutkWallFunction;
value $internalField;
}

Air_Surface
{
type nutkWallFunction;
value $internalField;
}
}


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

I would really appreciate if somebody could direct me towards the right direction.
thanks,
Michael

[harsha_kulkarni]

Try komega model or RSM model
Instead of kepsilon

[Robin.Kamenicky]

Hi Michael,

I have had a look at the case and if you are sure that the solution is wrong. Can you try a few things?

• Decrease relTol for U in fvSolution.Solver from 0.1 to 0.0. 0.1 is too big in my opinion. Eventually you can decrease also relTol for k and epsilon.

• Decrease the pressure and U residuals for SIMPLE algorithm in fvSolution.SIMPLE. Do you think that the solution converged?

• Change fvSolution.SIMPLE nNonOrthogonalCorrectors to 1

• In my opinion, the other problem can be your outlet BC where you use fixedPressure and zeroGradient for velocity. This is however not true since your flow is not fully developed. You might give a try to some other BC for U e.g. fluxCorrectedVelocity

Please let me know whether it works.

Cheers,
Robin

[MIchael_P]

Dear Robin,
thank you very much. I followed your advice and it looks much better now. Not yet fully matching the physical model but I am on the right track now. I am aware now that I myself have to set the rules for convergence and not, as I did, follow blindly the output ("solution has converged") of simplefoam.



However, I still have to do more analysis whether the solution has fully converged. After 2000 time steps the solution has not yet reached the residuals of 1e-6 that I had set for p and U.



Here would be my follow on questions. What are the best checks to see whether a solution has really converged? Are residuals of 1e-6 a good limit? What else would you check?



While I was waiting for responses from the forum I also started from the beginning and tried to check every step and found that my mesh, despite passing checkMesh might not be the best. I will post this separately. Here or in the meshing forum.

[MIchael_P]

Dear Harsha,
I will try your advice as soon as I can and let you know the results.
thanks,
Michael

[Robin.Kamenicky]

Originally Posted by MIchael_P

Here would be my follow on questions. What are the best checks to see whether a solution has really converged? Are residuals of 1e-6 a good limit? What else would you check?

Hi Michael,

Well, in general we want to evaluate particular flows based on the converged physical values. Hence, the most important is to be sure that the fields such as U, p, k, etc. converged. Indeed, it is very usefull to check residuals, but in the end you have the most interest in the fields. Those you can check visualy in Paraview, whether they change in subsequent iterations.

Surely, there are other ways to check it. In your case you can for instance plot U profile in the Trompete. Very often is also checked pressure drop across the computational domain. Sometime mass/volumetric flow at the inlet and outlet as they should match. You can also plot maximum velocity and pressure and see whether those changes with subsequent iterations or use probes. The choice can be, however, case dependent. You can use paraview or/and OpenFOAM Post-processing command line interface (CLI).

There is not anything as the right value of residuals which would tell you that the case converged. This is highly case dependent, sometime a case will converge with 1-e3, sometime 1-e6. Also for various fields you will be able to reach different residuals.

Regarding the mesh, surely mesh is crucial. OpenFOAM seems to be quite picky about mesh. CheckMesh values can be eventually setup by your preferences.

Hope this helps,
Let me know how your computation goes.
Robin