CorrectPhi.C
Go to the documentation of this file.
1 /*---------------------------------------------------------------------------*\
2  ========= |
3  \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
4  \\ / O peration |
5  \\ / A nd | Copyright (C) 2015-2017 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 "pimpleControl.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  pimpleControl& pimple
50 )
51 {
52  const fvMesh& mesh = U.mesh();
53  const Time& runTime = mesh.time();
54 
55  correctUphiBCs(U, phi);
56 
57  // Initialize 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
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("pcorr", p.dimensions(), 0.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 (pimple.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(0, 0);
105 
106  pcorrEqn.solve
107  (
108  mesh.solver(pcorr.select(pimple.finalNonOrthogonalIter()))
109  );
110 
111  if (pimple.finalNonOrthogonalIter())
112  {
113  phi -= pcorrEqn.flux();
114  }
115  }
116 }
117 
118 
119 template<class RAUfType, class DivRhoUType>
120 void Foam::CorrectPhi
121 (
122  volVectorField& U,
123  surfaceScalarField& phi,
124  const volScalarField& p,
125  const volScalarField& rho,
126  const volScalarField& psi,
127  const RAUfType& rAUf,
128  const DivRhoUType& divRhoU,
129  pimpleControl& pimple
130 )
131 {
132  const fvMesh& mesh = U.mesh();
133  const Time& runTime = mesh.time();
134 
135  correctUphiBCs(rho, U, phi);
136 
137  // Initialize BCs list for pcorr to zero-gradient
139  (
140  p.boundaryField().size(),
141  zeroGradientFvPatchScalarField::typeName
142  );
143 
144  // Set BCs of pcorr to fixed-value for patches at which p is fixed
146  {
147  if (p.boundaryField()[patchi].fixesValue())
148  {
149  pcorrTypes[patchi] = fixedValueFvPatchScalarField::typeName;
150  }
151  }
152 
153  volScalarField pcorr
154  (
155  IOobject
156  (
157  "pcorr",
158  runTime.timeName(),
159  mesh
160  ),
161  mesh,
162  dimensionedScalar("pcorr", p.dimensions(), 0.0),
163  pcorrTypes
164  );
165 
166  mesh.setFluxRequired(pcorr.name());
167 
168  while (pimple.correctNonOrthogonal())
169  {
170  // Solve for pcorr such that the divergence of the corrected flux
171  // matches the divRhoU provided (from previous iteration, time-step...)
172  fvScalarMatrix pcorrEqn
173  (
174  fvm::ddt(psi, pcorr)
175  + fvc::div(phi)
176  - fvm::laplacian(rAUf, pcorr)
177  ==
178  divRhoU
179  );
180 
181  pcorrEqn.solve
182  (
183  mesh.solver(pcorr.select(pimple.finalNonOrthogonalIter()))
184  );
185 
186  if (pimple.finalNonOrthogonalIter())
187  {
188  phi += pcorrEqn.flux();
189  }
190  }
191 }
192 
193 
194 // ************************************************************************* //
void CorrectPhi(volVectorField &U, surfaceScalarField &phi, const volScalarField &p, const RAUfType &rAUf, const DivUType &divU, pimpleControl &pimple)
Definition: CorrectPhi.C:43
Foam::surfaceFields.
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:428
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:508
const dictionary & solver(const word &name) const
Return the solver controls dictionary for the given field.
Definition: solution.C:353
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:644
const Time & time() const
Return the top-level database.
Definition: fvMesh.H:243
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
bool adjustPhi(surfaceScalarField &phi, const volVectorField &U, volScalarField &p)
Adjust the balance of fluxes to obey continuity.
Definition: adjustPhi.C:34
tmp< GeometricField< Type, fvPatchField, volMesh > > ddt(const dimensioned< Type > dt, const fvMesh &mesh)
Definition: fvcDdt.C:45
const dimensionSet & dimensions() const
Return dimensions.
dynamicFvMesh & mesh
bool correctNonOrthogonal()
Non-orthogonal corrector loop.
Calulate 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
bool finalNonOrthogonalIter() const
Helper function to identify final non-orthogonal iteration.
Calculate the divergence of the given field.
SolverPerformance< Type > solve(const dictionary &)
Solve segregated or coupled returning the solution statistics.
Definition: fvMatrixSolve.C:57
const Mesh & mesh() const
Return mesh.
Calculate the mesh motion flux and convert fluxes from absolute to relative and back.
correctUphiBCs(U, phi)
label patchi
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
A scalar instance of fvMatrix.
PIMPLE control class to supply convergence information/checks for the PIMPLE loop.
Definition: pimpleControl.H:52
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:876
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