facLaplacian.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) 2016-2017 Wikki Ltd
9 -------------------------------------------------------------------------------
10 License
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 "facLaplacian.H"
29 #include "faMesh.H"
30 #include "faLaplacianScheme.H"
31 
32 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
33 
34 namespace Foam
35 {
36 
37 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
38 
39 namespace fac
40 {
41 
42 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
43 
44 template<class Type>
45 tmp<GeometricField<Type, faPatchField, areaMesh>>
47 (
49  const word& name
50 )
51 {
53  (
54  vf.mesh(),
55  vf.mesh().laplacianScheme(name)
56  ).ref().facLaplacian(vf);
57 }
58 
59 template<class Type>
62 (
64  const word& name
65 )
66 {
68  (
69  fac::laplacian(tvf(), name)
70  );
71  tvf.clear();
72  return Laplacian;
73 }
74 
75 
76 template<class Type>
79 (
81 )
82 {
83  return fac::laplacian(vf, "laplacian(" + vf.name() + ')');
84 }
85 
86 template<class Type>
89 (
91 )
92 {
94  (
95  fac::laplacian(tvf())
96  );
97  tvf.clear();
98  return Laplacian;
99 }
100 
101 
102 template<class Type>
104 laplacian
105 (
106  const dimensionedScalar& gamma,
108  const word& name
109 )
110 {
111  return gamma*fac::laplacian(vf, name);
112 }
113 
114 template<class Type>
116 laplacian
117 (
118  const dimensionedScalar& gamma,
120  const word& name
121 )
122 {
124  (
125  fac::laplacian(gamma, tvf(), name)
126  );
127  tvf.clear();
128  return Laplacian;
129 }
130 
131 
132 template<class Type>
134 laplacian
135 (
136  const dimensionedScalar& gamma,
138 )
139 {
140  return gamma*fac::laplacian
141  (
142  vf, "laplacian(" + gamma.name() + ',' + vf.name() + ')'
143  );
144 }
145 
146 template<class Type>
148 laplacian
149 (
150  const dimensionedScalar& gamma,
152 )
153 {
155  (
156  fac::laplacian(gamma, tvf())
157  );
158  tvf.clear();
159  return Laplacian;
160 }
161 
162 
163 template<class Type>
165 laplacian
166 (
167  const areaScalarField& gamma,
169  const word& name
170 )
171 {
173  (
174  vf.mesh(),
175  vf.mesh().laplacianScheme(name)
176  ).ref().facLaplacian(gamma, vf);
177 }
178 
179 template<class Type>
181 laplacian
182 (
183  const tmp<areaScalarField>& tgamma,
185  const word& name
186 )
187 {
189  (
190  fac::laplacian(tgamma(), vf, name)
191  );
192  tgamma.clear();
193  return Laplacian;
194 }
195 
196 template<class Type>
198 laplacian
199 (
200  const areaScalarField& gamma,
202  const word& name
203 )
204 {
206  (
207  fac::laplacian(gamma, tvf(), name)
208  );
209  tvf.clear();
210  return Laplacian;
211 }
212 
213 template<class Type>
215 laplacian
216 (
217  const tmp<areaScalarField>& tgamma,
219  const word& name
220 )
221 {
223  (
224  fac::laplacian(tgamma(), tvf(), name)
225  );
226  tgamma.clear();
227  tvf.clear();
228  return Laplacian;
229 }
230 
231 
232 template<class Type>
234 laplacian
235 (
236  const areaScalarField& gamma,
238 )
239 {
240  return fac::laplacian
241  (
242  gamma,
243  vf,
244  "laplacian(" + gamma.name() + ',' + vf.name() + ')'
245  );
246 }
247 
248 template<class Type>
250 laplacian
251 (
252  const tmp<areaScalarField>& tgamma,
254 )
255 {
256  return fac::laplacian
257  (
258  tgamma,
259  vf,
260  "laplacian(" + tgamma().name() + ',' + vf.name() + ')'
261  );
262 }
263 
264 template<class Type>
266 laplacian
267 (
268  const areaScalarField& gamma,
270 )
271 {
272  return fac::laplacian
273  (
274  gamma,
275  tvf,
276  "laplacian(" + gamma.name() + ',' + tvf().name() + ')'
277  );
278 }
279 
280 template<class Type>
282 laplacian
283 (
284  const tmp<areaScalarField>& tgamma,
286 )
287 {
288  return fac::laplacian
289  (
290  tgamma,
291  tvf,
292  "laplacian(" + tgamma().name() + ',' + tvf().name() + ')'
293  );
294 }
295 
296 
297 template<class Type>
299 laplacian
300 (
301  const edgeScalarField& gamma,
303  const word& name
304 )
305 {
307  (
308  vf.mesh(),
309  vf.mesh().laplacianScheme(name)
310  ).ref().facLaplacian(gamma, vf);
311 }
312 
313 template<class Type>
315 laplacian
316 (
317  const tmp<edgeScalarField>& tgamma,
319  const word& name
320 )
321 {
323  (
324  fac::laplacian(tgamma(), vf, name)
325  );
326  tgamma.clear();
327  return Laplacian;
328 }
329 
330 template<class Type>
332 laplacian
333 (
334  const edgeScalarField& gamma,
336  const word& name
337 )
338 {
340  (
341  fac::laplacian(gamma, tvf(), name)
342  );
343  tvf.clear();
344  return Laplacian;
345 }
346 
347 template<class Type>
349 (
350  const tmp<edgeScalarField>& tgamma,
352  const word& name
353 )
354 {
356  (
357  fac::laplacian(tgamma(), tvf(), name)
358  );
359  tgamma.clear();
360  tvf.clear();
361  return Laplacian;
362 }
363 
364 
365 template<class Type>
367 laplacian
368 (
369  const edgeScalarField& gamma,
371 )
372 {
373  return fac::laplacian
374  (
375  gamma,
376  vf,
377  "laplacian(" + gamma.name() + ',' + vf.name() + ')'
378  );
379 }
380 
381 template<class Type>
383 laplacian
384 (
385  const tmp<edgeScalarField>& tgamma,
387 )
388 {
390  (
391  fac::laplacian(tgamma(), vf)
392  );
393  tgamma.clear();
394  return Laplacian;
395 }
396 
397 template<class Type>
399 laplacian
400 (
401  const edgeScalarField& gamma,
403 )
404 {
406  (
407  fac::laplacian(gamma, tvf())
408  );
409  tvf.clear();
410  return Laplacian;
411 }
412 
413 template<class Type>
415 (
416  const tmp<edgeScalarField>& tgamma,
418 )
419 {
421  (
422  fac::laplacian(tgamma(), tvf())
423  );
424  tgamma.clear();
425  tvf.clear();
426  return Laplacian;
427 }
428 
429 
430 /*
431 template<class Type>
432 tmp<GeometricField<Type, faPatchField, areaMesh>>
433 laplacian
434 (
435  const areaTensorField& gamma,
436  const GeometricField<Type, faPatchField, areaMesh>& vf
437 )
438 {
439  const faMesh& mesh = tvf().mesh();
440  return fac::laplacian
441  (
442  (
443  mesh.Sf() & fac::interpolate(tgamma) & mesh.Sf()
444  )/sqr(mesh.magSf()),
445  tvf
446  );
447 }
448 
449 template<class Type>
450 tmp<GeometricField<Type, faPatchField, areaMesh>>
451 laplacian
452 (
453  const areaTensorField& gamma,
454  const tmp<GeometricField<Type, faPatchField, areaMesh>>& tvf
455 )
456 {
457  tmp<GeometricField<Type, faPatchField, areaMesh>> Laplacian
458  (
459  fac::laplacian(gamma, tvf())
460  );
461  tvf.clear();
462  return Laplacian;
463 }
464 
465 template<class Type>
466 tmp<GeometricField<Type, faPatchField, areaMesh>>
467 laplacian
468 (
469  const tmp<areaTensorField>& tgamma,
470  const GeometricField<Type, faPatchField, areaMesh>& vf
471 )
472 {
473  tmp<GeometricField<Type, faPatchField, areaMesh>> Laplacian
474  (
475  fac::laplacian(tgamma(), vf)
476  );
477  tgamma.clear();
478  return Laplacian;
479 }
480 
481 template<class Type>
482 tmp<GeometricField<Type, faPatchField, areaMesh>>
483 laplacian
484 (
485  const tmp<areaTensorField>& tgamma,
486  const tmp<GeometricField<Type, faPatchField, areaMesh>>& tvf
487 )
488 {
489  tmp<GeometricField<Type, faPatchField, areaMesh>> Laplacian
490  (
491  fac::laplacian(tgamma(), tvf())
492  );
493  tgamma.clear();
494  tvf.clear();
495  return Laplacian;
496 }
497 
498 
499 template<class Type>
500 tmp<GeometricField<Type, faPatchField, areaMesh>> laplacian
501 (
502  const edgeTensorField& gamma,
503  const GeometricField<Type, faPatchField, areaMesh>& vf
504 )
505 {
506  const faMesh& mesh = tvf().mesh();
507 
508  return fac::laplacian
509  (
510  (
511  mesh.Sf() &
512  gamma
513  & mesh.Sf()
514  )/sqr(mesh.magSf()),
515  vf
516  );
517 }
518 
519 template<class Type>
520 tmp<GeometricField<Type, faPatchField, areaMesh>>
521 laplacian
522 (
523  const edgeTensorField& gamma,
524  const tmp<GeometricField<Type, faPatchField, areaMesh>>& tvf
525 )
526 {
527  tmp<GeometricField<Type, faPatchField, areaMesh>> Laplacian
528  (
529  fac::laplacian(gamma, tvf())
530  );
531  tvf.clear();
532  return Laplacian;
533 }
534 
535 template<class Type>
536 tmp<GeometricField<Type, faPatchField, areaMesh>>
537 laplacian
538 (
539  const tmp<edgeTensorField>& tgamma,
540  const GeometricField<Type, faPatchField, areaMesh>& vf
541 )
542 {
543  tmp<GeometricField<Type, faPatchField, areaMesh>> Laplacian
544  (
545  fac::laplacian(tgamma(), tvf())
546  );
547  tgamma.clear();
548  return Laplacian;
549 }
550 
551 template<class Type>
552 tmp<GeometricField<Type, faPatchField, areaMesh>>
553 laplacian
554 (
555  const tmp<edgeTensorField>& tgamma,
556  const tmp<GeometricField<Type, faPatchField, areaMesh>>& tvf
557 )
558 {
559  tmp<GeometricField<Type, faPatchField, areaMesh>> Laplacian
560  (
561  fac::laplacian(tgamma(), tvf())
562  );
563  tgamma.clear();
564  tvf.clear();
565  return Laplacian;
566 }
567 */
568 
569 
570 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
571 
572 } // End namespace fac
573 
574 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
575 
576 } // End namespace Foam
577 
578 // ************************************************************************* //
Foam::word
A class for handling words, derived from Foam::string.
Definition: word.H:65
Foam::tmp::clear
void clear() const noexcept
Definition: tmpI.H:287
Foam::tmp
A class for managing temporary objects.
Definition: PtrList.H:61
faMesh.H
ref
rDeltaT ref()
Foam::fa::laplacianScheme::New
static tmp< laplacianScheme< Type > > New(const faMesh &mesh, Istream &schemeData)
Return a pointer to a new laplacianScheme created on freestore.
Definition: faLaplacianScheme.C:47
facLaplacian.H
Calculate the laplacian of the given field.
faLaplacianScheme.H
Foam::dimensioned< scalar >
Foam
Namespace for OpenFOAM.
Definition: atmBoundaryLayer.C:33
fac
Calculate the second temporal derivative.
gamma
const scalar gamma
Definition: EEqn.H:9
Foam::name
word name(const expressions::valueTypeCode typeCode)
A word representation of a valueTypeCode. Empty for INVALID.
Definition: exprTraits.C:59
Foam::fac::laplacian
tmp< GeometricField< Type, faPatchField, areaMesh > > laplacian(const GeometricField< Type, faPatchField, areaMesh > &vf, const word &name)
Definition: facLaplacian.C:47
Foam::GeometricField
Generic GeometricField class.
Definition: areaFieldsFwd.H:53