
[Sponsors] 
February 17, 2011, 00:28 
Submarine with SimpleFoam

#1 
New Member
Alexandre Rubel
Join Date: Dec 2010
Location: Launceston, Tasmania AUSTRALIA
Posts: 28
Rep Power: 6 
Hi
I'm new with OpenFoam and I done a bit of CFD with CFX before. I doing a quite simple study about a submarine (without appendages) in OF 1.6. I'm using simpleFoam but for the moment I get quite bad results, and convergence problems. My inlet speed is 3.5 m/s and my mesh is around 6millions elements with a Y+ around 10. (this Mesh gives me very good results in CFX) The checkMesh is OK : Code:
Mesh stats points: 6238133 faces: 18602584 internal faces: 18491528 cells: 6182352 boundary patches: 4 point zones: 0 face zones: 0 cell zones: 0 Overall number of cells of each type: hexahedra: 6182352 prisms: 0 wedges: 0 pyramids: 0 tet wedges: 0 tetrahedra: 0 polyhedra: 0 Checking topology... Boundary definition 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 SURFACE_HULL 38368 38370 ok (closed singly connected) SURFACE_SIDE 46288 46464 ok (nonclosed singly connected) SURFACE_OUTLET 17600 17689 ok (nonclosed singly connected) SURFACE_INLET 8800 8889 ok (nonclosed singly connected) Checking geometry... Overall domain bounding box (7.696 8.712 8.712) (21.78 8.712 8.712) Mesh (nonempty, nonwedge) directions (1 1 1) Mesh (nonempty) directions (1 1 1) Boundary openness (5.04721e17 2.96585e16 3.75431e17) OK. Max cell openness = 7.73145e15 OK. Max aspect ratio = 496.544 OK. Minumum face area = 1.83581e08. Maximum face area = 0.376178. Face area magnitudes OK. Min volume = 6.63038e12. Max volume = 0.102174. Total volume = 6329.87. Cell volumes OK. Mesh nonorthogonality Max: 61.2197 average: 9.26162 Nonorthogonality check OK. Face pyramids OK. Max skewness = 2.01934 OK. Mesh OK. I use potentialFoam then lunch the simpleFoam without turbulences and finally set the turbulences on. I set the following relaxation factors Code:
relaxationFactors { p 0.3; U 0.7; k 0.5; omega 0.5; R 0.7; nuTilda 0.7; } I got a lot of question with no answers : 1/ Am I using a good scheme (potential,turb off, turb on) to solve my problem ? 2/ Should I use WallFunction or fixedValue (0) for my hull BC (omega,k,nut) ? 3/ Should I use nNonOrthogonalCorrectors ? 4/ Is my relaxation factors good ? 5/ Do you have any idea to help me ? Tanks a lot Alex 

February 17, 2011, 05:45 

#2 
Senior Member
Vesselin Krastev
Join Date: Jan 2010
Location: University of Tor Vergata, Rome
Posts: 361
Rep Power: 10 
Can you post your fvSchemes and fvSolution dictionaries?
V. 

February 17, 2011, 06:15 

#3 
Senior Member
maddalena
Join Date: Mar 2009
Posts: 436
Rep Power: 13 
Hi Alex,
a good way to enter in the OpenFOAM world of controlling any solver is this thread: pressure eq. "converges" after few time steps. I found it really inspirating (also because it helped me to solve quite a lot of problems ). V collaborate to that thread as well. mad 

February 17, 2011, 06:27 

#4 
New Member
Alexandre Rubel
Join Date: Dec 2010
Location: Launceston, Tasmania AUSTRALIA
Posts: 28
Rep Power: 6 
thanks for theese answers. I'm not a work (It's 9pm in Australia) so i can't try right know.
I gonna post my files tomorow and try your thread. Alex Last edited by alex_rubel; February 17, 2011 at 18:36. Reason: Add files 

February 17, 2011, 18:43 

#5 
New Member
Alexandre Rubel
Join Date: Dec 2010
Location: Launceston, Tasmania AUSTRALIA
Posts: 28
Rep Power: 6 
I posted my system file in the previous post
I red plenty of interesting threads this morning, but my understanding of english is not perfect so I really need your help for the fvScheme and fvSolution files. For the BC I found information which say that because my Y+ is between 1<Y+<30 and I'm using High Re model I should set : k> zeroGradient Omega > WF nut > nutSpalartAllamarasWallFunction What do you think about these settings ? Last edited by alex_rubel; February 17, 2011 at 20:25. 

February 18, 2011, 06:25 

#6 
Senior Member
Vesselin Krastev
Join Date: Jan 2010
Location: University of Tor Vergata, Rome
Posts: 361
Rep Power: 10 
1) There's no need to use potentialFoam or start the simulation with turbulence switched off: if the numerical settings are appropriate (and mesh is not a disaster) the run should start normally with turbulence switched on.
2) The BC's at the wall are substantially correct (for k you can also use kqRWallFunction, but is the same as zeroGradient). You can try also nutWallFunction for nut (it's a simplified wall function compared with the nutSpalartAllmaras one) 3) In the fvSchemes file there's an error (divergence and laplacian schemes are for epsilon instead of omega): I dont' think this is the file you employed for your runs, as the simulation should not have started if the fvSchemes entries were not consistent with the turbulent quantities involved in your case. However, if you change epsilon to omega this file is also quite ok as a first attempt (it is not the best choice if you are searching for accuracy, but it should work also on not so good meshes) 4) In the fvSolution you have to lower a lot the tolerance (NOT the relTol) parameter for all the quantities: something like 10^12 for p and 10^10 for the other quantities should be ok. Also (but this is not so fundamental) you can add 1 or 2 nonOrthogonalCorrectors. Finally, the PCG solver for the pressure is ok, but as your mesh is quite big probably you could have a faster solution using GAMG instead (have a look in some of the tutorials or here in the forum to see how to set properly this kind of solver: of course, tolerance and relTol parameters should be the same as for PCG) Good luck V. 

February 20, 2011, 18:55 

#7 
New Member
Alexandre Rubel
Join Date: Dec 2010
Location: Launceston, Tasmania AUSTRALIA
Posts: 28
Rep Power: 6 
Hi, i ran a simulation this week end with the following settings :
k>zeroGradient Omega>WF nut>SpalartAllmarasWF I ran the problem with the turbulence and the high speed and I get a quite good convergence using the fvSchemes and fvSolution attached. But I get 70% of error for my Cd, You said that my files are good for the convergence but might be improved to get a better accuracy. Can you help me for that Thanks Alex 

February 21, 2011, 01:41 

#8 
New Member
Alexandre Rubel
Join Date: Dec 2010
Location: Launceston, Tasmania AUSTRALIA
Posts: 28
Rep Power: 6 
i try the GAMG solver today and It works very well, after 1000 TStep I got the same Residuals "quality" and it reduce the ExecutionTime from 53674s (PCG) from 19778s (GAMG)
I red very interesting post about the nCellsInCoarsestLevel parameter and if I understand evrything I should set it up to sqrt(MeshNbrofElement) ? I curently using 10. If you can give me advice of the relTol parameter it would help me a lot. I may ask you too much but I hope you can help me Alex 

February 21, 2011, 07:02 

#9  
Senior Member
Vesselin Krastev
Join Date: Jan 2010
Location: University of Tor Vergata, Rome
Posts: 361
Rep Power: 10 
Quote:
1) In fvSchemes, try first to set div(phi,U) to Gauss linearUpwindV cellMDLimited Gauss linear 1; 2) If this change causes stability issues, put the limiter also on grad(U) (NOT on grad(p)), by setting: grad(U) cellMDLimited Gauss linear 1; 3) If you still have stability troubles (but I don't think this will be the case), put a limiter on the laplacianSchemes too, by setting: default Gauss linear limited 0.5; 4) If you will be successful in using the linearUpwind scheme for div(phi,U) (either with or without the additional limiting in points 2 and 3), you can try to extend it also on the convection terms for k and omega (but to my knowledge and experience the biggest improvement in accuracy tend to come from the momentum convection term dicretization choice), by setting: div(phi,k) Gauss linearUpwind cellMDLimited Gauss linear 1; and the same for omega. Hope this helps! V. 

February 21, 2011, 07:21 

#10  
Senior Member
Vesselin Krastev
Join Date: Jan 2010
Location: University of Tor Vergata, Rome
Posts: 361
Rep Power: 10 
Quote:
About the relTol parameters, It's my opinion that using too low values in a steadystate and underrelaxed solution practice (as is the SIMPLE algorithm in the simpleFoam solver) is not useful at all and produces only a big waste of time for each iteration...However, a minimum convergence criterion should be provided within the single iteration, thus I think that your values are quite ok (0.05 for p and 0.1 for the other quantities). If you want to make some tests, try to lower them to 0.001 for p and 0.01 for the other quantities, but I don't believe this will produce significant improvements in the convergence pattern of your residuals. Hope this helps, too! V. 

February 21, 2011, 07:23 

#11 
New Member
Alexandre Rubel
Join Date: Dec 2010
Location: Launceston, Tasmania AUSTRALIA
Posts: 28
Rep Power: 6 
Thanks, I gonna try with that.
Do you have information about the nCellsInCoarsestLevel parameter ? I red plenty of different explanation about this. Apparently it's set how many cellsare going to be solve by the same equation but I'm not sure at all Alex Edit Tanks again for your explanation it's a lot more clear in my head Last edited by alex_rubel; February 21, 2011 at 07:42. 

February 22, 2011, 19:47 

#12  
New Member
Alexandre Rubel
Join Date: Dec 2010
Location: Launceston, Tasmania AUSTRALIA
Posts: 28
Rep Power: 6 
Hi,
I get some good improvement for my simulation with a Y+ of 30 (KomegaSST with wallFunctions) using Code:
divSchemes { default none; div(phi,U) Gauss linearUpwindV Gauss linear; div(phi,k) Gauss linearUpwind Gauss linear; div(phi,omega) Gauss linearUpwind Gauss linear; div((nuEff*dev(grad(U).T()))) Gauss linear; I tried your settings but it leads to bounding problem for k and omega and an equivalent result. Quote:
g Code:
radSchemes { default Gauss linear; grad(p) Gauss linear; grad(U) Gauss linear; // upwind 1st order , linear = second } divSchemes { default none; div(phi,U) Gauss upwind; div(phi,k) Gauss upwind; div(phi,omega) Gauss upwind; div(phi,epsilon) Gauss upwind; div((nuEff*dev(grad(U).T()))) Gauss linear; } laplacianSchemes { default Gauss linear corrected; } For the WallFunctions I set zeroGradient for k, omegaWallFunction and nutSpalartAllmarasWallFunctions do I have to to set a very small fixedValue instead of zeroGradient, and what about nut ? I tried with zerGradient instead of omegaWallFunction but I didn't see any difference Thanks Alex ps : I used Gauss linearUpwind Gauss linear; instead of Gauss linearUpwind cellMDLimited Gauss linear 1; I know that cellMDLimited is a modifier for the gradient scheme Gauss linear but that's all is you have more information it will be great. pps : I hop I 'll be able to help you back one day Last edited by alex_rubel; February 23, 2011 at 00:13. 

February 23, 2011, 09:28 

#13  
Senior Member
Vesselin Krastev
Join Date: Jan 2010
Location: University of Tor Vergata, Rome
Posts: 361
Rep Power: 10 
Quote:
1) From your last post it is not totally clear what you want to do now to further improve your computed Cd value...However, if you want to pass from wall functions modeling to a direct solution of the boundary layer, you cannot definitely leave the first order upwind scheme on the convection terms, as accuracy is a milestone for the b.l. direct solving approach. 2) If you had a better behavior without the gradient limiters turned on, this is probably due to your hexamesh (hexameshes exhibits usually a more stable behavior then the tetradominant one in conjunction with higherorder interpolation schemes) 3) To directly solve the b. l. with the SST model you'll need to: Have the first nodes away from the wall inside the y+<1 limit Set k at the wall to fixedValue 10^10 (or, in general, to a very low value) Set omega at the wall to fixedValue <a very high value> (it should be something like 10^6, but there is a formula usually used in the literature which depends on 1/y^2, where y is the the absolute distance of the first nodes: I think you can find this formula somewhere here in the forum or at least in some paper related to the SST model). Set nut at the wall to calculated (an initial "value uniform 0" shoul be ok) Of course, lowering the y+ will probably lead to near wall cells with very high aspect ratios, but you should try to remain below the value of about 500 which you had in our y+=10 case (thus you have to reduce the cell's dimension also in the streamwise direction). In any case, as a first attempt you can try the same numerical setup which gave to you the best results in the WF configuration. About the cellMDLimited option, it is a muldidimensional (MD) gradient limiter (it means that the limiter is imposed separately in the three dimensions), but honestly still I don't know what exactly are the limiting criteria adopted (as in OF there are three other kind of limiters: cellLimited, faceMDLimited, and faceLimited). Good luck with your new runs! V. 

February 24, 2011, 01:51 

#14  
New Member
Alexandre Rubel
Join Date: Dec 2010
Location: Launceston, Tasmania AUSTRALIA
Posts: 28
Rep Power: 6 
Hi V,
Quote:
Quote:
And thanks a lot for the advices of your third point I gonna try with that and I'll post my results asap. I found something about a missing scalar of sqrt(2) in the KomegaSST model, I modified it (If I change the code, save it and make ./Allwmake it's ok ?), and for the moment my Cd is heading to 0.00125 instead of 0.0014 before Alex 

February 24, 2011, 05:57 

#15  
Senior Member
Vesselin Krastev
Join Date: Jan 2010
Location: University of Tor Vergata, Rome
Posts: 361
Rep Power: 10 
Quote:
Quote:
V. 

February 24, 2011, 19:24 

#16  
New Member
Alexandre Rubel
Join Date: Dec 2010
Location: Launceston, Tasmania AUSTRALIA
Posts: 28
Rep Power: 6 
Hi V.
Quote:
Code:
ddtSchemes { default steadyState; } gradSchemes { default Gauss linear; grad(p) Gauss linear; grad(U) Gauss linear; // grad(U) cellLimited Gauss linear 1; } divSchemes { default none; div(phi,U) Gauss linearUpwindV Gauss linear; div(phi,k) Gauss linearUpwind Gauss linear; div(phi,omega) Gauss linearUpwind Gauss linear; div((nuEff*dev(grad(U).T()))) Gauss linear; } laplacianSchemes { // default Gauss linear corrected; // default Gauss linear limited 0.5; default Gauss linear limited 0.333; } I ran a simulation for my low Y+ model with the following settings nut : Code:
SURFACE_HULL { type calculated; value uniform 0; } Code:
SURFACE_HULL { type fixedValue; value uniform 10^6; //To modify } Code:
SURFACE_HULL { type fixedValue; value uniform 10^10; } Code:
gradSchemes { default Gauss linear; grad(p) Gauss linear; // grad(U) Gauss linear; grad(U) cellMDLimited Gauss linear 1; } divSchemes { default none; div(phi,U) Gauss linearUpwindV cellMDLimited Gauss linear 1; div(phi,k) Gauss linearUpwind cellMDLimited Gauss linear 1; div(phi,omega) Gauss linearUpwind cellMDLimited Gauss linear 1; div((nuEff*dev(grad(U).T()))) Gauss linear; } laplacianSchemes { // default Gauss linear corrected; // default Gauss linear limited 0.5; default Gauss linear limited 0.333; } I used the following relaxation factors : Code:
relaxationFactors { p 0.1; U 0.5; k 0.3; omega 0.3; } Alex 

February 25, 2011, 06:31 

#17  
Senior Member
Vesselin Krastev
Join Date: Jan 2010
Location: University of Tor Vergata, Rome
Posts: 361
Rep Power: 10 
Quote:
About your settings for the highy+ mesh, I think they are pretty good (probably you will not be able to obtain better results in any case, as a consequence of the simplified wall modeling). Regards V. 

February 25, 2011, 06:57 

#18  
New Member
Alexandre Rubel
Join Date: Dec 2010
Location: Launceston, Tasmania AUSTRALIA
Posts: 28
Rep Power: 6 
Hi V
Quote:
There is a checkMesh in the simu7 zip file for the low Y+ model 

February 25, 2011, 07:12 

#19  
Senior Member
Vesselin Krastev
Join Date: Jan 2010
Location: University of Tor Vergata, Rome
Posts: 361
Rep Power: 10 
Quote:
<<Writing 34328 cells with high aspect ratio to set highAspectRatioCells Here it is the problem....I told you in one of my previous posts that you should keep a not extremely high aspect ratio for the nearwall cells, although this would constrain you to reduce the dimensions of the cells also in the streamwise and spanwise directions. At most you must keep the same max value of your high y+ case (which was around 500), but probably you'll need to reduce it further (I know that the number of cells will increase a lot, but this is a price you have to pay for the direct solution of the b. l.). So, reduce your max axpect ratio (at least to 500), set the correct b. c. for omega at the wall (see the reference in my previous post), keep all the other settings as before, and let's see what happens... Regards V. 

February 25, 2011, 07:17 

#20  
New Member
Alexandre Rubel
Join Date: Dec 2010
Location: Launceston, Tasmania AUSTRALIA
Posts: 28
Rep Power: 6 
Quote:
I'm doing a yaw angle study for my high Y+ model and the result are quite consistent I don't think I can improve a lot this model. I'm gonna try Kepsilon And Komega but my experiment on my case with CFX shows me that KomegaSST is the more accurate Alex 

Tags 
komegasst, simplefoam, submarine, wallfunction 
Thread Tools  
Display Modes  


Similar Threads  
Thread  Thread Starter  Forum  Replies  Last Post 
Laminar simpleFoam and inviscid simpleFoam  herenger  OpenFOAM Running, Solving & CFD  7  July 11, 2013 06:27 
MPI Error  simpleFoam  Floating Point Exception  scott  OpenFOAM Running, Solving & CFD  3  April 13, 2012 16:34 
simpleFoam ddt Euler ?  MoITB  OpenFOAM Running, Solving & CFD  2  June 12, 2010 13:36 
Naca0012 ke mpirun gives fpe whereas simpleFoam not  Pierpaolo  OpenFOAM  1  May 8, 2010 03:08 
Error running simpleFoam in parallel  skabilan  OpenFOAM Running, Solving & CFD  2  August 29, 2008 09:42 