code_saturne User's Forum

Hello.
Some time ago I noticed a problem with Saturne regarding any calculations with inflation mesh layers and k-omega turbulence model. First I encountered it in 5.x branch and recently I decided to check if it's still in 6.x versions so I compiled version 6.0.0 and gave it a try with 2 test cases (I read readmes and didn't find any changes that can affect the issue in current official version 6.0.1).

As I can see, the problem shows in any pure-aerodynamic isothermal case with inflation layers and k-omega model. I tried different settings in all cases but without success. The calculation flow is as follow. You start with default setting, add mesh, define inlet with given mass flow rate and an outlet. Set time stepping to spatial-variable with target Courant number of 5...10 as usual. Just standard basic things. Maybe change gradient reconstruction for stability. Then make first run with Upwing scheme for all variables and it's just fine. The you restart with SOLU with blend~0.8 for velocity and Upwind for turbulence. If there are no layers + k-epsilon it will be OK. But, if you have layers + k-omega, calculation looks great just after restart, then after some to tens of iteration pressure-velocity field start do diverge producing "cosmic" velocities. I didn't find any combination of settings that can solve this issue including pressure increase relaxation or pressure-velocity coupling option. If you let it calculate further it will never converge.
What can be the cause? How can I solve this? It's not a problem of particular geometry or velocity level or mesh. Even if the mesh is close to ideal and velocities are low (for example, 20 m/s in bulk) the issue is there.

As an example I attach the picture with diverged velocity field if GT snail duct. You can see that the mesh is not good in prism-tetra transition zone but it's near impossible to have this transition ideal in practical cases. In the snail itself velocities rise up to ~500 m/s that is correct (I compared to CFX that is our mainstream tool at work). But, as you can see, velocity becomes unrealistic high at the top part of the duct and it will never converge with "usual" simple settings.

I also prepared two test cases if someone will kindly check what may be wrong with setup.
GT snail case corresponding to the picture: https://www.dropbox.com/s/ew96d4shoxrty ... t.tgz?dl=0
Pipe bend with low velocities: https://www.dropbox.com/s/7obr8wfbt4f6i ... t.tgz?dl=0
There are complete cases without postprocessing dirs to reduce the size.

Will somebody, please, give some advice for this issue? Sometimes it's very important to have inflation layers so I need to make them work. Thanks for your attention.

Hello,

I'll try to take a look at this right after the holidays. There were a few changes in k-omega boundary condition handling, with an option to revert to the previous behavior and new options for an all-y+ variant, so it may be useful to check their influence (but if this is the cause of your issues, I would expect different behavior in v5 and v6).

Il'll certainly need Jean-Francois's help for the turbulence aspects, but I can take a first look. I did not download the case yet (answering from my phone), but am wondering whether these are incompressible or compressible flows (the latter may be less tested).

Best regards,

Yvan

Hello. Sorry for the late response.

Behaviour is similar in versions 5 and 6, but may be different in version 4. I ran a simple test case with layers and k-omega with success in v4 but it was years ago so I don't remember if I let it to iterate enough (anyway v4 is too old). It's important not to stop after 200-300 iterations (with numerical parameters from test cases above) because the divergence often starts indeed after 200...300 iterations (it's quite common for some kinds of solution divergensions in different software to start after initial, rough flow pattern has formed and more detail flow features start to develop). Regarding this issue in Saturne, I noticed that it's usually a near-wall region (but not completely on the wall) where velocity "runs away".

All cases are incompressible isothermal with constant density with different velocity levels so it's quite strange for me to see the issue. May be I set something wrong but I don't know what.

Thanks for your response, hope you will take a look after holidays.

Hello,

I will have a look at your test case.
If I understand well it works for k-eps without layers but does not work for k-w SST with layers.
First, could you tell me if it works with k-w SST without layers ?

Regards,

Jean-Francois

Hello, thanks for your response.

it works for k-eps without layers but does not work for k-w SST with layers
Yes.

But I found one strange thing with the pipe bend case. To answer your question about the mesh effect I ran the same case without layers, there was the same issue with k-omega. But, what is very strange for me, id didn't disappeared when I switched to k-epsilon (without layers). Velocities rise to unrealistic values near the outlet (standard zero pressure outlet). It's a simple pipe, nothing special, and there are no parameters for outlet BC in Saturne GUI, other options are as usual, Saturne version is 6.0.0 (8 threads), divergence becomes obvious around iteration 160.

I attached a picture for the bend case with velocity profile after solution has diverged (k-omega without layers), solver log for the k-epsilon (simplest) variant (run_solver.log) and an XML file for the k-epsilon variant (Bend-G000-NL-KE.xml).

Tetra mesh for the bend case without layers: https://www.dropbox.com/s/0al0634150lft ... .cgns?dl=0

If you want to reproduce an issue with any turbulence model and mesh without layers, download the mesh, open appropriate XML file (now you have exact test XMLs for k-omega and k-epsilon variants), replace the mesh in case properties, let it iterate for 50...100 iterations with upwind scheme, then switch to SOLU with blend 0.8...1.0 for velocity (upwind for turbulence) and wait for approximately 150...300 iteration depending on the model you chose. You will see that periodical maximum velocity spikes appear, then maximum pressure (and velocity) constantly rises up to unrealistic values

Hello,

Could you test adding an extrusion of the outlet over 1 or 2 layers, so that the outlet cells are "orthogonal" ? We know there are robustness issues the standard outlet BC when the cells are irregular at the outlet (and have some ideas to work around this, but need to test them).

Best regards,

Yvan

Ivan, I just tried to add the othogonal layer for the outlet using "Mesh extrusion" GUI feature for the OUT face group. Extrusion itself works but when the solver starts it fails because of incorrect boundary type: The workflow is very simple: use "Bend-G000-NL-KE.xml" case, add mesh extrusion, set OUT as selector. I attached the solver log. I think it's because OUT group is now includes interior faces: How can I work it around?

Hello,

Interior faces with boundary groups are ignored, as BC checks only loop on boundary faces.

By the coordinates of the face in error or the ERROR.ensight (or .med or .cgns) output, is the first face causing problems on the outlet ?

Best regards,

Yvan

Hello, Ivan.
I checked the error output closer. First, in ERROR.case, there is a couple of quad faces at the outlet part created by the Mesh extrusion feature (Mesh => Preprocessing => Other => Mesh extrusion) that have BC error. I attached the picture from ParaView (Bend-Extrusion-Error.png), faces with BC error form the blue stripe (actually there are two faces, I see it with the Multiblock Inspector if I select "Faces with B.C. error"). Also, in the listing, I see that the solver has added 106 faces to boundary_faces zone but only 104 faces to WallSide group (originally, in the running case, it's one of two pipe wall face groups: one wall group is located around the bend, second is at straight parts, for different roughness in our model case). It seems that the solver, although smart enough to modify boundaries, "forgets" about these two faces on extruded part so they has no parameters set and become errornous BCs. Looks like a bug... Could you, please, check this at your side to confirm? I provided the case XML and tetra mesh, you can test on any machine because the mesh is very small.

Yes, if you disable mesh extrusion, you will get divergence at the outlet, but the case will run without BC problems, so the mesh itself is correct. I checked it just now, the only thing I changed was disabled extrusion and the solver started fine.

Update
=====
I tried to play around with extrusion parameters (reduce cell size and number of extrusion layers from 2 to 1), the result is the same, you can see it on the picture №2 (Bend-Extrusion-Error-2.png).

Hello,

Is the case the same as the one you posted ? I can check.

Otherwise, in the mean time, you can always add "or no_group[]" to the wall faces selection criteria, or add a separate BC with 'no_group[]' so as to work around the undefined BC's error and see if at least the divergence is avoided.

Best regards,

Yvan