[tangd]

Hello,

I ran my case with both OpenFoam (sonicTurbFoam) and CFX5 and compared
results by generating graphs of temperature, pressure, mass and mass flow rate
as well as some contour screenshots (see attached images). In principle,
the axisymmetric-2D-system is a simple volume with an inlet. For OpenFoam
a 5° sector was used with wedges.

I found that temperature and pressure of OpenFoam are lower than the corresponding
results of CFX-5. Additionally, in OpenFoam the overall mass in the closed volume
is lower than the integral value of the inlet mass flow rate. Thus, it looks as if
mass is lost during the calculation. Note, the massflow input (and also the
total temperature) at inlet is nearly the same as in CFX5. The mesh was the same,
of course.

I have checked a lot of things: turbulence model, boundary conditions, thermal
properties, etc, but haven't found the reason why in OpenForm the mass balance (and
consequently the energy balance) is not fulfilled. The schemes I used are:

================================================== ============================
ddtSchemes
{
default CrankNicholson 0.5;
}

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,R) Gauss upwind;
div(R) Gauss linear;
div(phid,p) Gauss limitedLinear 1;
div(phiU,p) Gauss limitedLinear 1;
div(phi,h) Gauss limitedLinear 1;
div((muEff*dev2(grad(U).T()))) Gauss linear;
}

laplacianSchemes
{
default none;
laplacian(muEff,U) Gauss linear corrected;
laplacian(DkEff,k) Gauss linear corrected;
laplacian(DREff,R) Gauss linear corrected;
laplacian(DepsilonEff,epsilon) Gauss linear corrected;
laplacian((rho*1|A(U)),p) Gauss linear corrected;
laplacian(alphaEff,h) Gauss linear corrected;
}

interpolationSchemes
{
default linear;
}

snGradSchemes
{
default corrected;
}

fluxRequired
{
default no;
p;
}
================================================== ============================

I've also used Euler time discretisation, but nothing changed. So far, I didn't
change the spatial discretisation, it's still the default setting from prism.
The solution settings are:

================================================== ============================
solvers
{
p BICCG 1e-12 0;
rho ICCG 1e-12 0;
U BICCG 1e-12 0;
h BICCG 1e-12 0;
k BICCG 1e-12 0;
epsilon BICCG 1e-12 0;
R BICCG 1e-12 0;
}

PISO
{
momentumPredictor yes;
nCorrectors 2;
nNonOrthogonalCorrectors 0;
pRefCell 0;
pRefValue 0;
}
================================================== ============================

Thus, does anybody have an idea where the problem could be? The unconserevd
energy and mass balance make me uncomfortable ;-(

[lr103476]

Don't you need CrankNicholson 1.0 for fully second order accuracy?

Regards, Frank

[tangd]

I also used Euler here, it doesn't matter so much I think. Other ideas?

[david_h]

How did you calculate the system mass, m, the dotted curve in the t vs. (m0 - m) plot ?
Dave

[tangd]

The system mass is calculated by using pV=mRT, p and T are the mean values within the volume. Then I get the mass(m) value for each saved step, subtract the initial mass m0.

[david_h]

If I understood correctly you are calculating the mass as m = V/R {p}/{T}. Where { } denotes the volume average. Have you tried calculating the mass as m = V/R {p/T} ?

[tangd]

Yes, your idea turned out to be correct! I have checked the (m-m0) again by using your proposal m = V/R {p/T}, the mass is now conserved. Thanks a lot! But any idea about the pressure and Temperature?

[olwi]

If you really want to integrate the influx, you should look at *exactly* what goes in. That is exactly what's on the boundary patch. Not what's in the cell next to it... In addition, your own way of integarting plays a role. Is *your* integration 100% identical with the OpenFOAM time integration? For example using mid-point values matching the CN scheme?

It's relatively easy to assure mass conservation in a fv code like OpenFOAM, and from the quality in detail I've seen so far in OpenFOAM it would be profoundly surprising if you found a lack of mass conservation. "Surprising" is too weak a word, actually.

/Ola

[tangd]

You are right! My incorrect postprocessing calculation introduced the error. Now mass is conserved. But there are still problems during the simulation that make temperature and pressure smaller than that of CFX-5. any further suggestion?

[olwi]

Then remains to compare exactly how much influx of enthalpy you have in the two codes. If I understood correclty, there's no heat production inside? If the mass inventory is the same in both codes at time t, then you only need to look for the enthalpy influx. Same here: I doubt they've missed conservation in CFX5 or OpenFOAM... It's likley a matter of how sources or inlet properties are defined. Even if you think you do they same in both codes, the codes may have different conventions for how things are defined on the boundary.

[tangd]

The finding is that the enthalpy flow of OpenFOAM is very close to that of CFX5.

Here, we have checked the static part, mdot*h_static, as well as the kinetic part mdot*u²/2. Both are the same on the inlet for OpenFoam and CFX.

I inspected this by considering each cell at inlet then summing up them. Using the same method, I also checked the pressures and velocities at inlet from both tools. They are nearly the same. The inlet calculation seems to be no problem. But I still cannot figure out what kind of problem causes the different pressure and temperature values within the volume. Can you give me more suggestions?

[david_h]

Did you run CFX with viscous dissipation. The energy equation in sonicTurbFoam, "hEqn.H", does not appear to have a viscous dissipation term, which might explain a few degrees of temperature difference.

[tangd]

In the setting of heat transfer model "Total Energy" in CFX5, option "Incl. viscous Work Term" was _not_ activated.

[olwi]

Have you done the same check on energy/enthalpy conservation as you did for the mass? (total energy in the flow domain before/after in relation to the enthalpy influx.) If energy is conserved, then I guess there's only the consitutive relationship, etc., relating pressure to temperature and enthalpy...

[tangd]

Now the situation is that if I subtract the mean temperature <t> of OpenFoam with 293 and multiply with gamma (gamma=cp/cv), then I get almost the same temperature as in CFX5. So I guess the wrong equation is being used in calculating temperature in OpenFOAM. Any suggestion?

[olwi]

It's a bit hard to say OpenFOAM solves the "wrong" equation, isn't it? Sounds to me OF uses SI units for temperature (Kelvin), while CFX5 uses something else. The factor gamma can be explained if someone used cv instead of cp, or vice versa. Sounds to me it's more a question of reading in the manuals/case files to know what is actually calculated and ensuring you use the same physical properties and units in both codes. That's obviously the user's responsability. Should be easy enough in OF, just have a look in the source code.

[tangd]

[david_h]

One other thing you might take look at is the turbulence model. For example in:
/src/turbulenceModels/compressible/kEpsilon/kEpsilon.C,
the R() function includes a "(2/3)*k" term, however divRhoR() does not include this term.

The quickest fix is to add
-(2/3)*fvc::grad(rho*turbulence->k() )
after the fvc::grad(p) term when the momentum equation is solved.

Dave

[tangd]

Thank you so much. I have changed the code by using your suggestion. But unfortunately, the result is nearly the same as before.

I'm still thinking the "gamma" plays the key role here. So please give me more hints based on my previous post. Thanks again!

[chris1980]

Can someone comment if David is right because this would be a potential bug!