ops.H
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-2016 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 InClass
25  Foam::Pstream
26 
27 Description
28  Combination-Reduction operation for a parallel run.
29 
30  The information from all nodes is collected on the master node,
31  combined using the given combination function and the result is
32  broadcast to all nodes
33 
34 \*---------------------------------------------------------------------------*/
35 
36 #ifndef ops_H
37 #define ops_H
38 
39 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
40 
41 namespace Foam
42 {
43 
44 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
45 
46 #define EqOp(opName, op) \
47  \
48 template<class T1, class T2> \
49 class opName##Op2 \
50 { \
51 public: \
52  \
53  void operator()(T1& x, const T2& y) const \
54  { \
55  op; \
56  } \
57 }; \
58  \
59 template<class T> \
60 class opName##Op \
61 { \
62 public: \
63  \
64  void operator()(T& x, const T& y) const \
65  { \
66  op; \
67  } \
68 };
69 
70 EqOp(eq, x = y)
71 EqOp(plusEq, x += y)
72 EqOp(minusEq, x -= y)
73 EqOp(multiplyEq, x *= y)
74 EqOp(divideEq, x /= y)
75 EqOp(eqSqr, x = sqr(y))
76 EqOp(eqMag, x = mag(y))
77 EqOp(plusEqMagSqr, x += magSqr(y))
78 EqOp(maxEq, x = max(x, y))
79 EqOp(minEq, x = min(x, y))
80 EqOp(minMagSqrEq, x = (magSqr(x)<=magSqr(y) ? x : y))
81 EqOp(maxMagSqrEq, x = (magSqr(x)>=magSqr(y) ? x : y))
82 EqOp(andEq, x = (x && y))
83 EqOp(orEq, x = (x || y))
84 
85 EqOp(eqMinus, x = -y)
86 
87 EqOp(nopEq, (void)x)
88 
89 #undef EqOp
90 
91 
92 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
93 
94 #if __GNUC__
95 #define WARNRETURN __attribute__((warn_unused_result))
96 #else
97 #define WARNRETURN
98 #endif
99 
100 #define Op(opName, op) \
101  \
102  template<class T, class T1, class T2> \
103  class opName##Op3 \
104  { \
105  public: \
106  \
107  T operator()(const T1& x, const T2& y) const WARNRETURN \
108  { \
109  return op; \
110  } \
111  }; \
112  \
113  template<class T1, class T2> \
114  class opName##Op2 \
115  { \
116  public: \
117  \
118  T1 operator()(const T1& x, const T2& y) const WARNRETURN \
119  { \
120  return op; \
121  } \
122  }; \
123  \
124  template<class T> \
125  class opName##Op \
126  { \
127  public: \
128  \
129  T operator()(const T& x, const T& y) const WARNRETURN \
130  { \
131  return op; \
132  } \
133  };
134 
135 
136 #define weightedOp(opName, op) \
137  \
138  template<class Type, class CombineOp> \
139  class opName##WeightedOp \
140  { \
141  const CombineOp& cop_; \
142  \
143  public: \
144  \
145  opName##WeightedOp(const CombineOp& cop) \
146  : \
147  cop_(cop) \
148  {} \
149  \
150  void operator() \
151  ( \
152  Type& x, \
153  const label index, \
154  const Type& y, \
155  const scalar weight \
156  ) const \
157  { \
158  cop_(x, op); \
159  } \
160  }; \
161 
162 
163 Op(sum, x + y)
164 
165 Op(plus, x + y)
166 Op(minus, x - y)
167 Op(multiply, x * y)
168 Op(divide, x / y)
169 Op(cmptMultiply, cmptMultiply(x, y))
170 Op(cmptPow, cmptPow(x, y))
171 Op(cmptDivide, cmptDivide(x, y))
172 Op(stabilise, stabilise(x, y))
173 Op(max, max(x, y))
174 Op(min, min(x, y))
175 Op(minMagSqr, (magSqr(x)<=magSqr(y) ? x : y))
176 Op(maxMagSqr, (magSqr(x)>=magSqr(y) ? x : y))
177 Op(minMod, minMod(x, y))
178 Op(and, x && y)
179 Op(or, x || y)
180 Op(eqEq, x == y)
181 Op(less, x < y)
182 Op(lessEq, x <= y)
183 Op(greater, x > y)
184 Op(greaterEq, x >= y)
185 
186 weightedOp(multiply, (weight*y))
187 
188 #undef Op
189 #undef weightedOp
190 #undef WARNRETURN
191 
192 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
193 
194 } // End namespace Foam
195 
196 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
197 
198 #endif
199 
200 // ************************************************************************* //
void divide(FieldField< Field, Type > &f, const FieldField< Field, Type > &f1, const FieldField< Field, scalar > &f2)
void multiply(FieldField< Field, Type > &f, const FieldField< Field, Type > &f1, const FieldField< Field, scalar > &f2)
dimensioned< Type > max(const dimensioned< Type > &, const dimensioned< Type > &)
Type minMagSqr(const UList< Type > &f)
dimensionedSymmTensor sqr(const dimensionedVector &dv)
#define Op(opName, op)
Definition: ops.H:100
static bool less(const vector &x, const vector &y)
To compare normals.
dimensioned< Type > cmptDivide(const dimensioned< Type > &, const dimensioned< Type > &)
dimensioned< Type > sum(const DimensionedField< Type, GeoMesh > &df)
scalar y
#define EqOp(opName, op)
Definition: ops.H:46
Scalar cmptPow(const Scalar s1, const Scalar s2)
Definition: Scalar.H:293
dimensioned< Type > cmptMultiply(const dimensioned< Type > &, const dimensioned< Type > &)
dimensioned< scalar > magSqr(const dimensioned< Type > &)
dimensioned< Type > min(const dimensioned< Type > &, const dimensioned< Type > &)
Scalar minMod(const Scalar s1, const Scalar s2)
Definition: Scalar.H:227
tmp< DimensionedField< scalar, GeoMesh > > stabilise(const DimensionedField< scalar, GeoMesh > &dsf, const dimensioned< scalar > &ds)
Type maxMagSqr(const UList< Type > &f)
#define weightedOp(opName, op)
Definition: ops.H:136
dimensioned< scalar > mag(const dimensioned< Type > &)
Namespace for OpenFOAM.