pEqn.H
Go to the documentation of this file.
1 rho = thermo.rho();
2 
3 volScalarField rAU(1.0/UEqn.A());
6 
7 if (pimple.transonic())
8 {
10  (
11  "phid",
13  *(
14  (
15  fvc::flux(HbyA)
17  )
18  )
19  );
20 
21  fvc::makeRelative(phid, psi, U);
22  MRF.makeRelative(fvc::interpolate(psi), phid);
23 
24  while (pimple.correctNonOrthogonal())
25  {
26  fvScalarMatrix pEqn
27  (
28  fvm::ddt(psi, p)
29  + fvm::div(phid, p)
31  ==
32  fvOptions(psi, p, rho.name())
33  );
34 
35  pEqn.solve(mesh.solver(p.select(pimple.finalInnerIter())));
36 
37  if (pimple.finalNonOrthogonalIter())
38  {
39  phi == pEqn.flux();
40  }
41  }
42 }
43 else
44 {
46  (
47  "phiHbyA",
48  (
49  fvc::flux(rho*HbyA)
51  )
52  );
53 
54  fvc::makeRelative(phiHbyA, rho, U);
55  MRF.makeRelative(fvc::interpolate(rho), phiHbyA);
56 
57  // Update the pressure BCs to ensure flux consistency
58  constrainPressure(p, rho, U, phiHbyA, rhorAUf, MRF);
59 
60  while (pimple.correctNonOrthogonal())
61  {
62  fvScalarMatrix pEqn
63  (
64  fvm::ddt(psi, p)
65  + fvc::div(phiHbyA)
67  ==
68  fvOptions(psi, p, rho.name())
69  );
70 
71  pEqn.solve(mesh.solver(p.select(pimple.finalInnerIter())));
72 
73  if (pimple.finalNonOrthogonalIter())
74  {
75  phi = phiHbyA + pEqn.flux();
76  }
77  }
78 }
79 
80 #include "rhoEqn.H"
81 #include "compressibleContinuityErrs.H"
82 
83 U = HbyA - rAU*fvc::grad(p);
84 U.correctBoundaryConditions();
85 fvOptions.correct(U);
86 K = 0.5*magSqr(U);
87 
88 {
90  surfaceVectorField n(mesh.Sf()/mesh.magSf());
91  rhoUf += n*(fvc::absolute(phi, rho, U)/mesh.magSf() - (n & rhoUf));
92 }
93 
94 if (thermo.dpdt())
95 {
97 }
tmp< GeometricField< typename outerProduct< vector, Type >::type, fvPatchField, volMesh >> grad(const GeometricField< Type, fvsPatchField, surfaceMesh > &ssf)
Definition: fvcGrad.C:52
surfaceVectorField n(mesh.Sf()/mesh.magSf())
fvMatrix< scalar > fvScalarMatrix
Definition: fvMatricesFwd.H:42
U
Definition: pEqn.H:83
p
Definition: pEqn.H:50
phiHbyA
Definition: pEqn.H:20
tmp< GeometricField< Type, fvPatchField, volMesh > > div(const GeometricField< Type, fvsPatchField, surfaceMesh > &ssf)
Definition: fvcDiv.C:47
tmp< GeometricField< typename flux< Type >::type, fvsPatchField, surfaceMesh > > ddtCorr(const GeometricField< Type, fvPatchField, volMesh > &U, const GeometricField< Type, fvsPatchField, surfaceMesh > &Uf)
Definition: fvcDdt.C:170
GeometricField< vector, fvsPatchField, surfaceMesh > surfaceVectorField
const dictionary & pimple
GeometricField< vector, fvPatchField, volMesh > volVectorField
Definition: volFieldsFwd.H:55
tmp< GeometricField< Type, fvPatchField, volMesh > > laplacian(const GeometricField< Type, fvPatchField, volMesh > &vf, const word &name)
Definition: fvcLaplacian.C:45
CGAL::Exact_predicates_exact_constructions_kernel K
fv::options & fvOptions
surfaceScalarField rhorAUf("rhorAUf", fvc::interpolate(rho *rAU))
tmp< GeometricField< Type, fvPatchField, volMesh > > ddt(const dimensioned< Type > dt, const fvMesh &mesh)
Definition: fvcDdt.C:45
GeometricField< scalar, fvPatchField, volMesh > volScalarField
Definition: volFieldsFwd.H:52
psiReactionThermo & thermo
Definition: createFields.H:31
dynamicFvMesh & mesh
IOMRFZoneList & MRF
volScalarField & dpdt
tmp< volVectorField > constrainHbyA(const tmp< volVectorField > &tHbyA, const volVectorField &U, const volScalarField &p)
Definition: constrainHbyA.C:33
dimensioned< scalar > magSqr(const dimensioned< Type > &)
constrainPressure(p, rho, U, phiHbyA, rhorAUf, MRF)
static tmp< GeometricField< Type, fvsPatchField, surfaceMesh > > interpolate(const GeometricField< Type, fvPatchField, volMesh > &tvf, const surfaceScalarField &faceFlux, Istream &schemeData)
Interpolate field onto faces using scheme given by Istream.
surfaceScalarField phid("phid", fvc::interpolate(psi) *(fvc::flux(HbyA)+rhorAUf *fvc::ddtCorr(rho, U, phi)/fvc::interpolate(rho)))
tmp< surfaceScalarField > absolute(const tmp< surfaceScalarField > &tphi, const volVectorField &U)
Return the given relative flux in absolute form.
Definition: fvcMeshPhi.C:188
rhoUf
Definition: pEqn.H:91
volVectorField & HbyA
Definition: pEqn.H:13
fvVectorMatrix & UEqn
Definition: UEqn.H:13
tmp< surfaceScalarField > meshPhi(const volVectorField &U)
Definition: fvcMeshPhi.C:34
phi
Definition: pEqn.H:18
const volScalarField & psi
tmp< surfaceScalarField > flux(const volVectorField &vvf)
Return the face-flux field obtained from the given volVectorField.
Definition: fvcFlux.C:32
GeometricField< scalar, fvsPatchField, surfaceMesh > surfaceScalarField
rho
Definition: pEqn.H:1
volScalarField rAU(1.0/UEqn.A())
MRF makeRelative(fvc::interpolate(rho), phiHbyA)