I have created a moving mesh s
I have created a moving mesh solver by combining icoFoamAutoMotion and interFoam. The 2d motion of the cylinder is specified by movingFlapFvMesh, xampl=1.0, xfreq=0.1, yampl=1.0, yfreq=0.1 and thetaampl=0.0. The inlet velocity it 1 m/s and corresponds to Re=150.
The timestep is 0.005. Up to time=1.225 nice covergence is obtained, but at time=1.225 the solution begins to diverge, as can be seen from the Courant numbers shown below. Furthermore the number of iterations of solving motionUx and motionUy is increasing every iteration. Therefore I think that the motionSolver causes this problem, but I have no clue why. Any suggestions? Thanks. Greetings, Frank ================================================== Output ================================================== Mean and max Courant Numbers = 0.023978 0.948195 Time = 1.225 Correct mesh motion diffusion field. ICCG: Solving for motionUx, Initial residual = 0.00074787, Final residual = 9.94298e-06, No Iterations 280 ICCG: Solving for motionUy, Initial residual = 0.000497424, Final residual = 9.96133e-06, No Iterations 187 Correcting 2-D mesh motion ...done volume continuity errors : sum local = 4.06953e-15, global = -1.90209e-18 BICCG: Solving for Ux, Initial residual = 0.000139335, Final residual = 3.63493e-06, No Iterations 1 BICCG: Solving for Uy, Initial residual = 0.00080448, Final residual = 4.75049e-07, No Iterations 2 ICCG: Solving for p, Initial residual = 0.84611, Final residual = 9.43944e-07, No Iterations 290 time step continuity errors : sum local = 4.11744e-11, global = 1.12633e-14, cumulative = 6.75289e-10 ICCG: Solving for p, Initial residual = 0.859509, Final residual = 9.95797e-07, No Iterations 290 time step continuity errors : sum local = 4.79839e-11, global = 6.3742e-15, cumulative = 6.75296e-10 Wall patch = 3 Wall patch name = cyl_wall Uav = (0.282238 -0.717762 5.26235e-28) Aref = 1.44219 Lref = 0.564692 DragCoefficient = 5.82549 pressureDragCoefficient = 4.83947 viscDragCoefficient = 0.986017 LiftCoefficient = 2.19606 MomentCoefficient = 8.31774 ExecutionTime = 9948.97 s Mean and max Courant Numbers = 0.0239854 1.13739 Time = 1.23 Correct mesh motion diffusion field. ICCG: Solving for motionUx, Initial residual = 0.000752869, Final residual = 9.99251e-06, No Iterations 281 ICCG: Solving for motionUy, Initial residual = 0.000500919, Final residual = 9.95567e-06, No Iterations 187 Correcting 2-D mesh motion ...done volume continuity errors : sum local = 4.09456e-15, global = 2.59593e-18 BICCG: Solving for Ux, Initial residual = 0.000173378, Final residual = 4.67296e-06, No Iterations 1 BICCG: Solving for Uy, Initial residual = 0.00112283, Final residual = 7.11956e-07, No Iterations 2 ICCG: Solving for p, Initial residual = 0.890628, Final residual = 9.45922e-07, No Iterations 295 time step continuity errors : sum local = 5.84773e-11, global = -7.53318e-14, cumulative = 6.7522e-10 ICCG: Solving for p, Initial residual = 0.901723, Final residual = 9.32168e-07, No Iterations 296 time step continuity errors : sum local = 6.38908e-11, global = 1.57041e-13, cumulative = 6.75377e-10 Wall patch = 3 Wall patch name = cyl_wall Uav = (0.284429 -0.715571 5.79729e-29) Aref = 1.44221 Lref = 0.559905 DragCoefficient = 6.26111 pressureDragCoefficient = 5.26338 viscDragCoefficient = 0.997727 LiftCoefficient = 2.30624 MomentCoefficient = 10.1442 ExecutionTime = 10008.3 s Mean and max Courant Numbers = 0.0240167 1.89597 Time = 1.235 Correct mesh motion diffusion field. ICCG: Solving for motionUx, Initial residual = 0.000757996, Final residual = 9.99529e-06, No Iterations 281 ICCG: Solving for motionUy, Initial residual = 0.000504364, Final residual = 9.90425e-06, No Iterations 188 Correcting 2-D mesh motion ...done volume continuity errors : sum local = 4.10338e-15, global = 1.22539e-18 BICCG: Solving for Ux, Initial residual = 0.000229742, Final residual = 6.51542e-06, No Iterations 1 BICCG: Solving for Uy, Initial residual = 0.00160232, Final residual = 1.13808e-06, No Iterations 2 ICCG: Solving for p, Initial residual = 0.923514, Final residual = 9.85298e-07, No Iterations 299 time step continuity errors : sum local = 8.725e-11, global = -8.65098e-13, cumulative = 6.74512e-10 ICCG: Solving for p, Initial residual = 0.932281, Final residual = 9.66381e-07, No Iterations 301 time step continuity errors : sum local = 9.52416e-11, global = 8.27088e-13, cumulative = 6.75339e-10 Wall patch = 3 Wall patch name = cyl_wall Uav = (0.286627 -0.713373 -7.71628e-28) Aref = 1.44222 Lref = 0.55508 DragCoefficient = 7.31698 pressureDragCoefficient = 6.28853 viscDragCoefficient = 1.02845 LiftCoefficient = 2.55426 MomentCoefficient = 14.3745 ExecutionTime = 10067.6 s Mean and max Courant Numbers = 0.0240893 3.38602 Time = 1.24 Correct mesh motion diffusion field. ICCG: Solving for motionUx, Initial residual = 0.000763107, Final residual = 9.99048e-06, No Iterations 280 ICCG: Solving for motionUy, Initial residual = 0.000507756, Final residual = 9.93261e-06, No Iterations 188 Correcting 2-D mesh motion ...done volume continuity errors : sum local = 4.08207e-15, global = -1.90597e-18 BICCG: Solving for Ux, Initial residual = 0.000328665, Final residual = 1.3123e-06, No Iterations 2 BICCG: Solving for Uy, Initial residual = 0.00234532, Final residual = 2.49132e-06, No Iterations 2 ICCG: Solving for p, Initial residual = 0.945383, Final residual = 9.50596e-07, No Iterations 306 time step continuity errors : sum local = 1.20953e-10, global = -1.79658e-13, cumulative = 6.7516e-10 ICCG: Solving for p, Initial residual = 0.953924, Final residual = 9.49527e-07, No Iterations 308 time step continuity errors : sum local = 1.33918e-10, global = 5.25042e-14, cumulative = 6.75212e-10 Wall patch = 3 Wall patch name = cyl_wall Uav = (0.288832 -0.711168 5.26306e-28) Aref = 1.44221 Lref = 0.550215 DragCoefficient = 9.85785 pressureDragCoefficient = 8.74477 viscDragCoefficient = 1.11308 LiftCoefficient = 3.1429 MomentCoefficient = 24.2969 ExecutionTime = 10126.9 s Mean and max Courant Numbers = 0.0242311 6.68236 Time = 1.245 Correct mesh motion diffusion field. ICCG: Solving for motionUx, Initial residual = 0.000768247, Final residual = 9.8493e-06, No Iterations 282 ICCG: Solving for motionUy, Initial residual = 0.000511277, Final residual = 9.82936e-06, No Iterations 193 Correcting 2-D mesh motion ...done volume continuity errors : sum local = 4.09053e-15, global = 8.99978e-19 BICCG: Solving for Ux, Initial residual = 0.000517687, Final residual = 2.34858e-06, No Iterations 2 BICCG: Solving for Uy, Initial residual = 0.00348416, Final residual = 5.18035e-06, No Iterations 3 ICCG: Solving for p, Initial residual = 0.955581, Final residual = 9.54336e-07, No Iterations 312 time step continuity errors : sum local = 1.69908e-10, global = 1.25309e-12, cumulative = 6.76465e-10 ICCG: Solving for p, Initial residual = 0.96498, Final residual = 9.80104e-07, No Iterations 314 time step continuity errors : sum local = 1.91259e-10, global = -1.58831e-12, cumulative = 6.74877e-10 Wall patch = 3 Wall patch name = cyl_wall Uav = (0.291044 -0.708956 7.71984e-28) Aref = 1.44217 Lref = 0.545312 DragCoefficient = 15.0292 pressureDragCoefficient = 13.489 viscDragCoefficient = 1.54019 LiftCoefficient = 4.10684 MomentCoefficient = 48.2943 ExecutionTime = 10187 s Mean and max Courant Numbers = 0.0244966 14.8136 Time = 1.25 Correct mesh motion diffusion field. ICCG: Solving for motionUx, Initial residual = 0.000773288, Final residual = 9.89781e-06, No Iterations 282 ICCG: Solving for motionUy, Initial residual = 0.000514758, Final residual = 9.9967e-06, No Iterations 188 Correcting 2-D mesh motion ...done volume continuity errors : sum local = 4.04741e-15, global = 8.71871e-19 BICCG: Solving for Ux, Initial residual = 0.000889756, Final residual = 2.69979e-06, No Iterations 3 BICCG: Solving for Uy, Initial residual = 0.00519126, Final residual = 1.8881e-06, No Iterations 5 ICCG: Solving for p, Initial residual = 0.950479, Final residual = 9.4596e-07, No Iterations 318 time step continuity errors : sum local = 2.22139e-10, global = 1.45724e-12, cumulative = 6.76334e-10 ICCG: Solving for p, Initial residual = 0.957963, Final residual = 8.6076e-07, No Iterations 320 time step continuity errors : sum local = 2.21044e-10, global = -1.29554e-12, cumulative = 6.75039e-10 Wall patch = 3 Wall patch name = cyl_wall Uav = (0.293263 -0.706737 -1.80998e-28) Aref = 1.44214 Lref = 0.54037 DragCoefficient = 26.3516 pressureDragCoefficient = 23.0254 viscDragCoefficient = 3.32619 LiftCoefficient = 6.17307 MomentCoefficient = 104.607 ExecutionTime = 10246.6 s Mean and max Courant Numbers = 0.0249365 40.0673 Time = 1.255 Correct mesh motion diffusion field. ICCG: Solving for motionUx, Initial residual = 0.000778549, Final residual = 9.97842e-06, No Iterations 282 ICCG: Solving for motionUy, Initial residual = 0.000518411, Final residual = 9.89038e-06, No Iterations 193 Correcting 2-D mesh motion ...done volume continuity errors : sum local = 4.08214e-15, global = 1.5051e-18 BICCG: Solving for Ux, Initial residual = 0.00178422, Final residual = 3.98896e-06, No Iterations 7 BICCG: Solving for Uy, Initial residual = 0.00891935, Final residual = 7.66749e-06, No Iterations 8 ICCG: Solving for p, Initial residual = 0.924408, Final residual = 9.65802e-07, No Iterations 322 time step continuity errors : sum local = 2.86659e-10, global = 4.15199e-13, cumulative = 6.75454e-10 ICCG: Solving for p, Initial residual = 0.940882, Final residual = 9.8166e-07, No Iterations 322 time step continuity errors : sum local = 3.33021e-10, global = -5.05583e-13, cumulative = 6.74948e-10 Wall patch = 3 Wall patch name = cyl_wall Uav = (0.295489 -0.704511 2.74279e-28) Aref = 1.44219 Lref = 0.53539 DragCoefficient = 61.6947 pressureDragCoefficient = 52.343 viscDragCoefficient = 9.35167 LiftCoefficient = 18.4061 MomentCoefficient = 201.47 ExecutionTime = 10307 s Mean and max Courant Numbers = 0.0259063 127.598 Time = 1.26 Correct mesh motion diffusion field. ICCG: Solving for motionUx, Initial residual = 0.000783877, Final residual = 9.81292e-06, No Iterations 284 ICCG: Solving for motionUy, Initial residual = 0.000521931, Final residual = 9.95579e-06, No Iterations 194 Correcting 2-D mesh motion ...done volume continuity errors : sum local = 4.0859e-15, global = -2.14512e-18 BICCG: Solving for Ux, Initial residual = 0.00572169, Final residual = 7.68117e-06, No Iterations 12 BICCG: Solving for Uy, Initial residual = 0.0344701, Final residual = 7.918e-06, No Iterations 14 ICCG: Solving for p, Initial residual = 0.863708, Final residual = 9.84167e-07, No Iterations 324 time step continuity errors : sum local = 3.85863e-10, global = -8.52682e-13, cumulative = 6.74096e-10 ICCG: Solving for p, Initial residual = 0.89722, Final residual = 9.61545e-07, No Iterations 323 time step continuity errors : sum local = 4.96108e-10, global = -1.1644e-12, cumulative = 6.72931e-10 Wall patch = 3 Wall patch name = cyl_wall Uav = (0.297722 -0.702278 4.16906e-28) Aref = 1.44221 Lref = 0.53037 DragCoefficient = 194.264 pressureDragCoefficient = 178.762 viscDragCoefficient = 15.5021 LiftCoefficient = 68.4224 MomentCoefficient = 491.493 ExecutionTime = 10368 s Mean and max Courant Numbers = 0.0287663 420.85 Time = 1.265 Correct mesh motion diffusion field. ICCG: Solving for motionUx, Initial residual = 0.000788983, Final residual = 9.7713e-06, No Iterations 285 ICCG: Solving for motionUy, Initial residual = 0.00052554, Final residual = 9.99632e-06, No Iterations 194 Correcting 2-D mesh motion ...done volume continuity errors : sum local = 4.07999e-15, global = 4.6044e-19 BICCG: Solving for Ux, Initial residual = 0.0340294, Final residual = 3.50395e-06, No Iterations 20 BICCG: Solving for Uy, Initial residual = 0.290113, Final residual = 1.78799e-06, No Iterations 22 ICCG: Solving for p, Initial residual = 0.691357, Final residual = 9.72633e-07, No Iterations 380 time step continuity errors : sum local = 1.24537e-09, global = 6.79261e-12, cumulative = 6.79724e-10 ICCG: Solving for p, Initial residual = 0.997215, Final residual = 8.99401e-07, No Iterations 393 time step continuity errors : sum local = 1.00642e-07, global = 4.05372e-10, cumulative = 1.0851e-09 Wall patch = 3 Wall patch name = cyl_wall Uav = (0.299962 -0.700038 2.34445e-28) Aref = 1.44222 Lref = 0.525311 DragCoefficient = 3812.82 pressureDragCoefficient = 916.412 viscDragCoefficient = 2896.41 LiftCoefficient = 9426.75 MomentCoefficient = -102740 ExecutionTime = 10430 s Mean and max Courant Numbers = 0.571589 111497 Time = 1.27 Correct mesh motion diffusion field. ICCG: Solving for motionUx, Initial residual = 0.000794283, Final residual = 9.78199e-06, No Iterations 285 ICCG: Solving for motionUy, Initial residual = 0.000529197, Final residual = 9.90328e-06, No Iterations 196 Correcting 2-D mesh motion ...done volume continuity errors : sum local = 4.07254e-15, global = 1.28338e-19 BICCG: Solving for Ux, Initial residual = 0.999564, Final residual = 3.97533e-06, No Iterations 164 BICCG: Solving for Uy, Initial residual = 0.993602, Final residual = 7.81896e-06, No Iterations 160 ICCG: Solving for p, Initial residual = 0.45766, Final residual = 6.70541e-07, No Iterations 664 time step continuity errors : sum local = 1.8692e-08, global = -3.35061e-11, cumulative = 1.05159e-09 ICCG: Solving for p, Initial residual = 0.998543, Final residual = 9.70955e-07, No Iterations 750 time step continuity errors : sum local = 1.49225e-05, global = -3.6547e-08, cumulative = -3.54954e-08 Wall patch = 3 Wall patch name = cyl_wall Uav = (0.302209 -0.697791 -3.63624e-28) Aref = 1.4422 Lref = 0.520213 DragCoefficient = -2.70433e+07 pressureDragCoefficient = -2.69181e+07 viscDragCoefficient = -125286 LiftCoefficient = -6.50234e+06 MomentCoefficient = -1.03616e+08 ExecutionTime = 10500.2 s Mean and max Courant Numbers = 426.603 3.37259e+07 Time = 1.275 Correct mesh motion diffusion field. ICCG: Solving for motionUx, Initial residual = 0.000799647, Final residual = 9.80335e-06, No Iterations 285 ICCG: Solving for motionUy, Initial residual = 0.00053282, Final residual = 9.96991e-06, No Iterations 196 Correcting 2-D mesh motion ...done volume continuity errors : sum local = 4.08114e-15, global = -1.70846e-18 BICCG: Solving for Ux, Initial residual = 0.999992, Final residual = 3.38725, No Iterations 1001 BICCG: Solving for Uy, Initial residual = 0.999993, Final residual = 1.81754, No Iterations 1001 ICCG: Solving for p, Initial residual = 1, Final residual = 0.423881, No Iterations 5001 time step continuity errors : sum local = 5.22248e+06, global = -2.21042, cumulative = -2.21042 ICCG: Solving for p, Initial residual = 0.999492, Final residual = 10.9301, No Iterations 5001 time step continuity errors : sum local = 4.02462e+11, global = 272.815, cumulative = 270.604 Wall patch = 3 Wall patch name = cyl_wall Uav = (0.304463 -0.695537 7.76322e-28) Aref = 1.44217 Lref = 0.515075 DragCoefficient = -4.22157e+16 pressureDragCoefficient = -4.06296e+16 viscDragCoefficient = -1.58615e+15 LiftCoefficient = 1.55157e+16 MomentCoefficient = -5.22663e+17 |
In my experience, interFoam is
In my experience, interFoam is not ready for open channel flow, yet.
|
I use the icoFoamAutoMotion eq
I use the icoFoamAutoMotion equations and related algorithms, piso loop and stuff. The moving mesh part is copied from interFoam, since that solver works with the movingFvMesh class.
Regards, Frank |
Hoi Frank,
Je piso-loop con
Hoi Frank,
Je piso-loop convergeert niet. Als je kijkt naar de initial residual van p, moet deze minimaal 1 orde lager zijn bij de tweede iteratie in de loop, dit is niet het geval (Hier is die waarde zelfs iets hoger): ICCG: Solving for p, Initial residual = 0.84611, Final residual = 9.43944e-07, No Iterations 290 time step continuity errors : sum local = 4.11744e-11, global = 1.12633e-14, cumulative = 6.75289e-10 ICCG: Solving for p, Initial residual = 0.859509, Final residual = 9.95797e-07, No Iterations 290 Je kunt proberen om het aantal piso correctors op te voeren of een kleinere tijdsstap nemen. Succes Guido |
Ok, more piso correctors is no
Ok, more piso correctors is not improving the convergence.
On the other hand, a much smaller timestep is actually leading to convergence at the beginning of the motion. After a certain time (which is larger than the previous case) also this solution is diverges. When I turn off the mesh.move() function in my solver no problem occurs. So the problem lies in the piso loop when the mesh moves. Thanks, Frank |
Could someone please comment o
Could someone please comment on my solver, shown above.
When I move the mesh, the piso loop diverges, also when I decrease the timestep or increase the number of correctors. When I turn off the mesh motion, the flow solver works fine. In addition when I turn off the flow solver instead of the mesh motion, the moving mesh is calculated well. Do I have to include some extra terms in the piso loop to account for the mesh motion? Thanks in advance. Regards, Frank |
The solver looks fine (well, a
The solver looks fine (well, at least not obviously wrong). However, your residual behaviour indicates a serious problem and my guess is that you have messed up the boundary conditions.
Note that the wall boundary needs to be specified in terms of absolute velocity, so a moving wall with no flux through it will not have fixedValue (0 0 0). If this is the problem use the movingWallVelocity b.c., which corrects the normal component based on the mesh motion. As a test, have a look at the fluxes through your wall boundaries - if all is well, they should be exactly zero. Hrv |
As far as I know I already use
As far as I know I already used the movingWallVelocity b.c.
My boundary conditions in 0/U looks like (cyl_wall is the moving surface): dimensions [0 1 -1 0 0 0 0]; internalField uniform (1 0 0); boundaryField { far_field { type symmetryPlane; } outlet { type zeroGradient; } inlet { type fixedValue; value uniform (1 0 0); } cyl_wall { type movingWallVelocity; value uniform (0 0 0); } frontAndBackPlanes { type empty; } } =========================================== And 0/motionU looks like: dimensions [0 1 -1 0 0 0 0]; internalField uniform (0 0 0); boundaryField { far_field { type symmetryPlane; } outlet { type fixedValue; value uniform (0 0 0); } inlet { type fixedValue; value uniform (0 0 0); } cyl_wall { type fixedValue; value uniform (0 0 0); } frontAndBackPlanes { type empty; } } =========================================== When I used movingFlepFvMesh the phi of some surface elements after 1 sec looks like: 1/phi: cyl_wall { type calculated; value nonuniform List<scalar> 280 ( 0 1.73472e-18 0 1.73472e-18 -1.73472e-18 0 0 0 0 -1.73472e-18 0 0 0 1.73472e-18 0 0 8.67362e-19 -8.67362e-19 0 -8.67362e-19 . . . . . } What is wrong with this b.c. implementation. I think I mess thing up with U and motionU........The fixedValue of motionU is overruled by movingFlapFvMesh, right? Regards, Frank |
Summarising:
What b.c. do y
Summarising:
What b.c. do you generally have to apply to a moving body? * motionU: fixedValue on the moving surface (this is being overruled by the movingFvMesh class) * U: movingWallVelocity (uniform (0,0,0)) on the moving surface (according to dr Jasak this will correct the boundary flux) Is this correct? or must by movingWallVelocity be applied to both U and motionU? If these settings are correct, my problem of the diverging piso loop has its origins elsewhere. Thanks, Frank |
motionU: you can do whatever y
motionU: you can do whatever you like, this is only used to specify the motion. Typically, if you wish to specify the motion of an object, you would use the fixedValue b.c. on its surface (makes sense).
U: movingWallVelocity. I usually stick to this. It will correct the normal component of the velocity to achieve zero flux through the boundary based on the vertex motion (wherever it may come from). There's not connection between the motionU b.c. and the U b.c. - they are totally independent. I'm not sure what may be going wrong. Try setting very diffusive discretisation and if that still fails you've messed up the boundary conditions or something even more subtly. Have you tried looking at the mesh motion (motion only, no flow), to see if it behaves the way you think it should? This should work without any trouble, Hrv |
Thanks for your answers, but m
Thanks for your answers, but my set-up is still not working. I tried central, upwind and blended schemes for the convective flux, but the solution still diverges after 1 second or so.
When I read my ../1/phi file the values at the boundary are of the order 10e-18. Does this flux need to be exactly zero or is 10e-18 small enough. I also switched off the solver such that the mesh motion can be visualised. In that case the mesh looks very nice at all time up to t=20. No diverges occurs. The solution divergences only when the piso loop is solved. Maybe the movingFlapFvMesh coupling to the fixedValue of motionU is not done correctly.... Frank |
I'll have a look if you want.
I'll have a look if you want. Could you please pack the code (ready to run, if the mesh is made with blockMesh, please delete unnecessary files, with time-step, discretisation, initial fields etc set up) and the case (clean up the object files and similar garbage + please include everything I need to compile it up).
Hrv |
Thanks,
I just sent you th
Thanks,
I just sent you the case and solver files. Both are ready to run/compile. Frank |
Very easy: it does not work be
Very easy: it does not work because you've set the time step to give the Courant number of 3.5, which is jolly high. I have dropped the time step by a factor of 4 and with the Co number of around 1 it works beautifully for the first second or so and then the Co number goes higher again which makes it blow up. Nothing wrong with the solution, just the normal flow development on a newrrt fine mesh.
Then, I have added automatic time-step control and it works beautifully. Hrv |
Thanks again for your time :-)
Thanks again for your time :-).
In my first case (dt=0.02) also the Co of around 3.5 gives me a good solution for the first second. After that the Co always blows up which is rather independent of the chosen timestep. I really like a constant time step. Since I need to simulate for 100 seconds or so, this is easy to plan and estimate calculation times. Though, this seems to be impossible in my case. Frank |
I did the first 5 seconds and
I did the first 5 seconds and made a little animation for you (it is ugly but proves the point). The automatic time-step control is well-behaved and there's no problem with the simulation (I am using 1 non-orthogonal corrector). If you want the regular output dumps with automatic time-step control, e.g. every 1 second or similar, use adjustableRunTime writeControl option.
cylinder animation Hrv |
Allright, thanks. Is there act
Allright, thanks. Is there actually some sort of guideline how many non-orthogonal correctors are needed? Is the accuracy or the efficiency affected by this choice, or not?
Frank |
You have several options. I s
You have several options. I suspect the problem you are getting is that the mesh deformation may become quite sensitive to the convergence of the motion solver. The issue is that you've got a really big mesh qith a little bit of motion in the middle. Ideally, you would want to subset a piece of the mesh around the cylinder which is big enough to accomodate the motion but smaller than the complete domain. This would make yor simulation much faster (less work in mesh motion) and probably cause less problems. However, it is a bit more difficult to set up (you will be dealing with 2 meshes, where the motion mesh is a subset of the first one).
If you want to see exactly what is going on, try running the non-rothogonality check from time to time during the mesh motion. You will see that the non-orthogonality changes and this influences the convergence behaviour. You can tackle it in a number of different ways, but basically they fall into the two categories you mention - a less accurate but more robust discretisation for the laplacian: limited Gauss. You lose accuracy but it does not affect the speed - non-orthogonal correctors: (typically, 1 is plenty). Increases the cost of simulation but is more accurate than the limited snGrad stuff from above. Enjoy, Hrv |
Hi everybody,
In my simulat
Hi everybody,
In my simulations concerning moving meshes, I need a very small timestep in OpenFOAM compared to other CFD solvers (e.g. Fluent). In OpenFOAM I use second-order Crank-Nicholson, so in fact a smaller timestep than Fluent (which is first order in time for moving meshes) should work. The mesh motion is OK for larger timesteps, but for convergence of the icoDyMFoam solver a much much smaller timestep is necessary. Otherwise the Courant number increases till the solution diverges. Why is Crank-Nicholson in OpenFOAM not unconditionally stable, as proved in many textbooks (e.g Ferziger, Peric or Hirsch). Regards, Frank |
Why is Crank-Nicholson in Open
Quote:
I have shown in my PhD that Crank-Nicolson is not unconditionally unstable (based on numerical diffusion) so what you say is definitely news to me. Would you clarify please. Hrv |
mmm, its been a while that I r
mmm, its been a while that I read those books. I always thought that Von Neumann stability analysis shows this unconditional stability of Crank-Nicholson. I guess you're right, p.150 of Ferziger says that the above is true, but in some cases it may become unstable. It also says that the stability limits are problem dependent. (I borrowed Hirsch once, so I cannot review that book now).
In my simulations I need at least about 20000 time steps in one flapping (wing) period. This is 10 times more than in Fluent. Are there ways to tweak OpenFOAM such that a larger timestep can be used at the same accuracy? What are your ideas about this? Regards, Frank |
Hi Frank,
Please note that
Hi Frank,
Please note that the backward scheme is used for second-order temporal discretisation in Fluent and you are using the Crank-Nicolson scheme in OpenFoam. Regards, Zeljko |
Hi Zeljko,
Fluent 6.2 is on
Hi Zeljko,
Fluent 6.2 is only working with first-order temporal discretisation (Euler) when the dynamic mesh module is used. Therefore using second order Crank-Nicholson in OpenFOAM seems to be more efficient. But, up to now I need very small time-steps with which I am not happy. Regards, Frank |
I always use the backward temp
I always use the backward temporal scheme with the moving mesh in OpenFOAM.
Zeljko |
Hi Frank,
The answer is: go
Hi Frank,
The answer is: go back to basics. First, you know that running transient SIMPLE will have no Courant number limit. Start with that, and first order discretisation in time - if you like, you can now have a time step of half an hour. Of course, with the above the temporal accuracy is a disaster. The next step is to tune the time-step to get the right time resolution for your problem. When you do that, you can switch to second-order temporal scheme (backward it time, like Zeljko says) or Crank-Nicolson. Being aware that Crank Nicolson is borderline stable, you will now need to monitor the discretisation (dispersion) error to see if it fouls up your solution and again adjust the time-step accordingly. Keep in mind that you are perfectl welcome to configure the solver to run transient SIMPLE with Crank-Nicolson - I suspect you will need to write some top-level code to do that well. Finally, it is really important to compare like with like. I know exactly what commercial CFD is doing on cases like this (been there) and comments like " This is 10 times more than in Fluent" when you are comparing first order in a commercial code and full second order in space and time in OpenFOAM are completely meaningless. Hrv |
Thanks for your time again! Yo
Thanks for your time again! You are completely right that I should change one thing at a time to adjust and tune the OpenFOAM settings for my case.
To solve for the temporal resolution I need about at least 200 timesteps per plunging period. OpenFOAM (with Piso and backward or Crank-Nicholson) needs more than 40.000 timesteps to keep the max Courant number below 2. I use adjustableRunTime. It seems that I need to adjust icoDyMFoam to use SIMPLE pressure/velocity coupling. Why are all the transient solvers using PISO, by the way? Are you aware of any major problems w.r.t. implementing SIMPLE in the transient codes? Regards, Frank |
Well, my second PhD supervisor
Well, my second PhD supervisor at Imperial College was Raad Issa of the PISO fame. ;-) In any case, PISO is more efficient for cases where the time-step is determined by the need for temporal accuracy and the mesh is good-looking (ie. no small cells in regions of fast flow).
I have written several transient SIMPLE solvers (see development library) and there's never any problem. BTW, the "hot" thing at the moment is non-iterative time advancement for LES, which is shown to be second order and solves the pressure equation only once per time-step. Hrv |
Allright. I guess my cases are
Allright. I guess my cases are rather extreme, in the sense that large acceleration and velocities occur during plunging (and later rotation). Also the mesh is very nice in terms of skewness and orthogonality, but the cells near the wing's surface are very small.
Therefore, I will give SIMPLE a try since apparentrly the time-step not is dictated by the temporal accuracy. I have your development code, so I am going to try to combine transientSimpleFoam and icoDyMFoam. Regards and thanks again for your patience, Frank |
Dr. Jasak,
I read in your p
Dr. Jasak,
I read in your post above (from Sunday, November 19, 2006 - 06:37 am) that LES simulations may use non-iterative time advancement. Can this be done in OpenFoam 1.4.1? Which solvers and settings would one use to achieve this? Thank you, Alessandro |
Hello,
Sunday, November 19,
Hello,
Sunday, November 19, 2006 - 06:37 am: almost exactly a year ago... You got me confused for a moment. The reference I refer to came from Fluent News, 2004 and offers two options: References: 1. R.I. Issa, "Solution of the Implicitly Discretized Fluid Flow Equations by Operator-Splitting," Journal of Computational Physics 62, p. 40-65, 1985. 2. J.B. Perot, "An Analysis of the Fractional-Step Method," Journal of Computational Physics 108, p. 51-58, 1993. The first one is standard, ie. you will find it in most OpenFOAM solvers (among other things, the guy was my PhD supervisor so I knew about this practically from my first days in CFD). I haven't implemented the second one in OpenFOAM yet, but this should not be too hard. There also some modifications to the settings you can do to allow large time-steps when using PISO, but that's a different topic. Hope this is useful, Hrv |
Dr. Jasak,
Thank you for th
Dr. Jasak,
Thank you for the references. I will take a look at them and try to understand their implementation in OF. appreciate your help, Alessandro |
All times are GMT -4. The time now is 20:50. |