CorrectPhi.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) 2015-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 "CorrectPhi.H"
27 #include "volFields.H"
28 #include "surfaceFields.H"
29 #include "fvmDdt.H"
30 #include "fvmLaplacian.H"
31 #include "fvcDiv.H"
34 #include "adjustPhi.H"
35 #include "fvcMeshPhi.H"
36 #include "pressureReference.H"
38 
39 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
40 
41 template<class RAUfType>
43 (
44  surfaceScalarField& phi,
45  const volVectorField& U,
46  const volScalarField& p,
47  const RAUfType& rAU,
48  const autoPtr<volScalarField>& divU,
50  nonOrthogonalSolutionControl& pcorrControl
51 )
52 {
53  const fvMesh& mesh = phi.mesh();
54  const Time& runTime = mesh.time();
55 
56  // Initialise BCs list for pcorr to zero-gradient
57  wordList pcorrTypes
58  (
59  p.boundaryField().size(),
60  zeroGradientFvPatchScalarField::typeName
61  );
62 
63  // Set BCs of pcorr to fixed-value for patches at which p is fixed
64  forAll(p.boundaryField(), patchi)
65  {
66  if (p.boundaryField()[patchi].fixesValue())
67  {
68  pcorrTypes[patchi] = fixedValueFvPatchScalarField::typeName;
69  }
70  }
71 
72  volScalarField pcorr
73  (
74  IOobject
75  (
76  "pcorr",
77  runTime.name(),
78  mesh
79  ),
80  mesh,
81  dimensionedScalar(p.dimensions(), 0),
82  pcorrTypes
83  );
84 
85  if (pcorr.needReference())
86  {
87  fvc::makeRelative(phi, U);
88  adjustPhi(phi, U, pcorr);
89  fvc::makeAbsolute(phi, U);
90  }
91 
92  mesh.schemes().setFluxRequired(pcorr.name());
93 
94  while (pcorrControl.correctNonOrthogonal())
95  {
96  // Solve for pcorr such that the divergence of the corrected flux
97  // matches the divU provided (from previous iteration, time-step...)
98  fvScalarMatrix pcorrEqn
99  (
100  fvm::laplacian(rAU, pcorr)
101  ==
102  (
103  divU.valid()
104  ? fvc::div(phi) - divU()
105  : fvc::div(phi)
106  )
107  );
108 
109  pcorrEqn.setReference(pressureReference.refCell(), 0);
110 
111  pcorrEqn.solve();
112 
113  if (pcorrControl.finalNonOrthogonalIter())
114  {
115  phi -= pcorrEqn.flux();
116  }
117  }
118 }
119 
120 
121 template<class RAUfType>
123 (
124  surfaceScalarField& phi,
125  const volScalarField& p,
126  const volScalarField& psi,
127  const RAUfType& rAU,
128  const volScalarField& divRhoU,
129  nonOrthogonalSolutionControl& pcorrControl
130 )
131 {
132  const fvMesh& mesh = phi.mesh();
133  const Time& runTime = mesh.time();
134 
135  // Initialise BCs list for pcorr to zero-gradient
136  wordList pcorrTypes
137  (
138  p.boundaryField().size(),
139  zeroGradientFvPatchScalarField::typeName
140  );
141 
142  // Set BCs of pcorr to fixed-value for patches at which p is fixed
143  forAll(p.boundaryField(), patchi)
144  {
145  if (p.boundaryField()[patchi].fixesValue())
146  {
147  pcorrTypes[patchi] = fixedValueFvPatchScalarField::typeName;
148  }
149  }
150 
151  volScalarField pcorr
152  (
153  IOobject
154  (
155  "pcorr",
156  runTime.name(),
157  mesh
158  ),
159  mesh,
160  dimensionedScalar(p.dimensions(), 0),
161  pcorrTypes
162  );
163 
164  mesh.schemes().setFluxRequired(pcorr.name());
165 
166  while (pcorrControl.correctNonOrthogonal())
167  {
168  // Solve for pcorr such that the divergence of the corrected flux
169  // matches the divRhoU provided (from previous iteration, time-step...)
170  fvScalarMatrix pcorrEqn
171  (
172  fvm::ddt(psi, pcorr)
173  + fvc::div(phi)
174  - fvm::laplacian(rAU, pcorr)
175  ==
176  divRhoU
177  );
178 
179  pcorrEqn.solve();
180 
181  if (pcorrControl.finalNonOrthogonalIter())
182  {
183  phi += pcorrEqn.flux();
184  }
185  }
186 }
187 
188 
189 // ************************************************************************* //
Flux correction functions to ensure continuity.
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:434
For cases which do no have a pressure boundary adjust the balance of fluxes to obey continuity....
const Mesh & mesh() const
Return mesh.
Generic GeometricField class.
bool needReference() const
Does the field need a reference level for solution.
IOobject defines the attributes of an object for which implicit objectRegistry management is supporte...
Definition: IOobject.H:99
const word & name() const
Return name.
Definition: IOobject.H:310
Class to control time during OpenFOAM simulations that is also the top-level objectRegistry.
Definition: Time.H:76
An auto-pointer similar to the STL auto_ptr but with automatic casting to a reference to the type and...
Definition: autoPtr.H:51
bool valid() const
Return true if the autoPtr valid (ie, the pointer is set)
Definition: autoPtrI.H:83
const word & name() const
Return const reference to name.
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
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:99
const Time & time() const
Return the top-level database.
Definition: fvMesh.H:418
const fvSchemes & schemes() const
Return the fvSchemes.
Definition: fvMesh.C:1751
void setFluxRequired(const word &name) const
Definition: fvSchemes.C:503
Non-orthogonal solution control class. Provides non-orthogonal-loop control methods.
bool correctNonOrthogonal(const bool finalIter=false)
Non-orthogonal corrector loop.
bool finalNonOrthogonalIter() const
Flag to indicate the last non-orthogonal iteration.
Provides controls for the pressure reference in closed-volume simulations.
label refCell() const
Return the cell in which the reference pressure is set.
Calculate the divergence of the given field.
Calculate the mesh motion flux and convert fluxes from absolute to relative and back.
Calculate the matrix for the first temporal derivative.
Calculate the matrix for the laplacian of the field.
label patchi
const volScalarField & psi
volScalarField rAU(1.0/UEqn.A())
U
Definition: pEqn.H:72
adjustPhi(phiHbyA, U, p)
void CorrectPhi(surfaceScalarField &phi, const volVectorField &U, const volScalarField &p, const RAUfType &rAU, const autoPtr< volScalarField > &divU, const pressureReference &pressureReference, nonOrthogonalSolutionControl &pcorrControl)
Definition: CorrectPhi.C:43
tmp< VolField< Type > > ddt(const dimensioned< Type > dt, const fvMesh &mesh)
Definition: fvcDdt.C:45
tmp< VolField< Type > > laplacian(const VolField< Type > &vf, const word &name)
Definition: fvcLaplacian.C:45
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
volScalarField & p
Foam::surfaceFields.