pressureCorrector.C
Go to the documentation of this file.
1 /*---------------------------------------------------------------------------*\
2  ========= |
3  \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
4  \\ / O peration | Website: https://openfoam.org
5  \\ / A nd | Copyright (C) 2023 OpenFOAM Foundation
6  \\/ M anipulation |
7 -------------------------------------------------------------------------------
8 License
9  This file is part of OpenFOAM.
10 
11  OpenFOAM is free software: you can redistribute it and/or modify it
12  under the terms of the GNU General Public License as published by
13  the Free Software Foundation, either version 3 of the License, or
14  (at your option) any later version.
15 
16  OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
17  ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
18  FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
19  for more details.
20 
21  You should have received a copy of the GNU General Public License
22  along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
23 
24 \*---------------------------------------------------------------------------*/
25 
26 #include "twoPhaseSolver.H"
27 #include "constrainHbyA.H"
28 #include "constrainPressure.H"
29 #include "adjustPhi.H"
30 #include "findRefCell.H"
31 #include "fvcMeshPhi.H"
32 #include "fvcFlux.H"
33 #include "fvcDdt.H"
34 #include "fvcDiv.H"
35 #include "fvcSnGrad.H"
36 #include "fvcReconstruct.H"
37 #include "fvmLaplacian.H"
38 
39 // * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * * //
40 
42 (
44 )
45 {
46  volVectorField& U = U_;
48 
49  fvVectorMatrix& UEqn = tUEqn.ref();
50  setrAU(UEqn);
51 
52  const surfaceScalarField rAUf("rAUf", fvc::interpolate(rAU()));
53 
54  while (pimple.correct())
55  {
58  (
59  "phiHbyA",
62  );
63 
65 
66  if (p_rgh.needReference())
67  {
71  }
72 
74  (
75  (
78  )*rAUf*mesh.magSf()
79  );
80 
81  phiHbyA += phig;
82 
83  // Update the pressure BCs to ensure flux consistency
85 
86  // Cache any sources
87  fvScalarMatrix p_rghEqnSource
88  (
89  fvModels().sourceProxy(alpha1, p_rgh)
90  + fvModels().sourceProxy(alpha2, p_rgh)
91  );
92 
94  {
95  fvScalarMatrix p_rghEqn
96  (
98  == p_rghEqnSource
99  );
100 
101  p_rghEqn.setReference
102  (
103  pressureReference().refCell(),
105  );
106 
107  p_rghEqn.solve();
108 
110  {
111  phi = phiHbyA + p_rghEqn.flux();
112 
113  p_rgh.relax();
114 
115  U = HbyA
116  + rAU()*fvc::reconstruct((phig + p_rghEqn.flux())/rAUf);
119  }
120  }
121 
123 
124  // Correct Uf if the mesh is moving
125  fvc::correctUf(Uf, U, phi, MRF);
126 
127  // Make the fluxes relative to the mesh motion
129 
130  p == p_rgh + rho*buoyancy.gh;
131 
132  if (p_rgh.needReference())
133  {
135  (
136  "p",
137  p.dimensions(),
138  pressureReference().refValue()
139  - getRefCellValue(p, pressureReference().refCell())
140  );
141  p_rgh = p - rho*buoyancy.gh;
142  }
143  }
144 
145  clearrAU();
146  tUEqn.clear();
147 }
148 
149 
150 // ************************************************************************* //
For cases which do no have a pressure boundary adjust the balance of fluxes to obey continuity....
Generic GeometricField class.
void relax(const scalar alpha)
Relax field (for steady-state solution).
bool needReference() const
Does the field need a reference level for solution.
void correctBoundaryConditions()
Correct boundary field.
void makeRelative(volVectorField &U) const
Make the given absolute velocity relative within the MRF region.
Definition: MRFZoneList.C:157
bool constrain(fvMatrix< Type > &eqn) const
Apply constraints to an equation.
A special matrix type and solver, designed for finite volume solutions of scalar equations....
Definition: fvMatrix.H:118
SolverPerformance< Type > solve(const dictionary &)
Solve segregated or coupled returning the solution statistics.
Definition: fvMatrixSolve.C:58
tmp< SurfaceField< Type > > flux() const
Return the face-flux field from the matrix.
Definition: fvMatrix.C:964
void setReference(const label celli, const Type &value, const bool forceReference=false)
Set reference level for solution.
Definition: fvMatrix.C:563
const surfaceScalarField & magSf() const
Return cell face area magnitudes.
bool finalNonOrthogonalIter() const
Flag to indicate the last non-orthogonal iteration.
bool correct()
Piso loop within outer loop.
bool correctNonOrthogonal()
Non-orthogonal corrector loop.
Provides controls for the pressure reference in closed-volume simulations.
Foam::fvModels & fvModels() const
Return the fvModels that are created on demand.
Definition: solver.C:96
pimpleNoLoopControl pimple
PIMPLE inner-loop controls.
Definition: solver.H:107
Foam::fvConstraints & fvConstraints() const
Return the fvConstraints that are created on demand.
Definition: solver.C:107
const fvMesh & mesh
Region mesh.
Definition: solver.H:101
volScalarField & p_rgh
Reference to the buoyant pressure for buoyant cases.
Definition: VoFSolver.H:107
const surfaceScalarField & phi
Reference to the mass-flux field.
Definition: VoFSolver.H:212
volVectorField U_
Velocity field.
Definition: VoFSolver.H:94
tmp< fvVectorMatrix > tUEqn
Cached momentum matrix.
Definition: VoFSolver.H:145
autoPtr< surfaceVectorField > Uf
Pointer to the surface momentum field.
Definition: VoFSolver.H:133
autoPtr< volScalarField > rAU
Inverse momentum equation diagonal.
Definition: VoFSolver.H:129
const volVectorField & U
Reference to the velocity field.
Definition: VoFSolver.H:209
virtual const Foam::pressureReference & pressureReference() const =0
Return the pressure reference.
void clearrAU()
Clear the cached rAU is no longer needed.
Definition: VoFSolver.C:66
IOMRFZoneList MRF
MRF zone list.
Definition: VoFSolver.H:122
surfaceScalarField phi_
Volumetric flux field.
Definition: VoFSolver.H:97
const volScalarField & rho
Reference to the mixture continuity density field.
Definition: VoFSolver.H:110
void continuityErrors()
Calculate and print the continuity errors.
Definition: VoFSolver.C:45
void setrAU(const fvVectorMatrix &UEqn)
Set or update the cached rAU.
Definition: VoFSolver.C:53
Buoyancy related data for the Foam::solvers::isothermalFluid solver module when solving buoyant cases...
Definition: buoyancy.H:70
volScalarField gh
(g & h) - ghRef
Definition: buoyancy.H:95
surfaceScalarField ghf
(g & hf) - ghRef
Definition: buoyancy.H:98
virtual tmp< surfaceScalarField > surfaceTensionForce() const =0
Return the interface surface tension force for the momentum equation.
volScalarField & alpha1
Reference to the phase1-fraction.
volScalarField & alpha2
Reference to the phase2-fraction.
void incompressiblePressureCorrector(volScalarField &p)
Construct and solve the incompressible pressure equation.
fvVectorMatrix & UEqn
Definition: UEqn.H:11
Find the reference cell nearest (in index) to the given cell but which is not on a cyclic,...
Calculate the first temporal derivative.
Calculate the divergence of the given field.
Calculate the face-flux of the given field.
Calculate the mesh motion flux and convert fluxes from absolute to relative and back.
Reconstruct volField from a face flux field.
Calculate the snGrad of the given volField.
Calculate the matrix for the laplacian of the field.
volVectorField & HbyA
Definition: pEqn.H:13
surfaceScalarField phiHbyA("phiHbyA", fvc::interpolate(rho) *fvc::flux(HbyA))
tmp< SurfaceField< typename innerProduct< vector, Type >::type > > flux(const VolField< Type > &vf)
Return the face-flux field obtained from the given volVectorField.
static tmp< SurfaceField< Type > > interpolate(const VolField< Type > &tvf, const surfaceScalarField &faceFlux, Istream &schemeData)
Interpolate field onto faces using scheme given by Istream.
tmp< VolField< typename outerProduct< vector, Type >::type > > reconstruct(const SurfaceField< Type > &ssf)
void correctUf(autoPtr< surfaceVectorField > &Uf, const volVectorField &U, const surfaceScalarField &phi, const MRFType &MRF)
tmp< VolField< Type > > div(const SurfaceField< Type > &ssf)
Definition: fvcDiv.C:47
void makeAbsolute(surfaceScalarField &phi, const volVectorField &U)
Make the given flux absolute.
Definition: fvcMeshPhi.C:128
void makeRelative(surfaceScalarField &phi, const volVectorField &U)
Make the given flux relative.
Definition: fvcMeshPhi.C:89
tmp< SurfaceField< Type > > snGrad(const VolField< Type > &vf, const word &name)
Definition: fvcSnGrad.C:45
tmp< SurfaceField< typename Foam::flux< Type >::type > > ddtCorr(const VolField< Type > &U, const SurfaceField< Type > &Uf)
Definition: fvcDdt.C:196
tmp< fvMatrix< Type > > laplacian(const VolField< Type > &vf, const word &name)
Definition: fvmLaplacian.C:47
void constrainPressure(volScalarField &p, const RhoType &rho, const volVectorField &U, const surfaceScalarField &phiHbyA, const RAUType &rhorAU, const MRFType &MRF)
bool adjustPhi(surfaceScalarField &phi, const volVectorField &U, volScalarField &p)
Adjust the balance of fluxes to obey continuity.
Definition: adjustPhi.C:35
void constrainHbyA(volVectorField &HbyA, const volVectorField &U, const volScalarField &p)
Definition: constrainHbyA.C:34
scalar getRefCellValue(const volScalarField &field, const label refCelli)
Return the current value of field in the reference cell.
Definition: findRefCell.C:136
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
volScalarField & p