GeometricScalarField.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-2018 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 "GeometricScalarField.H"
27 
28 #define TEMPLATE template<template<class> class PatchField, class GeoMesh>
30 
31 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
32 
33 namespace Foam
34 {
35 
36 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
37 
38 template<template<class> class PatchField, class GeoMesh>
39 void stabilise
40 (
43  const dimensioned<scalar>& ds
44 )
45 {
46  stabilise(result.primitiveFieldRef(), gsf.primitiveField(), ds.value());
47  stabilise(result.boundaryFieldRef(), gsf.boundaryField(), ds.value());
48 }
49 
50 
51 template<template<class> class PatchField, class GeoMesh>
53 (
55  const dimensioned<scalar>& ds
56 )
57 {
59  (
61  (
62  "stabilise(" + gsf.name() + ',' + ds.name() + ')',
63  gsf.mesh(),
64  ds.dimensions() + gsf.dimensions()
65  )
66  );
67 
68  stabilise(tRes.ref(), gsf, ds);
69 
70  return tRes;
71 }
72 
73 
74 template<template<class> class PatchField, class GeoMesh>
76 (
78  const dimensioned<scalar>& ds
79 )
80 {
82 
84  (
85  New
86  (
87  tgsf,
88  "stabilise(" + gsf.name() + ',' + ds.name() + ')',
89  ds.dimensions() + gsf.dimensions()
90  )
91  );
92 
93  stabilise(tRes.ref(), gsf, ds);
94 
95  tgsf.clear();
96 
97  return tRes;
98 }
99 
100 
101 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
102 
103 BINARY_TYPE_OPERATOR(scalar, scalar, scalar, +, '+', add)
104 BINARY_TYPE_OPERATOR(scalar, scalar, scalar, -, '-', subtract)
105 
106 BINARY_OPERATOR(scalar, scalar, scalar, *, '*', multiply)
107 BINARY_OPERATOR(scalar, scalar, scalar, /, '|', divide)
108 
109 BINARY_TYPE_OPERATOR_SF(scalar, scalar, scalar, /, '|', divide)
110 
111 
112 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
113 
114 template<template<class> class PatchField, class GeoMesh>
115 void pow
116 (
120 )
121 {
122  pow(Pow.primitiveFieldRef(), gsf1.primitiveField(), gsf2.primitiveField());
123  pow(Pow.boundaryFieldRef(), gsf1.boundaryField(), gsf2.boundaryField());
124 }
125 
126 
127 template<template<class> class PatchField, class GeoMesh>
129 (
132 )
133 {
134  if (!gsf1.dimensions().dimensionless())
135  {
137  << "Base field is not dimensionless: " << gsf1.dimensions()
138  << exit(FatalError);
139  }
140 
141  if (!gsf2.dimensions().dimensionless())
142  {
144  << "Exponent field is not dimensionless: " << gsf2.dimensions()
145  << exit(FatalError);
146  }
147 
149  (
151  (
152  "pow(" + gsf1.name() + ',' + gsf2.name() + ')',
153  gsf1.mesh(),
154  dimless
155  )
156  );
157 
158  pow(tPow.ref(), gsf1, gsf2);
159 
160  return tPow;
161 }
162 
163 
164 template<template<class> class PatchField, class GeoMesh>
166 (
169 )
170 {
171  const GeometricField<scalar, PatchField, GeoMesh>& gsf1 = tgsf1();
172 
173  if (!gsf1.dimensions().dimensionless())
174  {
176  << "Base field is not dimensionless: " << gsf1.dimensions()
177  << exit(FatalError);
178  }
179 
180  if (!gsf2.dimensions().dimensionless())
181  {
183  << "Exponent field is not dimensionless: " << gsf2.dimensions()
184  << exit(FatalError);
185  }
186 
188  (
189  New
190  (
191  tgsf1,
192  "pow(" + gsf1.name() + ',' + gsf2.name() + ')',
193  dimless
194  )
195  );
196 
197  pow(tPow.ref(), gsf1, gsf2);
198 
199  tgsf1.clear();
200 
201  return tPow;
202 }
203 
204 
205 template<template<class> class PatchField, class GeoMesh>
207 (
210 )
211 {
212  const GeometricField<scalar, PatchField, GeoMesh>& gsf2 = tgsf2();
213 
214  if (!gsf1.dimensions().dimensionless())
215  {
217  << "Base field is not dimensionless: " << gsf1.dimensions()
218  << exit(FatalError);
219  }
220 
221  if (!gsf2.dimensions().dimensionless())
222  {
224  << "Exponent field is not dimensionless: " << gsf2.dimensions()
225  << exit(FatalError);
226  }
227 
229  (
230  New
231  (
232  tgsf2,
233  "pow(" + gsf1.name() + ',' + gsf2.name() + ')',
234  dimless
235  )
236  );
237 
238  pow(tPow.ref(), gsf1, gsf2);
239 
240  tgsf2.clear();
241 
242  return tPow;
243 }
244 
245 template<template<class> class PatchField, class GeoMesh>
247 (
250 )
251 {
252  const GeometricField<scalar, PatchField, GeoMesh>& gsf1 = tgsf1();
253  const GeometricField<scalar, PatchField, GeoMesh>& gsf2 = tgsf2();
254 
255  if (!gsf1.dimensions().dimensionless())
256  {
258  << "Base field is not dimensionless: " << gsf1.dimensions()
259  << exit(FatalError);
260  }
261 
262  if (!gsf2.dimensions().dimensionless())
263  {
265  << "Exponent field is not dimensionless: " << gsf2.dimensions()
266  << exit(FatalError);
267  }
268 
270  (
272  <scalar, scalar, scalar, scalar, PatchField, GeoMesh>::New
273  (
274  tgsf1,
275  tgsf2,
276  "pow(" + gsf1.name() + ',' + gsf2.name() + ')',
277  dimless
278  )
279  );
280 
281  pow(tPow.ref(), gsf1, gsf2);
282 
283  tgsf1.clear();
284  tgsf2.clear();
285 
286  return tPow;
287 }
288 
289 
290 template<template<class> class PatchField, class GeoMesh>
291 void pow
292 (
295  const dimensioned<scalar>& ds
296 )
297 {
298  pow(tPow.primitiveFieldRef(), gsf.primitiveField(), ds.value());
299  pow(tPow.boundaryFieldRef(), gsf.boundaryField(), ds.value());
300 }
301 
302 
303 template<template<class> class PatchField, class GeoMesh>
305 (
307  const dimensionedScalar& ds
308 )
309 {
310  if (!ds.dimensions().dimensionless())
311  {
313  << "Exponent is not dimensionless: " << ds.dimensions()
314  << exit(FatalError);
315  }
316 
318  (
320  (
321  "pow(" + gsf.name() + ',' + ds.name() + ')',
322  gsf.mesh(),
323  pow(gsf.dimensions(), ds)
324  )
325  );
326 
327  pow(tPow.ref(), gsf, ds);
328 
329  return tPow;
330 }
331 
332 template<template<class> class PatchField, class GeoMesh>
334 (
336  const dimensionedScalar& ds
337 )
338 {
339  if (!ds.dimensions().dimensionless())
340  {
342  << "Exponent is not dimensionless: " << ds.dimensions()
343  << exit(FatalError);
344  }
345 
347 
349  (
350  New
351  (
352  tgsf,
353  "pow(" + gsf.name() + ',' + ds.name() + ')',
354  pow(gsf.dimensions(), ds)
355  )
356  );
357 
358  pow(tPow.ref(), gsf, ds);
359 
360  tgsf.clear();
361 
362  return tPow;
363 }
364 
365 template<template<class> class PatchField, class GeoMesh>
367 (
369  const scalar& s
370 )
371 {
372  return pow(gsf, dimensionedScalar(s));
373 }
374 
375 template<template<class> class PatchField, class GeoMesh>
377 (
379  const scalar& s
380 )
381 {
382  return pow(tgsf, dimensionedScalar(s));
383 }
384 
385 
386 template<template<class> class PatchField, class GeoMesh>
387 void pow
388 (
390  const dimensioned<scalar>& ds,
392 )
393 {
394  pow(tPow.primitiveFieldRef(), ds.value(), gsf.primitiveField());
395  pow(tPow.boundaryFieldRef(), ds.value(), gsf.boundaryField());
396 }
397 
398 
399 template<template<class> class PatchField, class GeoMesh>
401 (
402  const dimensionedScalar& ds,
404 )
405 {
406  if (!ds.dimensions().dimensionless())
407  {
409  << "Base scalar is not dimensionless: " << ds.dimensions()
410  << exit(FatalError);
411  }
412 
413  if (!gsf.dimensions().dimensionless())
414  {
416  << "Exponent field is not dimensionless: " << gsf.dimensions()
417  << exit(FatalError);
418  }
419 
421  (
423  (
424  "pow(" + ds.name() + ',' + gsf.name() + ')',
425  gsf.mesh(),
426  dimless
427  )
428  );
429 
430  pow(tPow.ref(), ds, gsf);
431 
432  return tPow;
433 }
434 
435 
436 template<template<class> class PatchField, class GeoMesh>
438 (
439  const dimensionedScalar& ds,
441 )
442 {
444 
445  if (!ds.dimensions().dimensionless())
446  {
448  << "Base scalar is not dimensionless: " << ds.dimensions()
449  << exit(FatalError);
450  }
451 
452  if (!gsf.dimensions().dimensionless())
453  {
455  << "Exponent field is not dimensionless: " << gsf.dimensions()
456  << exit(FatalError);
457  }
458 
460  (
461  New
462  (
463  tgsf,
464  "pow(" + ds.name() + ',' + gsf.name() + ')',
465  dimless
466  )
467  );
468 
469  pow(tPow.ref(), ds, gsf);
470 
471  tgsf.clear();
472 
473  return tPow;
474 }
475 
476 template<template<class> class PatchField, class GeoMesh>
478 (
479  const scalar& s,
481 )
482 {
483  return pow(dimensionedScalar(s), gsf);
484 }
485 
486 template<template<class> class PatchField, class GeoMesh>
488 (
489  const scalar& s,
491 )
492 {
493  return pow(dimensionedScalar(s), tgsf);
494 }
495 
496 
497 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
498 
499 template<template<class> class PatchField, class GeoMesh>
500 void atan2
501 (
505 )
506 {
507  atan2
508  (
509  Atan2.primitiveFieldRef(),
510  gsf1.primitiveField(),
511  gsf2.primitiveField()
512  );
513  atan2
514  (
515  Atan2.boundaryFieldRef(),
516  gsf1.boundaryField(),
517  gsf2.boundaryField()
518  );
519 }
520 
521 
522 template<template<class> class PatchField, class GeoMesh>
524 (
527 )
528 {
530  (
532  (
533  "atan2(" + gsf1.name() + ',' + gsf2.name() + ')',
534  gsf1.mesh(),
535  atan2(gsf1.dimensions(), gsf2.dimensions())
536  )
537  );
538 
539  atan2(tAtan2.ref(), gsf1, gsf2);
540 
541  return tAtan2;
542 }
543 
544 
545 template<template<class> class PatchField, class GeoMesh>
547 (
550 )
551 {
552  const GeometricField<scalar, PatchField, GeoMesh>& gsf1 = tgsf1();
553 
555  (
556  New
557  (
558  tgsf1,
559  "atan2(" + gsf1.name() + ',' + gsf2.name() + ')',
560  atan2(gsf1.dimensions(), gsf2.dimensions())
561  )
562  );
563 
564  atan2(tAtan2.ref(), gsf1, gsf2);
565 
566  tgsf1.clear();
567 
568  return tAtan2;
569 }
570 
571 
572 template<template<class> class PatchField, class GeoMesh>
574 (
577 )
578 {
579  const GeometricField<scalar, PatchField, GeoMesh>& gsf2 = tgsf2();
580 
582  (
583  New
584  (
585  tgsf2,
586  "atan2(" + gsf1.name() + ',' + gsf2.name() + ')',
587  atan2( gsf1.dimensions(), gsf2.dimensions())
588  )
589  );
590 
591  atan2(tAtan2.ref(), gsf1, gsf2);
592 
593  tgsf2.clear();
594 
595  return tAtan2;
596 }
597 
598 template<template<class> class PatchField, class GeoMesh>
600 (
603 )
604 {
605  const GeometricField<scalar, PatchField, GeoMesh>& gsf1 = tgsf1();
606  const GeometricField<scalar, PatchField, GeoMesh>& gsf2 = tgsf2();
607 
609  (
611  <scalar, scalar, scalar, scalar, PatchField, GeoMesh>::New
612  (
613  tgsf1,
614  tgsf2,
615  "atan2(" + gsf1.name() + ',' + gsf2.name() + ')',
616  atan2(gsf1.dimensions(), gsf2.dimensions())
617  )
618  );
619 
620  atan2(tAtan2.ref(), gsf1, gsf2);
621 
622  tgsf1.clear();
623  tgsf2.clear();
624 
625  return tAtan2;
626 }
627 
628 
629 template<template<class> class PatchField, class GeoMesh>
630 void atan2
631 (
634  const dimensioned<scalar>& ds
635 )
636 {
637  atan2(tAtan2.primitiveFieldRef(), gsf.primitiveField(), ds.value());
638  atan2(tAtan2.boundaryFieldRef(), gsf.boundaryField(), ds.value());
639 }
640 
641 
642 template<template<class> class PatchField, class GeoMesh>
644 (
646  const dimensionedScalar& ds
647 )
648 {
650  (
652  (
653  "atan2(" + gsf.name() + ',' + ds.name() + ')',
654  gsf.mesh(),
655  atan2(gsf.dimensions(), ds)
656  )
657  );
658 
659  atan2(tAtan2.ref(), gsf, ds);
660 
661  return tAtan2;
662 }
663 
664 template<template<class> class PatchField, class GeoMesh>
666 (
668  const dimensionedScalar& ds
669 )
670 {
672 
674  (
675  New
676  (
677  tgsf,
678  "atan2(" + gsf.name() + ',' + ds.name() + ')',
679  atan2(gsf.dimensions(), ds)
680  )
681  );
682 
683  atan2(tAtan2.ref(), gsf, ds);
684 
685  tgsf.clear();
686 
687  return tAtan2;
688 }
689 
690 template<template<class> class PatchField, class GeoMesh>
692 (
694  const scalar& s
695 )
696 {
697  return atan2(gsf, dimensionedScalar(s));
698 }
699 
700 template<template<class> class PatchField, class GeoMesh>
702 (
704  const scalar& s
705 )
706 {
707  return atan2(tgsf, dimensionedScalar(s));
708 }
709 
710 
711 template<template<class> class PatchField, class GeoMesh>
712 void atan2
713 (
715  const dimensioned<scalar>& ds,
717 )
718 {
719  atan2(tAtan2.primitiveFieldRef(), ds.value(), gsf.primitiveField());
720  atan2(tAtan2.boundaryFieldRef(), ds.value(), gsf.boundaryField());
721 }
722 
723 
724 template<template<class> class PatchField, class GeoMesh>
726 (
727  const dimensionedScalar& ds,
729 )
730 {
732  (
734  (
735  "atan2(" + ds.name() + ',' + gsf.name() + ')',
736  gsf.mesh(),
737  atan2(ds, gsf.dimensions())
738  )
739  );
740 
741  atan2(tAtan2.ref(), ds, gsf);
742 
743  return tAtan2;
744 }
745 
746 
747 template<template<class> class PatchField, class GeoMesh>
749 (
750  const dimensionedScalar& ds,
752 )
753 {
755 
757  (
758  New
759  (
760  tgsf,
761  "atan2(" + ds.name() + ',' + gsf.name() + ')',
762  atan2(ds, gsf.dimensions())
763  )
764  );
765 
766  atan2(tAtan2.ref(), ds, gsf);
767 
768  tgsf.clear();
769 
770  return tAtan2;
771 }
772 
773 template<template<class> class PatchField, class GeoMesh>
775 (
776  const scalar& s,
778 )
779 {
780  return atan2(dimensionedScalar(s), gsf);
781 }
782 
783 template<template<class> class PatchField, class GeoMesh>
785 (
786  const scalar& s,
788 )
789 {
790  return atan2(dimensionedScalar(s), tgsf);
791 }
792 
793 
794 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
795 
796 UNARY_FUNCTION(scalar, scalar, pow3, pow3)
797 UNARY_FUNCTION(scalar, scalar, pow4, pow4)
798 UNARY_FUNCTION(scalar, scalar, pow5, pow5)
799 UNARY_FUNCTION(scalar, scalar, pow6, pow6)
800 UNARY_FUNCTION(scalar, scalar, pow025, pow025)
801 UNARY_FUNCTION(scalar, scalar, sqrt, sqrt)
802 UNARY_FUNCTION(scalar, scalar, cbrt, cbrt)
803 UNARY_FUNCTION(scalar, scalar, sign, sign)
804 UNARY_FUNCTION(scalar, scalar, pos, pos)
805 UNARY_FUNCTION(scalar, scalar, pos0, pos0)
806 UNARY_FUNCTION(scalar, scalar, neg, neg)
807 UNARY_FUNCTION(scalar, scalar, neg0, neg0)
808 UNARY_FUNCTION(scalar, scalar, posPart, posPart)
809 UNARY_FUNCTION(scalar, scalar, negPart, negPart)
810 
811 UNARY_FUNCTION(scalar, scalar, exp, trans)
812 UNARY_FUNCTION(scalar, scalar, log, trans)
813 UNARY_FUNCTION(scalar, scalar, log10, trans)
814 UNARY_FUNCTION(scalar, scalar, sin, trans)
815 UNARY_FUNCTION(scalar, scalar, cos, trans)
816 UNARY_FUNCTION(scalar, scalar, tan, trans)
817 UNARY_FUNCTION(scalar, scalar, asin, trans)
818 UNARY_FUNCTION(scalar, scalar, acos, trans)
819 UNARY_FUNCTION(scalar, scalar, atan, trans)
820 UNARY_FUNCTION(scalar, scalar, sinh, trans)
821 UNARY_FUNCTION(scalar, scalar, cosh, trans)
822 UNARY_FUNCTION(scalar, scalar, tanh, trans)
823 UNARY_FUNCTION(scalar, scalar, asinh, trans)
824 UNARY_FUNCTION(scalar, scalar, acosh, trans)
825 UNARY_FUNCTION(scalar, scalar, atanh, trans)
826 UNARY_FUNCTION(scalar, scalar, erf, trans)
827 UNARY_FUNCTION(scalar, scalar, erfc, trans)
828 UNARY_FUNCTION(scalar, scalar, lgamma, trans)
829 UNARY_FUNCTION(scalar, scalar, j0, trans)
830 UNARY_FUNCTION(scalar, scalar, j1, trans)
831 UNARY_FUNCTION(scalar, scalar, y0, trans)
832 UNARY_FUNCTION(scalar, scalar, y1, trans)
833 
834 
835 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
836 
837 #define BesselFunc(func) \
838  \
839 template<template<class> class PatchField, class GeoMesh> \
840 void func \
841 ( \
842  GeometricField<scalar, PatchField, GeoMesh>& gsf, \
843  const int n, \
844  const GeometricField<scalar, PatchField, GeoMesh>& gsf1 \
845 ) \
846 { \
847  func(gsf.primitiveFieldRef(), n, gsf1.primitiveField()); \
848  func(gsf.boundaryFieldRef(), n, gsf1.boundaryField()); \
849 } \
850  \
851 template<template<class> class PatchField, class GeoMesh> \
852 tmp<GeometricField<scalar, PatchField, GeoMesh>> func \
853 ( \
854  const int n, \
855  const GeometricField<scalar, PatchField, GeoMesh>& gsf \
856 ) \
857 { \
858  if (!gsf.dimensions().dimensionless()) \
859  { \
860  FatalErrorInFunction \
861  << "gsf not dimensionless" \
862  << abort(FatalError); \
863  } \
864  \
865  tmp<GeometricField<scalar, PatchField, GeoMesh>> tFunc \
866  ( \
867  GeometricField<scalar, PatchField, GeoMesh>::New \
868  ( \
869  #func "(" + gsf.name() + ')', \
870  gsf.mesh(), \
871  dimless \
872  ) \
873  ); \
874  \
875  func(tFunc.ref(), n, gsf); \
876  \
877  return tFunc; \
878 } \
879  \
880 template<template<class> class PatchField, class GeoMesh> \
881 tmp<GeometricField<scalar, PatchField, GeoMesh>> func \
882 ( \
883  const int n, \
884  const tmp<GeometricField<scalar, PatchField, GeoMesh>>& tgsf \
885 ) \
886 { \
887  const GeometricField<scalar, PatchField, GeoMesh>& gsf = tgsf(); \
888  \
889  if (!gsf.dimensions().dimensionless()) \
890  { \
891  FatalErrorInFunction \
892  << " : gsf not dimensionless" \
893  << abort(FatalError); \
894  } \
895  \
896  tmp<GeometricField<scalar, PatchField, GeoMesh>> tFunc \
897  ( \
898  New \
899  ( \
900  tgsf, \
901  #func "(" + gsf.name() + ')', \
902  dimless \
903  ) \
904  ); \
905  \
906  func(tFunc.ref(), n, gsf); \
907  \
908  tgsf.clear(); \
909  \
910  return tFunc; \
911 }
912 
915 
916 #undef BesselFunc
917 
918 
919 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
920 
921 } // End namespace Foam
922 
923 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
924 
925 #include "undefFieldFunctionsM.H"
926 
927 // ************************************************************************* //
void divide(FieldField< Field, Type > &f, const FieldField< Field, Type > &f1, const FieldField< Field, scalar > &f2)
dimensionedScalar sign(const dimensionedScalar &ds)
dimensionedScalar tanh(const dimensionedScalar &ds)
dimensionSet trans(const dimensionSet &)
Definition: dimensionSet.C:450
dimensionedScalar acos(const dimensionedScalar &ds)
#define BINARY_OPERATOR(ReturnType, Type1, Type2, Op, OpName, OpFunc)
void multiply(FieldField< Field, Type > &f, const FieldField< Field, Type > &f1, const FieldField< Field, scalar > &f2)
const word & name() const
Return name.
Definition: IOobject.H:295
dimensionedScalar log(const dimensionedScalar &ds)
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:124
error FatalError
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:319
#define UNARY_FUNCTION(ReturnType, Type1, Func, Dfunc)
const Boundary & boundaryField() const
Return const-reference to the boundary field.
T & ref() const
Return non-const reference or generate a fatal error.
Definition: tmpI.H:174
dimensionedScalar yn(const int n, const dimensionedScalar &ds)
#define BINARY_TYPE_OPERATOR_SF(ReturnType, Type1, Type2, Op, OpName, OpFunc)
dimensionedScalar sqrt(const dimensionedScalar &ds)
dimensionedScalar pow025(const dimensionedScalar &ds)
dimensionedScalar y0(const dimensionedScalar &ds)
const Internal::FieldType & primitiveField() const
Return a const-reference to the internal field.
Generic GeometricField class.
void subtract(FieldField< Field1, typename typeOfSum< Type1, Type2 >::type > &f, const FieldField< Field1, Type1 > &f1, const FieldField< Field2, Type2 > &f2)
dimensionedScalar posPart(const dimensionedScalar &ds)
dimensionedScalar neg(const dimensionedScalar &ds)
dimensionedScalar asin(const dimensionedScalar &ds)
dimensionedScalar pow5(const dimensionedScalar &ds)
dimensionedScalar j1(const dimensionedScalar &ds)
tmp< DimensionedField< TypeR, GeoMesh > > New(const tmp< DimensionedField< TypeR, GeoMesh >> &tdf1, const word &name, const dimensionSet &dimensions)
dimensionedScalar pos(const dimensionedScalar &ds)
const dimensionSet & dimensions() const
Return dimensions.
Scalar specific part of the implementation of GeometricField.
dimensionedScalar acosh(const dimensionedScalar &ds)
dimensionedScalar cos(const dimensionedScalar &ds)
#define BesselFunc(func)
dimensionedScalar exp(const dimensionedScalar &ds)
dimensionedScalar asinh(const dimensionedScalar &ds)
dimensionedScalar jn(const int n, const dimensionedScalar &ds)
dimensionedScalar cbrt(const dimensionedScalar &ds)
dimensionedScalar neg0(const dimensionedScalar &ds)
const Type & value() const
Return const reference to value.
dimensionedScalar atanh(const dimensionedScalar &ds)
dimensionedScalar y1(const dimensionedScalar &ds)
dimensionedScalar pos0(const dimensionedScalar &ds)
dimensionedScalar sin(const dimensionedScalar &ds)
Internal::FieldType & primitiveFieldRef()
Return a reference to the internal field.
void add(FieldField< Field1, typename typeOfSum< Type1, Type2 >::type > &f, const FieldField< Field1, Type1 > &f1, const FieldField< Field2, Type2 > &f2)
dimensionedScalar erf(const dimensionedScalar &ds)
dimensionedScalar atan2(const dimensionedScalar &x, const dimensionedScalar &y)
const Mesh & mesh() const
Return mesh.
const word & name() const
Return const reference to name.
Internal & ref()
Return a reference to the dimensioned internal field.
dimensionedScalar lgamma(const dimensionedScalar &ds)
dimensionedScalar pow(const dimensionedScalar &ds, const dimensionedScalar &expt)
dimensionedScalar pow3(const dimensionedScalar &ds)
tmp< DimensionedField< scalar, GeoMesh > > stabilise(const DimensionedField< scalar, GeoMesh > &dsf, const dimensioned< scalar > &ds)
dimensionedScalar erfc(const dimensionedScalar &ds)
const dimensionSet dimless(0, 0, 0, 0, 0, 0, 0)
Definition: dimensionSets.H:47
dimensionedScalar sinh(const dimensionedScalar &ds)
Boundary & boundaryFieldRef()
Return a reference to the boundary field.
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
dimensionedScalar atan(const dimensionedScalar &ds)
dimensionedScalar pow4(const dimensionedScalar &ds)
const dimensionSet & dimensions() const
Return const reference to dimensions.
dimensionedScalar pow6(const dimensionedScalar &ds)
dimensionedScalar cosh(const dimensionedScalar &ds)
Generic mesh wrapper used by volMesh, surfaceMesh, pointMesh etc.
Definition: GeoMesh.H:46
A class for managing temporary objects.
Definition: PtrList.H:53
dimensionedScalar tan(const dimensionedScalar &ds)
dimensionedScalar j0(const dimensionedScalar &ds)
#define BINARY_TYPE_OPERATOR(ReturnType, Type1, Type2, Op, OpName, OpFunc)
bool dimensionless() const
Return true if it is dimensionless.
Definition: dimensionSet.C:89
dimensionedScalar log10(const dimensionedScalar &ds)
Namespace for OpenFOAM.
dimensionedScalar negPart(const dimensionedScalar &ds)