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