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-2019 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 "fvScalarMatrix.H"
30 #include "fvmDdt.H"
31 #include "fvmLaplacian.H"
32 #include "fvcDiv.H"
35 #include "adjustPhi.H"
36 #include "fvcMeshPhi.H"
38 
39 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
40 
41 template<class RAUfType, class DivUType>
43 (
44  volVectorField& U,
45  surfaceScalarField& phi,
46  const volScalarField& p,
47  const RAUfType& rAUf,
48  const DivUType& divU,
49  nonOrthogonalSolutionControl& pcorrControl,
50  const bool evaluateUBCs
51 )
52 {
53  const fvMesh& mesh = U.mesh();
54  const Time& runTime = mesh.time();
55 
56  correctUphiBCs(U, phi, evaluateUBCs);
57 
58  // Initialize BCs list for pcorr to zero-gradient
60  (
61  p.boundaryField().size(),
62  zeroGradientFvPatchScalarField::typeName
63  );
64 
65  // Set BCs of pcorr to fixed-value for patches at which p is fixed
67  {
68  if (p.boundaryField()[patchi].fixesValue())
69  {
70  pcorrTypes[patchi] = fixedValueFvPatchScalarField::typeName;
71  }
72  }
73 
74  volScalarField pcorr
75  (
76  IOobject
77  (
78  "pcorr",
79  runTime.timeName(),
80  mesh
81  ),
82  mesh,
85  );
86 
87  if (pcorr.needReference())
88  {
89  fvc::makeRelative(phi, U);
90  adjustPhi(phi, U, pcorr);
91  fvc::makeAbsolute(phi, U);
92  }
93 
94  mesh.setFluxRequired(pcorr.name());
95 
96  while (pcorrControl.correctNonOrthogonal())
97  {
98  // Solve for pcorr such that the divergence of the corrected flux
99  // matches the divU provided (from previous iteration, time-step...)
100  fvScalarMatrix pcorrEqn
101  (
102  fvm::laplacian(rAUf, pcorr) == fvc::div(phi) - divU
103  );
104 
105  pcorrEqn.setReference(0, 0);
106 
107  pcorrEqn.solve();
108 
109  if (pcorrControl.finalNonOrthogonalIter())
110  {
111  phi -= pcorrEqn.flux();
112  }
113  }
114 }
115 
116 
117 template<class RAUfType, class DivRhoUType>
118 void Foam::CorrectPhi
119 (
120  volVectorField& U,
121  surfaceScalarField& phi,
122  const volScalarField& p,
123  const volScalarField& rho,
124  const volScalarField& psi,
125  const RAUfType& rAUf,
126  const DivRhoUType& divRhoU,
127  nonOrthogonalSolutionControl& pcorrControl,
128  const bool evaluateUBCs
129 )
130 {
131  const fvMesh& mesh = U.mesh();
132  const Time& runTime = mesh.time();
133 
134  correctUphiBCs(rho, U, phi, evaluateUBCs);
135 
136  // Initialize BCs list for pcorr to zero-gradient
138  (
139  p.boundaryField().size(),
140  zeroGradientFvPatchScalarField::typeName
141  );
142 
143  // Set BCs of pcorr to fixed-value for patches at which p is fixed
145  {
146  if (p.boundaryField()[patchi].fixesValue())
147  {
148  pcorrTypes[patchi] = fixedValueFvPatchScalarField::typeName;
149  }
150  }
151 
152  volScalarField pcorr
153  (
154  IOobject
155  (
156  "pcorr",
157  runTime.timeName(),
158  mesh
159  ),
160  mesh,
162  pcorrTypes
163  );
164 
165  mesh.setFluxRequired(pcorr.name());
166 
167  while (pcorrControl.correctNonOrthogonal())
168  {
169  // Solve for pcorr such that the divergence of the corrected flux
170  // matches the divRhoU provided (from previous iteration, time-step...)
171  fvScalarMatrix pcorrEqn
172  (
173  fvm::ddt(psi, pcorr)
174  + fvc::div(phi)
175  - fvm::laplacian(rAUf, pcorr)
176  ==
177  divRhoU
178  );
179 
180  pcorrEqn.solve();
181 
182  if (pcorrControl.finalNonOrthogonalIter())
183  {
184  phi += pcorrEqn.flux();
185  }
186  }
187 }
188 
189 
190 // ************************************************************************* //
Foam::surfaceFields.
bool finalNonOrthogonalIter() const
Flag to indicate the last non-orthogonal iteration.
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:434
const word & name() const
Return name.
Definition: IOobject.H:303
bool correctNonOrthogonal()
Non-orthogonal corrector loop.
const Boundary & boundaryField() const
Return const-reference to the boundary field.
tmp< GeometricField< Type, fvPatchField, volMesh > > div(const GeometricField< Type, fvsPatchField, surfaceMesh > &ssf)
Definition: fvcDiv.C:47
void setReference(const label celli, const Type &value, const bool forceReference=false)
Set reference level for solution.
Definition: fvMatrix.C:506
engineTime & runTime
Calculate the matrix for the laplacian of the field.
void setFluxRequired(const word &name) const
Definition: fvSchemes.C:497
static word timeName(const scalar, const int precision=precision_)
Return time name of given scalar time.
Definition: Time.C:622
const Time & time() const
Return the top-level database.
Definition: fvMesh.H:239
tmp< GeometricField< Type, fvPatchField, volMesh > > laplacian(const GeometricField< Type, fvPatchField, volMesh > &vf, const word &name)
Definition: fvcLaplacian.C:45
Class to control time during OpenFOAM simulations that is also the top-level objectRegistry.
Definition: Time.H:68
tmp< GeometricField< Type, fvPatchField, volMesh > > ddt(const dimensioned< Type > dt, const fvMesh &mesh)
Definition: fvcDdt.C:45
bool adjustPhi(surfaceScalarField &phi, const volVectorField &U, volScalarField &p)
Adjust the balance of fluxes to obey continuity.
bool needReference() const
Does the field need a reference level for solution.
const dimensionSet & dimensions() const
Return dimensions.
dynamicFvMesh & mesh
Calculate the matrix for the first temporal derivative.
A special matrix type and solver, designed for finite volume solutions of scalar equations. Face addressing is used to make all matrix assembly and solution loops vectorise.
Definition: fvPatchField.H:72
Calculate the divergence of the given field.
SolverPerformance< Type > solve(const dictionary &)
Solve segregated or coupled returning the solution statistics.
Definition: fvMatrixSolve.C:58
const Mesh & mesh() const
Return mesh.
void CorrectPhi(volVectorField &U, surfaceScalarField &phi, const volScalarField &p, const RAUfType &rAUf, const DivUType &divU, nonOrthogonalSolutionControl &pcorrControl, const bool evaluateUBCs)
Definition: CorrectPhi.C:43
Calculate the mesh motion flux and convert fluxes from absolute to relative and back.
label patchi
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
A scalar instance of fvMatrix.
wordList pcorrTypes(p.boundaryField().size(), zeroGradientFvPatchScalarField::typeName)
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:78
void makeAbsolute(surfaceScalarField &phi, const volVectorField &U)
Make the given flux absolute.
Definition: fvcMeshPhi.C:114
For cases which do no have a pressure boundary adjust the balance of fluxes to obey continuity...
tmp< GeometricField< Type, fvsPatchField, surfaceMesh > > flux() const
Return the face-flux field from the matrix.
Definition: fvMatrix.C:856
void makeRelative(surfaceScalarField &phi, const volVectorField &U)
Make the given flux relative.
Definition: fvcMeshPhi.C:75
IOobject defines the attributes of an object for which implicit objectRegistry management is supporte...
Definition: IOobject.H:92
Non-orthogonal solution control class. Provides non-orthogonal-loop control methods.
correctUphiBCs(U, phi, true)