[eugene]

Well at a glance, you shouldnt be using nonOrthogonal corrector steps in a steady state calculation.

Try to make your wall velocities tangential. You dont want flow coming in and/or going out since this will interact in unpredictable ways with the other variables (p,k,epsilon).

I am also trying to run one swirling flow case and unable to get converged solution. Is it recommended to run first potential Foam and after that start simpleFoam? I tried this trick to get good initial guess but velocity field after running potential foam was very very bad. Did I do any mistake? Any help?

Regards
GS

[hani]

Ghanshyam,

I get a good through-flow solution with potentialFoam in just about no time. Of course, I don't get the swirl from the potential solution.

To use potentialFoam with my present case I had to add to fvSchemes:
laplacian(1,p) Gauss linear corrected;

Regarding nonOrthogonalCorrectors, my experience with openFoam says that it should be used even with steady state case. With three nonOrthogonalCorrectors, it could run for some more iterations!!!! before diverging...Yes, so far I could not get converged solution for a swirling flow case.

Regards
GS

[henry]

nonOrthogonalCorrectors are not required for the SIMPLE p-U solution algorithm used in the steady-state codes but are available as an option but I have never found a case for which they are beneficial and often they are detrimental to the overall convergence rate.

[hjasak]

Sorry Eugene, I don't agree with the advice of not running non-orthogonal correctors. You should really strive to get a decent pressure field before you correct the velocity and if the mesh is bad, you're better off doing one non-orth corrector in SIMPLE. I have found a single corrector considerably improves the robustness and convergence in quite a lot of cases.

Also, I am pretty worried about not being able to get a decent potential flow solution. BTW, this one also requires non-orthogonal correctors and this is the first thing I would investigate. It will give you hints anout errors in the mesh, setup, boundary conditions etc.

Could you please tell us how bad this mesh is: non-orthogonality, skewness, nastiness near the wall, all-hex, all-tet or similar. Maybe some output from checkMesh would come in handy.

This should be bread-and-butter kind of problem for FOAM and the solver (in my experience) is fast, robust and accurate.

Hrv

I also added laplacian(1,p) Gauss linear corrected in fvScheme, still I dont get satisfactory solution. One more thing, potentialFoam over writes the "0" and it looks like in controlDict, for potentialFoam, endTime and writeInterval are pseudo. Am I correct?

Regards
GS

checkMesh gives following message:

----------------------------------
Number of non-orthogonality errors: 0. Number of severely non-orthogonal faces: 120.
Mesh non-orthogonality Max: 75.6717 average: 13.8299
Non-orthogonality check OK.

Writing 120 non-orthogonal faces to set nonOrthoFaces

Face pyramids OK.

Severe skewness for boundary face 2368028 skewness = 2.01081
--> FOAM Warning :
From function primitiveMesh::checkFaceSkewness(const bool report, labelHashSet* setPtr) const
in file meshes/primitiveMesh/primitiveMeshCheck.C at line 838
Large face skewness detected. Max skewness = 201.081 percent.
This may impair the quality of the result.
1 highly skew faces detected.
Writing 1 skew faces to set skewFaces

Minumum edge length = 1.30583e-05. Maximum edge length = 0.00244646.

All angles in faces are convex or less than 10 degrees concave.
----------------------------------

Thanks and regards
GS

[eugene]

To make sure your wall velocity is only tangential to the face, you can dot the boundary face unit normal vector with the velocity and subtract that from the specified velocity.

I.e.

const vectorField& faceN = mesh.boundary()[hubPatchID].nf()

U.boundaryField()[hubPatchID] -= faceN*(faceN & U.boundaryField()[hubPatchID]);

should do the trick.

Try to improve the quality of your mesh. In partcular, look out for distorted high aspect ratio cells (especially tets and pyramids if you have them). I have found these to be extremely detrimental to obtaining a steadyState solution.

Another trick is to take the difference between the fields at two timesteps to see where variations are occurring. Typically this will be at or near the problem, which will help with the diagnosis.

[hani]

Thanks for the help with the coding, I will try to implement it.

In this case I am unable to modify the grid since it is pre-determined by workshop organizers. I basically have to use this grid. There are however only hexahedral cells, and no tets and pyramids. I attach extracted results from checkMesh below, if you can help me interpret them.

Is there a way to visualize the distribution of local (cell) residuals directly (in paraFoam for instance)? Variations occur basically everywhere, but mostly in a region where I know that there is a risk for unsteadiness. One solution to this problem would be to say that the flow is by nature unsteady, and move on to turbFoam. However, I know many people who have been able to get a steady solution for this case with other codes. What settings should I use to 'force' the solution to be steady? However, at the end it is actually only healthy that a code doesn't generate a steady solution if there is none.

Here comes some extracted information from checkMesh (things that are not OK):
Severe non-orthogonality for face 982691 between cells 326306 and 351802: Angle
= 70.2071 deg.
Number of non-orthogonality errors: 0. Number of severely non-orthogonal faces:
117.
Mesh non-orthogonality Max: 77.9564 average: 19.0842
Non-orthogonality check OK.

Writing 117 non-orthogonal faces to set nonOrthoFaces

Severe skewness for boundary face 2947349 skewness = 2.09526
Severe skewness for boundary face 2947350 skewness = 2.23004
Severe skewness for boundary face 2947461 skewness = 2.09265
--> FOAM Warning :
From function primitiveMesh::checkFaceSkewness(const bool report, labelHashS
et* setPtr) const
in file meshes/primitiveMesh/primitiveMeshCheck.C at line 838
Large face skewness detected. Max skewness = 223.004 percent.
This may impair the quality of the result.
3 highly skew faces detected.
Writing 3 skew faces to set skewFaces

Failed 1 mesh checks.

Many thanks!

[hjasak]

Is this a publicly available geometry? If so, could you please tell me where to get it, I'd like to have a go. It would be even better if you've got the case set up to run :-)

Please let me know,

Hrv

[eugene]

You could try just pure upwind convection instead of linearLimited. Admittedly, I dont know the details of linearLimited, so I dont know whether it would improve matters.

[hani]

That is a great idea, Hrvoje! I am supposed to present my results for a number scientists and industrial people from all over the world, and I would really like to make some good marketing for OpenFoam. I think that OpenFoam has a huge potential, as soon as I learn how to use it. If you like, we could make a common contribution to the workshop.

You can find all the details about the case at www.turbine-99.org, and you can find the full OpenFoam case at the URL below. I have attached only the result from the last iteration, which you can consider as boundary conditions and a good initial guess for the internal field. Note that the problem with a small flow through the hub still remains. The code for including the boundary conditions is not included. If you need it I can give it to you, but I am working on improving it right now.

http://www.tfd.chalmers.se/~hani/out...rvojeJasak.tgz

About using pure upwind: I started off doing that, and I think that the flow was more stable (steady), but the levels of the residuals was similar to what I have described. Can I consider a computation to be converged if there are no more changes, although the residuals are not very small?

[amelia]

Hello,

I also get oscillations in my results when considering higher Reynolds, although the case is quite simple.

http://www.cfd-online.com/OpenFOAM_D...tml?1125888157

I am thinking of computing the average velocities over time and see if it converges. Do you think this is correct?

[hjasak]

Hi Hakan,

This looks like a decent mesh and you've got a recirculation at the outlet (from what I can see), so this should be easily fixed.

However, I would be happier to do the whole simulation myself so I would appreciate the code for the boundary conditions. Could you tell me what the inlet velocity and turbulence should be - I suspect I can do this pretty quickly :-)

Hrv

[hani]

Hrvoje,

I am really glad that you try to help me out! I also know other people who are interested in getting good results with OpenFoam for this case.

You may now download the boundary condition code at:
http://www.tfd.chalmers.se/~hani/outgoing/setT99BC.tgz

It uses the Uin and boundaryProperties files in constant to interpolate U, k and epsilon. It also adjusts for the correct volume flow. Feel free to improve it, I know that there is plenty of room for improvement. In the future I am planning to incorporate this with turbFoam to allow unsteady boundary conditions, but I have noticed when incorporating it with simpleFoam, that I cannot longer run in parallel and that I can no longer get 'clean' restarts without discontinuities in the residuals. These problems must be solved in the future.

I have tried to include the two lines by Eugene de Villiers (see the discussion) at the end, but at least when I compute the flow through the hub in Ensight, I still have the same amount of flow through that wall. This problem should also be solved in the future.

I just realized that checkMesh generates two files, nonOrthoFaces and skewFaces. Are these required for non-orthogonal correctors to work properly? I haven't yet had time to test if there is a difference.

Thanks a lot!
Håkan.

[hani]

Maybe you would also like the 0 time directory, where I set initial guesses for k and epsilon. Find it at:
http://www.tfd.chalmers.se/~hani/outgoing/0.tgz

The recirculation at the outlet that you mentioned is very small. I have seen results where maybe 1/4 of the upper right corner (looking in the streamwise direction, z-axis pointing up) of the outlet is recirculating. But let's start by letting OpenFoam decide what is the 'correct' solution.

Håkan.

[eugene]

Hmm, the lines I sent you should have worked. Try the following change:

(vectorField&) U.boundaryField()[hubPatchID] -= faceN*(faceN & U.boundaryField()[hubPatchID]);

And let me know if it makes any difference.

Note, unless your wall velocity changes with time, you only need to do this once.

[hani]

There is a difference, but not a good one. First of all, I got the following message while compiling:
warning: use of old-style cast

The flow trough the wall:
Without your modification: -2.1809e-7
First modification: -2.1809e-7
Second modification: -2.7649e-7

Note that I get these flow values from Ensight, so there might be a problem with post-processing. I would like to output it directly from my application, since then I get exactly what OpenFoam 'sees'. I could do something like:

scalar Qhub = -sum(phi.boundaryField()[hubPatchID]);
Info << "Volumetric flow through hub wall = " << Qhub << " m3/s " << endl;

But then I must be sure that phi is updated correctly at the hub wall after I have used your correction. Can you give me a one-liner for this purpose?

[eugene]

You are getting the warning because
(vectorField&) is depricated method of casting. If you want the warning to go away check the gcc-4's new casting rules:

static_cast<vectorfield&>

should to the trick I think.
You dont have to loop the faces, you can simply go
phi.boundaryField()[hubPatchID] = U.boundaryField()[hubPatchID] & mesh.boundaryField()[hubPatchID].Sf()

But then this shouldnt be necessary. Simply modify your U field before you construct phi and you should be set.