FreeFOAM The Cross-Platform CFD Toolkit
GeometricFieldFunctions.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) 1991-2010 OpenCFD Ltd.
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 
27 
28 #define TEMPLATE \
29  template<class Type, template<class> class PatchField, class GeoMesh>
31 
32 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
33 
34 namespace Foam
35 {
36 
37 // * * * * * * * * * * * * * * * Global functions * * * * * * * * * * * * * //
38 
39 template<class Type, template<class> class PatchField, class GeoMesh>
40 void component
41 (
43  <
45  PatchField,
46  GeoMesh
47  >& gcf,
49  const direction d
50 )
51 {
52  component(gcf.internalField(), gf.internalField(), d);
53  component(gcf.boundaryField(), gf.boundaryField(), d);
54 }
55 
56 
57 template<class Type, template<class> class PatchField, class GeoMesh>
58 void T
59 (
62 )
63 {
64  T(gf.internalField(), gf1.internalField());
65  T(gf.boundaryField(), gf1.boundaryField());
66 }
67 
68 
69 template<class Type, template<class> class PatchField, class GeoMesh, int r>
70 void pow
71 (
72  GeometricField<typename powProduct<Type, r>::type, PatchField, GeoMesh>& gf,
74 )
75 {
76  pow(gf.internalField(), gf1.internalField(), r);
77  pow(gf.boundaryField(), gf1.boundaryField(), r);
78 }
79 
80 template<class Type, template<class> class PatchField, class GeoMesh, int r>
82 pow
83 (
86 )
87 {
88  typedef typename powProduct<Type, r>::type powProductType;
89 
91  (
93  (
94  IOobject
95  (
96  "pow(" + gf.name() + ',' + name(r) + ')',
97  gf.instance(),
98  gf.db(),
101  ),
102  gf.mesh(),
103  pow(gf.dimensions(), r)
104  )
105  );
106 
107  pow<Type, r, PatchField, GeoMesh>(tPow(), gf);
108 
109  return tPow;
110 }
111 
112 
113 template<class Type, template<class> class PatchField, class GeoMesh, int r>
115 pow
116 (
119 )
120 {
121  typedef typename powProduct<Type, r>::type powProductType;
122 
124 
126  (
128  (
129  IOobject
130  (
131  "pow(" + gf.name() + ',' + name(r) + ')',
132  gf.instance(),
133  gf.db(),
136  ),
137  gf.mesh(),
138  pow(gf.dimensions(), r)
139  )
140  );
141 
142  pow<Type, r, PatchField, GeoMesh>(tPow(), gf);
143 
144  tgf.clear();
145 
146  return tPow;
147 }
148 
149 
150 template<class Type, template<class> class PatchField, class GeoMesh>
151 void sqr
152 (
154  <typename outerProduct<Type, Type>::type, PatchField, GeoMesh>& gf,
156 )
157 {
158  sqr(gf.internalField(), gf1.internalField());
159  sqr(gf.boundaryField(), gf1.boundaryField());
160 }
161 
162 template<class Type, template<class> class PatchField, class GeoMesh>
163 tmp
164 <
165  GeometricField
166  <
168  PatchField,
169  GeoMesh
170  >
171 >
173 {
174  typedef typename outerProduct<Type, Type>::type outerProductType;
175 
177  (
179  (
180  IOobject
181  (
182  "sqr(" + gf.name() + ')',
183  gf.instance(),
184  gf.db(),
187  ),
188  gf.mesh(),
189  sqr(gf.dimensions())
190  )
191  );
192 
193  sqr(tSqr(), gf);
194 
195  return tSqr;
196 }
197 
198 template<class Type, template<class> class PatchField, class GeoMesh>
199 tmp
200 <
201  GeometricField
202  <
204  PatchField,
205  GeoMesh
206  >
207 >
209 {
210  typedef typename outerProduct<Type, Type>::type outerProductType;
211 
213 
215  (
217  (
218  IOobject
219  (
220  "sqr(" + gf.name() + ')',
221  gf.instance(),
222  gf.db(),
225  ),
226  gf.mesh(),
227  sqr(gf.dimensions())
228  )
229  );
230 
231  sqr(tSqr(), gf);
232 
233  tgf.clear();
234 
235  return tSqr;
236 }
237 
238 
239 template<class Type, template<class> class PatchField, class GeoMesh>
240 void magSqr
241 (
244 )
245 {
246  magSqr(gsf.internalField(), gf.internalField());
247  magSqr(gsf.boundaryField(), gf.boundaryField());
248 }
249 
250 template<class Type, template<class> class PatchField, class GeoMesh>
251 tmp<GeometricField<scalar, PatchField, GeoMesh> > magSqr
252 (
254 )
255 {
257  (
259  (
260  IOobject
261  (
262  "magSqr(" + gf.name() + ')',
263  gf.instance(),
264  gf.db(),
267  ),
268  gf.mesh(),
269  sqr(gf.dimensions())
270  )
271  );
272 
273  magSqr(tMagSqr(), gf);
274 
275  return tMagSqr;
276 }
277 
278 template<class Type, template<class> class PatchField, class GeoMesh>
279 tmp<GeometricField<scalar, PatchField, GeoMesh> > magSqr
280 (
282 )
283 {
285 
287  (
289  (
290  IOobject
291  (
292  "magSqr(" + gf.name() + ')',
293  gf.instance(),
294  gf.db(),
297  ),
298  gf.mesh(),
299  sqr(gf.dimensions())
300  )
301  );
302 
303  magSqr(tMagSqr(), gf);
304 
305  tgf.clear();
306 
307  return tMagSqr;
308 }
309 
310 
311 template<class Type, template<class> class PatchField, class GeoMesh>
312 void mag
313 (
316 )
317 {
318  mag(gsf.internalField(), gf.internalField());
319  mag(gsf.boundaryField(), gf.boundaryField());
320 }
321 
322 template<class Type, template<class> class PatchField, class GeoMesh>
323 tmp<GeometricField<scalar, PatchField, GeoMesh> > mag
324 (
326 )
327 {
329  (
331  (
332  IOobject
333  (
334  "mag(" + gf.name() + ')',
335  gf.instance(),
336  gf.db(),
339  ),
340  gf.mesh(),
341  gf.dimensions()
342  )
343  );
344 
345  mag(tMag(), gf);
346 
347  return tMag;
348 }
349 
350 template<class Type, template<class> class PatchField, class GeoMesh>
351 tmp<GeometricField<scalar, PatchField, GeoMesh> > mag
352 (
354 )
355 {
357 
359  (
361  (
362  IOobject
363  (
364  "mag(" + gf.name() + ')',
365  gf.instance(),
366  gf.db(),
369  ),
370  gf.mesh(),
371  gf.dimensions()
372  )
373  );
374 
375  mag(tMag(), gf);
376 
377  tgf.clear();
378 
379  return tMag;
380 }
381 
382 
383 template<class Type, template<class> class PatchField, class GeoMesh>
384 void cmptAv
385 (
387  <
389  PatchField,
390  GeoMesh
391  >& gcf,
393 )
394 {
395  cmptAv(gcf.internalField(), gf.internalField());
396  cmptAv(gcf.boundaryField(), gf.boundaryField());
397 }
398 
399 template<class Type, template<class> class PatchField, class GeoMesh>
400 tmp
401 <
402  GeometricField
403  <
404  typename GeometricField<Type, PatchField, GeoMesh>::cmptType,
405  PatchField,
406  GeoMesh
407  >
408 >
410 {
412  cmptType;
413 
415  (
417  (
418  IOobject
419  (
420  "cmptAv(" + gf.name() + ')',
421  gf.instance(),
422  gf.db(),
425  ),
426  gf.mesh(),
427  gf.dimensions()
428  )
429  );
430 
431  cmptAv(CmptAv(), gf);
432 
433  return CmptAv;
434 }
435 
436 template<class Type, template<class> class PatchField, class GeoMesh>
437 tmp
438 <
439  GeometricField
440  <
441  typename GeometricField<Type, PatchField, GeoMesh>::cmptType,
442  PatchField,
443  GeoMesh
444  >
445 >
447 {
449  cmptType;
450 
452 
454  (
456  (
457  IOobject
458  (
459  "cmptAv(" + gf.name() + ')',
460  gf.instance(),
461  gf.db(),
464  ),
465  gf.mesh(),
466  gf.dimensions()
467  )
468  );
469 
470  cmptAv(CmptAv(), gf);
471 
472  tgf.clear();
473 
474  return CmptAv;
475 }
476 
477 
478 #define UNARY_REDUCTION_FUNCTION_WITH_BOUNDARY(returnType, func, gFunc) \
479  \
480 template<class Type, template<class> class PatchField, class GeoMesh> \
481 dimensioned<returnType> func \
482 ( \
483  const GeometricField<Type, PatchField, GeoMesh>& gf \
484 ) \
485 { \
486  return dimensioned<Type> \
487  ( \
488  #func "(" + gf.name() + ')', \
489  gf.dimensions(), \
490  Foam::func(gFunc(gf.internalField()), gFunc(gf.boundaryField())) \
491  ); \
492 } \
493  \
494 template<class Type, template<class> class PatchField, class GeoMesh> \
495 dimensioned<returnType> func \
496 ( \
497  const tmp<GeometricField<Type, PatchField, GeoMesh> >& tgf1 \
498 ) \
499 { \
500  dimensioned<returnType> res = func(tgf1()); \
501  tgf1.clear(); \
502  return res; \
503 }
504 
507 
508 #undef UNARY_REDUCTION_FUNCTION_WITH_BOUNDARY
509 
510 
511 #define UNARY_REDUCTION_FUNCTION(returnType, func, gFunc) \
512  \
513 template<class Type, template<class> class PatchField, class GeoMesh> \
514 dimensioned<returnType> func \
515 ( \
516  const GeometricField<Type, PatchField, GeoMesh>& gf \
517 ) \
518 { \
519  return dimensioned<Type> \
520  ( \
521  #func "(" + gf.name() + ')', \
522  gf.dimensions(), \
523  gFunc(gf.internalField()) \
524  ); \
525 } \
526  \
527 template<class Type, template<class> class PatchField, class GeoMesh> \
528 dimensioned<returnType> func \
529 ( \
530  const tmp<GeometricField<Type, PatchField, GeoMesh> >& tgf1 \
531 ) \
532 { \
533  dimensioned<returnType> res = func(tgf1()); \
534  tgf1.clear(); \
535  return res; \
536 }
537 
541 
542 #undef UNARY_REDUCTION_FUNCTION
543 
544 
545 BINARY_FUNCTION(Type, Type, Type, max)
546 BINARY_FUNCTION(Type, Type, Type, min)
547 BINARY_FUNCTION(Type, Type, Type, cmptMultiply)
548 BINARY_FUNCTION(Type, Type, Type, cmptDivide)
549 
550 BINARY_TYPE_FUNCTION(Type, Type, Type, max)
551 BINARY_TYPE_FUNCTION(Type, Type, Type, min)
552 BINARY_TYPE_FUNCTION(Type, Type, Type, cmptMultiply)
553 BINARY_TYPE_FUNCTION(Type, Type, Type, cmptDivide)
554 
555 
556 // * * * * * * * * * * * * * * * Global operators * * * * * * * * * * * * * //
557 
558 UNARY_OPERATOR(Type, Type, -, negate, transform)
559 
560 #ifndef __INTEL_COMPILER
561 BINARY_OPERATOR(Type, Type, scalar, *, '*', multiply)
562 BINARY_OPERATOR(Type, scalar, Type, *, '*', multiply)
563 #endif
564 BINARY_OPERATOR(Type, Type, scalar, /, '|', divide)
565 
566 BINARY_TYPE_OPERATOR_SF(Type, scalar, Type, *, '*', multiply)
567 BINARY_TYPE_OPERATOR_FS(Type, Type, scalar, *, '*', multiply)
568 
569 BINARY_TYPE_OPERATOR_FS(Type, Type, scalar, /, '|', divide)
570 
571 
572 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
573 
574 #define PRODUCT_OPERATOR(product, op, opFunc) \
575  \
576 template \
577 <class Type1, class Type2, template<class> class PatchField, class GeoMesh> \
578 void opFunc \
579 ( \
580  GeometricField \
581  <typename product<Type1, Type2>::type, PatchField, GeoMesh>& gf, \
582  const GeometricField<Type1, PatchField, GeoMesh>& gf1, \
583  const GeometricField<Type2, PatchField, GeoMesh>& gf2 \
584 ) \
585 { \
586  Foam::opFunc(gf.internalField(), gf1.internalField(), gf2.internalField());\
587  Foam::opFunc(gf.boundaryField(), gf1.boundaryField(), gf2.boundaryField());\
588 } \
589  \
590 template \
591 <class Type1, class Type2, template<class> class PatchField, class GeoMesh> \
592 tmp \
593 < \
594  GeometricField<typename product<Type1, Type2>::type, PatchField, GeoMesh> \
595 > \
596 operator op \
597 ( \
598  const GeometricField<Type1, PatchField, GeoMesh>& gf1, \
599  const GeometricField<Type2, PatchField, GeoMesh>& gf2 \
600 ) \
601 { \
602  typedef typename product<Type1, Type2>::type productType; \
603  tmp<GeometricField<productType, PatchField, GeoMesh> > tRes \
604  ( \
605  new GeometricField<productType, PatchField, GeoMesh> \
606  ( \
607  IOobject \
608  ( \
609  '(' + gf1.name() + #op + gf2.name() + ')', \
610  gf1.instance(), \
611  gf1.db(), \
612  IOobject::NO_READ, \
613  IOobject::NO_WRITE \
614  ), \
615  gf1.mesh(), \
616  gf1.dimensions() op gf2.dimensions() \
617  ) \
618  ); \
619  \
620  Foam::opFunc(tRes(), gf1, gf2); \
621  \
622  return tRes; \
623 } \
624  \
625 template \
626 <class Type1, class Type2, template<class> class PatchField, class GeoMesh> \
627 tmp \
628 < \
629  GeometricField<typename product<Type1, Type2>::type, PatchField, GeoMesh> \
630 > \
631 operator op \
632 ( \
633  const GeometricField<Type1, PatchField, GeoMesh>& gf1, \
634  const tmp<GeometricField<Type2, PatchField, GeoMesh> >& tgf2 \
635 ) \
636 { \
637  typedef typename product<Type1, Type2>::type productType; \
638  \
639  const GeometricField<Type2, PatchField, GeoMesh>& gf2 = tgf2(); \
640  \
641  tmp<GeometricField<productType, PatchField, GeoMesh> > tRes = \
642  reuseTmpGeometricField<productType, Type2, PatchField, GeoMesh>::New \
643  ( \
644  tgf2, \
645  '(' + gf1.name() + #op + gf2.name() + ')', \
646  gf1.dimensions() op gf2.dimensions() \
647  ); \
648  \
649  Foam::opFunc(tRes(), gf1, gf2); \
650  \
651  reuseTmpGeometricField<productType, Type2, PatchField, GeoMesh> \
652  ::clear(tgf2); \
653  \
654  return tRes; \
655 } \
656  \
657 template \
658 <class Type1, class Type2, template<class> class PatchField, class GeoMesh> \
659 tmp \
660 < \
661  GeometricField<typename product<Type1, Type2>::type, PatchField, GeoMesh> \
662 > \
663 operator op \
664 ( \
665  const tmp<GeometricField<Type1, PatchField, GeoMesh> >& tgf1, \
666  const GeometricField<Type2, PatchField, GeoMesh>& gf2 \
667 ) \
668 { \
669  typedef typename product<Type1, Type2>::type productType; \
670  \
671  const GeometricField<Type1, PatchField, GeoMesh>& gf1 = tgf1(); \
672  \
673  tmp<GeometricField<productType, PatchField, GeoMesh> > tRes = \
674  reuseTmpGeometricField<productType, Type1, PatchField, GeoMesh>::New \
675  ( \
676  tgf1, \
677  '(' + gf1.name() + #op + gf2.name() + ')', \
678  gf1.dimensions() op gf2.dimensions() \
679  ); \
680  \
681  Foam::opFunc(tRes(), gf1, gf2); \
682  \
683  reuseTmpGeometricField<productType, Type1, PatchField, GeoMesh> \
684  ::clear(tgf1); \
685  \
686  return tRes; \
687 } \
688  \
689 template \
690 <class Type1, class Type2, template<class> class PatchField, class GeoMesh> \
691 tmp \
692 < \
693  GeometricField<typename product<Type1, Type2>::type, PatchField, GeoMesh> \
694 > \
695 operator op \
696 ( \
697  const tmp<GeometricField<Type1, PatchField, GeoMesh> >& tgf1, \
698  const tmp<GeometricField<Type2, PatchField, GeoMesh> >& tgf2 \
699 ) \
700 { \
701  typedef typename product<Type1, Type2>::type productType; \
702  \
703  const GeometricField<Type1, PatchField, GeoMesh>& gf1 = tgf1(); \
704  const GeometricField<Type2, PatchField, GeoMesh>& gf2 = tgf2(); \
705  \
706  tmp<GeometricField<productType, PatchField, GeoMesh> > tRes = \
707  reuseTmpTmpGeometricField \
708  <productType, Type1, Type1, Type2, PatchField, GeoMesh>::New \
709  ( \
710  tgf1, \
711  tgf2, \
712  '(' + gf1.name() + #op + gf2.name() + ')', \
713  gf1.dimensions() op gf2.dimensions() \
714  ); \
715  \
716  Foam::opFunc(tRes(), gf1, gf2); \
717  \
718  reuseTmpTmpGeometricField \
719  <productType, Type1, Type1, Type2, PatchField, GeoMesh> \
720  ::clear(tgf1, tgf2); \
721  \
722  return tRes; \
723 } \
724  \
725 template \
726 <class Form, class Type, template<class> class PatchField, class GeoMesh> \
727 void opFunc \
728 ( \
729  GeometricField \
730  <typename product<Type, Form>::type, PatchField, GeoMesh>& gf, \
731  const GeometricField<Type, PatchField, GeoMesh>& gf1, \
732  const dimensioned<Form>& dvs \
733 ) \
734 { \
735  Foam::opFunc(gf.internalField(), gf1.internalField(), dvs.value()); \
736  Foam::opFunc(gf.boundaryField(), gf1.boundaryField(), dvs.value()); \
737 } \
738  \
739 template \
740 <class Form, class Type, template<class> class PatchField, class GeoMesh> \
741 tmp<GeometricField<typename product<Type, Form>::type, PatchField, GeoMesh> > \
742 operator op \
743 ( \
744  const GeometricField<Type, PatchField, GeoMesh>& gf1, \
745  const dimensioned<Form>& dvs \
746 ) \
747 { \
748  typedef typename product<Type, Form>::type productType; \
749  \
750  tmp<GeometricField<productType, PatchField, GeoMesh> > tRes \
751  ( \
752  new GeometricField<productType, PatchField, GeoMesh> \
753  ( \
754  IOobject \
755  ( \
756  '(' + gf1.name() + #op + dvs.name() + ')', \
757  gf1.instance(), \
758  gf1.db(), \
759  IOobject::NO_READ, \
760  IOobject::NO_WRITE \
761  ), \
762  gf1.mesh(), \
763  gf1.dimensions() op dvs.dimensions() \
764  ) \
765  ); \
766  \
767  Foam::opFunc(tRes(), gf1, dvs); \
768  \
769  return tRes; \
770 } \
771  \
772 template \
773 < \
774  class Form, \
775  class Cmpt, \
776  int nCmpt, \
777  class Type, template<class> class PatchField, \
778  class GeoMesh \
779 > \
780 tmp<GeometricField<typename product<Form, Type>::type, PatchField, GeoMesh> > \
781 operator op \
782 ( \
783  const GeometricField<Type, PatchField, GeoMesh>& gf1, \
784  const VectorSpace<Form,Cmpt,nCmpt>& vs \
785 ) \
786 { \
787  return gf1 op dimensioned<Form>(static_cast<const Form&>(vs)); \
788 } \
789  \
790  \
791 template \
792 <class Form, class Type, template<class> class PatchField, class GeoMesh> \
793 tmp<GeometricField<typename product<Type, Form>::type, PatchField, GeoMesh> > \
794 operator op \
795 ( \
796  const tmp<GeometricField<Type, PatchField, GeoMesh> >& tgf1, \
797  const dimensioned<Form>& dvs \
798 ) \
799 { \
800  typedef typename product<Type, Form>::type productType; \
801  \
802  const GeometricField<Type, PatchField, GeoMesh>& gf1 = tgf1(); \
803  \
804  tmp<GeometricField<productType, PatchField, GeoMesh> > tRes = \
805  reuseTmpGeometricField<productType, Type, PatchField, GeoMesh>::New \
806  ( \
807  tgf1, \
808  '(' + gf1.name() + #op + dvs.name() + ')', \
809  gf1.dimensions() op dvs.dimensions() \
810  ); \
811  \
812  Foam::opFunc(tRes(), gf1, dvs); \
813  \
814  reuseTmpGeometricField<productType, Type, PatchField, GeoMesh> \
815  ::clear(tgf1); \
816  \
817  return tRes; \
818 } \
819  \
820 template \
821 < \
822  class Form, \
823  class Cmpt, \
824  int nCmpt, \
825  class Type, template<class> class PatchField, \
826  class GeoMesh \
827 > \
828 tmp<GeometricField<typename product<Form, Type>::type, PatchField, GeoMesh> > \
829 operator op \
830 ( \
831  const tmp<GeometricField<Type, PatchField, GeoMesh> >& tgf1, \
832  const VectorSpace<Form,Cmpt,nCmpt>& vs \
833 ) \
834 { \
835  return tgf1 op dimensioned<Form>(static_cast<const Form&>(vs)); \
836 } \
837  \
838  \
839 template \
840 <class Form, class Type, template<class> class PatchField, class GeoMesh> \
841 void opFunc \
842 ( \
843  GeometricField \
844  <typename product<Form, Type>::type, PatchField, GeoMesh>& gf, \
845  const dimensioned<Form>& dvs, \
846  const GeometricField<Type, PatchField, GeoMesh>& gf1 \
847 ) \
848 { \
849  Foam::opFunc(gf.internalField(), dvs.value(), gf1.internalField()); \
850  Foam::opFunc(gf.boundaryField(), dvs.value(), gf1.boundaryField()); \
851 } \
852  \
853 template \
854 <class Form, class Type, template<class> class PatchField, class GeoMesh> \
855 tmp<GeometricField<typename product<Form, Type>::type, PatchField, GeoMesh> > \
856 operator op \
857 ( \
858  const dimensioned<Form>& dvs, \
859  const GeometricField<Type, PatchField, GeoMesh>& gf1 \
860 ) \
861 { \
862  typedef typename product<Form, Type>::type productType; \
863  tmp<GeometricField<productType, PatchField, GeoMesh> > tRes \
864  ( \
865  new GeometricField<productType, PatchField, GeoMesh> \
866  ( \
867  IOobject \
868  ( \
869  '(' + dvs.name() + #op + gf1.name() + ')', \
870  gf1.instance(), \
871  gf1.db(), \
872  IOobject::NO_READ, \
873  IOobject::NO_WRITE \
874  ), \
875  gf1.mesh(), \
876  dvs.dimensions() op gf1.dimensions() \
877  ) \
878  ); \
879  \
880  Foam::opFunc(tRes(), dvs, gf1); \
881  \
882  return tRes; \
883 } \
884  \
885 template \
886 < \
887  class Form, \
888  class Cmpt, \
889  int nCmpt, \
890  class Type, template<class> class PatchField, \
891  class GeoMesh \
892 > \
893 tmp<GeometricField<typename product<Form, Type>::type, PatchField, GeoMesh> > \
894 operator op \
895 ( \
896  const VectorSpace<Form,Cmpt,nCmpt>& vs, \
897  const GeometricField<Type, PatchField, GeoMesh>& gf1 \
898 ) \
899 { \
900  return dimensioned<Form>(static_cast<const Form&>(vs)) op gf1; \
901 } \
902  \
903 template \
904 <class Form, class Type, template<class> class PatchField, class GeoMesh> \
905 tmp<GeometricField<typename product<Form, Type>::type, PatchField, GeoMesh> > \
906 operator op \
907 ( \
908  const dimensioned<Form>& dvs, \
909  const tmp<GeometricField<Type, PatchField, GeoMesh> >& tgf1 \
910 ) \
911 { \
912  typedef typename product<Form, Type>::type productType; \
913  \
914  const GeometricField<Type, PatchField, GeoMesh>& gf1 = tgf1(); \
915  \
916  tmp<GeometricField<productType, PatchField, GeoMesh> > tRes = \
917  reuseTmpGeometricField<productType, Type, PatchField, GeoMesh>::New \
918  ( \
919  tgf1, \
920  '(' + dvs.name() + #op + gf1.name() + ')', \
921  dvs.dimensions() op gf1.dimensions() \
922  ); \
923  \
924  Foam::opFunc(tRes(), dvs, gf1); \
925  \
926  reuseTmpGeometricField<productType, Type, PatchField, GeoMesh> \
927  ::clear(tgf1); \
928  \
929  return tRes; \
930 } \
931  \
932 template \
933 < \
934  class Form, \
935  class Cmpt, \
936  int nCmpt, \
937  class Type, template<class> class PatchField, \
938  class GeoMesh \
939 > \
940 tmp<GeometricField<typename product<Form, Type>::type, PatchField, GeoMesh> > \
941 operator op \
942 ( \
943  const VectorSpace<Form,Cmpt,nCmpt>& vs, \
944  const tmp<GeometricField<Type, PatchField, GeoMesh> >& tgf1 \
945 ) \
946 { \
947  return dimensioned<Form>(static_cast<const Form&>(vs)) op tgf1; \
948 }
949 
952 
957 
958 #undef PRODUCT_OPERATOR
959 
960 
961 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
962 
963 } // End namespace Foam
964 
965 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
966 
968 
969 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //