
[Sponsors] 
Bounding epsilon error while running LRR/SSG using interFoam 

LinkBack  Thread Tools  Search this Thread  Display Modes 
February 26, 2021, 08:27 
Bounding epsilon error while running LRR/SSG using interFoam

#1 
New Member
Join Date: Jan 2021
Location: Norway
Posts: 22
Rep Power: 5 
Hello to every one,
I am trying to simulate supercritical multiphase flow (interFoam) using different RAS models to get the best solution. The simulation is running with ke, RNG, kwSST etc., but for LRR and SSG models, I am getting bounding epsilon message from the very beginning of the simulation. The epsilon blows up and the simulation stops after 12 seconds of simulation time. Further, in LRR, SSG Code:
div(((rho*nuEff)*dev2(T(grad(U))))) Gauss linear; was not accepted and div(((rho*nu)*dev2(T(grad(U))))) Gauss linear; was accepted. Code:
ExecutionTime = 89.92 s ClockTime = 98 s Courant Number mean: 0.298636 max: 0.895395 Interface Courant Number mean: 0.00270084 max: 0.895395 deltaT = 0.000471525 Time = 0.08411459 smoothSolver: Solving for alpha.water, Initial residual = 0.00299815, Final residual = 1.53669e09, No Iterations 3 Phase1 volume fraction = 0.97951 Min(alpha.water) = 2.46878e07 Max(alpha.water) = 1 MULES: Correcting alpha.water MULES: Correcting alpha.water Phase1 volume fraction = 0.97951 Min(alpha.water) = 0.000233734 Max(alpha.water) = 1.00009 GAMG: Solving for p_rgh, Initial residual = 0.00769725, Final residual = 3.59405e05, No Iterations 2 GAMG: Solving for p_rgh, Initial residual = 6.21739e05, Final residual = 2.90274e07, No Iterations 4 time step continuity errors : sum local = 6.90028e08, global = 3.90763e08, cumulative = 4.71314e06 GAMG: Solving for p_rgh, Initial residual = 6.95841e05, Final residual = 2.37694e07, No Iterations 3 GAMG: Solving for p_rgh, Initial residual = 8.27069e07, Final residual = 9.28369e09, No Iterations 3 time step continuity errors : sum local = 2.20722e09, global = 9.53455e10, cumulative = 4.71219e06 GAMG: Solving for p_rgh, Initial residual = 1.6805e05, Final residual = 4.85751e08, No Iterations 3 GAMG: Solving for p_rgh, Initial residual = 1.42639e07, Final residual = 6.06823e09, No Iterations 2 time step continuity errors : sum local = 1.44347e09, global = 8.06204e10, cumulative = 4.71138e06 DILUPBiCG: Solving for epsilon, Initial residual = 0.999478, Final residual = 3.1044e09, No Iterations 3 bounding epsilon, min: 1.43351e08 max: 0.00811 average: 4.5564e08 DILUPBiCG: Solving for Rxx, Initial residual = 0.967164, Final residual = 5.10025e10, No Iterations 1 DILUPBiCG: Solving for Rxy, Initial residual = 0.778008, Final residual = 1.6466e09, No Iterations 1 DILUPBiCG: Solving for Ryy, Initial residual = 0.98587, Final residual = 4.91203e10, No Iterations 1 DILUPBiCG: Solving for Rzz, Initial residual = 0.989121, Final residual = 5.00174e10, No Iterations 1 ExecutionTime = 90.44 s ClockTime = 99 s Subhojit Last edited by Tobi; March 3, 2021 at 04:45. 

March 2, 2021, 08:02 

#2  
Super Moderator
Tobias Holzmann
Join Date: Oct 2010
Location: Bad Wörishofen
Posts: 2,711
Blog Entries: 6
Rep Power: 51 
Quote:
Maybe a bug. To the other topic you mention. If you use the LLR turbulence model, are you underrelaxing the REqn? You can make it at least diagonal dominant by setting the matrix relaxation to 1 (same for epsilon). Same is valid for SSG. In addition, you did not provide any error message or other useful hints. fvSolution and the schemes you are using are of interest too.
__________________
Keep foaming, Tobias Holzmann 

March 2, 2021, 08:54 

#3 
New Member
Join Date: Jan 2021
Location: Norway
Posts: 22
Rep Power: 5 
Many thanks for your valuable comment Tobi. I have used relaxation of 1 for the abovementioned terms. The fvSolution and fvSchemes are:
Code:
ddtSchemes { default Euler; } gradSchemes { default faceLimited leastSquares 0.5; grad(U) faceLimited leastSquares 0.5; //grad(U) cellLimited leastSquares 0.5; //grad(k) cellLimited leastSquares 0.5; //grad(epsilon) cellLimited leastSquares 0.5; //default Gauss linear; } divSchemes { default none; div(rhoPhi,U) Gauss linearUpwind grad(U);//Gauss LUST grad(U); div(phi,alpha) Gauss vanLeer; div(phirb,alpha) Gauss linear; div(phi,k) Gauss upwind; div(phi,epsilon) Gauss upwind; div(phi,omega) Gauss upwind; //for kOmega, kOmegaSST div(((rho*nuEff)*dev2(T(grad(U))))) Gauss linear; div(phi,R) Gauss upwind; //For LRR, SSG div(R) Gauss linear; //For LRR, SSG div((rho*R)) Gauss linear; //For LRR, SSG div(((rho*nu)*dev2(T(grad(U))))) Gauss linear; //For LRR, SSG div(nonlinearStress) Gauss linear; //For LRR, SSG } laplacianSchemes { default Gauss linear limited corrected 0.33; //default Gauss linear corrected; } interpolationSchemes { default linear; } snGradSchemes { default limited corrected 0.33; //default corrected; } wallDist // for LRR, SSG { method meshWave; } // ************************************************************************* // Code:
solvers // water channel, dam break tutorials { "alpha.water.*" { nAlphaCorr 2; nAlphaSubCycles 1; cAlpha 1; MULESCorr yes; nLimiterIter 3; solver smoothSolver; smoother symGaussSeidel; tolerance 1e8; relTol 0; } "pcorr.*" { solver PCG; preconditioner { preconditioner GAMG; tolerance 1e05; relTol 0; smoother GaussSeidel; } tolerance 1e05; relTol 0; maxIter 50; } p_rgh { solver GAMG; tolerance 1e08; relTol 0.01; smoother GaussSeidel; } p_rghFinal { $p_rgh; relTol 0; } "(UkepsilonomegaR)" //add R for LRR and SSG, omega for kOmega, kOmegaSST { solver PBiCG; preconditioner DILU; tolerance 1e06; relTol 0; } "(UkepsilonomegaR)Final" //add R for LRR and SSG, omega for kOmega, kOmegaSST { $U; tolerance 1e08; relTol 0; } } PIMPLE { momentumPredictor off; // off for multiphase nCorrectors 3; nNonOrthogonalCorrectors 1; //residualControl //{ U 1e07; p_rgh 1e06; //} } relaxationFactors { equations { ".*" 1; "(UkepsilonomegaR)" 1; } } Code:
ExecutionTime = 804.44 s ClockTime = 827 s Courant Number mean: 1.97399e+30 max: 1.25118e+35 Interface Courant Number mean: 1.34134e+09 max: 1.12614e+13 deltaT = 2.07961e39 Time = 1.719277 smoothSolver: Solving for alpha.water, Initial residual = 1.14868e11, Final residual = 1.14868e11, No Iterations 0 Phase1 volume fraction = 0.741719 Min(alpha.water) = 1.10889e05 Max(alpha.water) = 1.00002 MULES: Correcting alpha.water MULES: Correcting alpha.water Phase1 volume fraction = 0.741719 Min(alpha.water) = 1.10889e05 Max(alpha.water) = 1.00002 GAMG: Solving for p_rgh, Initial residual = 1, Final residual = 0.0050873, No Iterations 6 GAMG: Solving for p_rgh, Initial residual = 0.00245794, Final residual = 1.94538e05, No Iterations 5 time step continuity errors : sum local = 7.6869e10, global = 5.68414e10, cumulative = 3.94798e+30 GAMG: Solving for p_rgh, Initial residual = 0.00318051, Final residual = 1.26821e05, No Iterations 4 GAMG: Solving for p_rgh, Initial residual = 7.93767e05, Final residual = 6.21374e07, No Iterations 4 time step continuity errors : sum local = 2.66209e11, global = 1.98525e11, cumulative = 3.94798e+30 GAMG: Solving for p_rgh, Initial residual = 0.000329261, Final residual = 2.90286e06, No Iterations 5 GAMG: Solving for p_rgh, Initial residual = 2.09426e05, Final residual = 7.86557e09, No Iterations 8 time step continuity errors : sum local = 3.3573e13, global = 2.25506e13, cumulative = 3.94798e+30 DILUPBiCG: Solving for epsilon, Initial residual = 0.999135, Final residual = 1.41005e14, No Iterations 1 bounding epsilon, min: 3.97239e+38 max: 3.03117e+52 average: 2.01431e+48 DILUPBiCG: Solving for Rxx, Initial residual = 0.979011, Final residual = 14.1644, No Iterations 1000 DILUPBiCG: Solving for Rxy, Initial residual = 0.000176047, Final residual = 0.0894335, No Iterations 1000 DILUPBiCG: Solving for Ryy, Initial residual = 0.999409, Final residual = 9.48781e09, No Iterations 5 DILUPBiCG: Solving for Rzz, Initial residual = 1, Final residual = 136864, No Iterations 1000 ExecutionTime = 813.45 s ClockTime = 836 s Courant Number mean: 3.99863e05 max: 0.212594 Interface Courant Number mean: 7.22945e09 max: 1.22022e05 deltaT = 2.49554e39 Time = 1.719277 smoothSolver: Solving for alpha.water, Initial residual = 6.81458e10, Final residual = 6.81458e10, No Iterations 0 Phase1 volume fraction = 0.741719 Min(alpha.water) = 1.10889e05 Max(alpha.water) = 1.00002 MULES: Correcting alpha.water MULES: Correcting alpha.water Phase1 volume fraction = 0.741719 Min(alpha.water) = 1.10889e05 Max(alpha.water) = 1.00002 GAMG: Solving for p_rgh, Initial residual = 0.26929, Final residual = 0.00255839, No Iterations 10 GAMG: Solving for p_rgh, Initial residual = 3.55552e73, Final residual = 3.55552e73, No Iterations 0 time step continuity errors : sum local = 2.61007e07, global = 1.32084e07, cumulative = 3.94798e+30 GAMG: Solving for p_rgh, Initial residual = 4.16704e73, Final residual = 4.16704e73, No Iterations 0 GAMG: Solving for p_rgh, Initial residual = 4.16704e73, Final residual = 4.16704e73, No Iterations 0 time step continuity errors : sum local = 3.05898e07, global = 1.08856e07, cumulative = 3.94798e+30 GAMG: Solving for p_rgh, Initial residual = 4.22678e73, Final residual = 4.22678e73, No Iterations 0 GAMG: Solving for p_rgh, Initial residual = 4.22678e73, Final residual = 4.22678e73, No Iterations 0 time step continuity errors : sum local = 3.10284e07, global = 1.06426e07, cumulative = 3.94798e+30 DILUPBiCG: Solving for epsilon, Initial residual = 5.59504e21, Final residual = 5.59504e21, No Iterations 0 DILUPBiCG: Solving for Rxx, Initial residual = 1, Final residual = 1, No Iterations 1000 DILUPBiCG: Solving for Rxy, Initial residual = 0.999986, Final residual = 0.999986, No Iterations 1000 DILUPBiCG: Solving for Ryy, Initial residual = 1, Final residual = 1, No Iterations 1000 DILUPBiCG: Solving for Rzz, Initial residual = 1, Final residual = 1, No Iterations 1000 ExecutionTime = 825.17 s ClockTime = 848 s Primary job terminated normally, but 1 process returned a nonzero exit code.. Per userdirection, the job has been aborted.   mpirun detected that one or more processes exited with nonzero status, thus causing the job to be terminated. The first process to do so was: Process name: [[21610,1],0] Exit code: 145 In ReynoldsStress.c file at line 243 Code:
 fvc::div(this>alpha_*rho*this>nu()*dev2(T(fvc::grad(U)))) Code:
const scalar minVsf = min(vsf).value(); if (minVsf < lowerBound.value()) { Info<< "bounding " << vsf.name() << ", min: " << minVsf << " max: " << max(vsf).value() << " average: " << gAverage(vsf.primitiveField()) << endl; vsf.primitiveFieldRef() = max ( max ( vsf.primitiveField(), fvc::average(max(vsf, lowerBound))().primitiveField() * pos0(vsf.primitiveField()) ), lowerBound.value() ); vsf.boundaryFieldRef() = max(vsf.boundaryField(), lowerBound.value()); } return vsf; } Thanking you. Subhojit Last edited by Tobi; March 3, 2021 at 04:46. Reason: Tobi: Cleaned the mess 

March 2, 2021, 11:22 

#4 
Super Moderator
Tobias Holzmann
Join Date: Oct 2010
Location: Bad Wörishofen
Posts: 2,711
Blog Entries: 6
Rep Power: 51 
First of all, please clean your posts. Its just a mess. We have the code tag to include dictionaries, code and others in order to keep the format.
Probably your problem starts somewhere before: Code:
Courant Number mean: 1.97399e+30 max: 1.25118e+35 Interface Courant Number mean: 1.34134e+09 max: 1.12614e+13 deltaT = 2.07961e39 Time = 1.719277 Anything went wrong here previously. Check out your states in paraview, e.g., at 1.70 s, 1,71 s 1,715 s to check where the problem starts and to analyze why your simulation crashes. As you can see, the Co Number kills everything. Probably somewhere is a high velocity (which can be limited by the fvOptions to make your simulation more stable).
__________________
Keep foaming, Tobias Holzmann 

March 2, 2021, 11:57 

#5 
New Member
Join Date: Jan 2021
Location: Norway
Posts: 22
Rep Power: 5 
Thanks. I tried several fvSchemes to limit the alpha between 0 and 1, but still some ve values are there which is possibly generating the ve epsilon. The problem is starting with bounding epsilon and after some time getting a very high value of epsilon (both +ve and ve).
Subhojit 

March 2, 2021, 12:16 

#6 
Super Moderator
Tobias Holzmann
Join Date: Oct 2010
Location: Bad Wörishofen
Posts: 2,711
Blog Entries: 6
Rep Power: 51 
Bounding problems occur from time to time and therefore, we have the bounding functions that change the values to physical ones. The alpha field can have small negative values (this is a numerical topic) but I had such values in my calculations too. However, I never used the turbulence models you are using and hence, no idea what can cause the problem. In addition, we (the community) do not have insight into your case so well ... a bit tricky
__________________
Keep foaming, Tobias Holzmann 

March 3, 2021, 04:32 

#7  
New Member
Join Date: Jan 2021
Location: Norway
Posts: 22
Rep Power: 5 
Quote:


March 3, 2021, 04:47 

#8 
Super Moderator
Tobias Holzmann
Join Date: Oct 2010
Location: Bad Wörishofen
Posts: 2,711
Blog Entries: 6
Rep Power: 51 
Very important > https://letmegooglethat.com/?q=fvOptions+limit+velocity
Including limiters in fvOptions
__________________
Keep foaming, Tobias Holzmann 

March 3, 2021, 08:01 

#9  
New Member
Join Date: Jan 2021
Location: Norway
Posts: 22
Rep Power: 5 
Quote:
Subhojit 

March 3, 2021, 08:05 

#10 
Super Moderator
Tobias Holzmann
Join Date: Oct 2010
Location: Bad Wörishofen
Posts: 2,711
Blog Entries: 6
Rep Power: 51 
For epsilon you don't need anything and definitely no source term . I guess you meant if you need some settings in fvOptions for epsilon. Well, the limitVelocity fvOptions will try to limit the velocity if you exceed the maximum number you specified. However, that does not mean that your simulation will run fine/correct/without errors. Its just a try to stabilize the simulation.
__________________
Keep foaming, Tobias Holzmann 

March 3, 2021, 08:45 

#11  
New Member
Join Date: Jan 2021
Location: Norway
Posts: 22
Rep Power: 5 
Quote:
What I noticed is that at most of the cells in the domain, the epsilon is close to the minimum value of 2.2e16, and the calculated k (from R) is close to the minimum value of 3.3e16. And at some cells close to the interface, the epsilon is increasing rapidly. The maximum velocity is below 10 m/s which is okay for my case. Subhojit 

March 3, 2021, 13:13 

#13 
New Member
Join Date: Jan 2021
Location: Norway
Posts: 22
Rep Power: 5 
Many thanks Tobi. Lets see how it goes.


March 8, 2022, 04:50 

#14 
New Member
SunTime
Join Date: Nov 2020
Posts: 14
Rep Power: 5 
Hi, have you solved your problem? Recently， I meet the same problem，do you have some experience to share with me?


Tags 
bounding error, interfoam diverging, multiphase flow, reynolds stress model 
Thread Tools  Search this Thread 
Display Modes  


Similar Threads  
Thread  Thread Starter  Forum  Replies  Last Post 
interfoam parallel running exit code 145  fatlem  OpenFOAM Running, Solving & CFD  6  February 25, 2023 07:58 
SimpleFoam k and epsilon bounded  nedved  OpenFOAM Running, Solving & CFD  16  March 4, 2017 08:30 
How to write k and epsilon before the abnormal end  xiuying  OpenFOAM Running, Solving & CFD  8  August 27, 2013 15:33 
SimpleFoam k and epsilon bounded  nedved  OpenFOAM Running, Solving & CFD  1  November 25, 2008 20:21 
Bounding epsilon and K with rasInterFoam  openfoam_user  OpenFOAM Running, Solving & CFD  0  October 23, 2008 08:48 