pEqnPimpleDyM.H
Go to the documentation of this file.
1 volScalarField rAU(1.0 / UEqn.A());
2 //***************** NOTE *******************
3 // constrainHbyA has been used since OpenFOAM-v1606; however, it may degrade the accuracy of derivatives
4 // because constraining variables will create discontinuity. Here we have a option to use the old
5 // implementation in OpenFOAM-3.0+ and before (no constraint for HbyA)
6 autoPtr<volVectorField> HbyAPtr = nullptr;
7 label useConstrainHbyA = daOptionPtr_->getOption<label>("useConstrainHbyA");
9 {
10  HbyAPtr.reset(new volVectorField(constrainHbyA(rAU * UEqn.H(), U, p)));
11 }
12 else
13 {
14  HbyAPtr.reset(new volVectorField("HbyA", U));
15  HbyAPtr() = rAU * UEqn.H();
16 }
17 volVectorField& HbyA = HbyAPtr();
18 surfaceScalarField phiHbyA(
19  "phiHbyA",
20  fvc::flux(HbyA));
21 
22 tmp<volScalarField> rAtU(rAU);
23 
24 if (pimple.nCorrPISO() <= 1)
25 {
26  tUEqn.clear();
27 }
28 
29 // Non-orthogonal pressure corrector loop
30 while (pimple.correctNonOrthogonal())
31 {
32  fvScalarMatrix pEqn(
33  fvm::laplacian(rAtU(), p) == fvc::div(phiHbyA));
34 
35  pEqn.setReference(pRefCell, pRefValue);
36 
37  // get the solver performance info such as initial
38  // and final residuals
39  SolverPerformance<scalar> solverP = pEqn.solve(mesh.solver(p.select(pimple.finalInnerIter())));
40 
41  DAUtility::primalResidualControl(solverP, pimplePrintToScreen, "p", daGlobalVarPtr_->primalMaxRes);
42 
43  if (pimple.finalNonOrthogonalIter())
44  {
45  phi = phiHbyA - pEqn.flux();
46  }
47 }
48 
49 if (pimplePrintToScreen)
50 {
51 #include "continuityErrsPython.H"
52 }
53 
54 // Explicitly relax pressure for momentum corrector
55 p.relax();
56 
57 // bound p
58 DAUtility::boundVar(allOptions, p, pimplePrintToScreen);
59 
60 U = HbyA - rAtU * fvc::grad(p);
61 // bound U
62 DAUtility::boundVar(allOptions, U, pimplePrintToScreen);
63 U.correctBoundaryConditions();
64 
65 // Correct Uf if the mesh is moving
66 this->correctUfPimpleDyM(Uf, U, phi);
67 
68 // Make the fluxes relative to the mesh motion
pRefValue
scalar & pRefValue
Definition: createRefsPimpleDyM.H:12
U
U
Definition: pEqnPimpleDyM.H:60
pimple
pimpleControlDF & pimple
Definition: createRefsPimpleDyM.H:5
useConstrainHbyA
label useConstrainHbyA
Definition: pEqnPimpleDyM.H:7
UEqn
fvVectorMatrix & UEqn
Definition: UEqnPimpleDyM.H:15
continuityErrsPython.H
pRefCell
label & pRefCell
Definition: createRefsPimpleDyM.H:11
HbyA
volVectorField & HbyA
Definition: pEqnPimpleDyM.H:17
mesh
fvMesh & mesh
Definition: createRefsHeatTransfer.H:4
p
volScalarField & p
Definition: createRefsPimpleDyM.H:6
phiHbyA
surfaceScalarField phiHbyA("phiHbyA", fvc::flux(HbyA))
phi
surfaceScalarField & phi
Definition: createRefsPimpleDyM.H:8
tUEqn
tmp< fvVectorMatrix > tUEqn(fvm::ddt(U)+fvm::div(phi, U)+turbulencePtr_->divDevReff(U) - fvSource)
rAtU
tmp< volScalarField > rAtU(rAU)
makeRelative
MRF makeRelative(fvc::interpolate(rho), phiHbyA)
allOptions
const dictionary & allOptions
Definition: createRefsPimpleDyM.H:14
Uf
surfaceVectorField & Uf
Definition: createRefsPimpleDyM.H:15
rAU
volScalarField rAU(1.0/UEqn.A())
HbyAPtr
autoPtr< volVectorField > HbyAPtr
Definition: pEqnPimpleDyM.H:6
correctUfPimpleDyM
this correctUfPimpleDyM(Uf, U, phi)