[clef1]

Hello everybody,

I am trying to run simulations on 2D pseudo mesh in OpenFoam for the calculation of drag coefficient around a cylinder for various Reynolds numbers. I want to validate my data with experimental data.

I created my mesh through SALOME and I created a template that contains the mesh in unv and constant/polyMesh folder (along with transport and turbulence properties files), 0/U,p,k,omega,nut and my system folder. Then I copy-paste this template folder and make the nscescary changes (in U,k, omega and controlDict, in drag coefficent calculation) for new simulations with different Reynolds number.

HOWEVER, when I do this, I keep getting the SAME results for drag coefficient Cd with my template case. I came to realize that the solution is the same for each mesh I created (i tried finer and coarser meshes).

Any ideas why this happen?

Thank you

[Tobermory]

Have a think about the physics, i.e. why does Cd vary as the Re number varies? What are the physical processes that are causing the changes in the wake and the drag? This should get you thinking about laminar/tubulent transition, flow separation etc. etc.

Now think about the simulation approach you are applying and ask youself: which of the above physical processes is it capturing? Clearly, the CFD code cannot reproduce real life if it is missing the physics. Good luck.

[clef1]

Originally Posted by Tobermory

Have a think about the physics, i.e. why does Cd vary as the Re number varies? What are the physical processes that are causing the changes in the wake and the drag? This should get you thinking about laminar/tubulent transition, flow separation etc. etc.

Now think about the simulation approach you are applying and ask youself: which of the above physical processes is it capturing? Clearly, the CFD code cannot reproduce real life if it is missing the physics. Good luck.

I might haven't explained clearly my problem. My problem is that the solution is fully mesh dependent. For example, if I use Mesh_1 for Re=10000, I will get cd = 0.70. If I now use the same Mesh_1, for Re=100000 (with exact same boundary and initial conditions for U,k,omega as the case in Re =10000. Also same fvschemes and fvsolution), I get the same cd=0.70.

Another approach of my problem is that for two Meshes (for example Mesh_1 and Mesh_2 for same geometry) and for same Re number, same files and initial/boundary conditions, I get different results. For example for case_1 with Re =10000, for Mesh_1 I will get cd=0.65 and for Mesh_2 I will get cd = 1.00.


I use simpleFoam for a range 10000 < Re < 500000. I have seen studies that use steady state solvers and simpleFoam for this range.

As for the physics behind the problem, can you please be more precise in what I miss?

Thank you,

[Tobermory]

No - you provided enough information in your first post, although the extra info is interesting. Note that this "simple case" is actually a pretty complicated flow to simulate accurately, so take heart.

The short answer is that the code is working fine, i.e. is doing what you are telling it to do. From your posts, I am assuming that you are running a k-omega closure model with a wall function on the cylinder surface. The wall function has no Re number dependence, and so you'd expect no variation in Cd with Re. You would expect variation in Cd with wall mesh resolution, of course - research wall functions to understand why.

Finally - on the question of physics, let me throw a question back at you - why does the Cd vary with Re in reality? What is causing this? Once you understand this, then you can start asking yourself, am I using the right modelling approach (turbulence model, RANS/URANS/LES, transition model, etc etc etc)? However, step 1 is always to get a suitably fine mesh that matches your wall treatment method.

[clef1]

In fact, when I try various variations of my mesh in the cylinder wall boundary (on inflation layers mainly), I get different results for same Re and initial conditions. So the chances are my problem can be solved based on your assumption for mesh refinement, especially on cylinder wall, and correct use of wall function. I will do a more analytical research on wall functions and revert.

For your reference, I use k-omega SST turbulence model with wall functions (kqRWallFunction) on cylinder wall patches.

As for the physics behind the problem, I have an Naval Engineering background and understand the basics (at least) of the problem. However, as I am relatively new in CFD, I try to figure out how those physics are translated in CFD terminology/methods.

Thank you for your reply

[clef1]

I have made some simulations since my last reply. Here are my comments:

From your posts, I am assuming that you are running a k-omega closure model with a wall function on the cylinder surface. The wall function has no Re number dependence, and so you'd expect no variation in Cd with Re. You would expect variation in Cd with wall mesh resolution, of course - research wall functions to understand why.

After some research, for k-omega SST model, i tried to keep y+ below 5 in order to avoid log-law wall functions. For k, on my cylinder wall I used fixedValue equal to 0 (as the kinetic energy is zero on the walls), for omega I used again a fixed small value and for nut for cylinder I tried two approaches. 1st approach is to use calculated for small y+ and the other is to use nutUBlendedWallFunction, which can be used along a big range of y+ values (including viscous sublayer, buffer and log law regions). Both seem to work well.

Finally - on the question of physics, let me throw a question back at you - why does the Cd vary with Re in reality? What is causing this? Once you understand this, then you can start asking yourself, am I using the right modelling approach (turbulence model, RANS/URANS/LES, transition model, etc etc etc)? However, step 1 is always to get a suitably fine mesh that matches your wall treatment method.

My simulations give much more better results in higher Reynolds numbers (for Re >= 500000) than in lower Reynolds number 10000 > Re >100000. After some research, this must be due to the transitient region for fully turbulent wake.
As I read in some papers, URANS is better in this case than RANS, but RANS also can give satisfactory results.

Is that right? And If so, what model do you suggest for simulations in that region (10000 > Re >100000)?

[Tobermory]

Yes - exactly, you are on the right track now. Note that the kW SST model assumes fully turbulent flow, and so if there is laminar->turbulent transition in your boundary layer (which occurs at lower Re#s), then you will not be able to pick that up (without including a transition model), and the wake width / Cd will be wrong.

As I said initially - at lo Re, this is a really difficult flow to model! Luckily, for most applications, the Re is high enough that the flow is fully turbulent and Cd is approx constant.

[clef1]

Another question/comformation regarding mesh. As Re gets bigger, y+ gets bigger too. So for simulations that have remarkable difference in Re, I cannot use the same mesh and same boundary conditions, as y+ will differ significantly.
That means, If I want to use the same boundary conditions (and same wall treatment), I need more refined mesh as Re grows, in order to have a y+ close to the one used in lower Re.

Is that correct?

And finally, do you know how wall functions implement regarding U? As I have seen so far, wall functions only exist for k,omega,nut for the k omega SST case.

[Tobermory]

Originally Posted by clef1

That means, If I want to use the same boundary conditions (and same wall treatment), I need more refined mesh as Re grows, in order to have a y+ close to the one used in lower Re.

Is that correct?

Yes.

And finally, do you know how wall functions implement regarding U? As I have seen so far, wall functions only exist for k,omega,nut for the k omega SST case.

The nut wall function sets the wall shear stress, which acts as a wall function for velocity.

[clef1]

Thank you for your help and fast reply,
Your previous guidelines have helped me a lot in my simulatiouns.