Simulation seems to converge but crashes suddenly
2 Attachment(s)
Hello,
I'm trying to run a simulation on an unstructured mesh. The simulation seems to converge but, after a thousand of iterations, it suddenly crashes. The strange things is that, if I restart it from just some iterations before the crash and without modify anything, it starts to converge again and crashes further beyond the previous crash point. The simulation is steady state, MRFSimpleFoam in OF2.1.1. The case is decomposed in 32 processors. I use cyclicAMI for periodic faces. RASmodel is kEpsilon. BCs are :  inlet : p=zeroGradient U=flowRateInletVelocity  outlet : p=fixedValue  all walls : U = fixedValue ( 0 0 0); I declare inlet, outlet and cyclicAMI as non rotating patches in MRFmodel. The checkMesh result is Mesh stats points: 24362 faces: 232715 internal faces: 211841 cells: 111139 boundary patches: 9 point zones: 0 face zones: 0 cell zones: 1 Overall number of cells of each type: hexahedra: 0 prisms: 0 wedges: 0 pyramids: 0 tet wedges: 0 tetrahedra: 111139 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: 1 (OK). Checking patch topology for multiply connected surfaces ... Patch Faces Points Surface topology inlet 353 206 ok (nonclosed singly connected) albero 434 250 ok (nonclosed singly connected) wall_inlet 542 305 ok (nonclosed singly connected) ciclico_ps 1878 1059 ok (nonclosed singly connected) dischi 8215 4469 ok (nonclosed singly connected) blade 2916 1606 ok (nonclosed singly connected) ciclico_ss 1928 1084 ok (nonclosed singly connected) top_bott_out 4036 2158 ok (nonclosed singly connected) outlet 572 356 ok (nonclosed singly connected) Checking geometry... Overall domain bounding box (0.00825 0.0229597 0.0123) (0.0594161 0.055631 0.0032) Mesh (nonempty, nonwedge) directions (1 1 1) Mesh (nonempty) directions (1 1 1) Boundary openness (9.90405e17 1.02286e16 2.27575e17) OK. Max cell openness = 2.85274e16 OK. Max aspect ratio = 28.7209 OK. Minumum face area = 6.8635e09. Maximum face area = 1.76476e06. Face area magnitudes OK. Min volume = 6.18733e13. Max volume = 6.85671e10. Total volume = 7.09264e06. Cell volumes OK. Mesh nonorthogonality Max: 86.2371 average: 17.0449 *Number of severely nonorthogonal faces: 203. Nonorthogonality check OK. <<Writing 203 nonorthogonal faces to set nonOrthoFaces Face pyramids OK. Max skewness = 1.07316 OK. Coupled point location match (average 0) OK. Mesh OK. I attached some pictures of residuals. You can see that residuals keeps steady or falling until they rapidly increase. When I restart the case from the latest "useful" (=not crashed) time, the simulation proceeds without problems for another thousand of iterations.By the first steps after restart, global results are much different from the ones before the crash. After each crash, results get closer to that calculated for the same geometry with structured mesh (and without any problem!). After the latest crash, I got almost the same results that I reached with structured mesh. But after some further iteration, the case crashed again. I would like to ask you if you have any ideas about why the simulation crashes. I also would like to know if I have to consider correct the results I get in a similar case. I mean that, in my opinion, if simulation crashes it means that results are not reliable. But at the same time, they look so close to the ones by structured mesh, witch I consider to be reliable ( and close to the ones by starCD ). Please, can someone help me? Thank in advance. 
Interesting....
Hello!
One thing that does not look good is the nonorthogonality of the mesh. Can you post your log file (atleast the residual log of the first 10 iterations) ? Use gradient limiters ...'cellLimited' for your grad scheme. This should make your simulation more stable. Akshay 
Thank you Akshay for the advise.
Now i cannot post the residual log but I will do it as soon as i can. Do you think the problem is the mesh? Actually I do think so but many other times meshes with a high number of nonorthogonal cells didn't give me such problems. Actually, i often find problems when i try to use unstructured mesh. I use Salome for meshing, but sometimes i would like to have a more powerful tool. May i ask you what you use? 
OpenFOAM is highly sensitive to the nonorthogonality of meshes. I use snappyHexMesh to create my meshes. It does a great job. Other meshes work as well(even if they are unstructured) as long as orthogonality is maintained.

Thank you again.
So you think the problem could be orthogonality. Sorry if I bother you, but how can I control nonorthogonality, especially with an unstructured mesh? Recently I tried to use snappyHexMesh too. Do you use it also for internal fluid dynamics? Because it seems to me that it works better for external fluid dynamics, but I admit I dont know it well. I have tried helixOS interface to simplify the use of snappyHexMesh. Thank in advance. 
Hello Akshay
I posted the residual log file. I tried to use cellLimited as grad scheme but OF 211 does not recognize it. The only limited scheme are: limitWith limitedCubic limitedCubicV limitedLinear limitedLinearV >>>>>>>>> Sorry, I'm not able to have the LOG file uploaded. I post you the first lines: Create time Create mesh for time = 34380 Reading field p Reading field U Reading/calculating face flux field phi AMI: Creating addressing and weights between 748 source faces and 758 target faces AMI: Patch source weights min/max/average = 0.999877, 1.00032, 1.00001 AMI: Patch target weights min/max/average = 0.999917, 1.0003, 1.00001 Selecting incompressible transport model Newtonian Selecting RAS turbulence model kEpsilon No field sources present SIMPLE: convergence criteria field p tolerance 1e07 field Urel tolerance 1e06 field nuTilda tolerance 1e06 field omega tolerance 1e06 field k tolerance 1e06 field epsilon tolerance 1e06 Starting time loop Time = 34381 DILUPBiCG: Solving for Ux, Initial residual = 0.000116484, Final residual = 1.95806e07, No Iterations 2 DILUPBiCG: Solving for Uy, Initial residual = 0.000190872, Final residual = 4.10964e07, No Iterations 2 DILUPBiCG: Solving for Uz, Initial residual = 0.00104798, Final residual = 1.49606e06, No Iterations 2 DICPCG: Solving for p, Initial residual = 0.00671232, Final residual = 6.54732e06, No Iterations 62 time step continuity errors : sum local = 0.00145536, global = 0.000191862, cumulative = 0.000191862 DILUPBiCG: Solving for epsilon, Initial residual = 3.7014e06, Final residual = 5.76453e08, No Iterations 1 DILUPBiCG: Solving for k, Initial residual = 2.98994e05, Final residual = 2.74708e08, No Iterations 2 ExecutionTime = 0.49 s ClockTime = 0 s Time = 34382 DILUPBiCG: Solving for Ux, Initial residual = 6.917e05, Final residual = 8.73941e08, No Iterations 2 DILUPBiCG: Solving for Uy, Initial residual = 0.000102618, Final residual = 2.82595e07, No Iterations 2 DILUPBiCG: Solving for Uz, Initial residual = 0.000612589, Final residual = 8.99979e07, No Iterations 2 DICPCG: Solving for p, Initial residual = 0.00271528, Final residual = 2.61552e06, No Iterations 123 time step continuity errors : sum local = 0.000581223, global = 3.52181e05, cumulative = 0.00022708 DILUPBiCG: Solving for epsilon, Initial residual = 3.77728e06, Final residual = 5.73571e08, No Iterations 1 DILUPBiCG: Solving for k, Initial residual = 2.99839e05, Final residual = 1.0709e08, No Iterations 2 ExecutionTime = 0.73 s ClockTime = 0 s Time = 34383 DILUPBiCG: Solving for Ux, Initial residual = 3.61533e05, Final residual = 2.63406e07, No Iterations 2 DILUPBiCG: Solving for Uy, Initial residual = 5.50816e05, Final residual = 1.01818e07, No Iterations 2 DILUPBiCG: Solving for Uz, Initial residual = 0.000340453, Final residual = 5.88818e07, No Iterations 2 DICPCG: Solving for p, Initial residual = 0.00283995, Final residual = 2.71022e06, No Iterations 30 time step continuity errors : sum local = 0.000601727, global = 8.04178e05, cumulative = 0.000146662 DILUPBiCG: Solving for epsilon, Initial residual = 4.42084e06, Final residual = 6.33828e08, No Iterations 1 bounding epsilon, min: 0.00650573 max: 1060.77 average: 10.5298 DILUPBiCG: Solving for k, Initial residual = 3.0096e05, Final residual = 5.41345e09, No Iterations 2 ExecutionTime = 0.87 s ClockTime = 0 s Time = 34384 DILUPBiCG: Solving for Ux, Initial residual = 2.26295e05, Final residual = 6.13416e08, No Iterations 2 DILUPBiCG: Solving for Uy, Initial residual = 3.24352e05, Final residual = 4.73366e08, No Iterations 2 DILUPBiCG: Solving for Uz, Initial residual = 0.00019446, Final residual = 1.83417e08, No Iterations 3 DICPCG: Solving for p, Initial residual = 0.00185662, Final residual = 1.73255e06, No Iterations 118 time step continuity errors : sum local = 0.000384432, global = 3.08246e05, cumulative = 0.000115837 DILUPBiCG: Solving for epsilon, Initial residual = 3.71761e06, Final residual = 2.17537e07, No Iterations 1 DILUPBiCG: Solving for k, Initial residual = 2.99508e05, Final residual = 7.43986e09, No Iterations 2 ExecutionTime = 1.09 s ClockTime = 1 s Time = 34385 DILUPBiCG: Solving for Ux, Initial residual = 1.64447e05, Final residual = 4.3641e08, No Iterations 2 DILUPBiCG: Solving for Uy, Initial residual = 2.18851e05, Final residual = 5.3845e08, No Iterations 2 DILUPBiCG: Solving for Uz, Initial residual = 0.000130193, Final residual = 1.85792e07, No Iterations 2 DICPCG: Solving for p, Initial residual = 0.0010382, Final residual = 9.90235e07, No Iterations 114 time step continuity errors : sum local = 0.000219681, global = 3.99351e05, cumulative = 7.59024e05 DILUPBiCG: Solving for epsilon, Initial residual = 3.57012e06, Final residual = 8.03344e08, No Iterations 1 DILUPBiCG: Solving for k, Initial residual = 2.98805e05, Final residual = 6.51274e09, No Iterations 2 ExecutionTime = 1.31 s ClockTime = 1 s Time = 34386 DILUPBiCG: Solving for Ux, Initial residual = 1.18153e05, Final residual = 3.32574e08, No Iterations 2 DILUPBiCG: Solving for Uy, Initial residual = 1.49857e05, Final residual = 8.53877e08, No Iterations 2 DILUPBiCG: Solving for Uz, Initial residual = 8.9257e05, Final residual = 2.14161e07, No Iterations 2 DICPCG: Solving for p, Initial residual = 0.00057172, Final residual = 5.15406e07, No Iterations 189 time step continuity errors : sum local = 0.00011434, global = 7.86249e07, cumulative = 7.51161e05 DILUPBiCG: Solving for epsilon, Initial residual = 3.54648e06, Final residual = 1.1617e07, No Iterations 1 DILUPBiCG: Solving for k, Initial residual = 2.98313e05, Final residual = 6.89284e09, No Iterations 2 ExecutionTime = 1.61 s ClockTime = 1 s Time = 34387 DILUPBiCG: Solving for Ux, Initial residual = 8.47096e06, Final residual = 2.02699e08, No Iterations 2 DILUPBiCG: Solving for Uy, Initial residual = 1.05874e05, Final residual = 7.0427e08, No Iterations 2 DILUPBiCG: Solving for Uz, Initial residual = 6.13598e05, Final residual = 3.56812e09, No Iterations 3 DICPCG: Solving for p, Initial residual = 0.000366482, Final residual = 3.65723e07, No Iterations 138 time step continuity errors : sum local = 8.11399e05, global = 7.36762e06, cumulative = 6.77485e05 DILUPBiCG: Solving for epsilon, Initial residual = 3.60755e06, Final residual = 7.63014e08, No Iterations 1 DILUPBiCG: Solving for k, Initial residual = 2.97936e05, Final residual = 8.57472e09, No Iterations 2 ExecutionTime = 1.87 s ClockTime = 1 s Time = 34388 DILUPBiCG: Solving for Ux, Initial residual = 6.20941e06, Final residual = 1.95689e08, No Iterations 2 DILUPBiCG: Solving for Uy, Initial residual = 7.69862e06, Final residual = 5.75046e08, No Iterations 2 DILUPBiCG: Solving for Uz, Initial residual = 4.28168e05, Final residual = 1.70379e07, No Iterations 2 DICPCG: Solving for p, Initial residual = 0.000244954, Final residual = 2.29733e07, No Iterations 180 time step continuity errors : sum local = 5.09679e05, global = 6.05531e07, cumulative = 6.8354e05 DILUPBiCG: Solving for epsilon, Initial residual = 3.87334e06, Final residual = 6.84854e08, No Iterations 1 DILUPBiCG: Solving for k, Initial residual = 2.97638e05, Final residual = 9.06926e09, No Iterations 2 ExecutionTime = 2.15 s ClockTime = 2 s Time = 34389 DILUPBiCG: Solving for Ux, Initial residual = 4.77713e06, Final residual = 1.43383e08, No Iterations 2 DILUPBiCG: Solving for Uy, Initial residual = 5.88675e06, Final residual = 1.4134e08, No Iterations 2 DILUPBiCG: Solving for Uz, Initial residual = 3.13716e05, Final residual = 1.14339e07, No Iterations 2 DICPCG: Solving for p, Initial residual = 0.000170138, Final residual = 1.67529e07, No Iterations 38 time step continuity errors : sum local = 3.71664e05, global = 1.77426e06, cumulative = 7.01283e05 DILUPBiCG: Solving for epsilon, Initial residual = 6.4676e06, Final residual = 8.48483e08, No Iterations 1 bounding epsilon, min: 0.0854313 max: 1060.77 average: 10.5295 DILUPBiCG: Solving for k, Initial residual = 2.981e05, Final residual = 1.25845e07, No Iterations 2 ExecutionTime = 2.3 s ClockTime = 2 s Time = 34390 DILUPBiCG: Solving for Ux, Initial residual = 3.80104e06, Final residual = 1.01539e08, No Iterations 2 DILUPBiCG: Solving for Uy, Initial residual = 4.65179e06, Final residual = 5.97925e09, No Iterations 2 DILUPBiCG: Solving for Uz, Initial residual = 2.38548e05, Final residual = 4.96177e08, No Iterations 2 DICPCG: Solving for p, Initial residual = 0.000118929, Final residual = 1.15022e07, No Iterations 108 time step continuity errors : sum local = 2.55176e05, global = 3.39004e06, cumulative = 6.67383e05 DILUPBiCG: Solving for epsilon, Initial residual = 4.35129e06, Final residual = 6.37832e08, No Iterations 1 DILUPBiCG: Solving for k, Initial residual = 2.97125e05, Final residual = 6.79241e09, No Iterations 2 ExecutionTime = 2.51 s ClockTime = 2 s Time = 34391 DILUPBiCG: Solving for Ux, Initial residual = 3.12327e06, Final residual = 8.26615e09, No Iterations 2 DILUPBiCG: Solving for Uy, Initial residual = 3.79e06, Final residual = 5.2095e09, No Iterations 2 DILUPBiCG: Solving for Uz, Initial residual = 1.87191e05, Final residual = 2.82144e08, No Iterations 2 DICPCG: Solving for p, Initial residual = 8.62629e05, Final residual = 9.76075e08, No Iterations 113 time step continuity errors : sum local = 2.16542e05, global = 1.94692e06, cumulative = 6.47914e05 DILUPBiCG: Solving for epsilon, Initial residual = 3.68029e06, Final residual = 6.11616e08, No Iterations 1 DILUPBiCG: Solving for k, Initial residual = 2.96539e05, Final residual = 7.78769e09, No Iterations 2 ExecutionTime = 2.73 s ClockTime = 2 s Time = 34392 DILUPBiCG: Solving for Ux, Initial residual = 2.63464e06, Final residual = 6.09912e09, No Iterations 2 DILUPBiCG: Solving for Uy, Initial residual = 3.18131e06, Final residual = 4.0811e09, No Iterations 2 DILUPBiCG: Solving for Uz, Initial residual = 1.50655e05, Final residual = 5.82442e08, No Iterations 2 DICPCG: Solving for p, Initial residual = 6.69605e05, Final residual = 9.75085e08, No Iterations 17 time step continuity errors : sum local = 2.16325e05, global = 5.92714e07, cumulative = 6.41986e05 DILUPBiCG: Solving for epsilon, Initial residual = 3.46206e06, Final residual = 5.47954e08, No Iterations 1 DILUPBiCG: Solving for k, Initial residual = 2.96205e05, Final residual = 8.97014e09, No Iterations 2 ExecutionTime = 2.86 s ClockTime = 2 s Time = 34393 DILUPBiCG: Solving for Ux, Initial residual = 2.259e06, Final residual = 5.01807e09, No Iterations 2 DILUPBiCG: Solving for Uy, Initial residual = 2.72394e06, Final residual = 5.76879e09, No Iterations 2 DILUPBiCG: Solving for Uz, Initial residual = 1.23467e05, Final residual = 4.15094e08, No Iterations 2 DICPCG: Solving for p, Initial residual = 5.2557e05, Final residual = 9.71596e08, No Iterations 178 time step continuity errors : sum local = 2.15551e05, global = 7.09067e07, cumulative = 6.49077e05 DILUPBiCG: Solving for epsilon, Initial residual = 3.35479e06, Final residual = 5.28543e08, No Iterations 1 DILUPBiCG: Solving for k, Initial residual = 2.95863e05, Final residual = 8.85696e09, No Iterations 2 ExecutionTime = 3.14 s ClockTime = 3 s Time = 34394 DILUPBiCG: Solving for Ux, Initial residual = 1.96087e06, Final residual = 3.74087e09, No Iterations 2 DILUPBiCG: Solving for Uy, Initial residual = 2.36285e06, Final residual = 3.33676e09, No Iterations 2 DILUPBiCG: Solving for Uz, Initial residual = 1.0291e05, Final residual = 2.50191e08, No Iterations 2 DICPCG: Solving for p, Initial residual = 4.31965e05, Final residual = 9.70765e08, No Iterations 19 time step continuity errors : sum local = 2.15365e05, global = 3.91339e06, cumulative = 6.88211e05 DILUPBiCG: Solving for epsilon, Initial residual = 3.29362e06, Final residual = 5.36464e08, No Iterations 1 DILUPBiCG: Solving for k, Initial residual = 2.95555e05, Final residual = 8.68175e09, No Iterations 2 ExecutionTime = 3.27 s ClockTime = 3 s limiterBlended What scheme did you mean? Is it in OF 220? 
Hey
try this... gradSchemes { default cellLimited leastSquares 1.0; } And... You did mention you were using an unstructured grid right? 
Quote:
1. For the fvScheme file in your case check the following post: http://www.cfdonline.com/Forums/ope...orrectors.html I hope that will solve your problem regarding scheme implementation. If you are not getting proper convergence try this one too. 2. Since, your mesh is "severely nonorthogonal faces: 203". I would rather run a transient run instead of SIMPLE and with transient case try increasing the nCorrectors and nNonOrthogonalCorrectors in fvSolution file. Do, try these steps and do let us know. 
xxxx: You've got all Tets and a total count of 111139. You are decomposing this size into 32 partitions!! that is an overkill!!
Tushar: I don't think he needs to bring in more complexities through a transient run. MRFSimpleFoam should do the job. NonOrtho of 86 can be handled through some stable solver settings. I've run cases with nonortho of 87 successfully. And... xxxxx: Please use another name here and not 'xxxxx' ;) 
Hello,
One more reason why the solver had convergence problem. If you see the checkMesh output the max skewness is: Code:
Max skewness = 1.07316 OK. :) 
Tushar: Skewness of 1.07 is as low as it can get. Above 5 and i'll call it a MINOR problem ;)

Hello everybody.
Akshay: thank you. I will try cellLimited soon. The case was partitioned between 32 processors to see how it works ( I have just started to use multiprocessors and I still don't know exactly how it works).. You can call me 'zzzzz' if you prefer :) . In the meantime I talked with a person more expert than me in CFD, but not in OF. He supposes that, since crash occurs in only one step, this could be caused by errors in reading/writing in RAM. Tushar: I'm not an expert but skewness as low as 1.07 should not be a problem. 
Well,
Akshay Quote:
Quote:
http://www.tchpc.tcd.ie/fluent/Unpac...ug/node318.htm http://www.cfdonline.com/Forums/flu...an098a.html http://en.wikipedia.org/wiki/Types_o...uiangular_skew You, will realize the actual problem comes when skewness exceeds 1 (>1). Also, when mesh is distorted then the error induced by convection term increases significantly which causes unboundedness. Now, if suitable scheme is not used then the solver crashes. I hope you guys got my point. :) 
Tushar:
OpenFOAM has its very own definition of skewness. What it shows in checkMesh is based on its definition of skewness. In OF skewness above 4 is usually called high skewness. Please don't confuse fluent with openfoam. They are very different! 
Quote:
I know one can solve with higher skewness in OpenFOAM. But, for that you need to find the correct scheme. May be, some more efforts to skewcorrect the solution. That I think was a reason behind the divergence of the mentioned solution. What do you think solving is more important ? or, solution accuracy? with respect to the actual results. When the user uses higher skewness it incorporates more chances of errors. 
2 Attachment(s)
Akshay, may I ask you another thing?
Following your advise I started to learn how to use snappyHexMesh. I still have some doubts. First of all, I saw that when I use snappyHexMesh (with snappyHex disctionary taken from simpleFoam/turbine tutorial and modified to adapt to a simple geometry) OF creates to folders, named 1 and 2, each containig a polyMesh folder. I think this are, let me say, a first trial mesh and a refined one, rigth? Than only the mesh in folder "2" is copied to constant. If I perform checkMesh I see that mesh at time 1 seems to be better that at time 2 ( smaller skewness, and so on..). If I copy polyMesh folder from time 1 to constant and open paraview, the mesh doesn't look much different from the one at time 2. Am I wrong? In this case, should I use the mesh of folder 1? Another question is, can you please explain briefly what the two numbers in " levels ((2 4)); " in snappyHexMeshDict mean? I know they are related to mesh levels and using them one can control the resulting mesh, but how do they actually work? The last question is about mesh quality. It seems to me that snappyHexMesh produces some holes in the surface. I highlighted some of these in the attached picture.You can also see there's a sort of bubble in the foreground (pointed by an arrow). The colours are cell normals, so different colours mean in neighbour cells mean steeply variations in orientation. Can this problem cause the crash of the case? How can I remove it? Thank in advance. 
Quote:
Please, excuse me for my earlier comment it was too rude. For few cases skewness of 14 are considered OK with OpenFOAM. Quote:
 Best Regards! 
All times are GMT 4. The time now is 17:32. 