两类pEqn.H求解可压缩流动计算时稳定性问题



  • 有对pimple算法求解pEqn.H有研究的老哥吗,关于这个算法有一些困惑,(1)https://github.com/OpenFOAM/OpenFOAM-2.3.x/blob/master/applications/solvers/combustion/reactingFoam/reactingFoam.C
    中密度方程为什么要求解两次,一次出现在runtime++后面,一次包含在pEqn.H里面。
    (2)reactingFoam和rhoreactingFoam中的pEqn.H求解的压力方程并不相同,而且reactingFoam的pEqn.H和rhoPimpleFoam里面的pEqn.H基本上是一致的。有使用过rhoreactingFoam的大佬吗?rhoreactingFoam在求解可压缩流动是是否比reactingFoam有优势?我现在使用的是reactingFoam这一类的pEqn.H,计算过程中一个时间步内压力变化很大,导致温度出现了震荡现象,导致物性也发生了震荡,无法收敛。
    下面贴出这两中pEqn.H
    rhoPimpleFoam/reactingFoam

    if (pimple.transonic())
    {
        surfaceScalarField phid
        (
            "phid",
            fvc::interpolate(psi)
           *(
                (fvc::interpolate(HbyA) & mesh.Sf())
              + rhorAUf*fvc::ddtCorr(rho, U, phi)/fvc::interpolate(rho)
            )
        );
    
        fvOptions.makeRelative(fvc::interpolate(psi), phid);
    
        while (pimple.correctNonOrthogonal())
        {
            fvScalarMatrix pEqn
            (
                fvm::ddt(psi, p)
              + fvm::div(phid, p)
              - fvm::laplacian(rhorAUf, p)
              ==
                fvOptions(psi, p, rho.name())
            );
    
            fvOptions.constrain(pEqn);
    
            pEqn.solve(mesh.solver(p.select(pimple.finalInnerIter())));
    
            if (pimple.finalNonOrthogonalIter())
            {
                phi == pEqn.flux();
            }
        }
    }
    

    rhoreactingFoam

     thermo.rho() -= psi*p;
    
        volScalarField rAU(1.0/UEqn.A());
        surfaceScalarField rhorAUf("rhorAUf", fvc::interpolate(rho*rAU));
    
        volVectorField HbyA("HbyA", U);
        HbyA = rAU*UEqn.H();
    
        if (pimple.transonic())
        {
            surfaceScalarField phiHbyA
            (
                "phiHbyA",
                (
                    (fvc::interpolate(rho*HbyA) & mesh.Sf())
                  + rhorAUf*fvc::ddtCorr(rho, U, phi)
                )/fvc::interpolate(rho)
            );
    
            fvOptions.makeRelative(phiHbyA);
    
            surfaceScalarField phid("phid", fvc::interpolate(thermo.psi())*phiHbyA);
    
            phiHbyA *= fvc::interpolate(rho);
    
            fvScalarMatrix pDDtEqn
            (
                fvc::ddt(rho) + fvc::div(phiHbyA)
              + correction(psi*fvm::ddt(p) + fvm::div(phid, p))
            );
    
            while (pimple.correctNonOrthogonal())
            {
                fvScalarMatrix pEqn
                (
                    pDDtEqn
                  - fvm::laplacian(rhorAUf, p)
                 ==
                    fvOptions(psi, p, rho.name())
                );
    
                fvOptions.constrain(pEqn);
    
                pEqn.solve(mesh.solver(p.select(pimple.finalInnerIter())));
    
                if (pimple.finalNonOrthogonalIter())
                {
                    phi = phiHbyA + pEqn.flux();
                }
            }
        }
        thermo.rho() += psi*p;
    

    只贴了transonic部分的


 

Forest
Mountains