[sebastian]

Hi everybody!

Trying to compare a Fluent calculation with OF.

With Fluent, everything worked okay without any problems.

Unfortunatelly the OF case is not running so far.

I have generated a first starting solution with potentialFoam. Then I want to switch to simpleFoam, to run the case with the k-epsilon model.
But after two or three iterations, the solution is blowing up.
Its working with switching off the turbulence model (at least its still running with low relaxation factors and now I have about 200 iterations). But turning on the k-e model again, the continuity is exploding immediatelly, as well as k and epsilon and then also momentum.

The mesh seems to be okay. Its a tet mesh and the cell quality is decent.

The discretisation scheme is upwind for div and laplacian. Grad scheme is linear and time is steady state.
I have tried several initalizations for k and epsilon but the result is always the same..


I know there can be a lot of reasons for the prob.. but I would be happy about any comment.


Thanks in advance!

Sebastian

[kiran]

can u paste the last error message lines from log file. so we could see where exactly the problem is.

[sebastian]

Hey. Sure, to copy the error message here might help



Time = 4

GAMG: Solving for Ux, Initial residual = 0.0306223, Final residual = 1.9577e-07, No Iterations 2
GAMG: Solving for Uy, Initial residual = 0.19433, Final residual = 1.27689e-06, No Iterations 2
GAMG: Solving for Uz, Initial residual = 0.152602, Final residual = 8.53599e-07, No Iterations 2
GAMG: Solving for p, Initial residual = 0.000902567, Final residual = 9.04656e-07, No Iterations 8
GAMG: Solving for p, Initial residual = 2.08343e-06, Final residual = 3.03159e-07, No Iterations 1
GAMG: Solving for p, Initial residual = 8.41737e-07, Final residual = 8.41737e-07, No Iterations 0
GAMG: Solving for p, Initial residual = 8.41737e-07, Final residual = 8.41737e-07, No Iterations 0
time step continuity errors : sum local = 5.95441, global = 0.0633777, cumulative = 0.0676273
GAMG: Solving for epsilon, Initial residual = 0.852492, Final residual = 2.86155e-07, No Iterations 2
bounding epsilon, min: -3.11834e+06 max: 1.85789e+12 average: 4.50946e+07
GAMG: Solving for k, Initial residual = 0.936029, Final residual = 1.83872e-07, No Iterations 2
bounding k, min: -1618.5 max: 7.23768e+09 average: 186799
ExecutionTime = 532.15 s ClockTime = 532 s

Time = 5

GAMG: Solving for Ux, Initial residual = 0.390574, Final residual = 2.22292e-06, No Iterations 2
GAMG: Solving for Uy, Initial residual = 0.337754, Final residual = 9.8492e-07, No Iterations 2
GAMG: Solving for Uz, Initial residual = 0.853074, Final residual = 2.59431e-06, No Iterations 2
GAMG: Solving for p, Initial residual = 0.47455, Final residual = 9.87328e-07, No Iterations 100
GAMG: Solving for p, Initial residual = 3.54943e-06, Final residual = 5.0544e-07, No Iterations 2
GAMG: Solving for p, Initial residual = 2.20942e-06, Final residual = 3.87525e-07, No Iterations 1
GAMG: Solving for p, Initial residual = 9.97362e-07, Final residual = 9.97362e-07, No Iterations 0
time step continuity errors : sum local = 1.79565e+10, global = -1.20617e+06, cumulative = -1.20617e+06
GAMG: Solving for epsilon, Initial residual = 1, Final residual = 2.89978e-07, No Iterations 3
bounding epsilon, min: -9.32203e+21 max: 3.1884e+31 average: 2.41736e+25
GAMG: Solving for k, Initial residual = 1, Final residual = 7.01714e-17, No Iterations 1
ExecutionTime = 824.21 s ClockTime = 824 s

Time = 6

GAMG: Solving for Ux, Initial residual = 0.0527437, Final residual = 2.19226e-09, No Iterations 2
GAMG: Solving for Uy, Initial residual = 0.0422177, Final residual = 2.22655e-08, No Iterations 2
GAMG: Solving for Uz, Initial residual = 0.0738953, Final residual = 3.43556e-08, No Iterations 2
GAMG: Solving for p, Initial residual = 1, Final residual = 9.71148e-07, No Iterations 405
GAMG: Solving for p, Initial residual = 1.63209e-07, Final residual = 1.63209e-07, No Iterations 0
GAMG: Solving for p, Initial residual = 1.63209e-07, Final residual = 1.63209e-07, No Iterations 0
GAMG: Solving for p, Initial residual = 1.63209e-07, Final residual = 1.63209e-07, No Iterations 0
time step continuity errors : sum local = 3.9722e+48, global = 1.65518e+31, cumulative = 1.65518e+31
[0] #0 Foam::error:rintStack(Foam::Ostream&) in "/home/sebastian/OpenFOAM/OpenFOAM-1.6/lib/linuxGccDPOpt/libOpenFOAM.so"
[0] #1 Foam::sigFpe::sigFpeHandler(int) in "/home/sebastian/OpenFOAM/OpenFOAM-1.6/lib/linuxGccDPOpt/libOpenFOAM.so"
[0] #2 Uninterpreted:
[0] #3 Foam::GaussSeidelSmoother::smooth(Foam::Field<doub le>&, Foam::Field<double> const&, unsigned char, int) const in "/home/sebastian/OpenFOAM/OpenFOAM-1.6/lib/linuxGccDPOpt/libOpenFOAM.so"
[0] #4 Foam::GAMGSolver::Vcycle(Foam::PtrList<Foam::lduMa trix::smoother> const&, Foam::Field<double>&, Foam::Field<double> const&, Foam::Field<double>&, Foam::Field<double>&, Foam::Field<double>&, Foam::PtrList<Foam::Field<double> >&, Foam::PtrList<Foam::Field<double> >&, unsigned char) const in "/home/sebastian/OpenFOAM/OpenFOAM-1.6/lib/linuxGccDPOpt/libOpenFOAM.so"
[0] #5 Foam::GAMGSolver::solve(Foam::Field<double>&, Foam::Field<double> const&, unsigned char) const in "/home/sebastian/OpenFOAM/OpenFOAM-1.6/lib/linuxGccDPOpt/libOpenFOAM.so"
[0] #6 Foam::fvMatrix<double>::solve(Foam::dictionary const&) in "/home/sebastian/OpenFOAM/OpenFOAM-1.6/lib/linuxGccDPOpt/libfiniteVolume.so"
[0] #7 Foam::lduMatrix::solverPerformance Foam::solve<double>(Foam::tmp<Foam::fvMatrix<doubl e> > const&) in "/home/sebastian/OpenFOAM/OpenFOAM-1.6/lib/linuxGccDPOpt/libincompressibleRASModels.so"
[0] #8 Foam::incompressible::RASModels::kEpsilon::correct () in "/home/sebastian/OpenFOAM/OpenFOAM-1.6/lib/linuxGccDPOpt/libincompressibleRASModels.so"
[0] #9 main in "/home/sebastian/OpenFOAM/OpenFOAM-1.6/applications/bin/linuxGccDPOpt/simpleFoam"
[0] #10 __libc_start_main in "/lib/tls/libc.so.6"
[0] #11 _start at ../sysdeps/i386/elf/start.S:122
[node07:10323] *** Process received signal ***
[node07:10323] Signal: Floating point exception (8)
[node07:10323] Signal code: (-6)
[node07:10323] Failing at address: 0x2853
[node07:10323] [ 0] [0xffffe440]
[node07:10323] [ 1] /home/sebastian/OpenFOAM/OpenFOAM-1.6/lib/linuxGccDPOpt/libOpenFOAM.so(_ZN4Foam6sigFpe13sigFpeHandlerEi+0x 61) [0x411873a1]
[node07:10323] [ 2] [0xffffe420]
[node07:10323] [ 3] /home/sebastian/OpenFOAM/OpenFOAM-1.6/lib/linuxGccDPOpt/libOpenFOAM.so(_ZNK4Foam19GaussSeidelSmoother6smoo thERNS_5FieldIdEERKS2_hi+0x52) [0x40fe0c52]
[node07:10323] [ 4] /home/sebastian/OpenFOAM/OpenFOAM-1.6/lib/linuxGccDPOpt/libOpenFOAM.so(_ZNK4Foam10GAMGSolver6VcycleERKNS_7 PtrListINS_9lduMatrix8smootherEEERNS_5FieldIdEERKS 8_S9_S9_S9_RNS1_IS8_EESD_h+0xeb7) [0x40ff1577]
[node07:10323] [ 5] /home/sebastian/OpenFOAM/OpenFOAM-1.6/lib/linuxGccDPOpt/libOpenFOAM.so(_ZNK4Foam10GAMGSolver5solveERNS_5Fi eldIdEERKS2_h+0x3b9) [0x40ff28a9]
[node07:10323] [ 6] /home/sebastian/OpenFOAM/OpenFOAM-1.6/lib/linuxGccDPOpt/libfiniteVolume.so(_ZN4Foam8fvMatrixIdE5solveERKNS _10dictionaryE+0x178) [0x407b3bd8]
[node07:10323] [ 7] /home/sebastian/OpenFOAM/OpenFOAM-1.6/lib/linuxGccDPOpt/libincompressibleRASModels.so(_ZN4Foam5solveIdEENS _9lduMatrix17solverPerformanceERKNS_3tmpINS_8fvMat rixIT_EEEE+0x58) [0x400af6d8]
[node07:10323] [ 8] /home/sebastian/OpenFOAM/OpenFOAM-1.6/lib/linuxGccDPOpt/libincompressibleRASModels.so(_ZN4Foam14incompress ible9RASModels8kEpsilon7correctEv+0xdc0) [0x400a9c20]
[node07:10323] [ 9] /home/sebastian/OpenFOAM/OpenFOAM-1.6/applications/bin/linuxGccDPOpt/simpleFoam [0x805b950]
[node07:10323] [10] /lib/tls/libc.so.6(__libc_start_main+0xe0) [0x4136b500]
[node07:10323] [11] /home/sebastian/OpenFOAM/OpenFOAM-1.6/applications/bin/linuxGccDPOpt/simpleFoam [0x8059c81]
[node07:10323] *** End of error message ***


Sebastian

[bastil]

Originally Posted by sebastian

The mesh seems to be okay. Its a tet mesh and the cell quality is decent.

• What is decent? What does checkMesh say?
• In general from you log: Using GAMG for velocity and turbulence is a waste of time - use DILUPBiCG.
• Turn non-orthogonal correctors to maximum 1 (from my experience 0 runs or it will never run). another waste of time
• What version of OpenFOAM do you use?
• Does the potentialFoam solution look "reasonable" - at least as reasonable as it can?
• Can you share the case?
• Is it a Hi-Re-Mesh which you run with HI-Re k-epsilon model?

Regards Bastian

[Michel_HB]

Hi Sebastian,

i had the same problem.
make sure that all files in the 0-directory have the same patches and boundary conditions...

Michel

[sebastian]

Hi!

- checkMesh says everything is okay. I also checked the mesh quality in ICEM and the worst element has a quality of 0.18 (if that says anything to you)
- I use OpenFoam 1.6
- Well youre right, the potential foam solution is very unphysical. But starting without using the potentialFoam generated solution, the problem is the same.
- Yes, I expect high Reynolds numbers in my domain. What do you mean by Hi-Re-k-epsilon model?

Thanks for the speeding up tips!

I also checked the boundary conditions again. The patches should be okay. For a first try I am using fixed values for k and epsilon at the inlets and outlets of the domain (estimated with formulas), and zero gradients at the walls.

[bastil]

Originally Posted by sebastian

- Yes, I expect high Reynolds numbers in my domain. What do you mean by Hi-Re-k-epsilon model?

There are two kinds of k-epsilon turbulence models that are more or less expensive an accurate: Hi- and Low-Re-k-eps models. None of them have to do with the question if you expect high reynolds numbers in your model. Furthermore, Low-Re-models are an extension of the models that are also valid at regions with Low-Re numbers (typically present at walls). The requirements for your mesh are different for Hi- and Low reynolds models (typically estimated by an that called yplus-number). I was asking this since I wanted to know what kind of turbulence model you exactly use in FLUENT and OpenFOAM and what kind of mesh you have?
If you can send me the model I can take a look at it.

Regards Bastian

[sebastian]

Thanks Bastian! But I can not send you the case.

Maybe we can solve the prob here in the forum.

The mesh I try to calculate is an coarser version of a finer mesh I later want to work with. In the fine mesh, y+ is around 30. In this mesh its about 100. In Fluent I used standard k-e model with standard wall functions.


Sebastian

[maddalena]

Hello... Mmm... I do not now too much of your case but... are you sure of this?

Originally Posted by sebastian

For a first try I am using fixed values for k and epsilon at the inlets and outlets of the domain (estimated with formulas), and zero gradients at the walls.

Shouldn't be k and epsilon fixedValue at the inlets and zeroGradient at the outlets? A solution that is blowing up after only some iterations sounds to me or a bc problems or a fvSchemes that needs tuning...
Hope this help.

mad

[bastil]

Originally Posted by maddalena

Shouldn't be k and epsilon fixedValue at the inlets and zeroGradient at the outlets? A solution that is blowing up after only some iterations sounds to me or a bc problems or a fvSchemes that needs tuning...

I agree....

[Chris Lucas]

Hi,

you said that you used a zeroGradient BC at the wall for k and epsilon. You should use the wall function BC's "kqRWallFunction" and "epsilonWallFunction".

I would also recommend to use "turbulentIntensityKineticEnergyInlet" and "turbulentMixingLengthDissipationRateInlet" at the Inlet

Regards,
Christian

[sebastian]

Hi!

Thanks a lot! Really seems the fixed outlet conditions have been the prob. At least now its running longer and its still looking okay

Also thanks for the suggested bc at the inlets!


Sebastian

[sebastian]

Hi again!

Sorry, I am asking again. But this is the only way for me to get some help..
I still have the problems with my calculation, that its blowing up.

My case: its a nozzle flow with a core stream and fan stream entering in a farfield.

I think its blowing up, because of bounding epsilon (I know this has been discussed a lot, but it didnt really help me so far)

Please could anybody have a look at my case, I tried to present it as detailed as possible. Please ask for additional informations.

After about 600 iterations, the simulation stops:


Time = 606

DILUPBiCG: Solving for Ux, Initial residual = 0.999913, Final residual = 7.55902e-08, No Iterations 2
DILUPBiCG: Solving for Uy, Initial residual = 0.999633, Final residual = 6.0743e-08, No Iterations 2
DILUPBiCG: Solving for Uz, Initial residual = 0.999938, Final residual = 6.09759e-08, No Iterations 2
GAMG: Solving for p, Initial residual = 1, Final residual = 9.34594e-07, No Iterations 191
GAMG: Solving for p, Initial residual = 3.59854e-08, Final residual = 3.59854e-08, No Iterations 0
time step continuity errors : sum local = 2.89543e+37, global = 2.81514e+22, cumulative = 2.81514e+22
DILUPBiCG: Solving for epsilon, Initial residual = 1, Final residual = 1.65232e-08, No Iterations 4
bounding epsilon, min: -1.91717e+65 max: 8.52119e+80 average: 1.40408e+75
DILUPBiCG: Solving for k, Initial residual = 0.999998, Final residual = 5.26354e-16, No Iterations 1
bounding k, min: -2.04901e+54 max: 7.78681e+72 average: 8.76528e+66
ExecutionTime = 19859.3 s ClockTime = 19876 s


Only a few iterations before, everything seems to be okay! The residuals are low, even for epsilon. Only the message of bounding epsilon appears a few times during the calculation. Compared with a Fluent calculation so far it looks good!

Here the schemes I used:

ddtSchemes
{
default steadyState;
}
gradSchemes
{
default Gauss linear;
grad(p) Gauss linear;
grad(U) Gauss linear;
}
divSchemes
{
default none;
div(phi,U) Gauss upwind;
div(phi,k) Gauss upwind;
div(phi,epsilon) Gauss upwind;
div(U) Gauss upwind;
div((nuEff*dev(grad(U).T()))) Gauss linear;
}
laplacianSchemes
{
default none
laplacian(nu,U) Gauss upwind phi corrected;
laplacian((1|A(U)),p) Gauss upwind phi corrected;
laplacian(nuEff,U) Gauss upwind phi corrected;
laplacian(DkEff,k) Gauss upwind phi corrected;
laplacian(DepsilonEff,epsilon) Gauss upwind phi corrected;
}
interpolationSchemes
{
default upwind phi;
}
snGradSchemes
{
default corrected;
}
fluxRequired
{
default no;
p ;
}

And the solver:

solvers
{
U
{
solver PBiCG;
preconditioner DILU;
tolerance 1e-05;
relTol 0;
}
p
{
solver GAMG;
preconditioner DIC;
tolerance 1e-06;
relTol 0;
smoother GaussSeidel;
nCellsInCoarsestLevel 20;
agglomerator faceAreaPair;
mergeLevels 1;
}
k
{
solver PBiCG;
preconditioner DILU;
tolerance 1e-05;
relTol 0;
}
epsilon
{
solver PBiCG;
preconditioner DILU;
tolerance 1e-05;
relTol 0;
}
}
SIMPLE
{
nNonOrthogonalCorrectors 1;
}

relaxationFactors
{
p 0.1;
U 0.2;
k 0.2;
epsilon 0.2;
}


I use fixed BC at the inlets for k and epsilon. I have calculated them according my initialisation values. I also tried to update them according to a Fluent simulation after around 500 iterations to better values, but it didnt help. I know that variable BCs like a turbulent lenght scale based one or based on the eddy viscosity ratio would be better. But using them blows my calculation up after only a few iterations.

Please help and dont hesitate to ask for more information.

Thanks a lot!

Sebastian

[Chris Lucas]

Hi,

could you please upload your case file (zero folder, system folder). It is difficult to give you a good answer without all the information.

Regards,
Christian

[sebastian]

Hi Chris,

how do I upload it exactly?

Sebastian

[Chris Lucas]

Hi,

there is a button with a paper clip (next to the simley button or the Undo button). However, each file must be smaller than 100kB.

If you initials your calculation with potentialFoam, please, only include the boundary conditions.

Regards,
Christian

[sebastian]

Okay, here is my case (/0 and /system folder)

Thanks to everybody who can take a look at it!

Sebastian

[Chris Lucas]

Hi,

I have a few suggestions. Firstly, at the outlet, I would recommend a inletOutet BC for U and a fixedMeanValue for p.

fixedMeanValue :
http://www.cfd-online.com/Forums/ope...condition.html

Have you tried to use the turbulent inlet BC I mentioned above?

Have you tried to use a different interpolation scheme like linear or QUICK. In a recent simulation of a compressible nozzle flow, I had some problems with the upwind scheme

Have you look at the solution before the simulation crashed?

Regards,
Christian

[sebastian]

Hi Christian,

first of all thanks for your effort!

Yes, I have tried to run it with the turbulent BC you mentioned. You can find the way I used them in the folder /0 named epsilon.functions an k.functions.
But the calculation crushed after a few iterations, so I guess there went something wrong with them... Maybe you want to have a look on them.

No, I only used upwind for interpolation so far. I thought that would be the most stable variant.

Yes, I looked at the solution right before it crashed. So far it looked quite reasonable!


Best wishes,
Sebastian

[maddalena]

Hi, I will also try to limit divergence... this is similar to what Fluent has:

Originally Posted by sebastian

divSchemes
{
default none;
div(phi,U) Gauss linearUpwindV cellLimited Gauss linear 1;
div(phi,k) Gauss linearUpwind cellLimited Gauss linear 1;
div(phi,epsilon) Gauss linearUpwind cellLimited Gauss linear 1;
div(U) Gauss linearUpwindV cellLimited Gauss linear 1;
div((nuEff*dev(grad(U).T()))) Gauss linear;
}

and maybe increase the relative tolerance of the solvers: for example:

Originally Posted by sebastian

solvers
{
U
{
solver PBiCG;
preconditioner DILU;
tolerance 1e-05;
relTol 1e-2;
}

hope this help,

mad