
[Sponsors] 
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 pressurevelocity 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(); } } 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" } //  Pressurevelocity 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(); } } 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" Code:
thermo.correctRho(psi*p  psip0, rhoMin, rhoMax) ; rho = thermo.rho(); 

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. 

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 rewrite 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) Adding all the terms, we get: Code:
fvc::ddt(rho) + psi*correction(fvm::ddt(p)) + fvc::div(phiHbyA)  fvm::laplacian(rhorAUf, p) 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. 

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 

September 18, 2019, 01:50 

#6  
Member

Quote:
__________________
Kai 

Tags 
compressible, openfoam, peqn, pressure correction 
Thread Tools  Search this Thread 
Display Modes  


Similar Threads  
Thread  Thread Starter  Forum  Replies  Last Post 
How rho_ and psi_ are calculated in compressible solvers of OpenFOAM?  JBUNSW  OpenFOAM Programming & Development  17  February 12, 2021 08:59 
Flux discretization of 2D Compressible NavierStokes  selig5576  Main CFD Forum  7  March 29, 2017 22:25 
Simulation of high pressure diesel injector  all phases compressible with cavitation  fivos  CFX  4  July 30, 2015 06:48 
One phase is compressible, can I use compressibleinterfoam?  sharonyue  OpenFOAM Running, Solving & CFD  2  March 31, 2014 09:40 
momentum underrelaxation for compressible flow with SIMPLE  Mihai ARGHIR  Main CFD Forum  0  April 7, 2000 04:58 