[dowlee]

Issues on the simulation of high-speed compressible flow within turbomachinery

Greeting to All Foamers,

Since I have employed Openfoam (including openfoam-1.6-ext, foam-extend-3.0/3.1 and openfoam-2.3.0) to run the simulation of high-speed compressible and turbulent flows (Usually the flow is transonic) of turbomachinery (including compressor and turbine) for nearly seven months, here I want to post some related issues on this topic, and really need some discussions with Foamers on this forum. And any critical comments are welcome and appreciated. If there are any errors, please help me to correct them, thanks very much in advance. In the following, the flows only and absolutely represent the high-speed compressible and turbulent flows inside turbomachinery.

According to my investigation, I find that at present OF is not very widely used to simulate such kind of flows inside turbomachinery. Except some applications for steady-RANS calculations which are reported on journal/conference papers (typical work done by e.g. Oliver Borm and Luca Mangani) and on this forum, it is rare to use OF to run the unsteady-RANS/DES/LES calculations for turbomachinery (I believe somebody does this, only just the results are not published). I feel that to get some good results by using OF, it is not very direct, compared with commercial software and in-house solvers. And as I know some Foamers (including me) have been puzzled with the accuracy of turbomachinery flow simulation by using OF for a very long time. I really want to know why OF is not widely used in turbomachinery community, and want Foamers guide me to figure out what should be further done inside the code of OF to promote the massive application of OF into such kind of flows inside turbomachinery.

I have no experience to modify the code inside OF, but only just read some codes and run some simulations for axial and centrifugal compressors. From my viewpoint as a OF user, I want to talk something from my experience based on my limited knowledge.

(1)Mesh
I use ICEM to generate the mesh of single blade passage. Usually the mesh is checked by ICEM without any problems. And based on this mesh, the results are very good by using CFX. But for the same mesh, the checkMesh utility always complains the cell determinant is too small. For example, based on same mesh, ICEM reports the minimum cell determinant is 0.3 while OF reports the minimum cell determinant is 1.0E-07. The small-determinant cells locate near the wall. This is due to that the normal distance between the first cell center away from the wall and the wall is small enough to make sure acceptable yPlus distribution. Some Foamers suggest that the mesh should always make checkMesh happy. The only way to diminish the warning on the minimum cell determinant is to increase the grid density near the wall along streamwise and spanwise direction, or increase the normal distance which resulting in bigger yPLus value. Do you think the cell determinant is a serious problem to impact the accuracy of flow simulation in your case?

(2)Flow Fields Initialization
The initialization of flow fields in OF is not direct. For the initialization of compressor flow fields, from the start it is better to set a uniform flow with a lower shaft speed and a lower outlet static pressure. Then gradually increase the shaft speed and the outlet pressure. In one of my cases, the design shaft speed is 58000 rpm and the outlet pressure at design point is nearly 7 times of the standard-atmosphere pressure while the inlet flow condition is from standard atmosphere. To get a steady-RANS result at design point with mesh nodes of 3 million, it took me around 2 weeks to increase the shaft speed from 0 rpm to 58000 rpm and the outlet pressure from 101325 pa to the design-point pressure. Compared with the commercial software and in-house solvers which just take one or two days to get steady-RANS results, this sounds too time-consuming.

Though the potentialFoam solver has the ability to provide an initial field, but after my several tests for turbomachinery, I never succeed to use it to provide a good initial field.

This also makes some indirect practice for grid-dependency comparison (e.g. i. check the sensitivity of overall aerodynamic performance to different grid levels based on the same geometry; ii. get a satisfied yPlus distribution on the wall). The mapField utility can be used to map the field of source mesh onto the field of target mesh, but sometimes it doesn’t work very well. For example, this utility may complain some cells of target mesh are outside of the source mesh, so those cells of target mesh fail to be mapped and result in the numerical iteration process based on the target mesh collapsed. To alleviate or avoid this collapse, the flow parameters (e.g. pressure/velocity/temperature) of those unmapped cells can be modified manually, but when the number of those cells is more, to do this is not realistic any more.

Are there any other robust and efficient tools or methods to initialize the flow fields of turbomachinery in OF? I mean the initialization can be finished in one step, but not by gradually adjusting the shaft speed and the outlet pressure.

(2)Boundary conditions
For typical turbomachinery flow simulation, the common used boundary conditions are : i) specification of total pressure, total temperature and flow direction at inlet; ii) specification of outlet static pressure; iii) zero velocity, zero-gradient temperature and zero-gradient pressure at wall, e.g. hub/shroud/blade; iv) cyclic boundary (cyclicGgi/cyclic/cyclicAMI).

i)Inlet boundary
The available boundary conditions provided by OF for inlet type in turbomachinery are: totalPressure, totalTemperature and pressureDirectedInletVelocity (or pressureInletVelocity, and so on for velocity boundary).
The totalPressure boundary calculates the static pressure from the velocity while the pressureDirectedInletVelocity boundary calculates the velocity from the flux with the specified inlet direction.

According to the propagation direction of flow character information for subsonic inlet flow, only the static pressure at inlet is calculated from the domain. And once the static pressure is updated, the velocity can be calculated from total pressure and inlet flow direction while the static temperature can be calculated from the isentropic equation. The question here is whether the totalPressure and pressureDirectedInletVelocity boundary conditions are identical with general-means numerical treatment of turbomachinery subsonic inlet.

In addition to the above available boundary conditions, Oliver Borm has developed another kind of boundary conditions for turbomachinery inlet (see openfoam-extend-DensityBasedTurbo): isentropicTotalTemperature for temperature boundary, temperatureDirectedInletVelocity for velocity boundary. And the pressure boundary is given with zeroGradient. If the distance between the inlet and first blade leading edge is short, I am not sure whether the zeroGradient pressure boundary is still suitable. I have compared the results calculated by OF and CFX based on the same mesh. It shows that for the results by using CFX, the inlet pressure gradient is very big. Though the numerical convergence by using OF with zeroGradient pressure at inlet is good, but the mass flow rate is too much lower.

ii)Outlet
At turbomachinery outlet, static pressure is specified while zeroGradient is employed for temperature and velocity (it should be note that fluxCorrectedVelocity boundary can also be used for velocity).

If the distance between the blade trailing edge and outlet is not short and the flow is uniform near the outlet, the flow streamwise gradient is small and so zeroGradient boundary is acceptable for the outlet boundary.

In my centrifugal impeller case (without diffusor), the distance is short, and at the trailing edge big separation and typical jet-wake flow exist. So the outlet flow is highly non-uniform and the flow streamwise gradient at outlet should not be neglected, again, is the zeroGradient boundary for velocity and temperature still acceptable at the outlet? I have checked my cases for many times, and I have been stuck for several months, the mass-flow –rate difference between inlet and outlet calculated by using OF is always very big.

And I have thoroughly read one professional code for compressor and turbine widely used in industry and academy. Within the code, for the inlet and outlet boundaries, instead of zeroGradient, extrapolation from neighboring cells inside the domain is adopted. As OF is a kind of general flow simulator, I don`t know whether such extrapolation is appreciated. And what is your opinion about zeroGradient boundary condition?

(4)Sovler
The available solvers for compressible and turbulent flow include steady and transient options.

i)Steady-state solver
Luca Mangani in his paper states that since the standard OpenFoam library does not contain a solver properly working for compressible flow fields. So he develops a new solver by using SIMPLE-C algorithm for steady calculation. Oliver Borm developed transonicMRFDyMFoam and transonicSteadySRFFoam solvers. Recently Foam community has released steadyCompressibleMRFFoam solver.

ii)transient solver
In OpenFoam-1.6-ext and Foam-extend-3.0/3.1, the rhoPorousMRFPimpleSolver is available for turbomachinery. Oliver Borm developed transonicUnsteadyMRFDyMFoam. In OpenFOAM-2.3.0, fvOptions makes several solvers available for turbomachinery.

Does anybody test the performance of these above solvers? I think the above solvers are general to solve a kind of flow problems. For turbomachinery, is it necessary to include some special treatment in the solver, e.g. inlet and outlet boundary, robustness with flow separation, more suitable turbulent wall function and so on. I have been told that sometimes a solver maybe work well with not very complex flow inside compressor, but if the blade loading is very high and big separation exists, the solver may be no longer suitable and more advanced physical /numerical model is needed. For my centrifugal impeller case, big separation exists and the total pressure ratio exceeds 10. For this case, I have tested different schemes (e.g. upwind/central/limitedLinear/filterLinear schemes) based on three kind of meshes (the cell nodes are 380 k with maximum yPlus of 3, 3 million with maximum yPlus of 50, 3 million with maximum yPlus of 3), but mass-flow-rate difference is never within the general reference value of 0.5%.

From the above, I don’t want to leave you such an impression that I complains the feasibility and functionality of OF for turbomachinery flow simulations. But actually I am very interested to OF with its superiority. I want to know whether there are problems that stop the massive application in turbomachinery community and if there are, what should be done in the next future. And at most I want to solve my problems at present. Welcome Fomers give me some hints.

[dowlee]

Hi, does anybody give me some advice?

[jensi_t]

Hi,

I use dbnsTurbFoam (foam ext 3.1 based on Borms solver) with LRR and ke modelling (wall functions and y+ about 30). It is a transient solver, but i just managed to run it with local time-stepping (EulerLocal). This should help you getting a steady state solution more quickly. I don't have turning parts, but i think you can use ggi interfaces.

Greetings,
Jens

P.S.: I use fluentDataToFoam to create an initial solution

[dowlee]

Dear Jens
1) Before fortunately, I used rhoPorousMRFPimpleFoam with LES model, and got very good results. In the future, I will investigate this dbnsTurbFoam solver.
2) Yes, using Fluent results to create an initial solution is very convenient. Thanks very much.
Kindly regards
Dowlee

Originally Posted by jensi_t

Hi,

I use dbnsTurbFoam (foam ext 3.1 based on Borms solver) with LRR and ke modelling (wall functions and y+ about 30). It is a transient solver, but i just managed to run it with local time-stepping (EulerLocal). This should help you getting a steady state solution more quickly. I don't have turning parts, but i think you can use ggi interfaces.

Greetings,
Jens

P.S.: I use fluentDataToFoam to create an initial solution

[normy]

Originally Posted by jensi_t

Hi,

I use dbnsTurbFoam (foam ext 3.1 based on Borms solver) with LRR and ke modelling (wall functions and y+ about 30). It is a transient solver, but i just managed to run it with local time-stepping (EulerLocal). This should help you getting a steady state solution more quickly.

Hi Jens, I'm new in Openfoam as you'll see
I wanted to know how you managed to run dbnsTurbFoam with local time-stepping. Did you build a new solver or change only the fvSchemes?
Because I'm trying to modify the source code but I had some problems...
Could you help me, please?
Thanks,
Norman

[jensi_t]

Hey Norman,

Sorry I totally forgot to answer. You have to create the physikalDealtaT field in createFields.H:

    Info<< "Reading thermophysical properties\n" << endl;

    autoPtr<basicPsiThermo> thermo
    (
        basicPsiThermo::New(mesh)
    );

    // Primitive variables

    volScalarField& h = thermo->h();
    volScalarField& p = thermo->p();
    const volScalarField& T = thermo->T();

    Info<< "Reading field rho\n" << endl;
    volScalarField rho
    (
        IOobject
        (
            "rho",
            runTime.timeName(),
            mesh,
            IOobject::NO_READ,
            IOobject::AUTO_WRITE
        ),
        thermo->rho()
    );

    Info<< "Reading field U\n" << endl;
    volVectorField U
    (
        IOobject
        (
            "U",
            runTime.timeName(),
            mesh,
            IOobject::MUST_READ,
            IOobject::AUTO_WRITE
        ),
        mesh
    );


    // Conservative variables

    volVectorField rhoU
    (
        IOobject
        (
            "rhoU",
            runTime.timeName(),
            mesh,
            IOobject::NO_READ,
            IOobject::NO_WRITE
        ),
        rho*U
    );

    volScalarField rhoE
    (
        IOobject
        (
            "rhoE",
            runTime.timeName(),
            mesh,
            IOobject::NO_READ,
            IOobject::NO_WRITE
        ),
        rho*(h + 0.5*magSqr(U)) - p
    );


    // Create numeric flux
    //numericFlux<roeFlux, firstOrderLimiter> dbnsFlux    
    numericFlux<roeFlux, VenkatakrishnanLimiter> dbnsFlux
    //numericFlux<roeFlux, BarthJespersenLimiter> dbnsFlux
    //numericFlux<rusanovFlux, BarthJespersenLimiter> dbnsFlux
    (
        p,
        U,
        T,
        thermo()
    );

    // Create mass flux alias for easier coupling with other code components
    const surfaceScalarField& phi = dbnsFlux.rhoFlux();

    Info<< "Creating turbulence model\n" << endl;
    autoPtr<compressible::turbulenceModel> turbulence
    (
        compressible::turbulenceModel::New
        (
            rho,
            U,
            phi,
            thermo
        )
    );

    
    
    //Test Jens
     localTimeStep localTimeStep(U, thermo, turbulence());
    
     IOField<scalar> physDeltaT
     (
      IOobject
      (
         "physDeltaT",
          runTime.timeName(),
          mesh,
          IOobject::NO_READ,
          IOobject::NO_WRITE
      ),
      1.0
      );

and then I think the Runge Kutta time integration is already implemented in the localTimeStep.H, so you also have to modify the solver

#include "fvCFD.H"
#include "basicPsiThermo.H"
#include "turbulenceModel.H"
#include "bound.H"
#include "hllcFlux.H"
#include "roeFlux.H"
#include "rusanovFlux.H"
#include "betaFlux.H"
#include "MDLimiter.H"
#include "firstOrderLimiter.H"
#include "BarthJespersenLimiter.H"
#include "VenkatakrishnanLimiter.H"
#include "numericFlux.H"
#include "localTimeStep.H"


// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

int main(int argc, char *argv[])
{
#   include "setRootCase.H"
#   include "createTime.H"
#   include "createMesh.H"
#   include "createFields.H"

// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

    Info<< "\nStarting time loop\n" << endl;

    // Runge-Kutta coefficient
    scalarList beta(4);
    beta[0] = 0.1100;
    beta[1] = 0.2766;
    beta[2] = 0.5000;
    beta[3] = 1.0000;

    while (runTime.run())
    {
#       include "readTimeControls.H"
#       include "readFieldBounds.H"
//#       include "compressibleCourantNo.H"
//#       include "setDeltaT.H"

        runTime++;

        Info<< "\n Time = " << runTime.value() << endl;

        // Switch off solver messages for diagonal solver RK
        lduMatrix::debug = 0;

    
    localTimeStep.update(maxCo, adjustTimeStep);
    
        // Low storage Runge-Kutta time integration
        forAll (beta, i)
        {
            // Solve the approximate Riemann problem 
            dbnsFlux.computeFlux();
        

        physDeltaT[0]=beta[i];

            // Time integration
            solve
            (
                fvm::ddt(rho)
              + fvc::div(dbnsFlux.rhoFlux())
            );

            solve
            (
                fvm::ddt(rhoU)
              + fvc::div(dbnsFlux.rhoUFlux())
              + fvc::div(turbulence->devRhoReff())
            );

            solve
            (
                fvm::ddt(rhoE)
              + fvc::div(dbnsFlux.rhoEFlux())
              + fvc::div(turbulence->devRhoReff() & U)
              - fvc::laplacian(turbulence->alphaEff(), h)

            );

#           include "updateFields.H"
        }

        // Switch on solver messages for turbulence
        lduMatrix::debug = 1;

        turbulence->correct();

        runTime.write();

        Info<< "    ExecutionTime = "
            << runTime.elapsedCpuTime()
            << " s\n" << endl;
        

    Info<< "\n end \n";

    return(0);
}

Good Luck!

[hmbo2]

Hi,
I am new to use openFoam/2.3.0 and my case study is turbomachinery as one and half stage(stator-rotor-stator) and I create my mesh using blockMeshDG and the mesh is very good. Now I want to start my solver, so which solver in openFoam/2.3.0 is suitable for turbomachinery (axial turbine and axial compressor) and how can I change the fvoptions to work suitable for the solver?
Thanks In advanced

Hayder

[hmbo2]

Hello, Anybody help?

I am still using the openFoam/2.3.0 with rhoSimplecFoam as a solver for turbine stator with total pressure at inlet and static pressure at the outlet. I have got a good results in the Yaw angle but the values of meridonal, circumeferntial and absolute components of velocity is higher than the experimental data. please can any body help me?

Thanks in advanced

[hakolekar]

Hi Dowlee,

I am working on turbomachinery flows - more specifically turbines with small pressure drops.

It will be great to hear how you resolved some of the issues mentioned in your discussion and whether OpenFoam 2.3.0 has the capability to provide a good solution,

Thank you very much.
Regards,

Harshal

[anthony761]

Hello Dowlee, hello dear foamers,

Did you succeed to perform any compressible flow simulation for turbomachinery application? What is you feedback? I have also encoutered some problems you were talking about. This is my case:

OpenFOAM-v2006
Geometry: Centrifugal compressor, only one blade passage with cyclic/cyclicAMI (with diffuser vanedless and without volute)
BCs: flowRateInletVelocity with fixedValue temperature inlet and fixedValue pressure outlet
MRF model
Turbulence: kOmegaSST
Thermophysical: hePsiThermo, pureMixture, transport const, hConst, perfectGas, sensibleInternalEnergy

Tried Solvers: rhoSimpleFoam, rhoCentralFoam, sonicFoam

Results are strange, the pressure increased as expected at the exit but not the temperature which seems remain almost constant. Do you have any idea of the problem?

[Fouch]

Hello Anthony,


I want to discuss about turbomachinery and rhoSimpleFoam.


I am calculating a centrifugal fan so very close to a compressor.


In my case, I cannot reach the expected pressure but CFX can. Could you share your fvSolution and fvSchemes ?


You can find mine on this topic :
rhoSimpleFoam accuracy


Thank you in advance.

[anthony761]

Originally Posted by Fouch

Hello Anthony,


I want to discuss about turbomachinery and rhoSimpleFoam.


I am calculating a centrifugal fan so very close to a compressor.


In my case, I cannot reach the expected pressure but CFX can. Could you share your fvSolution and fvSchemes ?


You can find mine on this topic :
rhoSimpleFoam accuracy


Thank you in advance.

Ok, let's discuss in your post