zhyang-dev/openfoam
1
0
Fork 0
mirror of https://develop.openfoam.com/Development/openfoam.git synced 2025-11-28 03:28:01 +00:00
Code Issues Actions Packages Projects Releases Wiki Activity

Corrected the transonic formulation.

Browse Source
This commit is contained in:
henry
2008-08-25 11:18:13 +01:00
parent 4164b8efe8
commit e612b8e4da
2 changed files with 75 additions and 27 deletions
Download Patch File Download Diff File Expand all files Collapse all files

2
applications/solvers/compressible/rhoSimpleFoam/hEqn.H
View File

@ -5,7 +5,7 @@
- fvm::Sp(fvc::div(phi), h) - fvm::Sp(fvc::div(phi), h)
- fvm::laplacian(turbulence->alphaEff(), h) - fvm::laplacian(turbulence->alphaEff(), h)
== ==
fvc::div(phi/fvc::interpolate(rho)*fvc::interpolate(p, "div(U,p)")) fvc::div(phi/fvc::interpolate(rho), p, "div(U,p)")
- p*fvc::div(phi/fvc::interpolate(rho)) - p*fvc::div(phi/fvc::interpolate(rho))
); );

72
applications/solvers/compressible/rhoSimpleFoam/pEqn.H
View File

@ -1,17 +1,63 @@
volScalarField AU = UEqn().A(); rho = thermo->rho();
U = UEqn().H()/AU;
UEqn.clear();
phi = fvc::interpolate(rho*U) & mesh.Sf();
bool closedVolume = adjustPhi(phi, U, p);
for (int nonOrth=0; nonOrth<=nNonOrthCorr; nonOrth++) volScalarField rUA = 1.0/UEqn().A();
U = rUA*UEqn().H();
UEqn.clear();
bool closedVolume = false;
if (transonic)
{ {
surfaceScalarField phid
(
"phid",
fvc::interpolate(thermo->psi())*(fvc::interpolate(U) & mesh.Sf())
);
for (int nonOrth=0; nonOrth<=nNonOrthCorr; nonOrth++)
{
fvScalarMatrix pEqn0
(
fvm::div(phid, p)
- fvm::laplacian(rho*rUA, p)
);
fvScalarMatrix pEqn = pEqn0;
pEqn.relax(mesh.relaxationFactor("pEqn"));
pEqn.setReference(pRefCell, pRefValue);
// retain the residual from the first iteration
if (nonOrth == 0)
{
eqnResidual = pEqn.solve().initialResidual();
maxResidual = max(eqnResidual, maxResidual);
}
else
{
pEqn.solve();
}
if (nonOrth == nNonOrthCorr)
{
phi == pEqn0.flux();
}
}
}
else
{
phi = fvc::interpolate(rho)*(fvc::interpolate(U) & mesh.Sf());
//phi = fvc::interpolate(rho*U) & mesh.Sf();
closedVolume = adjustPhi(phi, U, p);
for (int nonOrth=0; nonOrth<=nNonOrthCorr; nonOrth++)
{
fvScalarMatrix pEqn fvScalarMatrix pEqn
( (
fvm::laplacian(rho/AU, p) == fvc::div(phi) fvm::laplacian(rho*rUA, p) == fvc::div(phi)
); );
pEqn.setReference(pRefCell, pRefValue); pEqn.setReference(pRefCell, pRefValue);
// retain the residual from the first iteration // retain the residual from the first iteration
if (nonOrth == 0) if (nonOrth == 0)
{ {
@ -27,14 +73,20 @@ for (int nonOrth=0; nonOrth<=nNonOrthCorr; nonOrth++)
{ {
phi -= pEqn.flux(); phi -= pEqn.flux();
} }
}
} }
#include "incompressible/continuityErrs.H" #include "incompressible/continuityErrs.H"
// Explicitly relax pressure for momentum corrector // Explicitly relax pressure for momentum corrector
p.relax(); p.relax();
U -= fvc::grad(p)/AU; rho = thermo->rho();
rho.relax();
Info<< "rho max/min : " << max(rho).value() << " " << min(rho).value() << endl;
U -= rUA*fvc::grad(p);
U.correctBoundaryConditions(); U.correctBoundaryConditions();
bound(p, pMin); bound(p, pMin);
@ -46,7 +98,3 @@ if (closedVolume)
p += (initialMass - fvc::domainIntegrate(thermo->psi()*p)) p += (initialMass - fvc::domainIntegrate(thermo->psi()*p))
/fvc::domainIntegrate(thermo->psi()); /fvc::domainIntegrate(thermo->psi());
} }
rho = thermo->rho();
rho.relax();
Info<< "rho max/min : " << max(rho).value() << " " << min(rho).value() << endl;