# Understanding Terms in Compressible pEqn

 Register Blogs Members List Search Today's Posts Mark Forums Read

 November 29, 2017, 08:47 Understanding Terms in Compressible pEqn #1 New Member   Join Date: Oct 2017 Posts: 12 Rep Power: 8 I am trying to understand some of the terms in the compressible version of the pEqn (rhoPimpleFoam) and had some questions that I have not been able to figure out on my own. My first question has to do with the pressure equation itself. In "Computational Methods for Fluid Dynamics" chapter 10, the pressure-velocity correction is explained as inputting the predicted velocity and density into the continuity equation with a small correction parameter (, ): Where (in OpenFOAM notation) and according to Ferziger In OpenFOAM programming, this is implemented as: Code: // Calculate flux of HbyA and add contribution from gravity term surfaceScalarField phiHbyA ( "phiHbyA", fvc::interpolate(rho)*fvc::flux(HbyA) + rhorAUf*fvc::ddtCorr(rho, U, phi) + rhorAUf*fvc::interpolate(rho)*(g & mesh.Sf()) ); Code: // Setup pressure correction equation fvScalarMatrix pDDtEqn ( fvc::ddt(rho) + psi*correction(fvm::ddt(p)) + fvc::div(phiHbyA) == fvOptions(psi, p, rho.name()) ); while (pimple.correctNonOrthogonal()) { fvScalarMatrix pEqn(pDDtEqn - fvm::laplacian(rhorAUf, p)); pEqn.solve(mesh.solver(p.select(pimple.finalInnerIter()))); // Correct velocity flux by (subtracting?) flux from P Equation if (pimple.finalNonOrthogonalIter()) { phi = phiHbyA + pEqn.flux(); } } What I don't understand is the following: 1. Why is "fvc::div(phiHbyA) - fvm::laplacian(rhorAUf, P)" is included instead of "phiHbyA - fvm::grad(rhorAUf,P)"? Is there a particular reason the divergence of those terms taken? "(phiHbyA - fvm::grad(rhoAUf,P))*S" should be the correct velocity flux, so taking the divergence seems unnecessary to me. And even if we take the divergence, shouldn't we have to take the divergence of every term in that equation including the 'ddt' terms? 2. Is psi*correction(fvm::ddt(p)) the equivalent term? Or is this simply the density () correction? My second question has to do with how the density is calculated. In rhoPimpleFoam, the density is initially predicted by "rhoEqn.H" prior to calculating the velocity terms: Code: if (pimple.nCorrPIMPLE() <= 1) { #include "rhoEqn.H" } // --- Pressure-velocity PIMPLE corrector loop while (pimple.loop()) { #include "UEqn.H" #include "EEqn.H" // --- Pressure corrector loop while (pimple.correct()) { if (pimple.consistent()) { #include "pcEqn.H" } else { #include "pEqn.H" } } if (pimple.turbCorr()) { turbulence->correct(); } } But the code then solves again the "rhoEqn.H" after the PHI term is updated: Code: while (pimple.correctNonOrthogonal()) { fvScalarMatrix pEqn(pDDtEqn - fvm::laplacian(rhorAUf, p)); pEqn.solve(mesh.solver(p.select(pimple.finalInnerIter()))); // Correct velocity flux by (subtracting?) flux from P Equation if (pimple.finalNonOrthogonalIter()) { phi = phiHbyA + pEqn.flux(); } } #include "rhoEqn.H" #include "compressibleContinuityErrs.H" And then it is corrected again afterwards by: Code: thermo.correctRho(psi*p - psip0, rhoMin, rhoMax) ; rho = thermo.rho(); Is there a reason why the density is corrected so many times? And shouldn't the equation state act as the constraint on the density instead of solving the continuity equation again? kaifu, massive_turbulence, lukasf and 2 others like this.

 November 29, 2017, 11:54 #2 New Member   Join Date: Oct 2017 Posts: 12 Rep Power: 8 In regard to my first question, you can always write the continuity equation in differential form instead of integral form: If we break it into the (.)* and (.)' terms, we get for the time derivative: and for the divergence term The first term in the time derivative should be "fvc::ddt(rho)" while the second term would then be approximated by "psi*fvm::ddt(p)", where p here is the correction pressure. For the divergence terms, the first term should be "fvc::div(phiHbyA)-fvc::div(laplacian(rhorAUf, p))", but in rhoPimpleFoam, it is written as "fvc::div(phiHbyA)-fvm::div(laplacian(rhorAUf, p))". Why is that? And what about the remaining terms in the divergence term (i.e. all the u' and rho' terms)? I think I am missing something, as I am not sure how the pressure equation is derived in rhoPimpleFoam. Kummi likes this.

 December 2, 2017, 12:28 #3 New Member   Join Date: Oct 2017 Posts: 12 Rep Power: 8 I think I have (semi) figured out my first question. If we re-write the continuity equation as: Where (*) is the predicted values and (') are the corrected values. We can then break each term down to: But from the definition of , we can rewrite the above equation as: This simply equals (in OpenFOAM notation): Code: fvc::div(phiHbyA) - fvm::laplacian(rhorAUf,P) And for the density correction: Adding all the terms, we get: Code: fvc::ddt(rho) + psi*correction(fvm::ddt(p)) + fvc::div(phiHbyA) - fvm::laplacian(rhorAUf, p) Which is exactly what's given in rhoPimpleFoam solver. The only remaining term that is left that I can't figure out is the following: There doesn't seem to be any additional terms in the OpenFOAM pressure equation to account for this. Is this somehow hidden inside psi*correction(fvm::ddt(p)) and I just don't know? Any help is appreciated. mickbatti, bharat_aero, peyman.havaej and 4 others like this.

 February 19, 2018, 11:17 Have you been able to solve for the missing term ? #4 New Member   Ahmed Shoeb Join Date: Sep 2017 Location: Germany Posts: 3 Rep Power: 8 Hello thank you for your detailed explanation, have you been able to find a solution to account for the missing term ? Also could you please explain the role of the the fvc::ddtCorr(rho, U, phi) term ?

 October 16, 2018, 09:59 #5 New Member   Jan Gaertner Join Date: Nov 2017 Posts: 20 Rep Power: 8 Hello, the term is found when you use the transonic option. As I understand it they do not use the correction term here and instead the fvm description as there is no Code:  fvc::div(phi,rho) Zhiheng Wang and lukasf like this.

September 18, 2019, 01:50
#6
Member

Kai
Join Date: May 2010
Location: Shanghai
Posts: 61
Blog Entries: 1
Rep Power: 16
Quote:
 Originally Posted by AMK53 #include "rhoEqn.H" #include "compressibleContinuityErrs.H" [/CODE]And then it is corrected again afterwards by: Code: thermo.correctRho(psi*p - psip0, rhoMin, rhoMax) ; rho = thermo.rho(); Is there a reason why the density is corrected so many times? And shouldn't the equation state act as the constraint on the density instead of solving the continuity equation again?
When it goes to convergece, both of density update by thermo model and rhoEqn.H are all the same. I guess during calculaltion, this method could bring robustness to the solver.
__________________
Kai

 Tags compressible, openfoam, peqn, pressure correction