pEqnRhoSimple.H
Go to the documentation of this file.
1 volScalarField rAU(1.0 / UEqn.A());
2 surfaceScalarField rhorAUf("rhorAUf", fvc::interpolate(rho* rAU));
3 //***************** NOTE *******************
4 // constrainHbyA has been used since OpenFOAM-v1606; however, it may degrade the accuracy of derivatives
5 // because constraining variables will create discontinuity. Here we have a option to use the old
6 // implementation in OpenFOAM-3.0+ and before (no constraint for HbyA)
7 autoPtr<volVectorField> HbyAPtr = nullptr;
8 label useConstrainHbyA = daOptionPtr_->getOption<label>("useConstrainHbyA");
10 {
11  HbyAPtr.reset(new volVectorField(constrainHbyA(rAU * UEqn.H(), U, p)));
12 }
13 else
14 {
15  HbyAPtr.reset(new volVectorField("HbyA", U));
16  HbyAPtr() = rAU * UEqn.H();
17 }
18 volVectorField& HbyA = HbyAPtr();
19 
20 tUEqn.clear();
21 
22 bool closedVolume = false;
23 
24 surfaceScalarField phiHbyA("phiHbyA", fvc::interpolate(rho) * fvc::flux(HbyA));
25 
26 MRF.makeRelative(fvc::interpolate(rho), phiHbyA);
27 
28 // Update the pressure BCs to ensure flux consistency
30 
31 // NOTE: we don't support transonic = true
32 
34 
35 while (simple.correctNonOrthogonal())
36 {
37  fvScalarMatrix pEqn(
38  fvc::div(phiHbyA)
39  - fvm::laplacian(rhorAUf, p)
40  - fvOptions(psi, p, rho.name()));
41 
42  pEqn.setReference(
43  pressureControl.refCell(),
44  pressureControl.refValue());
45 
46  // get the solver performance info such as initial
47  // and final residuals
48  SolverPerformance<scalar> solverP = pEqn.solve();
49 
50  DAUtility::primalResidualControl(solverP, printToScreen_, "p", daGlobalVarPtr_->primalMaxRes);
51 
52  if (simple.finalNonOrthogonalIter())
53  {
54  phi = phiHbyA + pEqn.flux();
55  }
56 }
57 
58 if (printToScreen_)
59 {
60 #include "continuityErrsPython.H"
61 }
62 
63 // Explicitly relax pressure for momentum corrector
64 p.relax();
65 
66 // bound p
67 DAUtility::boundVar(allOptions, p, printToScreen_);
68 
69 U = HbyA - rAU * fvc::grad(p);
70 // bound U
71 DAUtility::boundVar(allOptions, U, printToScreen_);
72 U.correctBoundaryConditions();
73 fvOptions.correct(U);
74 
75 bool pLimited = pressureControl.limit(p);
76 
77 // For closed-volume cases adjust the pressure and density levels
78 // to obey overall mass continuity
80 {
81  p += (initialMass_ - fvc::domainIntegrate(psi * p))
82  / fvc::domainIntegrate(psi);
83 }
84 
86 {
87  p.correctBoundaryConditions();
88 }
89 
90 rho = thermo.rho();
91 
92 rho.relax();
93 
94 // bound rho
95 DAUtility::boundVar(allOptions, rho, printToScreen_);
closedVolume
bool closedVolume
Definition: pEqnRhoSimple.H:22
pressureControl
pressureControl & pressureControl
Definition: createRefsRhoSimpleC.H:12
simple
simpleControl & simple
Definition: createRefsRhoSimpleC.H:5
initialMass_
initialMass_
Definition: createRefsRhoSimpleC.H:13
rAU
volScalarField rAU(1.0/UEqn.A())
UEqn
fvVectorMatrix & UEqn
Definition: UEqnPimpleDyM.H:15
MRF
IOMRFZoneListDF & MRF
Definition: createRefsRhoSimple.H:18
rho
rho
Definition: pEqnRhoSimple.H:90
thermo
fluidThermo & thermo
Definition: createRefsRhoPimple.H:6
HbyA
volVectorField & HbyA
Definition: pEqnRhoSimple.H:18
useConstrainHbyA
label useConstrainHbyA
Definition: pEqnRhoSimple.H:8
continuityErrsPython.H
pLimited
bool pLimited
Definition: pEqnRhoSimple.H:75
phiHbyA
surfaceScalarField phiHbyA("phiHbyA", fvc::interpolate(rho) *fvc::flux(HbyA))
p
volScalarField & p
Definition: createRefsPimpleDyM.H:6
U
U
Definition: pEqnRhoSimple.H:69
phi
surfaceScalarField & phi
Definition: createRefsPimpleDyM.H:8
rhorAUf
surfaceScalarField rhorAUf("rhorAUf", fvc::interpolate(rho *rAU))
tUEqn
tmp< fvVectorMatrix > tUEqn(fvm::ddt(U)+fvm::div(phi, U)+turbulencePtr_->divDevReff(U) - fvSource)
HbyAPtr
autoPtr< volVectorField > HbyAPtr
Definition: pEqnRhoSimple.H:7
allOptions
const dictionary & allOptions
Definition: createRefsPimpleDyM.H:14
adjustPhi
adjustPhi(phiHbyA, U, p)
constrainPressure
constrainPressure(p, rho, U, phiHbyA, rhorAUf, MRF)
psi
const volScalarField & psi
Definition: createRefsRhoPimple.H:14