
[Sponsors] 
Simulation of mixing tank (RSTM) diverges suddenly 

LinkBack  Thread Tools  Search this Thread  Display Modes 
August 19, 2019, 09:22 
Simulation of mixing tank (RSTM) diverges suddenly

#1 
New Member
Mazen Draw
Join Date: Sep 2018
Posts: 21
Rep Power: 8 
Dear all,
I am running a simulation of a mixing tank on OF 7 where the only boundary condition that I have is Wall. No Inlet nor Outlet. I have a rotating zone for the impeller and stationary zone for the tank with the baffles. The geometry contains only half of the tank, due to the symmetry of the flow. The patches, where the tank is cut in half, communicate with each other via rotational cyclicAMI BC. The patches, where the rotating zone is connected to the stationary zone, communicate via a noOrdering cyclicAMI BC. here is a checkMesh output: Code:
/**\ =========  \\ / F ield  OpenFOAM: The Open Source CFD Toolbox \\ / O peration  Website: https://openfoam.org \\ / A nd  Version: 7 \\/ M anipulation  \**/ Build : 7109ba3c8d53a Exec : checkMesh Date : Aug 19 2019 Time : 13:46:04 Host : "fwd250" PID : 12968 I/O : uncollated Case : /mnt/sshfs/bigdata/openfoam/draw33/tutorials/multiphase/reactingTwoPhaseEulerFoam/laminar/SSG/2007_Montante_onePhase_steadySSG nProcs : 1 sigFpe : Enabling floating point exception trapping (FOAM_SIGFPE). fileModificationChecking : Monitoring runtime modified files using timeStampMaster (fileModificationSkew 10) allowSystemOperations : Allowing usersupplied system call operations // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // Create time Create polyMesh for time = 0 Time = 0 Mesh stats points: 1134003 faces: 3298814 internal faces: 3197032 cells: 1082641 faces per cell: 6 boundary patches: 18 point zones: 0 face zones: 2 cell zones: 2 Overall number of cells of each type: hexahedra: 1082641 prisms: 0 wedges: 0 pyramids: 0 tet wedges: 0 tetrahedra: 0 polyhedra: 0 Checking topology... Boundary definition OK. Cell to face addressing OK. Point usage OK. Upper triangular ordering OK. Face vertices OK. *Number of regions: 2 The mesh has multiple regions which are not connected by any face. <<Writing region information to "0/cellToRegion" <<Writing region 0 with 403920 cells to cellSet region0 <<Writing region 1 with 678721 cells to cellSet region1 Checking patch topology for multiply connected surfaces... Patch Faces Points Surface topology IMPELLERSIDE_2 11880 12103 ok (nonclosed singly connected) BAFFLES 8712 9044 ok (nonclosed singly connected) TANKWALL 11088 11571 ok (nonclosed singly connected) TANKPERIODIC2 4620 4788 ok (nonclosed singly connected) TANKPERIODIC1 4620 4788 ok (nonclosed singly connected) TANKBOTTOM 3060 3216 ok (nonclosed singly connected) TANKTOP 3060 3216 ok (nonclosed singly connected) IMPELLERSIDE_1 11880 12103 ok (nonclosed singly connected) IMPELLERBOTTOM 5108 5243 ok (nonclosed singly connected) IMPELLERTOP 5108 5243 ok (nonclosed singly connected) IMPELLERPERIODIC2 9435 9690 ok (nonclosed singly connected) IMPELLERPERIODIC1 9435 9690 ok (nonclosed singly connected) IMPELLER3 1359 1392 ok (nonclosed singly connected) IMPELLER2 1359 1392 ok (nonclosed singly connected) IMPELLER1 1359 1392 ok (nonclosed singly connected) SHAFT 8000 8270 ok (nonclosed singly connected) SPARGERSURFACE 1200 1274 ok (nonclosed singly connected) INLET 499 548 ok (nonclosed singly connected) Checking geometry... Overall domain bounding box (0.118 0.118 0) (0.118 1.44504e17 0.236) Mesh has 3 geometric (nonempty/wedge) directions (1 1 1) Mesh has 3 solution (nonempty) directions (1 1 1) Boundary openness (1.3176e16 1.39255e15 4.29463e17) OK. Max cell openness = 4.42129e16 OK. Max aspect ratio = 29.1252 OK. Minimum face area = 6.71041e08. Maximum face area = 7.10408e05. Face area magnitudes OK. Min volume = 3.7373e11. Max volume = 1.51012e07. Total volume = 0.00510649. Cell volumes OK. Mesh nonorthogonality Max: 45.2912 average: 10.76 Nonorthogonality check OK. Face pyramids OK. Max skewness = 0.920013 OK. Coupled point location match (average 0) OK. Mesh OK. End First I ran the simulation with kOmegaSST. I did not get the best convergence for the pressure but at least the simulation was stable. When I am running the SSG model, the simulation first runs normally then suddenly diverges and crashes. I noticed a bounding epsilon message which then leads to a very big time step continuity error. This leads afterwards to very small residual in velocity (1e19 e.g.) and big one in pressure (1 e.g.). These residuals then oscillates then the whole simulation stops. I started then reading about these kind of issues in the forum. All changes of numerical schemes in fvSchemes or the solvers in fvSolution did not really help too much. fvSchemes: Code:
/** C++ **\  =========    \\ / F ield  cfMesh: A library for mesh generation   \\ / O peration    \\ / A nd  Author: Franjo Juretic   \\/ M anipulation  Email: franjo.juretic@cfields.com  \**/ FoamFile { version 2.0; format ascii; class dictionary; object fvSchemes; } // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // ddtSchemes { default steadyState; } gradSchemes { default Gauss linear; } divSchemes { default none; div(phi,U) bounded Gauss limitedLinearV 1; div(phi,k) Gauss upwind; //bounded Gauss limitedLinear 1; div(phi,epsilon) Gauss upwind; //bounded Gauss limitedLinear 1; div(phi,omega) bounded Gauss limitedLinear 1; div(phi,R) bounded Gauss limitedLinear 1; //upwind; div((nuEff*dev2(T(grad(U))))) Gauss linear; div((nu*dev2(T(grad(U))))) Gauss linear; div(R) Gauss linear; } laplacianSchemes { default Gauss linear uncorrected; } interpolationSchemes { default linear; } snGradSchemes { default uncorrected; } fluxRequired { default no; p; } wallDist { method Poisson; nRequired true; } // ************************************************************************* // Code:
/** C++ **\  =========    \\ / F ield  cfMesh: A library for mesh generation   \\ / O peration    \\ / A nd  Author: Franjo Juretic   \\/ M anipulation  Email: franjo.juretic@cfields.com  \**/ FoamFile { version 2.0; format ascii; class dictionary; object fvSolution; } // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // solvers { "(kepsilonomegaR).*" { solver PBiCG; //solver smoothSolver; preconditioner DILU; tolerance 1e8; relTol 0.001; minIter 1; } yPsi { solver PCG; preconditioner none; tolerance 1e10; relTol 0; } p { solver GAMG; tolerance 1e8; relTol 1e8; minIter 5; maxIter 100; smoother GaussSeidel; nPreSweeps 1; nPostSweeps 3; nFinestSweeps 3; scaleCorrection true; directSolveCoarsest false; cacheAgglomeration on; nCellsInCoarsestLevel 1000; agglomerator faceAreaPair; mergeLevels 1; minIter 2; }; U { solver PBiCG; preconditioner DILU; tolerance 1e09; relTol 0; }; } SIMPLE { nNonOrthogonalCorrectors 4; pRefCell 0; pRefValue 0; residualControl { p 1e5; U 1e5; "(kepsilonomegaR)" 1e5; } } relaxationFactors { fields { p 0.3; } equations { U 0.5; k 0.3; epsilon 0.3; omega 0.3; R 0.3; } } PISO { nCorrectors 5; nNonOrthogonalCorrectors 2; pRefCell 0; pRefValue 0; } // ************************************************************************* // So I started playing around with the relaxation factors. When setting the under relaxation factor of velocity and epsilon equal to 0.001 and reynolds stresses to 0.0001, the simulation is much more stable. The residual of reynolds stresses still increase in value but much more slower than before. After running the simulation with these relaxation factors for a while, the pressure residual starts to oscillates around 0.0002 and continues like this. Code:
relaxationFactors { fields { p 0.3; } equations { U 0.001; k 0.3; epsilon 0.001; omega 0.3; R 0.0001; } } I want to map the best converged solution of the steadystate case to the transient case and start again. When I run the transient case now, I am getting oscillations faster and the simulation crashes also quicker than in steadystate case. 

August 19, 2019, 09:57 

#2 
New Member
Mazen Draw
Join Date: Sep 2018
Posts: 21
Rep Power: 8 
BTW I am using simpleFoam solver and this is a one phase simulation. The transient case will be then solved with pimpleFoam solver. Afterwards, the whole model will be extended to 2 phase and then 3 phase simulation using reactingEulerFoam.


July 20, 2022, 14:45 

#3 
Senior Member
Join Date: Dec 2019
Location: Cologne, Germany
Posts: 367
Rep Power: 8 
Hi,
Did you solve your problem? I am facing the same kind of problem, I am simulating the whole tank with MRF, Kepsilon works fine, although the power number is 20% lower, i can achieve higher values for komegasst but the presidual stalls at the order 1e2, which is really bad. 

August 22, 2022, 06:52 

#4 
New Member
Mazen Draw
Join Date: Sep 2018
Posts: 21
Rep Power: 8 
Hello Geth,
This has been a while ago, but AFAIR simpleFoam wasn't able to obtain the steadystate solution with "small" residuals. Therefore, I switched to pimpleFoam using Local Time Stepping. You can do so by selecting localEuler as the time derivative scheme in fvSchemes. I believe the reason why simpleFoam doesn't perform well is that there is no steadystate solution for the flow per se. Instead the flow variables are oscillating around a statistical average value. Therefore, using a transient solver like pimpleFoam even in the pseudotransient mode (when using localEuler) will yield smaller residuals. Honestly if I have to do it all over again, I would probably just use pimpleFoam in Piso Mode (1 nOuterCorrectors with 2 or 3 nCorrectors) until the statistical steadystate solution is approached. Then I would switch to Pimple Mode (3 nOuterCorrectors with 1 or 2 nCorrectors) and time average the solution over a reasonable time interval. This should give you a steadystate solution with "small" residuals. 

February 4, 2024, 06:28 

#5 
Member
ASR
Join Date: Jan 2023
Location: Leeds, UK
Posts: 52
Rep Power: 3 
I've also simulated a mixing tank running a single phase simulation for SST KOmega. I found the residuals for p stalled at 10^2. It was due to BC. The solution extremely sensitive to BC.
For highRe use nutkWallFunction (y+ > 30); For scalable wall functions (1< y+ < 300) use nutUWallFunction or nutUSpaldingWallFunction; For lowRe use nutLowReWallFunction (y+ < 1). 

Tags 
bounding epsilon, bounding error, divergence ssg, mixing tank, ssg 
Thread Tools  Search this Thread 
Display Modes  


Similar Threads  
Thread  Thread Starter  Forum  Replies  Last Post 
Continuous mixing flow simulation in a Tank  Meca_Ali  FLUENT  7  March 14, 2018 06:50 
3D mixing tank using MRF/SRF  cfd_confuse  Main CFD Forum  1  September 13, 2016 03:43 
Simulation of thermocline in a tank  fkhan7  CFX  25  October 22, 2015 12:23 
Meshing industrial scale mixing tank  mlbontbs87  FLUENT  0  April 26, 2011 16:18 
3D mixing tank using MRF/SRF  cfd_confuse  FLUENT  0  October 2, 2010 04:46 