[gquixley]

First of all, thank you for visiting this thread. This is my first post to these forums so I will be as thorough as possible.

Objective
I am currently running a series of URANS 3D Cylinder simulations inside OpenFOAM V17.12 using the k-omega SST and SSTLM models to observe the transitional effects on the separation point at Re=33,000.

What has been done already
2D simulations have been completed to verify the SST and SSTLM models are working as expected. A steady 3D simulation has been completed with a reported Cp value at the stagnation point equal to 1.

What is the problem?
The problem I am seeing is the mean extracted pressure coefficient(Cp) around the cylinder is higher than 1 for an unsteady SST simulation.

This Cp value above one, as many of us were taught, should not be possible in incompressible flows however Issa wrote in his technical note, "Rise of Total Pressure in Frictional Flow", this assumption is in fact erroneous. The rise in pressure at the stagnation point can be attributed to the stress gradients locally contributing at the expense of depleting the total pressure elsewhere.

Now you may be thinking at this point, you've just answered your own question - you fool. However, Issa conducted two flow configurations and reported Cp values up to 4.8% higher than 1. I am observing Cp values in the range of 10-20% higher. Has anyone seen or experienced this before?

After conducting a literature review around the subject I stumbled across several papers with similar issues to my own:

• Benim - Modeling turbulent flow past a circular cylinder by RANS, URANS, LES and DES
• Wenjun - Numerical investigation of flow around one finite circular cylinder with two free ends

There are several other papers but these two show the point I am trying to illustrate. Both papers report Cp values above 1 at the stagnation point however neither report the reasons attributed to it.

The setup
controlDict: pisoFoam

fvSchemes:

• ddtSchemes: backward;
• gradSchemes: default -> cellLimited Gauss linear 1; grad(U) -> Gauss linear;
• divSchemes: Gauss linear;
• laplacianSchemes: default -> Gauss linear limited 0.5;
• interpolationSchemes: deafult -> linear;
• snGradSchemes: default -> limited 0.5;
• wallDist: method -> meshWave;

fvSolution solvers:

• p -> PCG with DIC preconditioner
• U, k, omega -> PBiCGStab with DILU preconditioner
• PISO -> nCorrectors = 2, nNonOrthogonalCorrectors = 1

I have run these cases from an initial condition using a potential flow to initialize the simulation then a steady formulation for 200 or so iterations before switching to an unsteady formulation. This was done for a coarse mesh resolution then using the mapFields function, mapped to finer meshes. The resulting Cp>1 appears in all solutions regardless of mesh density.

At this point I am reaching out because I cannot find any literature addressing this issue and am looking for peoples experience and if you notice anything wrong in the set up. I would be happy to provide follow up information for anyone interested and finally thank you for reading through this whole post.

Kind regards,
GQ

[flotus1]

My intuition tells me that cp>1 might be possible in unsteady flows.
cp is normalized to return a value of 1 in the case where all kinetic energy of the free stream is transferred into pressure.
With periodic vortex shedding at a cylinder, we have a situation where some rotational component can get added to the free stream velocity.
Imagine the flow is at a stage where the wake directly behind the cylinder is pushed downwards. Now the vortex "snaps back" into the opposite direction, which causes some rotation of the fluid around the cylinder. This rotational component is added to the free stream velocity and might temporarily lead to cases of cp>1. Not sure if this also means that the time-averaged cp could be higher than 1. Or how large that effect might even be.
That's what I would draw on a napkin, please don't quote me on that

[gquixley]

This Cp > 1 is occurring at, what would be, the stagnation point. At this Re, I'm running a cylinder with an effective diameter of 0.5m and velocity = 1m/s. Do you believe the rotation of fluid is strong enough to contribute to the free stream velocity and in turn increase Cp?

The complexity of this simple case is impressive. I ran another 3D Cylinder simulation without free-ends by imposing a periodic boundary condition. Once I added free-ends and replaced the periodic boundary condition with a farfield condition I started to see Cp > 1...

[FMDenaro]

Just a question: how do you assume a “stagnation point” in viscous (and dissipative) flows?

[gquixley]

Originally Posted by FMDenaro

Just a question: how do you assume a “stagnation point” in viscous (and dissipative) flows?

A stagnation point is the point in the vicinity of a solid body, in this case the cylinder, where the streamlines split. I suppose the assumption is made for this incompressible case that there will be a stagnation point of near zero velocity in direction of flow, that would typically correspond to a Cp value equal to 1.

Could the value above 1 mean there is no stagnation point? Or as noted in Issa's technical note could this be a localized effect?

[gquixley]

Not sure if this is a trick question or not... I'm assuming there will be a stagnation point as with most scenarios. The stagnation point being the region where streamlines split. This is a generic answer but I don't know if this changes for a viscous and dissipative flow.

Are you suggesting that there won't be a stagnation point thus the Cp value exceeds 1?

[FMDenaro]

Originally Posted by gquixley

Not sure if this is a trick question or not... I'm assuming there will be a stagnation point as with most scenarios. The stagnation point being the region where streamlines split. This is a generic answer but I don't know if this changes for a viscous and dissipative flow.

Are you suggesting that there won't be a stagnation point thus the Cp value exceeds 1?

The assumptions for having the Cp<=1 are: incompressible and inviscid steady flows. Just for a compressible flow you can have Cp>1.

On the other hand, defining one "stagnation point" for viscous flows is not well posed. It is not sufficient to see a split in the streamlines to have a stagnation point.

[LuckyTran]

If this is incompressible pisoFoam then we are forced to rule out all attempts to explain it using compressible phenomenon.

I think it's worth some time to discuss how you actually calculate the total pressure or get the Cp, even if you just use whatever the software gives you. These formulae involved almost always assume constant specific heats and don't actually do the proper integration as per theory. So I can see Cp and Total pressure being a little wonky if you have temperature dependent specific heats in your thermophysical properties for example... It's very common that there is an inconsistency involved.

But yeah, total pressure doesn't increase unless there's boundary work being done. This is numerical CFD, so I'm always willing to tolerate a little % error (whether + or -).


Furthermore, due to the way algebra works, there is a slight difference between calculating cp at every time-step and averaging in time, or just averaging the pressure field and calculating the cp from this already averaged pressure and velocity.

[gquixley]

Thank you for your responses.

To clarify, I agree with you LuckyTran with the fact this is an incompressible pisoFoam simulation so there should be no compressible effects. To my knowledge there are on temperature components (for this simulation, in actuality I'm sure there would be).

The way in which I am calculating Cp is according to the equation:
Cp = (p-p0)/0.5*(rho*U^2)

Due to rho being included in the pressure correction step of the Piso algorithm I am not diviving by rho, also p0 is set to 0 throughout the simulation and the velocity used is 1 m/s. I have conducted this equation on both instanteous values of p and averaged over multiple seconds of simulation. The instantaneous values always provided Cp much higher than 1 whereas the averaged results were closer to one but still much above the expected 1. Thus the above equation boils down to:
Cp=2*p

I have been running a few tests to observe what happens to Cp when certain factors change. I first altered the Reynolds numbers to see how that would effect the results and there was virtually no change in the Cp, ~1.2, although the turbulent kinetic energy changed at the same magnitude as with the Reynolds. Next I ran just K-Omega to see if this was just a localised issue and the Cp came back in significantly over 3, so there is something happening there too. I have tried with different mesh densities and even using pimplefoam to little success... Currently testing different turbulent intensity values and changing the ddtScheme from backward to CrankNicolson (spoiler: the simulation crashed)

All in all that's where I'm at for the moment. If anyone can try to reproduce this then that would be great! I simply made a 1.5 aspect ratio cylinder using the above mentioned settings and am reporting Cp values around 1.2.

[LuckyTran]

I hope you remembered that p and U in pisoFOAM is p/rho and U/rho.

We can appreciate the fact that you can easily calculate the total pressure field, plot it, look at it, and see where the total pressure is increasing and confirm whether this is consistent with your knowledge of how the pressure coefficient is behaving.

Nobody wants to hear about a simulation crashing at this stage.... because now we have to question whether anything was ever converged, ever.