fvcMeshPhi.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) 2011-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 "fvcMeshPhi.H"
27 #include "fvMesh.H"
28 #include "ddtScheme.H"
29 #include "surfaceInterpolate.H"
30 
31 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
32 
34 (
35  const volVectorField& vf
36 )
37 {
39  (
40  vf.mesh(),
41  vf.mesh().ddtScheme("ddt(" + vf.name() + ')')
42  ).ref().meshPhi(vf);
43 }
44 
45 
47 (
48  const dimensionedScalar& rho,
49  const volVectorField& vf
50 )
51 {
53  (
54  vf.mesh(),
55  vf.mesh().ddtScheme("ddt(" + rho.name() + ',' + vf.name() + ')')
56  ).ref().meshPhi(vf);
57 }
58 
59 
61 (
62  const volScalarField& rho,
63  const volVectorField& vf
64 )
65 {
67  (
68  vf.mesh(),
69  vf.mesh().ddtScheme("ddt(" + rho.name() + ',' + vf.name() + ')')
70  ).ref().meshPhi(vf);
71 }
72 
73 
75 (
76  surfaceScalarField& phi,
77  const volVectorField& U
78 )
79 {
80  if (phi.mesh().moving())
81  {
82  phi -= fvc::meshPhi(U);
83  }
84 }
85 
87 (
88  surfaceScalarField& phi,
89  const dimensionedScalar& rho,
90  const volVectorField& U
91 )
92 {
93  if (phi.mesh().moving())
94  {
95  phi -= rho*fvc::meshPhi(rho, U);
96  }
97 }
98 
100 (
101  surfaceScalarField& phi,
102  const volScalarField& rho,
103  const volVectorField& U
104 )
105 {
106  if (phi.mesh().moving())
107  {
108  phi -= fvc::interpolate(rho)*fvc::meshPhi(rho, U);
109  }
110 }
111 
112 
114 (
115  surfaceScalarField& phi,
116  const volVectorField& U
117 )
118 {
119  if (phi.mesh().moving())
120  {
121  phi += fvc::meshPhi(U);
122  }
123 }
124 
126 (
127  surfaceScalarField& phi,
128  const dimensionedScalar& rho,
129  const volVectorField& U
130 )
131 {
132  if (phi.mesh().moving())
133  {
134  phi += rho*fvc::meshPhi(rho, U);
135  }
136 }
137 
139 (
140  surfaceScalarField& phi,
141  const volScalarField& rho,
142  const volVectorField& U
143 )
144 {
145  if (phi.mesh().moving())
146  {
147  phi += fvc::interpolate(rho)*fvc::meshPhi(rho, U);
148  }
149 }
150 
151 
153 (
154  const tmp<surfaceScalarField>& tphi,
155  const volVectorField& U
156 )
157 {
158  if (tphi().mesh().moving())
159  {
160  return tphi - fvc::meshPhi(U);
161  }
162  else
163  {
164  return tmp<surfaceScalarField>(tphi, true);
165  }
166 }
167 
168 
170 (
171  const tmp<surfaceScalarField>& tphi,
172  const volScalarField& rho,
173  const volVectorField& U
174 )
175 {
176  if (tphi().mesh().moving())
177  {
178  return tphi - fvc::interpolate(rho)*fvc::meshPhi(rho, U);
179  }
180  else
181  {
182  return tmp<surfaceScalarField>(tphi, true);
183  }
184 }
185 
186 
188 (
189  const tmp<surfaceScalarField>& tphi,
190  const volVectorField& U
191 )
192 {
193  if (tphi().mesh().moving())
194  {
195  const word phiName(tphi().name());
196 
198  (
199  phiName,
200  tphi + fvc::meshPhi(U)
201  );
202  }
203  else
204  {
205  return tmp<surfaceScalarField>(tphi, true);
206  }
207 }
208 
209 
211 (
212  const tmp<surfaceScalarField>& tphi,
213  const volScalarField& rho,
214  const volVectorField& U
215 )
216 {
217  if (tphi().mesh().moving())
218  {
219  const word phiName(tphi().name());
220 
222  (
223  phiName,
224  tphi + fvc::interpolate(rho)*fvc::meshPhi(rho, U)
225  );
226  }
227  else
228  {
229  return tmp<surfaceScalarField>(tphi, true);
230  }
231 }
232 
233 
235 (
236  const tmp<surfaceScalarField>& tphi,
237  const volScalarField& alpha,
238  const volScalarField& rho,
239  const volVectorField& U
240 )
241 {
242  if (tphi().mesh().moving())
243  {
244  const word phiName(tphi().name());
245 
247  (
248  phiName,
249  tphi
250  + fvc::interpolate(alpha)*fvc::interpolate(rho)*fvc::meshPhi(rho, U)
251  );
252  }
253  else
254  {
255  return tmp<surfaceScalarField>(tphi, true);
256  }
257 }
258 
259 
261 (
263  const volVectorField& U,
264  const surfaceScalarField& phi
265 )
266 {
267  const fvMesh& mesh = U.mesh();
268 
269  if (mesh.dynamic())
270  {
271  Uf() = fvc::interpolate(U);
272  surfaceVectorField n(mesh.Sf()/mesh.magSf());
273  Uf() += n*(phi/mesh.magSf() - (n & Uf()));
274  }
275 }
276 
277 
279 (
281  const volScalarField& rho,
282  const volVectorField& U,
283  const surfaceScalarField& phi
284 )
285 {
286  const fvMesh& mesh = U.mesh();
287 
288  if (mesh.dynamic())
289  {
290  rhoUf() = fvc::interpolate(rho*U);
291  surfaceVectorField n(mesh.Sf()/mesh.magSf());
292  rhoUf() += n*(fvc::absolute(phi, rho, U)/mesh.magSf() - (n & rhoUf()));
293  }
294 }
295 
296 
297 // ************************************************************************* //
const surfaceVectorField & Sf() const
Return cell face area vectors.
const word & name() const
Return name.
Definition: IOobject.H:303
dynamicFvMesh & mesh
autoPtr< BasicCompressibleMomentumTransportModel > New(const volScalarField &rho, const volVectorField &U, const surfaceScalarField &phi, const typename BasicCompressibleMomentumTransportModel::transportModel &transport)
A class for handling words, derived from string.
Definition: word.H:59
void correctUf(autoPtr< surfaceVectorField > &Uf, const volVectorField &U, const surfaceScalarField &phi)
Definition: fvcMeshPhi.C:261
autoPtr< surfaceVectorField > Uf
const surfaceScalarField & magSf() const
Return cell face area magnitudes.
const Mesh & mesh() const
Return mesh.
const word & name() const
Return const reference to name.
word name(const complex &)
Return a string representation of a complex.
Definition: complex.C:47
static tmp< GeometricField< Type, fvsPatchField, surfaceMesh > > interpolate(const GeometricField< Type, fvPatchField, volMesh > &tvf, const surfaceScalarField &faceFlux, Istream &schemeData)
Interpolate field onto faces using scheme given by Istream.
Calculate the mesh motion flux and convert fluxes from absolute to relative and back.
tmp< surfaceScalarField > absolute(const tmp< surfaceScalarField > &tphi, const volVectorField &U)
Return the given relative flux in absolute form.
Definition: fvcMeshPhi.C:188
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
tmp< surfaceScalarField > meshPhi(const volVectorField &U)
Definition: fvcMeshPhi.C:34
virtual bool dynamic() const
Is mesh dynamic.
Definition: polyMesh.H:506
void correctRhoUf(autoPtr< surfaceVectorField > &rhoUf, const volScalarField &rho, const volVectorField &U, const surfaceScalarField &phi)
Definition: fvcMeshPhi.C:279
label n
An auto-pointer similar to the STL auto_ptr but with automatic casting to a reference to the type and...
Definition: PtrList.H:52
A class for managing temporary objects.
Definition: PtrList.H:53
void makeRelative(surfaceScalarField &phi, const volVectorField &U)
Make the given flux relative.
Definition: fvcMeshPhi.C:75
autoPtr< surfaceVectorField > rhoUf
tmp< surfaceScalarField > relative(const tmp< surfaceScalarField > &tphi, const volVectorField &U)
Return the given absolute flux in relative form.
Definition: fvcMeshPhi.C:153