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 | www.openfoam.com
6 \\/ M anipulation |
7-------------------------------------------------------------------------------
8 Copyright (C) 2015-2017 OpenFOAM Foundation
9-------------------------------------------------------------------------------
10License
11 This file is part of OpenFOAM.
12
13 OpenFOAM is free software: you can redistribute it and/or modify it
14 under the terms of the GNU General Public License as published by
15 the Free Software Foundation, either version 3 of the License, or
16 (at your option) any later version.
17
18 OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
19 ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
20 FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
21 for more details.
22
23 You should have received a copy of the GNU General Public License
24 along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
25
26\*---------------------------------------------------------------------------*/
27
28#include "CorrectPhi.H"
29#include "volFields.H"
30#include "surfaceFields.H"
31#include "fvScalarMatrix.H"
32#include "fvmDdt.H"
33#include "fvmLaplacian.H"
34#include "fvcDiv.H"
37#include "adjustPhi.H"
38#include "fvcMeshPhi.H"
39#include "pimpleControl.H"
40
41// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
42
43template<class RAUfType, class DivUType>
45(
48 const volScalarField& p,
49 const RAUfType& rAUf,
50 const DivUType& divU,
52)
53{
54 const fvMesh& mesh = U.mesh();
55 const Time& runTime = mesh.time();
56
58
59 // Initialize BCs list for pcorr to zero-gradient
61 (
62 p.boundaryField().size(),
64 );
65
66 // Set BCs of pcorr to fixed-value for patches at which p is fixed
67 forAll(p.boundaryField(), patchi)
68 {
69 if (p.boundaryField()[patchi].fixesValue())
70 {
72 }
73 }
74
76 (
78 (
79 "pcorr",
80 runTime.timeName(),
81 mesh
82 ),
83 mesh,
84 dimensionedScalar(p.dimensions(), Zero),
86 );
87
88 if (pcorr.needReference())
89 {
90 fvc::makeRelative(phi, U);
92 fvc::makeAbsolute(phi, U);
93 }
94
95 mesh.setFluxRequired(pcorr.name());
96
97 while (pimple.correctNonOrthogonal())
98 {
99 // Solve for pcorr such that the divergence of the corrected flux
100 // matches the divU provided (from previous iteration, time-step...)
101 fvScalarMatrix pcorrEqn
102 (
103 fvm::laplacian(rAUf, pcorr) == fvc::div(phi) - divU
104 );
105
106 pcorrEqn.setReference(0, 0);
107
108 pcorrEqn.solve
109 (
110 mesh.solver(pcorr.select(pimple.finalNonOrthogonalIter()))
111 );
112
113 if (pimple.finalNonOrthogonalIter())
114 {
115 phi -= pcorrEqn.flux();
116 }
117 }
118}
119
120
121template<class RAUfType, class DivRhoUType>
123(
126 const volScalarField& p,
127 const volScalarField& rho,
128 const volScalarField& psi,
129 const RAUfType& rAUf,
130 const DivRhoUType& divRhoU,
132)
133{
134 const fvMesh& mesh = U.mesh();
135 const Time& runTime = mesh.time();
136
138
139 // Initialize BCs list for pcorr to zero-gradient
141 (
142 p.boundaryField().size(),
144 );
145
146 // Set BCs of pcorr to fixed-value for patches at which p is fixed
147 forAll(p.boundaryField(), patchi)
148 {
149 if (p.boundaryField()[patchi].fixesValue())
150 {
152 }
153 }
154
156 (
158 (
159 "pcorr",
160 runTime.timeName(),
161 mesh
162 ),
163 mesh,
164 dimensionedScalar(p.dimensions(), Zero),
166 );
167
168 mesh.setFluxRequired(pcorr.name());
169
170 while (pimple.correctNonOrthogonal())
171 {
172 // Solve for pcorr such that the divergence of the corrected flux
173 // matches the divRhoU provided (from previous iteration, time-step...)
174 fvScalarMatrix pcorrEqn
175 (
176 fvm::ddt(psi, pcorr)
177 + fvc::div(phi)
178 - fvm::laplacian(rAUf, pcorr)
179 ==
180 divRhoU
181 );
182
183 pcorrEqn.solve
184 (
185 mesh.solver(pcorr.select(pimple.finalNonOrthogonalIter()))
186 );
187
188 if (pimple.finalNonOrthogonalIter())
189 {
190 phi += pcorrEqn.flux();
191 }
192 }
193}
194
195
196// ************************************************************************* //
For cases which do no have a pressure boundary adjust the balance of fluxes to obey continuity....
surfaceScalarField & phi
pimpleControl & pimple
Defines the attributes of an object for which implicit objectRegistry management is supported,...
Definition: IOobject.H:170
Class to control time during OpenFOAM simulations that is also the top-level objectRegistry.
Definition: Time.H:80
A special matrix type and solver, designed for finite volume solutions of scalar equations....
Definition: fvMatrix.H:121
tmp< GeometricField< Type, fvsPatchField, surfaceMesh > > flux() const
Return the face-flux field from the matrix.
Definition: fvMatrix.C:1443
SolverPerformance< Type > solve(const dictionary &)
Solve returning the solution statistics.
void setReference(const label celli, const Type &value, const bool forceReference=false)
Set reference level for solution.
Definition: fvMatrix.C:1002
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:91
PIMPLE control class to supply convergence information/checks for the PIMPLE loop.
Definition: pimpleControl.H:58
U
Definition: pEqn.H:72
volScalarField & p
const volScalarField & psi
volScalarField pcorr(IOobject("pcorr", runTime.timeName(), mesh, IOobject::NO_READ, IOobject::NO_WRITE), mesh, dimensionedScalar(p.dimensions(), Zero), pcorrTypes)
wordList pcorrTypes(p.boundaryField().size(), zeroGradientFvPatchScalarField::typeName)
dynamicFvMesh & mesh
engineTime & runTime
A scalar instance of fvMatrix.
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.
surfaceScalarField rAUf("rAUf", fvc::interpolate(rAU))
adjustPhi(phiHbyA, U, p_rgh)
zeroField divU
Definition: alphaSuSp.H:3
correctUphiBCs(U, phi)
void CorrectPhi(volVectorField &U, surfaceScalarField &phi, const volScalarField &p, const RAUfType &rAUf, const DivUType &divU, pimpleControl &pimple)
Definition: CorrectPhi.C:45
#define forAll(list, i)
Loop across all elements in list.
Definition: stdFoam.H:333
static const char *const typeName
The type name used in ensight case files.
Foam::surfaceFields.