GeometricScalarField.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 | Copyright (C) 2011-2017 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  IOobject
63  (
64  "stabilise(" + gsf.name() + ',' + ds.name() + ')',
65  gsf.instance(),
66  gsf.db(),
69  ),
70  gsf.mesh(),
71  ds.dimensions() + gsf.dimensions()
72  )
73  );
74 
75  stabilise(tRes.ref(), gsf, ds);
76 
77  return tRes;
78 }
79 
80 
81 template<template<class> class PatchField, class GeoMesh>
83 (
85  const dimensioned<scalar>& ds
86 )
87 {
89 
91  (
92  New
93  (
94  tgsf,
95  "stabilise(" + gsf.name() + ',' + ds.name() + ')',
96  ds.dimensions() + gsf.dimensions()
97  )
98  );
99 
100  stabilise(tRes.ref(), gsf, ds);
101 
102  tgsf.clear();
103 
104  return tRes;
105 }
106 
107 
108 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
109 
110 BINARY_TYPE_OPERATOR(scalar, scalar, scalar, +, '+', add)
111 BINARY_TYPE_OPERATOR(scalar, scalar, scalar, -, '-', subtract)
112 
113 BINARY_OPERATOR(scalar, scalar, scalar, *, '*', multiply)
114 BINARY_OPERATOR(scalar, scalar, scalar, /, '|', divide)
115 
116 BINARY_TYPE_OPERATOR_SF(scalar, scalar, scalar, /, '|', divide)
117 
118 
119 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
120 
121 template<template<class> class PatchField, class GeoMesh>
122 void pow
123 (
127 )
128 {
129  pow(Pow.primitiveFieldRef(), gsf1.primitiveField(), gsf2.primitiveField());
130  pow(Pow.boundaryFieldRef(), gsf1.boundaryField(), gsf2.boundaryField());
131 }
132 
133 
134 template<template<class> class PatchField, class GeoMesh>
136 (
139 )
140 {
141  if (!gsf1.dimensions().dimensionless())
142  {
144  << "Base field is not dimensionless: " << gsf1.dimensions()
145  << exit(FatalError);
146  }
147 
148  if (!gsf2.dimensions().dimensionless())
149  {
151  << "Exponent field is not dimensionless: " << gsf2.dimensions()
152  << exit(FatalError);
153  }
154 
156  (
158  (
159  IOobject
160  (
161  "pow(" + gsf1.name() + ',' + gsf2.name() + ')',
162  gsf1.instance(),
163  gsf1.db(),
166  ),
167  gsf1.mesh(),
168  dimless
169  )
170  );
171 
172  pow(tPow.ref(), gsf1, gsf2);
173 
174  return tPow;
175 }
176 
177 
178 template<template<class> class PatchField, class GeoMesh>
180 (
183 )
184 {
185  const GeometricField<scalar, PatchField, GeoMesh>& gsf1 = tgsf1();
186 
187  if (!gsf1.dimensions().dimensionless())
188  {
190  << "Base field is not dimensionless: " << gsf1.dimensions()
191  << exit(FatalError);
192  }
193 
194  if (!gsf2.dimensions().dimensionless())
195  {
197  << "Exponent field is not dimensionless: " << gsf2.dimensions()
198  << exit(FatalError);
199  }
200 
202  (
203  New
204  (
205  tgsf1,
206  "pow(" + gsf1.name() + ',' + gsf2.name() + ')',
207  dimless
208  )
209  );
210 
211  pow(tPow.ref(), gsf1, gsf2);
212 
213  tgsf1.clear();
214 
215  return tPow;
216 }
217 
218 
219 template<template<class> class PatchField, class GeoMesh>
221 (
224 )
225 {
226  const GeometricField<scalar, PatchField, GeoMesh>& gsf2 = tgsf2();
227 
228  if (!gsf1.dimensions().dimensionless())
229  {
231  << "Base field is not dimensionless: " << gsf1.dimensions()
232  << exit(FatalError);
233  }
234 
235  if (!gsf2.dimensions().dimensionless())
236  {
238  << "Exponent field is not dimensionless: " << gsf2.dimensions()
239  << exit(FatalError);
240  }
241 
243  (
244  New
245  (
246  tgsf2,
247  "pow(" + gsf1.name() + ',' + gsf2.name() + ')',
248  dimless
249  )
250  );
251 
252  pow(tPow.ref(), gsf1, gsf2);
253 
254  tgsf2.clear();
255 
256  return tPow;
257 }
258 
259 template<template<class> class PatchField, class GeoMesh>
261 (
264 )
265 {
266  const GeometricField<scalar, PatchField, GeoMesh>& gsf1 = tgsf1();
267  const GeometricField<scalar, PatchField, GeoMesh>& gsf2 = tgsf2();
268 
269  if (!gsf1.dimensions().dimensionless())
270  {
272  << "Base field is not dimensionless: " << gsf1.dimensions()
273  << exit(FatalError);
274  }
275 
276  if (!gsf2.dimensions().dimensionless())
277  {
279  << "Exponent field is not dimensionless: " << gsf2.dimensions()
280  << exit(FatalError);
281  }
282 
284  (
286  <scalar, scalar, scalar, scalar, PatchField, GeoMesh>::New
287  (
288  tgsf1,
289  tgsf2,
290  "pow(" + gsf1.name() + ',' + gsf2.name() + ')',
291  dimless
292  )
293  );
294 
295  pow(tPow.ref(), gsf1, gsf2);
296 
297  tgsf1.clear();
298  tgsf2.clear();
299 
300  return tPow;
301 }
302 
303 
304 template<template<class> class PatchField, class GeoMesh>
305 void pow
306 (
309  const dimensioned<scalar>& ds
310 )
311 {
312  pow(tPow.primitiveFieldRef(), gsf.primitiveField(), ds.value());
313  pow(tPow.boundaryFieldRef(), gsf.boundaryField(), ds.value());
314 }
315 
316 
317 template<template<class> class PatchField, class GeoMesh>
319 (
321  const dimensionedScalar& ds
322 )
323 {
324  if (!ds.dimensions().dimensionless())
325  {
327  << "Exponent is not dimensionless: " << ds.dimensions()
328  << exit(FatalError);
329  }
330 
332  (
334  (
335  IOobject
336  (
337  "pow(" + gsf.name() + ',' + ds.name() + ')',
338  gsf.instance(),
339  gsf.db(),
342  ),
343  gsf.mesh(),
344  pow(gsf.dimensions(), ds)
345  )
346  );
347 
348  pow(tPow.ref(), gsf, ds);
349 
350  return tPow;
351 }
352 
353 template<template<class> class PatchField, class GeoMesh>
355 (
357  const dimensionedScalar& ds
358 )
359 {
360  if (!ds.dimensions().dimensionless())
361  {
363  << "Exponent is not dimensionless: " << ds.dimensions()
364  << exit(FatalError);
365  }
366 
368 
370  (
371  New
372  (
373  tgsf,
374  "pow(" + gsf.name() + ',' + ds.name() + ')',
375  pow(gsf.dimensions(), ds)
376  )
377  );
378 
379  pow(tPow.ref(), gsf, ds);
380 
381  tgsf.clear();
382 
383  return tPow;
384 }
385 
386 template<template<class> class PatchField, class GeoMesh>
388 (
390  const scalar& s
391 )
392 {
393  return pow(gsf, dimensionedScalar(s));
394 }
395 
396 template<template<class> class PatchField, class GeoMesh>
398 (
400  const scalar& s
401 )
402 {
403  return pow(tgsf, dimensionedScalar(s));
404 }
405 
406 
407 template<template<class> class PatchField, class GeoMesh>
408 void pow
409 (
411  const dimensioned<scalar>& ds,
413 )
414 {
415  pow(tPow.primitiveFieldRef(), ds.value(), gsf.primitiveField());
416  pow(tPow.boundaryFieldRef(), ds.value(), gsf.boundaryField());
417 }
418 
419 
420 template<template<class> class PatchField, class GeoMesh>
422 (
423  const dimensionedScalar& ds,
425 )
426 {
427  if (!ds.dimensions().dimensionless())
428  {
430  << "Base scalar is not dimensionless: " << ds.dimensions()
431  << exit(FatalError);
432  }
433 
434  if (!gsf.dimensions().dimensionless())
435  {
437  << "Exponent field is not dimensionless: " << gsf.dimensions()
438  << exit(FatalError);
439  }
440 
442  (
444  (
445  IOobject
446  (
447  "pow(" + ds.name() + ',' + gsf.name() + ')',
448  gsf.instance(),
449  gsf.db(),
452  ),
453  gsf.mesh(),
454  dimless
455  )
456  );
457 
458  pow(tPow.ref(), ds, gsf);
459 
460  return tPow;
461 }
462 
463 
464 template<template<class> class PatchField, class GeoMesh>
466 (
467  const dimensionedScalar& ds,
469 )
470 {
472 
473  if (!ds.dimensions().dimensionless())
474  {
476  << "Base scalar is not dimensionless: " << ds.dimensions()
477  << exit(FatalError);
478  }
479 
480  if (!gsf.dimensions().dimensionless())
481  {
483  << "Exponent field is not dimensionless: " << gsf.dimensions()
484  << exit(FatalError);
485  }
486 
488  (
489  New
490  (
491  tgsf,
492  "pow(" + ds.name() + ',' + gsf.name() + ')',
493  dimless
494  )
495  );
496 
497  pow(tPow.ref(), ds, gsf);
498 
499  tgsf.clear();
500 
501  return tPow;
502 }
503 
504 template<template<class> class PatchField, class GeoMesh>
506 (
507  const scalar& s,
509 )
510 {
511  return pow(dimensionedScalar(s), gsf);
512 }
513 
514 template<template<class> class PatchField, class GeoMesh>
516 (
517  const scalar& s,
519 )
520 {
521  return pow(dimensionedScalar(s), tgsf);
522 }
523 
524 
525 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
526 
527 template<template<class> class PatchField, class GeoMesh>
528 void atan2
529 (
533 )
534 {
535  atan2
536  (
537  Atan2.primitiveFieldRef(),
538  gsf1.primitiveField(),
539  gsf2.primitiveField()
540  );
541  atan2
542  (
543  Atan2.boundaryFieldRef(),
544  gsf1.boundaryField(),
545  gsf2.boundaryField()
546  );
547 }
548 
549 
550 template<template<class> class PatchField, class GeoMesh>
552 (
555 )
556 {
558  (
560  (
561  IOobject
562  (
563  "atan2(" + gsf1.name() + ',' + gsf2.name() + ')',
564  gsf1.instance(),
565  gsf1.db(),
568  ),
569  gsf1.mesh(),
570  atan2(gsf1.dimensions(), gsf2.dimensions())
571  )
572  );
573 
574  atan2(tAtan2.ref(), gsf1, gsf2);
575 
576  return tAtan2;
577 }
578 
579 
580 template<template<class> class PatchField, class GeoMesh>
582 (
585 )
586 {
587  const GeometricField<scalar, PatchField, GeoMesh>& gsf1 = tgsf1();
588 
590  (
591  New
592  (
593  tgsf1,
594  "atan2(" + gsf1.name() + ',' + gsf2.name() + ')',
595  atan2(gsf1.dimensions(), gsf2.dimensions())
596  )
597  );
598 
599  atan2(tAtan2.ref(), gsf1, gsf2);
600 
601  tgsf1.clear();
602 
603  return tAtan2;
604 }
605 
606 
607 template<template<class> class PatchField, class GeoMesh>
609 (
612 )
613 {
614  const GeometricField<scalar, PatchField, GeoMesh>& gsf2 = tgsf2();
615 
617  (
618  New
619  (
620  tgsf2,
621  "atan2(" + gsf1.name() + ',' + gsf2.name() + ')',
622  atan2( gsf1.dimensions(), gsf2.dimensions())
623  )
624  );
625 
626  atan2(tAtan2.ref(), gsf1, gsf2);
627 
628  tgsf2.clear();
629 
630  return tAtan2;
631 }
632 
633 template<template<class> class PatchField, class GeoMesh>
635 (
638 )
639 {
640  const GeometricField<scalar, PatchField, GeoMesh>& gsf1 = tgsf1();
641  const GeometricField<scalar, PatchField, GeoMesh>& gsf2 = tgsf2();
642 
644  (
646  <scalar, scalar, scalar, scalar, PatchField, GeoMesh>::New
647  (
648  tgsf1,
649  tgsf2,
650  "atan2(" + gsf1.name() + ',' + gsf2.name() + ')',
651  atan2(gsf1.dimensions(), gsf2.dimensions())
652  )
653  );
654 
655  atan2(tAtan2.ref(), gsf1, gsf2);
656 
657  tgsf1.clear();
658  tgsf2.clear();
659 
660  return tAtan2;
661 }
662 
663 
664 template<template<class> class PatchField, class GeoMesh>
665 void atan2
666 (
669  const dimensioned<scalar>& ds
670 )
671 {
672  atan2(tAtan2.primitiveFieldRef(), gsf.primitiveField(), ds.value());
673  atan2(tAtan2.boundaryFieldRef(), gsf.boundaryField(), ds.value());
674 }
675 
676 
677 template<template<class> class PatchField, class GeoMesh>
679 (
681  const dimensionedScalar& ds
682 )
683 {
685  (
687  (
688  IOobject
689  (
690  "atan2(" + gsf.name() + ',' + ds.name() + ')',
691  gsf.instance(),
692  gsf.db(),
695  ),
696  gsf.mesh(),
697  atan2(gsf.dimensions(), ds)
698  )
699  );
700 
701  atan2(tAtan2.ref(), gsf, ds);
702 
703  return tAtan2;
704 }
705 
706 template<template<class> class PatchField, class GeoMesh>
708 (
710  const dimensionedScalar& ds
711 )
712 {
714 
716  (
717  New
718  (
719  tgsf,
720  "atan2(" + gsf.name() + ',' + ds.name() + ')',
721  atan2(gsf.dimensions(), ds)
722  )
723  );
724 
725  atan2(tAtan2.ref(), gsf, ds);
726 
727  tgsf.clear();
728 
729  return tAtan2;
730 }
731 
732 template<template<class> class PatchField, class GeoMesh>
734 (
736  const scalar& s
737 )
738 {
739  return atan2(gsf, dimensionedScalar(s));
740 }
741 
742 template<template<class> class PatchField, class GeoMesh>
744 (
746  const scalar& s
747 )
748 {
749  return atan2(tgsf, dimensionedScalar(s));
750 }
751 
752 
753 template<template<class> class PatchField, class GeoMesh>
754 void atan2
755 (
757  const dimensioned<scalar>& ds,
759 )
760 {
761  atan2(tAtan2.primitiveFieldRef(), ds.value(), gsf.primitiveField());
762  atan2(tAtan2.boundaryFieldRef(), ds.value(), gsf.boundaryField());
763 }
764 
765 
766 template<template<class> class PatchField, class GeoMesh>
768 (
769  const dimensionedScalar& ds,
771 )
772 {
774  (
776  (
777  IOobject
778  (
779  "atan2(" + ds.name() + ',' + gsf.name() + ')',
780  gsf.instance(),
781  gsf.db(),
784  ),
785  gsf.mesh(),
786  atan2(ds, gsf.dimensions())
787  )
788  );
789 
790  atan2(tAtan2.ref(), ds, gsf);
791 
792  return tAtan2;
793 }
794 
795 
796 template<template<class> class PatchField, class GeoMesh>
798 (
799  const dimensionedScalar& ds,
801 )
802 {
804 
806  (
807  New
808  (
809  tgsf,
810  "atan2(" + ds.name() + ',' + gsf.name() + ')',
811  atan2(ds, gsf.dimensions())
812  )
813  );
814 
815  atan2(tAtan2.ref(), ds, gsf);
816 
817  tgsf.clear();
818 
819  return tAtan2;
820 }
821 
822 template<template<class> class PatchField, class GeoMesh>
824 (
825  const scalar& s,
827 )
828 {
829  return atan2(dimensionedScalar(s), gsf);
830 }
831 
832 template<template<class> class PatchField, class GeoMesh>
834 (
835  const scalar& s,
837 )
838 {
839  return atan2(dimensionedScalar(s), tgsf);
840 }
841 
842 
843 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
844 
845 UNARY_FUNCTION(scalar, scalar, pow3, pow3)
846 UNARY_FUNCTION(scalar, scalar, pow4, pow4)
847 UNARY_FUNCTION(scalar, scalar, pow5, pow5)
848 UNARY_FUNCTION(scalar, scalar, pow6, pow6)
849 UNARY_FUNCTION(scalar, scalar, pow025, pow025)
850 UNARY_FUNCTION(scalar, scalar, sqrt, sqrt)
851 UNARY_FUNCTION(scalar, scalar, cbrt, cbrt)
852 UNARY_FUNCTION(scalar, scalar, sign, sign)
853 UNARY_FUNCTION(scalar, scalar, pos, pos)
854 UNARY_FUNCTION(scalar, scalar, pos0, pos0)
855 UNARY_FUNCTION(scalar, scalar, neg, neg)
856 UNARY_FUNCTION(scalar, scalar, neg0, neg0)
857 UNARY_FUNCTION(scalar, scalar, posPart, posPart)
858 UNARY_FUNCTION(scalar, scalar, negPart, negPart)
859 
860 UNARY_FUNCTION(scalar, scalar, exp, trans)
861 UNARY_FUNCTION(scalar, scalar, log, trans)
862 UNARY_FUNCTION(scalar, scalar, log10, trans)
863 UNARY_FUNCTION(scalar, scalar, sin, trans)
864 UNARY_FUNCTION(scalar, scalar, cos, trans)
865 UNARY_FUNCTION(scalar, scalar, tan, trans)
866 UNARY_FUNCTION(scalar, scalar, asin, trans)
867 UNARY_FUNCTION(scalar, scalar, acos, trans)
868 UNARY_FUNCTION(scalar, scalar, atan, trans)
869 UNARY_FUNCTION(scalar, scalar, sinh, trans)
870 UNARY_FUNCTION(scalar, scalar, cosh, trans)
871 UNARY_FUNCTION(scalar, scalar, tanh, trans)
872 UNARY_FUNCTION(scalar, scalar, asinh, trans)
873 UNARY_FUNCTION(scalar, scalar, acosh, trans)
874 UNARY_FUNCTION(scalar, scalar, atanh, trans)
875 UNARY_FUNCTION(scalar, scalar, erf, trans)
876 UNARY_FUNCTION(scalar, scalar, erfc, trans)
877 UNARY_FUNCTION(scalar, scalar, lgamma, trans)
878 UNARY_FUNCTION(scalar, scalar, j0, trans)
879 UNARY_FUNCTION(scalar, scalar, j1, trans)
880 UNARY_FUNCTION(scalar, scalar, y0, trans)
881 UNARY_FUNCTION(scalar, scalar, y1, trans)
882 
883 
884 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
885 
886 #define BesselFunc(func) \
887  \
888 template<template<class> class PatchField, class GeoMesh> \
889 void func \
890 ( \
891  GeometricField<scalar, PatchField, GeoMesh>& gsf, \
892  const int n, \
893  const GeometricField<scalar, PatchField, GeoMesh>& gsf1 \
894 ) \
895 { \
896  func(gsf.primitiveFieldRef(), n, gsf1.primitiveField()); \
897  func(gsf.boundaryFieldRef(), n, gsf1.boundaryField()); \
898 } \
899  \
900 template<template<class> class PatchField, class GeoMesh> \
901 tmp<GeometricField<scalar, PatchField, GeoMesh>> func \
902 ( \
903  const int n, \
904  const GeometricField<scalar, PatchField, GeoMesh>& gsf \
905 ) \
906 { \
907  if (!gsf.dimensions().dimensionless()) \
908  { \
909  FatalErrorInFunction \
910  << "gsf not dimensionless" \
911  << abort(FatalError); \
912  } \
913  \
914  tmp<GeometricField<scalar, PatchField, GeoMesh>> tFunc \
915  ( \
916  new GeometricField<scalar, PatchField, GeoMesh> \
917  ( \
918  IOobject \
919  ( \
920  #func "(" + gsf.name() + ')', \
921  gsf.instance(), \
922  gsf.db(), \
923  IOobject::NO_READ, \
924  IOobject::NO_WRITE \
925  ), \
926  gsf.mesh(), \
927  dimless \
928  ) \
929  ); \
930  \
931  func(tFunc.ref(), n, gsf); \
932  \
933  return tFunc; \
934 } \
935  \
936 template<template<class> class PatchField, class GeoMesh> \
937 tmp<GeometricField<scalar, PatchField, GeoMesh>> func \
938 ( \
939  const int n, \
940  const tmp<GeometricField<scalar, PatchField, GeoMesh>>& tgsf \
941 ) \
942 { \
943  const GeometricField<scalar, PatchField, GeoMesh>& gsf = tgsf(); \
944  \
945  if (!gsf.dimensions().dimensionless()) \
946  { \
947  FatalErrorInFunction \
948  << " : gsf not dimensionless" \
949  << abort(FatalError); \
950  } \
951  \
952  tmp<GeometricField<scalar, PatchField, GeoMesh>> tFunc \
953  ( \
954  New \
955  ( \
956  tgsf, \
957  #func "(" + gsf.name() + ')', \
958  dimless \
959  ) \
960  ); \
961  \
962  func(tFunc.ref(), n, gsf); \
963  \
964  tgsf.clear(); \
965  \
966  return tFunc; \
967 }
968 
971 
972 #undef BesselFunc
973 
974 
975 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
976 
977 } // End namespace Foam
978 
979 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
980 
981 #include "undefFieldFunctionsM.H"
982 
983 // ************************************************************************* //
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:291
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)
const fileName & instance() const
Definition: IOobject.H:386
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
const objectRegistry & db() const
Return the local objectRegistry.
Definition: IOobject.C:331
dimensionedScalar tan(const dimensionedScalar &ds)
dimensionedScalar j0(const dimensionedScalar &ds)
#define BINARY_TYPE_OPERATOR(ReturnType, Type1, Type2, Op, OpName, OpFunc)
IOobject defines the attributes of an object for which implicit objectRegistry management is supporte...
Definition: IOobject.H:92
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)