All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
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-2021 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"
37 #include "pressureReference.H"
39 
40 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
41 
42 template<class RAUfType, class DivUType>
44 (
45  surfaceScalarField& phi,
46  const volVectorField& U,
47  const volScalarField& p,
48  const RAUfType& rAUf,
49  const DivUType& divU,
50  const pressureReference& pressureReference,
51  nonOrthogonalSolutionControl& pcorrControl
52 )
53 {
54  const fvMesh& mesh = phi.mesh();
55  const Time& runTime = mesh.time();
56 
57  // Initialise BCs list for pcorr to zero-gradient
59  (
60  p.boundaryField().size(),
61  zeroGradientFvPatchScalarField::typeName
62  );
63 
64  // Set BCs of pcorr to fixed-value for patches at which p is fixed
65  forAll(p.boundaryField(), patchi)
66  {
67  if (p.boundaryField()[patchi].fixesValue())
68  {
69  pcorrTypes[patchi] = fixedValueFvPatchScalarField::typeName;
70  }
71  }
72 
73  volScalarField pcorr
74  (
75  IOobject
76  (
77  "pcorr",
78  runTime.timeName(),
79  mesh
80  ),
81  mesh,
82  dimensionedScalar(p.dimensions(), 0),
84  );
85 
86  if (pcorr.needReference())
87  {
88  fvc::makeRelative(phi, U);
89  adjustPhi(phi, U, pcorr);
90  fvc::makeAbsolute(phi, U);
91  }
92 
93  mesh.setFluxRequired(pcorr.name());
94 
95  while (pcorrControl.correctNonOrthogonal())
96  {
97  // Solve for pcorr such that the divergence of the corrected flux
98  // matches the divU provided (from previous iteration, time-step...)
99  fvScalarMatrix pcorrEqn
100  (
101  fvm::laplacian(rAUf, pcorr) == fvc::div(phi) - divU
102  );
103 
104  pcorrEqn.setReference(pressureReference.refCell(), 0);
105 
106  pcorrEqn.solve();
107 
108  if (pcorrControl.finalNonOrthogonalIter())
109  {
110  phi -= pcorrEqn.flux();
111  }
112  }
113 }
114 
115 
116 template<class RAUfType, class DivRhoUType>
117 void Foam::CorrectPhi
118 (
119  surfaceScalarField& phi,
120  const volScalarField& p,
121  const volScalarField& rho,
122  const volScalarField& psi,
123  const RAUfType& rAUf,
124  const DivRhoUType& divRhoU,
125  nonOrthogonalSolutionControl& pcorrControl
126 )
127 {
128  const fvMesh& mesh = phi.mesh();
129  const Time& runTime = mesh.time();
130 
131  // Initialise BCs list for pcorr to zero-gradient
133  (
134  p.boundaryField().size(),
135  zeroGradientFvPatchScalarField::typeName
136  );
137 
138  // Set BCs of pcorr to fixed-value for patches at which p is fixed
139  forAll(p.boundaryField(), patchi)
140  {
141  if (p.boundaryField()[patchi].fixesValue())
142  {
143  pcorrTypes[patchi] = fixedValueFvPatchScalarField::typeName;
144  }
145  }
146 
147  volScalarField pcorr
148  (
149  IOobject
150  (
151  "pcorr",
152  runTime.timeName(),
153  mesh
154  ),
155  mesh,
156  dimensionedScalar(p.dimensions(), 0),
157  pcorrTypes
158  );
159 
160  mesh.setFluxRequired(pcorr.name());
161 
162  while (pcorrControl.correctNonOrthogonal())
163  {
164  // Solve for pcorr such that the divergence of the corrected flux
165  // matches the divRhoU provided (from previous iteration, time-step...)
166  fvScalarMatrix pcorrEqn
167  (
168  fvm::ddt(psi, pcorr)
169  + fvc::div(phi)
170  - fvm::laplacian(rAUf, pcorr)
171  ==
172  divRhoU
173  );
174 
175  pcorrEqn.solve();
176 
177  if (pcorrControl.finalNonOrthogonalIter())
178  {
179  phi += pcorrEqn.flux();
180  }
181  }
182 }
183 
184 
185 // ************************************************************************* //
Foam::surfaceFields.
fvMatrix< scalar > fvScalarMatrix
Definition: fvMatricesFwd.H:42
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:434
tmp< GeometricField< Type, fvPatchField, volMesh > > div(const GeometricField< Type, fvsPatchField, surfaceMesh > &ssf)
Definition: fvcDiv.C:47
Calculate the matrix for the laplacian of the field.
GeometricField< vector, fvPatchField, volMesh > volVectorField
Definition: volFieldsFwd.H:58
tmp< GeometricField< Type, fvPatchField, volMesh > > laplacian(const GeometricField< Type, fvPatchField, volMesh > &vf, const word &name)
Definition: fvcLaplacian.C:45
bool adjustPhi(surfaceScalarField &phi, const volVectorField &U, volScalarField &p)
Adjust the balance of fluxes to obey continuity.
Definition: adjustPhi.C:34
void CorrectPhi(surfaceScalarField &phi, const volVectorField &U, const volScalarField &p, const RAUfType &rAUf, const DivUType &divU, const pressureReference &pressureReference, nonOrthogonalSolutionControl &pcorrControl)
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:57
dynamicFvMesh & mesh
Calculate the matrix for the first temporal derivative.
Calculate the divergence of the given field.
Calculate the mesh motion flux and convert fluxes from absolute to relative and back.
List< word > wordList
A List of words.
Definition: fileName.H:54
label patchi
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
A scalar instance of fvMatrix.
wordList pcorrTypes(p.boundaryField().size(), zeroGradientFvPatchScalarField::typeName)
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...
GeometricField< scalar, fvsPatchField, surfaceMesh > surfaceScalarField
void makeRelative(surfaceScalarField &phi, const volVectorField &U)
Make the given flux relative.
Definition: fvcMeshPhi.C:75