
[Sponsors] 
October 19, 2020, 09:43 
PIMPLE algorithm ignores residual control

#1 
New Member
Dhruv
Join Date: Jan 2020
Location: Blacksburg, Virginia
Posts: 6
Rep Power: 6 
Hello all
I am trying to simulate a cavitating flow inside a venturisection using InterPhaseChangeFoam an unsteady, multiphase solver. Here is my fvSolution Code:
/** C++ **\  =========    \\ / F ield  OpenFOAM: The Open Source CFD Toolbox   \\ / O peration  Version: 2.3.0   \\ / A nd  Web: www.OpenFOAM.org   \\/ M anipulation   \**/ FoamFile { version 2.0; format ascii; class dictionary; object fvSolution; } // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // solvers { "alpha.water.*" { cAlpha 0; nAlphaCorr 2; nAlphaSubCycles 1; MULESCorr yes; nLimiterIter 5; solver PBiCG; preconditioner DILU; tolerance 1e8;//1e8 relTol 0; }; "U.*" { solver PBiCG; preconditioner DILU; tolerance 1e6; //1e6 relTol 0; }; "(komega)" { solver PBiCG; preconditioner DILU; tolerance 1e6; //1e6 relTol 0; }; "(Ukomega)Final" { solver PBiCG; preconditioner DILU; tolerance 1e6; //1e6 relTol 0; } p_rgh { solver GAMG; tolerance 1e8; //1e8 relTol 0.1; smoother DICGaussSeidel; nPreSweeps 0; nPostSweeps 2; cacheAgglomeration true; nCellsInCoarsestLevel 10; agglomerator faceAreaPair; mergeLevels 1; maxIter 50; }; nuTilda { solver PBiCG; preconditioner DILU; smoother DICGaussSeidel; nSweeps 2; tolerance 1e6; relTol 0; }; p_rghFinal { solver PCG; preconditioner { preconditioner GAMG; tolerance 1e6; //1e6 relTol 0; nVcycles 2; smoother DICGaussSeidel; nPreSweeps 0; nPostSweeps 2; nFinestSweeps 2; cacheAgglomeration true; nCellsInCoarsestLevel 10; agglomerator faceAreaPair; mergeLevels 1; }; tolerance 1e7; //1e7 relTol 0; maxIter 50; }; } PIMPLE { momentumPredictor yes; nOuterCorrectors 100; nCorrectors 3; nNonOrthogonalCorrectors 0; residualControl //outerCorrectorResidualControl { "(p_rghU)" { tolerance 5e5;//5e5 relTol 0;//0.1 } "(komegaepsilon)" { tolerance 1e4; relTol 0; } } } relaxationFactors { fields { p 0.8; } equations { U 0.8; // 0.9 is more stable but 0.95 more convergent ".*" 0.8; // 0.9 is more stable but 0.95 more convergent k 0.8; omega 0.8; } } // ************************************************************************* // However, looking at the residual files, the iterations continue even after dropping to the required level Code:
Courant Number mean: 0.0528581 max: 2.48597 Time = 1.39901 PIMPLE: iteration 1 DILUPBiCG: Solving for alpha.water, Initial residual = 3.68011e07, Final residual = 1.36401e09, No Iterations 1 Phase1 volume fraction = 0.991636 Min(alpha.water) = 0.0154223 Max(alpha.water) = 1 MULES: Correcting alpha.water MULES: Correcting alpha.water Liquid phase volume fraction = 0.991616 Min(alpha.water) = 0.0154223 Max(alpha.water) = 1 DILUPBiCG: Solving for Ux, Initial residual = 0.00127737, Final residual = 7.84762e07, No Iterations 1 DILUPBiCG: Solving for Uy, Initial residual = 0.00227295, Final residual = 3.61585e09, No Iterations 2 GAMG: Solving for p_rgh, Initial residual = 6.75951e06, Final residual = 5.29966e10, No Iterations 1 GAMG: Solving for p_rgh, Initial residual = 1.63286e07, Final residual = 3.04435e10, No Iterations 1 GAMGPCG: Solving for p_rgh, Initial residual = 1.90892e08, Final residual = 1.90892e08, No Iterations 0 This is kOmegaSSTF DILUPBiCG: Solving for omega, Initial residual = 0.000120964, Final residual = 7.63802e08, No Iterations 1 bounding omega, min: 29.5739 max: 2.83886e+06 average: 82685.8 DILUPBiCG: Solving for k, Initial residual = 0.00256498, Final residual = 2.6164e09, No Iterations 2 bounding k, min: 0.100634 max: 22.1979 average: 0.703296 PIMPLE: iteration 2 DILUPBiCG: Solving for alpha.water, Initial residual = 1.51834e07, Final residual = 3.8706e11, No Iterations 1 Phase1 volume fraction = 0.991644 Min(alpha.water) = 0.0154281 Max(alpha.water) = 1 MULES: Correcting alpha.water MULES: Correcting alpha.water Liquid phase volume fraction = 0.991624 Min(alpha.water) = 0.0154281 Max(alpha.water) = 1 DILUPBiCG: Solving for Ux, Initial residual = 0.000288568, Final residual = 1.69436e07, No Iterations 1 DILUPBiCG: Solving for Uy, Initial residual = 0.000502123, Final residual = 2.93426e07, No Iterations 1 GAMG: Solving for p_rgh, Initial residual = 6.27618e07, Final residual = 6.02198e10, No Iterations 1 GAMG: Solving for p_rgh, Initial residual = 2.07927e07, Final residual = 3.39953e10, No Iterations 1 GAMGPCG: Solving for p_rgh, Initial residual = 1.89539e08, Final residual = 1.89539e08, No Iterations 0 This is kOmegaSSTF DILUPBiCG: Solving for omega, Initial residual = 3.04491e05, Final residual = 1.35276e08, No Iterations 1 bounding omega, min: 22.4256 max: 2.83873e+06 average: 82690.8 DILUPBiCG: Solving for k, Initial residual = 0.000716519, Final residual = 5.31525e07, No Iterations 1 bounding k, min: 0.0272711 max: 22.274 average: 0.703536 PIMPLE: iteration 3 DILUPBiCG: Solving for alpha.water, Initial residual = 1.46881e07, Final residual = 3.71256e11, No Iterations 1 Phase1 volume fraction = 0.991645 Min(alpha.water) = 0.0154305 Max(alpha.water) = 1 MULES: Correcting alpha.water MULES: Correcting alpha.water Liquid phase volume fraction = 0.991625 Min(alpha.water) = 0.0154305 Max(alpha.water) = 1 DILUPBiCG: Solving for Ux, Initial residual = 6.57397e05, Final residual = 3.9333e08, No Iterations 1 DILUPBiCG: Solving for Uy, Initial residual = 0.000114833, Final residual = 7.40035e08, No Iterations 1 GAMG: Solving for p_rgh, Initial residual = 1.85128e07, Final residual = 3.77726e10, No Iterations 1 GAMG: Solving for p_rgh, Initial residual = 2.61274e08, Final residual = 2.8472e10, No Iterations 1 GAMGPCG: Solving for p_rgh, Initial residual = 1.89681e08, Final residual = 1.89681e08, No Iterations 0 This is kOmegaSSTF DILUPBiCG: Solving for omega, Initial residual = 8.63384e06, Final residual = 3.5632e09, No Iterations 1 bounding omega, min: 39.1196 max: 2.83868e+06 average: 82692 DILUPBiCG: Solving for k, Initial residual = 0.000213868, Final residual = 1.30345e07, No Iterations 1 bounding k, min: 0.0417151 max: 22.2951 average: 0.703467 PIMPLE: iteration 4 DILUPBiCG: Solving for alpha.water, Initial residual = 1.45561e07, Final residual = 3.72721e11, No Iterations 1 Phase1 volume fraction = 0.991645 Min(alpha.water) = 0.0154314 Max(alpha.water) = 1 MULES: Correcting alpha.water MULES: Correcting alpha.water Liquid phase volume fraction = 0.991625 Min(alpha.water) = 0.0154314 Max(alpha.water) = 1 DILUPBiCG: Solving for Ux, Initial residual = 1.53741e05, Final residual = 1.0046e08, No Iterations 1 DILUPBiCG: Solving for Uy, Initial residual = 2.80869e05, Final residual = 2.3609e08, No Iterations 1 GAMG: Solving for p_rgh, Initial residual = 7.86704e08, Final residual = 2.95063e10, No Iterations 1 GAMG: Solving for p_rgh, Initial residual = 2.12339e08, Final residual = 2.34074e10, No Iterations 1 GAMGPCG: Solving for p_rgh, Initial residual = 1.88666e08, Final residual = 1.88666e08, No Iterations 0 This is kOmegaSSTF DILUPBiCG: Solving for omega, Initial residual = 4.02313e06, Final residual = 1.46929e09, No Iterations 1 bounding omega, min: 22.9509 max: 2.83865e+06 average: 82692.7 DILUPBiCG: Solving for k, Initial residual = 0.000139817, Final residual = 4.32277e08, No Iterations 1 bounding k, min: 0.0272965 max: 22.301 average: 0.703634 PIMPLE: iteration 5 DILUPBiCG: Solving for alpha.water, Initial residual = 1.44821e07, Final residual = 3.73634e11, No Iterations 1 Phase1 volume fraction = 0.991646 Min(alpha.water) = 0.0154317 Max(alpha.water) = 1 MULES: Correcting alpha.water MULES: Correcting alpha.water Liquid phase volume fraction = 0.991625 Min(alpha.water) = 0.0154317 Max(alpha.water) = 1 DILUPBiCG: Solving for Ux, Initial residual = 3.97181e06, Final residual = 3.51636e09, No Iterations 1 DILUPBiCG: Solving for Uy, Initial residual = 9.34745e06, Final residual = 1.30711e08, No Iterations 1 GAMG: Solving for p_rgh, Initial residual = 5.11634e08, Final residual = 2.30832e10, No Iterations 1 GAMG: Solving for p_rgh, Initial residual = 2.08343e08, Final residual = 2.0801e10, No Iterations 1 GAMGPCG: Solving for p_rgh, Initial residual = 1.8864e08, Final residual = 1.8864e08, No Iterations 0 This is kOmegaSSTF DILUPBiCG: Solving for omega, Initial residual = 2.51697e06, Final residual = 1.89748e09, No Iterations 1 bounding omega, min: 50.392 max: 2.83865e+06 average: 82692.4 DILUPBiCG: Solving for k, Initial residual = 8.56674e05, Final residual = 6.59347e07, No Iterations 1 bounding k, min: 0.0407519 max: 22.3026 average: 0.703485 PIMPLE: iteration 6 DILUPBiCG: Solving for alpha.water, Initial residual = 1.44401e07, Final residual = 3.73939e11, No Iterations 1 Phase1 volume fraction = 0.991646 Min(alpha.water) = 0.0154318 Max(alpha.water) = 1 MULES: Correcting alpha.water MULES: Correcting alpha.water Liquid phase volume fraction = 0.991625 Min(alpha.water) = 0.0154318 Max(alpha.water) = 1 DILUPBiCG: Solving for Ux, Initial residual = 1.59784e06, Final residual = 2.50779e09, No Iterations 1 DILUPBiCG: Solving for Uy, Initial residual = 5.68396e06, Final residual = 3.10249e08, No Iterations 1 GAMG: Solving for p_rgh, Initial residual = 6.23817e08, Final residual = 2.25845e10, No Iterations 1 GAMG: Solving for p_rgh, Initial residual = 2.04297e07, Final residual = 1.95536e10, No Iterations 1 GAMGPCG: Solving for p_rgh, Initial residual = 1.87125e08, Final residual = 1.87125e08, No Iterations 0 This is kOmegaSSTF DILUPBiCG: Solving for omega, Initial residual = 2.74728e06, Final residual = 7.57324e10, No Iterations 1 bounding omega, min: 28.6162 max: 2.83864e+06 average: 82692 DILUPBiCG: Solving for k, Initial residual = 0.000144201, Final residual = 4.12669e08, No Iterations 1 bounding k, min: 0.0272926 max: 22.3035 average: 0.703574 PIMPLE: converged in 6 iterations ExecutionTime = 0.6 s ClockTime = 1 s For reference, here is my fvSchemes Code:
/** C++ **\  =========    \\ / F ield  OpenFOAM: The Open Source CFD Toolbox   \\ / O peration  Version: 2.3.0   \\ / A nd  Web: www.OpenFOAM.org   \\/ M anipulation   \**/ FoamFile { version 2.0; format ascii; class dictionary; object fvSchemes; } // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // ddtSchemes { default Euler; } interpolationSchemes { default linear; } divSchemes { default none; div(rhoPhi,U) Gauss linearUpwind grad(U); div(phi,omega) Gauss linearUpwind grad(omega); div(phi,k) Gauss linearUpwind grad(k); div(phi,omega) Gauss linearUpwind grad(omega); div(phi,alpha) Gauss vanLeer; div(phirb,alpha) Gauss linear; div(phi,nuTilda) Gauss linearUpwind grad(nuTilda); // div((muEff*dev(T(grad(U))))) Gauss linear; div(((rho*nuEff)*dev2(T(grad(U))))) Gauss linear; } gradSchemes { default Gauss linear; } laplacianSchemes { default Gauss linear limited corrected 0.5; } snGradSchemes { default limited corrected 0.5; } fluxRequired { default none; p_rgh; pcorr; alpha.water; } wallDist { method meshWave; } // ************************************************************************* // Best TurbulentGuy 

October 21, 2020, 03:55 

#2 
New Member
Join Date: Feb 2020
Location: Germany
Posts: 17
Rep Power: 6 
Hello TurbulentGuy,
this might arise from your underrelaxation. You put the factor 0.8 for basically all quantities, which is then used between the PIMPLE correctors. However, as it is not specifically stated, the last corrector will not use underrelaxation. So to speak, "UFinal", "kFinal", ... are set to 1. The solver probably even recognises the time step is converged after 5 iterations and calculates one last corrector step (#6) with no underrelaxation. This is just a guess and worth a try, if it still appears with no underrelaxation at all. I have to admit, that I am also confused about the "residualControl" and "outerResidualControl" stuff. One last thing for the turbulence, as I have similar problems right now. The bounding of omega looks bad to be honest. If the values differ that much, there is probably something wrong going on in the calculated physics. Just as a side note. Good luck and let me know, if that was the source of error. Best regards Doerschlott 

October 21, 2020, 07:45 

#3 
New Member
Join Date: Mar 2020
Posts: 18
Rep Power: 6 
Hi TurbulentGuy,
I am not sure if it's the same with PIMPLE, but for SimpleFoam and in OpenFoam version 7 the syntax is like this SIMPLE { nNonOrthogonalCorrectors 10; nCorrectors 5, nOuterCorrectors 5; consistent yes; residualControl { p 1e4; U 1e4; } } You can see if it worked right after you started the solver. It'll say "Create time Create mesh for time = 0 SIMPLE: Convergence criteria found p: tolerance 0.0001 U: tolerance 0.0001 ......" Best Maike 

October 22, 2020, 05:10 

#4 
Senior Member
Carlo_P
Join Date: May 2019
Location: Italy
Posts: 176
Rep Power: 7 
Hey, on which version are you running?
In the openfoam7, it was a bug. In opefoam 8 is working. Next trick... the residual control should be written like this: residualControl { p 1e4; U 1e4; ".*" 1e4; } so you are sure that you will include all the variables 

October 23, 2020, 21:57 

#5  
New Member
Dhruv
Join Date: Jan 2020
Location: Blacksburg, Virginia
Posts: 6
Rep Power: 6 
Quote:
First of all thank you for your insights. I did not know the final underelaxation terms. However, even after I write them, the calculation continues even after the residual control limits are reached. I thought maybe it's because of k and omega bounding limits(as you mentioned) but they don't change much per iteration. Do you have an idea how to take care of that? Best TurbulentGuy 

October 23, 2020, 22:05 

#6  
New Member
Dhruv
Join Date: Jan 2020
Location: Blacksburg, Virginia
Posts: 6
Rep Power: 6 
Quote:
I am using OpenFOAMv1912. And it does work as you guys say(the instructions are printed when I run the program). Best TurbulentGuy 

October 27, 2020, 07:58 

#7 
New Member
Join Date: Feb 2020
Location: Germany
Posts: 17
Rep Power: 6 
Hey TurbulentGuy,
the thing is, that as soon as the convergence criteria is met, a final corrector loop (iteration #6) is conducted, regardless of the stated underrelaxation. Though I'm glad I could teach you about the final terms. The residual control definitely worked, as you put 100 for the outerCorrectors and it clearly finishes the time step before that. For the bounding problem: I had the problem when I increased the time steps too much. I was calculating a steadystate with a transient solver and as soon as the flow field is developed, you can go well beyond CFL>1. But if you overstretch, things start to oscillate like the turbulence fields, in which case OpenFOAM starts to bound the solutions for them. As far as I know, this happens only to turbulence fields and not for "normal" fields like velocity, temperature and so on. As the CFL number is in the range of ~1 in your case, I can only assume the boundary conditions of turbulence fields as the problem. You might check and compare them with tutorial cases of SST. Although I never simulated multiphase with OpenFOAM, so my knowledge is very limited there. Sorry. Best regards doerschlott 

October 27, 2020, 08:07 

#8 
New Member
Dhruv
Join Date: Jan 2020
Location: Blacksburg, Virginia
Posts: 6
Rep Power: 6 
Hello doerschlott
Thank you for helping to clear my confusion. Regarding the bounding problem, I agree with you. After your comment, I ran a small calculation on adjustable deltaT which the calculation showed around 5e6. I am currently running the calculation at deltaT of 1e6 and I do not see any bounding limits. I also think that having a high bound limit is fine to avoid a dividebyzero situation. Please correct me if I am wrong. Warm regards TurbulentGuy 

December 19, 2020, 06:50 

#9 
New Member
Join Date: Feb 2020
Location: Germany
Posts: 17
Rep Power: 6 
Hey TurbulentGuy, sorry for the late response
This is true to some degree, although the real problem is having omega at near zerovalues, as you divide by omega. Have a look into "kOmegaSSTBase.C", function "correctNut": Code:
// Correct the turbulence viscosity this>nut_ = a1_*k_/max(a1_*omega_, b1_*F23()*sqrt(S2)); Best regards doerschlott 

Thread Tools  Search this Thread 
Display Modes  


Similar Threads  
Thread  Thread Starter  Forum  Replies  Last Post 
Segmentation fault when using reactingFOAM for Fluids  Tommy Floessner  OpenFOAM Running, Solving & CFD  4  April 22, 2018 12:30 
chtMultiRegionSimpleFoam turbulent case  Aditya Patil  OpenFOAM Running, Solving & CFD  6  April 24, 2017 22:13 
Floating point exception error  lpz_michele  OpenFOAM Running, Solving & CFD  53  October 19, 2015 02:50 
Moving mesh  Niklas Wikstrom (Wikstrom)  OpenFOAM Running, Solving & CFD  122  June 15, 2014 06:20 
Orifice Plate with a fully developed flow  Problems with convergence  jonmec  OpenFOAM Running, Solving & CFD  3  July 28, 2011 05:24 