[daniel_mills]

We have a need to run a pressure\density based compressible steady state solver. I have been fighting with rhoSimpleFoam for a few weeks now.

We would like to be able to specify:
2 inlets with TotalPressure using a fixed pressure.
2 outlets with possibilty of reverse flow and fixed static pressure.

I am currently using totalPressure with fixed values for the inlets and outletInlet for the outlets in the p file. I am using pressureInletVelocity for the inlets and inletOutlet for the outlets with the values set to 0 for the U file.

The flow field will develop a swirl as this runs.

The geometry checks fine using checkMesh.

This also runs fine in simpleFoam developing a fully converged solution.

For some reason when we run it with the pressure\density solver like rhoSimpleFoam we get an instability in pressure\density after a few hundred steps.

I have been able to change the instability by playing with the "phi -= pEqn.flux()" statement in the the code, but never actually get rid of it. I think the problem may reside pressure\density equation for the SIMPLE solvers.

Is OpenFoam capable of solving this type of problem?

HAS ANYONE BEEN ABLE TO SOLVE A PROBLEM SIMILAR TO WHAT I HAVE DESCRIBED WITH A STEADY STATE SOLVER?

[daniel_mills]

I have noticed when someone asks a question about stability problems in rhoSimpleFoam or density correctors in steady state solvers they are never answered. Here are two other posts posing similar questions from other users that never got further response. If there is a problem with rhoSimpleFoam thats alright, but it would be nice to get some sort of response?

Just tell me there is a problem and no one currently has a solution.

http://www.cfd-online.com/OpenFOAM_D...es/1/2068.html

http://www.cfd-online.com/OpenFOAM_D...es/1/2718.html

[hjasak]

There's no problem in the solver. What you are seeing is that compressible flows are much more sensitive to solver setting, initial field, time-step and boundary conditions. Unless you get all the components right together, the solver will blow up in your case.

All in all, this just requires lots of care.

Hrv

[daniel_mills]

Alright I can easily believe that the problem is in my setup somewhere. Maybe someone would be so kind as to spot where my potential problems are. I will post my checkMesh, fvSolution, fvSchemes, and my 0 files. Please anyone with any ideas let me know what else I should try. No matter what I try I seem to get what looks very much like a low wave mode pressure instability after about 500-600 steps.

Remember this mesh converges in simpleFoam. I also get much more stability by decreasing the pressure drop or increasing viscosity.



checkMesh:

Number of points: 124105
edges: 649440
faces: 979095
internal faces: 920158
cells: 453762
boundary patches: 19
point zones: 0
face zones: 0
cell zones: 0

Checking topology and geometry ...

Point usage check OK.

Upper triangular ordering OK.

Topological cell zip-up check OK.

Face vertices OK.

Face-face connectivity OK.

Basic topo ok ...

Checking patch topology for multiply connected surfaces ...

Patch Faces Points Surface

Periodic-Center 816 875 ok (not multiply connected)
Periodic-Combustor 16178 10957 ok (not multiply connected)
Periodic-Holes 328 204 ok (not multiply connected)
Periodic-Plenum-2 2792 1522 ok (not multiply connected)
Periodic-Plenum-1 2046 1095 ok (not multiply connected)
SmallTubes-2 2538 1324 ok (not multiply connected)
LargeTubes-2 3107 1707 ok (not multiply connected)
SmallTubes-1 2078 1092 ok (not multiply connected)
LargeTubes-1 6047 3191 ok (not multiply connected)
Wall-6 1032 687 ok (not multiply connected)
WallSet-5 7102 3751 ok (not multiply connected)
WallSet-4 3742 2016 ok (not multiply connected)
WallSet-2 4455 2328 ok (not multiply connected)
Wall-3 860 957 ok (not multiply connected)
Wall-1 742 409 ok (not multiply connected)
BC-PressOut-center 95 113 ok (not multiply connected)
BC-PressOut-1 3088 2180 ok (not multiply connected)
BC-PressInlet-2 1048 568 ok (not multiply connected)
BC-PressInlet-1 843 461 ok (not multiply connected)


Patch topo ok ...

Topology check done.

Domain bounding box: min = (-1.42109e-17 -2.75474e-51 0) max = (0.0757305 0.0245345 0.0245345) meters.

Checking geometry...

Boundary openness in x-direction = -4.32336e-18
Boundary openness in y-direction = -8.77453e-18
Boundary openness in z-direction = -3.35318e-18

Boundary closed (OK).

Max cell openness = 4.27239e-22 Max aspect ratio = 2.58574. All cells OK.

Minumum face area = 1.86851e-09. Maximum face area = 1.79996e-06. Face area magnitudes OK.

Min volume = 3.02777e-14. Max volume = 1.42021e-09. Total volume = 1.34424e-05. Cell volumes OK.

Mesh non-orthogonality Max: 65.6633 average: 17.5979

Non-orthogonality check OK.

Face pyramids OK.

Max skewness = 119.41 percent. Face skewness OK.

Minumum edge length = 5.08262e-05. Maximum edge length = 0.00214438.

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

Face flatness (1 = flat, 0 = butterfly) : average = 0.999994 min = 0.999846

All faces are flat in that the ratio between projected and actual area is > 0.8

Geometry check done.

Number of cells by type:

hexahedra: 19560
prisms: 44800
wedges: 0
pyramids: 285
tet wedges: 0
tetrahedra: 389117
polyhedra: 0

Number of regions: 1 (OK).

Mesh OK.


fvSolution: (note I update the relaxation values to .2 on p and .6 on everything else after about 150 steps. It continues to run fine for another 400 steps)


solvers

{
p AMG 1e-06 0 100;
U BICCG 1e-06 0;
h BICCG 1e-06 0;
k BICCG 1e-08 0;
epsilon BICCG 1e-08 0;
}

SIMPLE

{
nNonOrthogonalCorrectors 0;
pRefCell 0;
pRefValue 0;
}

relaxationFactors

{
p 0.05;
rho 0.05;
U 0.3;
h 0.1;
k 0.05;
epsilon 0.05;
}




fvSchemes:

ddtSchemes

{
default steadyState;
}

gradSchemes

{
default Gauss linear;
}

divSchemes

{
default none;
div(phi,U) Gauss limitedLinearV 1;
div(phi,k) Gauss upwind;
div(phi,epsilon) Gauss upwind;
div(phi,h) Gauss limitedLinear 1;
div((muEff*dev2(grad(U).T()))) Gauss linear;
}

laplacianSchemes

{
default none;
laplacian(muEff,U) Gauss linear limited 0.5;
laplacian((rho|A(U)),p) Gauss linear limited 0.5;
laplacian(alphaEff,h) Gauss linear limited 0.5;
laplacian(DkEff,k) Gauss linear limited 0.5;
laplacian(DepsilonEff,epsilon) Gauss linear limited 0.5;
laplacian(mut,U) Gauss linear limited 0.5;
}

interpolationSchemes

{
default linear;
interpolate(U) linear;
}

snGradSchemes

{
default corrected;
}

fluxRequired

{
default no;
p;
}


Epsilon: (note I crash almost immediately when I run laminar. I also get increased stability by fiddling with viscosity.)

dimensions [0 2 -3 0 0 0 0];

internalField uniform 1.5;

boundaryField

{
Periodic-Center
{
type cyclic;
}
Periodic-Combustor
{
type cyclic;
}
Periodic-Holes
{
type cyclic;
}
Periodic-Plenum-2
{
type cyclic;
}
Periodic-Plenum-1
{
type cyclic;
}
SmallTubes-2
{
type zeroGradient;
}
LargeTubes-2
{
type zeroGradient;
}
SmallTubes-1
{
type zeroGradient;
}
LargeTubes-1
{
type zeroGradient;
}
Wall-6
{
type zeroGradient;
}
WallSet-5
{
type zeroGradient;
}
WallSet-4
{
type zeroGradient;
}
WallSet-2
{
type zeroGradient;
}
Wall-3
{
type zeroGradient;
}
Wall-1
{
type zeroGradient;
}
BC-PressOut-center
{
type inletOutlet;
inletValue uniform 1.5;
}
BC-PressOut-1
{
type inletOutlet;
inletValue uniform 1.5;
}
BC-PressInlet-2
{
type fixedValue;
value uniform 1.5;
}
BC-PressInlet-1
{
type fixedValue;
value uniform 1.5;
}


k:

dimensions [0 2 -2 0 0 0 0];

internalField uniform 0.055;

boundaryField

{
Periodic-Center
{
type cyclic;
}
Periodic-Combustor
{
type cyclic;
}
Periodic-Holes
{
type cyclic;
}
Periodic-Plenum-2
{
type cyclic;
}
Periodic-Plenum-1
{
type cyclic;
}
SmallTubes-2
{
type zeroGradient;
}
LargeTubes-2
{
type zeroGradient;
}
SmallTubes-1
}
type zeroGradient;
}
LargeTubes-1
{
type zeroGradient;
}
Wall-6
{
type zeroGradient;
}
WallSet-5
{
type zeroGradient;
}
WallSet-4
{
type zeroGradient;
}
WallSet-2
{
type zeroGradient;
}
Wall-3
{
type zeroGradient;
}
Wall-1
{
type zeroGradient;
}
BC-PressOut-center
{
type inletOutlet;
inletValue uniform .055;
}
BC-PressOut-1
}
type inletOutlet;
inletValue uniform .055;
}
BC-PressInlet-2
{
type fixedValue;
value uniform .055;
}
BC-PressInlet-1
{
type fixedValue;
value uniform .055;
}


p: (note I have also tried totalPressure and fixed value for the outlet)

dimensions [1 -1 -2 0 0 0 0];

internalField uniform 102928;

boundaryField
{
Periodic-Center
{
type cyclic;
}
Periodic-Combustor
{
type cyclic;
}
Periodic-Holes
{
type cyclic;
}
Periodic-Plenum-2
{
type cyclic;
}
Periodic-Plenum-1
{
type cyclic;
}
SmallTubes-2
{
type zeroGradient;
}
LargeTubes-2
{
type zeroGradient;
}
SmallTubes-1
{
type zeroGradient;
}
LargeTubes-1
{
type zeroGradient;
}
Wall-6
{
type zeroGradient;
}
WallSet-5
{
type zeroGradient;
}
WallSet-4
{
type zeroGradient;
}
WallSet-2
{
type zeroGradient;
}
Wall-3
{
type zeroGradient;
}
Wall-1
{
type zeroGradient;
}
BC-PressOut-center
{
type fixedValue;
value uniform 101325;
}
BC-PressOut-1
{
type fixedValue;
value uniform 101325;
}
BC-PressInlet-2
{
type totalPressure;
p0 uniform 104531;
value uniform 104531;
}
BC-PressInlet-1
{
type totalPressure;
p0 uniform 104531;
value uniform 104531;
}


T:

dimensions [0 0 0 1 0 0 0];

internalField uniform 300;

boundaryField
{
Periodic-Center
{
type cyclic;
}
Periodic-Combustor
{
type cyclic;
}
Periodic-Holes
{
type cyclic;
}
Periodic-Plenum-2
{
type cyclic;
}
Periodic-Plenum-1
{
type cyclic;
}
SmallTubes-2
{
type zeroGradient;
}
LargeTubes-2
{
type zeroGradient;
}
SmallTubes-1
{
type zeroGradient;
}
LargeTubes-1
{
type zeroGradient;
}
Wall-6
{
type zeroGradient;
}
WallSet-5
{
type zeroGradient;
}
WallSet-4
{
type zeroGradient;
}
WallSet-2
{
type zeroGradient;
}
Wall-3
{
type zeroGradient;
}
Wall-1
{
type zeroGradient;
}
BC-PressOut-center
{
type inletOutlet;
inletValue uniform 300;
}
BC-PressOut-1
{
type inletOutlet;
inletValue uniform 300;
}
BC-PressInlet-2
{
type fixedValue;
value uniform 300;
}
BC-PressInlet-1
{
type fixedValue;
value uniform 300;
}


U: (note I have also tried pressureInletOutletVelocity for the outlet)

dimensions [0 1 -1 0 0 0 0];

internalField uniform (1 0.1 0.1);

boundaryField
{
Periodic-Center
{
type cyclic;
}
Periodic-Combustor
{
type cyclic;
}
Periodic-Holes
{
type cyclic;
}
Periodic-Plenum-2
{
type cyclic;
}
Periodic-Plenum-1
{
type cyclic;
}
SmallTubes-2
{
type fixedValue;
value uniform (0 0 0);
}
LargeTubes-2
{
type fixedValue;
value uniform (0 0 0);
}
SmallTubes-1
{
type fixedValue;
value uniform (0 0 0);
}
LargeTubes-1
{
type fixedValue;
value uniform (0 0 0);
}
Wall-6
{
type fixedValue;
value uniform (0 0 0);
}
WallSet-5
{
type fixedValue;
value uniform (0 0 0);
}
WallSet-4
{
type fixedValue;
value uniform (0 0 0);
}
WallSet-2
{
type fixedValue;
value uniform (0 0 0);
}
Wall-3
{
type fixedValue;
value uniform (0 0 0);
}
Wall-1
{
type fixedValue;
value uniform (0 0 0);
}
BC-PressOut-center
{
type inletOutlet;
inletValue uniform (0 0 0);
}
BC-PressOut-1
{
type inletOutlet;
inletValue uniform (0 0 0);
}
BC-PressInlet-2
{
type pressureInletVelocity;
phi phi;
rho rho;
value uniform (0 0 0);
}
BC-PressInlet-1
{
type pressureInletVelocity;
phi phi;
rho rho;
value uniform (0 0 0);
}

[daniel_mills]

Even though this should be a low mach number(<0.3) would it still make sense to experiment with the pressureTransmissive outlet BC since it seems to be a growing oscillating pressure wave in a cylinder that seems to be causing the problems???

[pierre]

Hi Daniel,

Your set-up seems fine. Apart from that "low wave instabiltiy" what happens to the other fields, do they converge towards reasonable values, do they oscillate consequently too?
You can always try to take the immoral and cowardly step of using upwind for all your variables. Maybe you might also try to tighten your pressure solver tolerance.

Take it easy

Pierre

[hartinger]

probably try one of the transient codes like rhoSonicFoam and rhopSonicFoam. Using different codes helps spotting errors in the setup.

The total pressure boundary condition should be alright.

markus

[daniel_mills]

Thank you very much for your replies,

We get the same thing when upwinding all values and I have tried playing with the tolerances of pressure, it still develops a unstable pressure wave.

The problem with pressure is usually preluded by a sharp drop in the minimum Ux. A drop from around -33 to around -60. After a couple steps velocity is back to normal, but pressure develops a growing oscillating wave. Rho of course ocsillates with pressure. Temperature and velocity look pretty good for a few steps until variance in pressure begins driving them.

We have played with rhoTurbFoam a little but haven't spent much time with the transient codes yet.

[pierre]

As Markus suggested maybe rhopSonicFoam might help, as it uses a slightly different way of biulding fluxes. Quite a stable code I hear.

Pierre

[daniel_mills]

Does anyone have anything else we can try on this? We get the low wave pressure mode on both rhoTurbFoam and rhoSimpleFoam. It always develops in the channel(cylinder). It always develops after the swirling flow field begins to setup. We have tried more generic cases and the problem persists. The swirl develops then pressure oscillates and becomes unstable.

We have began trying rhopSonicFoam, but have not gotten very far.

Any other suggestions would be a great help.

Thanks

[aerothermal]

Originally Posted by daniel_mills

Does anyone have anything else we can try on this? We get the low wave pressure mode on both rhoTurbFoam and rhoSimpleFoam. It always develops in the channel(cylinder). It always develops after the swirling flow field begins to setup. We have tried more generic cases and the problem persists. The swirl develops then pressure oscillates and becomes unstable.

We have began trying rhopSonicFoam, but have not gotten very far.

Any other suggestions would be a great help.

Thanks

Hi Daniel,

Did you solve your issue?
If not, I will post my issue here. If affirmatively, I will try your solution.

Regards,

Guilherme

[alberto]

Hi, the last comment before yours is from 2006. It is probably more convenient if you open a new thread, and maybe provide a sample case that reproduces the problem.

Best,

[aerothermal]

Hi Alberto,

The last post is from 2006 but the user joined in 2009!
Certainly I will open a new thread. Thanks.


Regards,

Guilherme

[alberto]

Oh, I guess this is due to the forum migration.

[rock.senthilkumar]

Hai all,
I am new to OpenFoam.With my experience in Fluent ( simulated supersonic flow at Mach number 2.2 over Gemini Reentry capsule) thought to try OF. Simulation condition of my Fluent work is steady, Density based, M=2.2, inviscid both inviscid and turbulent flow. So far , I converted my fluent mesh file (.msh) as OpenFOAM mesh file. For steady , compressible flow rhoSimpleFoam seems to be a equivalent solver, but in OpenFOAM 1.7.1 (which is the one I installed in my machine) only rhoPorousSimpleFoam solver is available. Kindly explain what does this porous in rhoPorousSimpleFoam means for, and suggest me with correct solver.. Thanks in advance.

Originally Posted by alberto

Oh, I guess this is due to the forum migration.

[alberto]

rhoSimpleFoam is part of OpenFOAM 1.7.1

[rock.senthilkumar]

Dear Alberto,
I would like to thank you for your reply . And I feel very Sorry for my mistake., can you please suggest me any other solver other than rhoSimpleFoam for my case..

[alberto]

rhoSimpleFoam is the only compressible steady solver available.

There is rhoCentralFoam, which is suitable for viscid/inviscid flows with high Mach number, but it is unsteady. Similarly, sonicFoam is fine for transonic/supersonic turbulent unsteady flows.

Best,

[rock.senthilkumar]

Dear Alberto,
Thanks a lot for your reply.Now I decided the rhoSimpleFoam solver for my work. Sorry for bothering you again, Kindly send any test case for rhoSimpleFoam to my mail id rock.senthilkumar@gmail.com ., so that I can do any case study also.

Originally Posted by alberto

rhoSimpleFoam is the only compressible steady solver available.

There is rhoCentralFoam, which is suitable for viscid/inviscid flows with high Mach number, but it is unsteady. Similarly, sonicFoam is fine for transonic/supersonic turbulent unsteady flows.

Best,

[alberto]

Take a look at the tutorials :-)