
[Sponsors] 
January 24, 2006, 05:57 
I have OpenFOAM converging a c

#1 
Senior Member
Håkan Nilsson
Join Date: Mar 2009
Location: Gothenburg, Sweden
Posts: 193
Rep Power: 9 
I have OpenFOAM converging a cyclic water turbine runner computation. I am however not pleased with the convergence rate, and I think that my setup of the case could be improved (and I could of course use your ideas)
My case is shown in the following figure. It shows the cyclic 1/5 of the runner which is cyclic about the yaxis. The geometry is colored gray. The cyclic boundaries cut through the center of the blades, and they are colored by the turbulent kinetic energy. The flow comes in at the top, where I set steady highly swirling profiles for U,k,epsilon and dp/dn=0. At the walls I use wall functions, and I use the standard kepsilon model. At the outlet I have tried different alternatives, with the following common settings: p=0, dU/dn=0, and the following alternatives: dk/dn=0, depsilon/dn=0 or inletOutlet on k and epsilon, with an inlet value similar to the average at the true inlet. The computation converges nicely in the upper half of the computational domain, but below the runner (in the diffusing part) the convergence is extremely slow, particularly close to the rotational axis (in the region colored red below the hub). In part of this region there is a recirculation. There are several possible resons for this behaviour: 1) the flow in this region is not steady in reality (however, there does not seem to be any oscillations in the solution but the flow seems to be extremely slowly converging towards a solution), 2) the outlet boundary condition is not good enough (this is likely, I don't like the concept of constant pressure for outlets with swirling flow), 3) there is something wrong with my cyclics (in another thread I have described my problems using this mesh in OpenFOAM, and I had to change the too strict geometrical tolerance to make it work. I have however checked that the cyclic coupling looks alright), 4) there is a problem with the prism cells at the axis (they are generated due to the collapse of one of the faces in the hexcells, the rest of the grid cells are hexas), 5) etc. etc. (please give me your ideas) The residuals at the end of the computation are shown in the following figure: The flow is computed with a modified simpleFoam, where the Coriolis and centrifugal forces are added as: dimensionedVector Omega ( "Omega", dimensionSet(0, 0, 1, 0, 0), vector(0, 2*M_PI/60*(595), 0) ); volVectorField Fcent = (Omega ^ (Omega ^ mesh.C()));  tmp<fvvectormatrix> UEqn ( fvm::div(phi, U) + turbulence>divR(U) + (2*Omega ^ U) + Fcent ); I hope that someone out there has any idea on what might cause this slow convergence below the hub. Below you can see my fvSchemes, fvSolution Håkan. fvSchemes: ddtSchemes { default steadyState; } gradSchemes { default Gauss linear; grad(p) Gauss linear; grad(U) Gauss linear; } divSchemes { default none; div(phi,U) Gauss GammaV 1; div(phi,k) Gauss Gamma 1; div(phi,epsilon) Gauss Gamma 1; // div(phi,U) Gauss upwind; // div(phi,k) Gauss upwind; // div(phi,epsilon) Gauss upwind; div(phi,R) Gauss upwind; div(R) Gauss linear; div(phi,nuTilda) Gauss upwind; div((nuEff*dev(grad(U).T()))) Gauss linear; } laplacianSchemes { default none; laplacian(nuEff,U) Gauss linear corrected; laplacian(1A(U),p) Gauss linear corrected; laplacian(DkEff,k) Gauss linear corrected; laplacian(DepsilonEff,epsilon) Gauss linear corrected; laplacian(DREff,R) Gauss linear corrected; laplacian(DnuTildaEff,nuTilda) Gauss linear corrected; } interpolationSchemes { default linear; interpolate(U) linear; } snGradSchemes { default corrected; } fluxRequired { default no; p; }  fvSolution: solvers { // p ICCG 1e06 0.01; p AMG 1e07 0.01 100; U BICCG 1e06 0.1; k BICCG 1e06 0.1; epsilon BICCG 1e06 0.1; R BICCG 1e06 0.1; nuTilda BICCG 1e06 0.1; } SIMPLE { nNonOrthogonalCorrectors 0; pRefCell 28; pRefValue 0; } relaxationFactors { p 0.2; U 0.7; k 0.5; epsilon 0.5; R 0.5; nuTilda 0.5; } 

January 24, 2006, 08:49 
Unless you have high nonortho

#2 
Senior Member
Eugene de Villiers
Join Date: Mar 2009
Posts: 725
Rep Power: 13 
Unless you have high nonorthogonality, your settings look fine. Might want to play with your convection schemes.
At a guess I would say your problem is due the erronous constant outlet pressure. At least that is the first thing I would check. 

January 24, 2006, 09:26 
Thanks for your response, Euge

#3 
Senior Member
Håkan Nilsson
Join Date: Mar 2009
Location: Gothenburg, Sweden
Posts: 193
Rep Power: 9 
Thanks for your response, Eugene.
Hexmeshes in these kinds of geometries tend to have both nonorthogonal and skew cells, which is the case also in the present mesh. Are you proposing that I should increase nNonOrthogonalCorrectors? I haven't yet compared the reults quantitatively with experiments and previous computations, but the results look fine in the upper part of the domain, which also has nonorthogonal and skew cells. The convergence problem is only below the hub. I have tried both upwind and Gamma convection schemes with the same problem, but I actually don't think that the convergence problem is due to the convection schemes. What I do know is that upwind is likely to give bad results for this kind of flow. I agree that there is probably a problem with the outlet boundary condition. What boundary conditions would you recommend at the outlet for this case with highly swirling flow? In previous research with another code I have used Neumann for all variables at the outlet. I saw the newly discussed comparison between Fluent and OpenFOAM for the cyclinder in crossflow and that there was also a problem with the outlet boundary condition for that case in OpenFOAM. It is really important to have good outlet boundary conditions unless you can afford to extend your computational domain, so if anyone has some good suggestion I would be very happy. Håkan. 

January 24, 2006, 10:05 
If you run checkmesh and your

#4 
Senior Member
Eugene de Villiers
Join Date: Mar 2009
Posts: 725
Rep Power: 13 
If you run checkmesh and your nonorthogonality is above say 70 degrees, you will want to put a limiter on the corrector of the laplacian schemes like:
laplacian(1A(U),p) Gauss linear limited 0.3; But as you say, this is unlikely to be your problem. As mentioned in the other thread, you can simulate an avergae pressure outlet boundary by updating the outlet each iteration in the top level code provided the outlet is on a single processor. Basically something like: //before solver starts scalar pMeanOut = 0; //check if this processor has an outlet boundary label outletPatchID = mesh.boundary().patch().findPatchID("yourOutletBou ndaryNameHere") if(outletPatchID != 1) { pMeanOut = mean(p.boundaryfield()[outletPatchID]); } //inside the time loop if(outletPatchID != 1) { scalarField newP = p.boundaryField()[outletPatchID]; const labelList& outletPatchCells = mesh.boundary()[outletPatchID] .faceCells(); forAll(newP, faceI) { newP[faceI] = p[outletPatchCells[faceI]]; } newP /= pMeanOut/mean(newP); p.boundaryField()[outletPatchID] == newP; } This hasnt been tested, so take it only as a suggestion. And before you try this just check if zero gradient on outlet pressure doesnt produce a stable solution. 

January 24, 2006, 11:04 
I will try the limiter on the

#5 
Senior Member
Håkan Nilsson
Join Date: Mar 2009
Location: Gothenburg, Sweden
Posts: 193
Rep Power: 9 
I will try the limiter on the corrector on the laplacian scheme in the future. I have many cells with nonorthogonality above 70 degrees. I think that this will have impact on the quality on the results, but to start with the issue is to get rid of the convergence problem, which I don't think will be dependent on this limiter.
About setting a mean pressure outlet boundary condition as you suggest, what is the difference numerically between that and a pure Neumann bc? As I understand it the only difference is that the mean pressure at the outlet is set to the specified reference pressure instead of specifying the reference pressure in a point when using pure Neumann. This would only result in a difference in the level of the pressure of the final result. I think that the way the reference pressure is set is mostly important when doing unsteady simulations, where it should be set in a way that it does not affect the overall pressure level, i.e. don't set the reference pressure in a point where the pressure is likely to oscillate. If you do your whole pressure field will start oscillating and you will probably have slow computations due to the imposed global fluctuating pressure level. In my case the computation is steady, so it should not matter how I specify the reference pressure (as long as I don't say that the pressure should be constant at the outlet). My conclusion is that a Neumann outlet boundary condition for the pressure (together with a reference pressure point) should be sufficient, and I'm testing that right now. The reason why I didn't test this already is that I had problems with this for another case. Maybe it is too unstable to use from scratch? Now I use it when restarting the solution I show above, and it is still converging after some 700 iterations. I really appreciate your help Eugene! Håkan. 

January 25, 2006, 04:11 
Håkan,
Normally I will look

#6 
New Member
ShiuhHwa Shyu (Steven)
Join Date: Mar 2009
Location: Taiwan
Posts: 17
Rep Power: 9 
Håkan,
Normally I will look at the residual distribution and how it changes during the iterations. I guess that will give you more feeling about your problem. 

Thread Tools  
Display Modes  


Similar Threads  
Thread  Thread Starter  Forum  Replies  Last Post 
Recommended PC Spec for STARCCM+  Alex  CDadapco  12  November 8, 2013 02:41 
Recommended setup for flow over cylinder  rswbroers  OpenFOAM Running, Solving & CFD  19  September 27, 2007 07:03 
Recommended setup for the Singly rotated frame  jaswi  OpenFOAM Running, Solving & CFD  1  June 13, 2007 05:12 
Recommended PC for CFD  Paul  Main CFD Forum  10  March 13, 2006 16:16 
Recommended values on turbulence B.C.  Jihwan  CDadapco  2  August 10, 2005 07:45 