fvmLaplacian.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-2022 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 "volFields.H"
27 #include "surfaceFields.H"
28 #include "fvMatrix.H"
29 #include "laplacianScheme.H"
30 #include "gaussLaplacianScheme.H"
31 
32 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
33 
34 namespace Foam
35 {
36 
37 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
38 
39 namespace fvm
40 {
41 
42 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
43 
44 template<class Type>
45 tmp<fvMatrix<Type>>
47 (
48  const VolField<Type>& vf,
49  const word& name
50 )
51 {
52  surfaceScalarField Gamma
53  (
54  IOobject
55  (
56  "1",
57  vf.time().constant(),
58  vf.mesh(),
60  ),
61  vf.mesh(),
63  );
64 
65  return fvm::laplacian(Gamma, vf, name);
66 }
67 
68 
69 template<class Type>
72 (
73  const VolField<Type>& vf
74 )
75 {
76  surfaceScalarField Gamma
77  (
78  IOobject
79  (
80  "1",
81  vf.time().constant(),
82  vf.mesh(),
84  ),
85  vf.mesh(),
87  );
88 
89  return fvm::laplacian
90  (
91  Gamma,
92  vf,
93  "laplacian(" + vf.name() + ')'
94  );
95 }
96 
97 
98 template<class Type>
101 (
102  const zero&,
103  const VolField<Type>& vf,
104  const word& name
105 )
106 {
107  return tmp<fvMatrix<Type>>
108  (
109  new fvMatrix<Type>(vf, dimensionSet(0, 0, -2, 0, 0))
110  );
111 }
112 
113 
114 template<class Type>
117 (
118  const zero&,
119  const VolField<Type>& vf
120 )
121 {
122  return tmp<fvMatrix<Type>>
123  (
124  new fvMatrix<Type>(vf, dimensionSet(0, 0, -2, 0, 0))
125  );
126 }
127 
128 
129 template<class Type>
132 (
133  const one&,
134  const VolField<Type>& vf,
135  const word& name
136 )
137 {
138  return fvm::laplacian(vf, name);
139 }
140 
141 
142 template<class Type>
145 (
146  const one&,
147  const VolField<Type>& vf
148 )
149 {
150  return fvm::laplacian(vf);
151 }
152 
153 
154 template<class Type, class GType>
157 (
158  const dimensioned<GType>& gamma,
159  const VolField<Type>& vf,
160  const word& name
161 )
162 {
163  const SurfaceField<GType> Gamma
164  (
165  IOobject
166  (
167  gamma.name(),
168  vf.instance(),
169  vf.mesh(),
171  ),
172  vf.mesh(),
173  gamma
174  );
175 
176  return fvm::laplacian(Gamma, vf, name);
177 }
178 
179 
180 template<class Type, class GType>
183 (
184  const dimensioned<GType>& gamma,
185  const VolField<Type>& vf
186 )
187 {
188  const SurfaceField<GType> Gamma
189  (
190  IOobject
191  (
192  gamma.name(),
193  vf.instance(),
194  vf.mesh(),
196  ),
197  vf.mesh(),
198  gamma
199  );
200 
201  return fvm::laplacian(Gamma, vf);
202 }
203 
204 
205 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
206 
207 template<class Type, class GType>
210 (
211  const VolField<GType>& gamma,
212  const VolField<Type>& vf,
213  const word& name
214 )
215 {
217  (
218  vf.mesh(),
219  vf.mesh().schemes().laplacian(name)
220  ).ref().fvmLaplacian(gamma, vf);
221 }
222 
223 
224 template<class Type, class GType>
227 (
228  const tmp<VolField<GType>>& tgamma,
229  const VolField<Type>& vf,
230  const word& name
231 )
232 {
233  tmp<fvMatrix<Type>> tLaplacian(fvm::laplacian(tgamma(), vf, name));
234  tgamma.clear();
235  return tLaplacian;
236 }
237 
238 
239 template<class Type, class GType>
242 (
243  const VolField<GType>& gamma,
244  const VolField<Type>& vf
245 )
246 {
247  return fvm::laplacian
248  (
249  gamma,
250  vf,
251  "laplacian(" + gamma.name() + ',' + vf.name() + ')'
252  );
253 }
254 
255 
256 template<class Type, class GType>
259 (
260  const tmp<VolField<GType>>& tgamma,
261  const VolField<Type>& vf
262 )
263 {
264  tmp<fvMatrix<Type>> tLaplacian(fvm::laplacian(tgamma(), vf));
265  tgamma.clear();
266  return tLaplacian;
267 }
268 
269 
270 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
271 
272 template<class Type, class GType>
275 (
276  const SurfaceField<GType>& gamma,
277  const VolField<Type>& vf,
278  const word& name
279 )
280 {
282  (
283  vf.mesh(),
284  vf.mesh().schemes().laplacian(name)
285  ).ref().fvmLaplacian(gamma, vf);
286 }
287 
288 
289 template<class Type, class GType>
292 (
293  const tmp<SurfaceField<GType>>& tgamma,
294  const VolField<Type>& vf,
295  const word& name
296 )
297 {
298  tmp<fvMatrix<Type>> tLaplacian = fvm::laplacian(tgamma(), vf, name);
299  tgamma.clear();
300  return tLaplacian;
301 }
302 
303 
304 template<class Type, class GType>
307 (
308  const SurfaceField<GType>& gamma,
309  const VolField<Type>& vf
310 )
311 {
312  return fvm::laplacian
313  (
314  gamma,
315  vf,
316  "laplacian(" + gamma.name() + ',' + vf.name() + ')'
317  );
318 }
319 
320 
321 template<class Type, class GType>
324 (
325  const tmp<SurfaceField<GType>>& tGamma,
326  const VolField<Type>& vf
327 )
328 {
329  tmp<fvMatrix<Type>> tLaplacian(fvm::laplacian(tGamma(), vf));
330  tGamma.clear();
331  return tLaplacian;
332 }
333 
334 
335 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
336 
337 template<class Type>
340 (
341  const VolField<scalar>& gamma,
342  const VolField<Type>& vf
343 )
344 {
345  return fvm::laplacianCorrection(fvc::interpolate(gamma), vf);
346 }
347 
348 
349 template<class Type>
352 (
353  const tmp<VolField<scalar>>& tgamma,
354  const VolField<Type>& vf
355 )
356 {
357  tmp<fvMatrix<Type>> tLaplacianCorrection
358  (
359  fvm::laplacianCorrection(tgamma(), vf)
360  );
361  tgamma.clear();
362  return tLaplacianCorrection;
363 }
364 
365 
366 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
367 
368 template<class Type>
371 (
372  const SurfaceField<scalar>& gamma,
373  const VolField<Type>& vf
374 )
375 {
376  return correction
377  (
379  (
380  gamma*vf.mesh().magSf(),
381  vf.mesh().nonOrthDeltaCoeffs(),
382  vf
383  )
384  );
385 }
386 
387 
388 template<class Type>
391 (
392  const tmp<SurfaceField<scalar>>& tGamma,
393  const VolField<Type>& vf
394 )
395 {
396  tmp<fvMatrix<Type>> tLaplacianCorrection
397  (
398  fvm::laplacianCorrection(tGamma(), vf)
399  );
400  tGamma.clear();
401  return tLaplacianCorrection;
402 }
403 
404 
405 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
406 
407 } // End namespace fvm
408 
409 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
410 
411 } // End namespace Foam
412 
413 // ************************************************************************* //
const Mesh & mesh() const
Return mesh.
Generic GeometricField class.
IOobject defines the attributes of an object for which implicit objectRegistry management is supporte...
Definition: IOobject.H:99
const Time & time() const
Return time.
Definition: IOobject.C:318
fileName & instance() const
Return the instance directory, constant, system, <time> etc.
Definition: IOobject.C:355
const word & name() const
Return name.
Definition: IOobject.H:310
static const word & constant()
Return constant name.
Definition: TimePaths.H:122
Dimension set for the base types.
Definition: dimensionSet.H:125
const word & name() const
Return const reference to name.
A special matrix type and solver, designed for finite volume solutions of scalar equations....
Definition: fvMatrix.H:118
Basic second-order laplacian using face-gradients and Gauss' theorem.
static tmp< laplacianScheme< Type, GType > > New(const fvMesh &mesh, Istream &schemeData)
Return a pointer to a new laplacianScheme created on freestore.
A class representing the concept of 1 (scalar(1)) used to avoid unnecessary manipulations for objects...
Definition: one.H:51
A class for managing temporary objects.
Definition: tmp.H:55
A class for handling words, derived from string.
Definition: word.H:62
A class representing the concept of 0 used to avoid unnecessary manipulations for objects that are kn...
Definition: zero.H:50
static tmp< SurfaceField< Type > > interpolate(const VolField< Type > &tvf, const surfaceScalarField &faceFlux, Istream &schemeData)
Interpolate field onto faces using scheme given by Istream.
tmp< fvMatrix< Type > > laplacian(const VolField< Type > &vf, const word &name)
Definition: fvmLaplacian.C:47
tmp< fvMatrix< Type > > laplacianCorrection(const VolField< scalar > &gamma, const VolField< Type > &vf)
Definition: fvmLaplacian.C:340
Namespace for OpenFOAM.
word name(const bool)
Return a word representation of a bool.
Definition: boolIO.C:39
const dimensionSet dimless
tmp< fvMatrix< Type > > correction(const fvMatrix< Type > &)
Return the correction form of the given matrix.
Foam::surfaceFields.