dimensionedScalar.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) 2004-2010 OpenCFD Ltd.
6  \\/ M anipulation |
7 -------------------------------------------------------------------------------
8  | Copyright (C) 2011-2017 OpenFOAM Foundation
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 "dimensionedScalar.H"
29 
30 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
31 
32 namespace Foam
33 {
34 
35 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
36 
37 dimensionedScalar operator+(const dimensionedScalar& ds1, const scalar s2)
38 {
39  return ds1 + dimensionedScalar(s2);
40 }
41 
42 
43 dimensionedScalar operator+(const scalar s1, const dimensionedScalar& ds2)
44 {
45  return dimensionedScalar(s1) + ds2;
46 }
47 
48 
49 dimensionedScalar operator-(const dimensionedScalar& ds1, const scalar s2)
50 {
51  return ds1 - dimensionedScalar(s2);
52 }
53 
54 
55 dimensionedScalar operator-(const scalar s1, const dimensionedScalar& ds2)
56 {
57  return dimensionedScalar(s1) - ds2;
58 }
59 
60 
61 dimensionedScalar operator*(const dimensionedScalar& ds1, const scalar s2)
62 {
63  return ds1 * dimensionedScalar(s2);
64 }
65 
66 
67 dimensionedScalar operator/(const scalar s1, const dimensionedScalar& ds2)
68 {
69  return dimensionedScalar(s1)/ds2;
70 }
71 
72 
74 (
75  const dimensionedScalar& ds,
76  const dimensionedScalar& expt
77 )
78 {
79  return dimensionedScalar
80  (
81  "pow(" + ds.name() + ',' + expt.name() + ')',
82  pow(ds.dimensions(), expt),
83  ::pow(ds.value(), expt.value())
84  );
85 }
86 
87 
89 {
90  return dimensionedScalar
91  (
92  "pow3(" + ds.name() + ')',
93  pow3(ds.dimensions()),
94  pow3(ds.value())
95  );
96 }
97 
98 
100 {
101  return dimensionedScalar
102  (
103  "pow4(" + ds.name() + ')',
104  pow4(ds.dimensions()),
105  pow4(ds.value())
106  );
107 }
108 
109 
111 {
112  return dimensionedScalar
113  (
114  "pow5(" + ds.name() + ')',
115  pow5(ds.dimensions()),
116  pow5(ds.value())
117  );
118 }
119 
120 
122 {
123  return dimensionedScalar
124  (
125  "pow6(" + ds.name() + ')',
126  pow6(ds.dimensions()),
127  pow6(ds.value())
128  );
129 }
130 
131 
133 {
134  return dimensionedScalar
135  (
136  "pow025(" + ds.name() + ')',
137  pow025(ds.dimensions()),
138  pow025(ds.value())
139  );
140 }
141 
142 
144 {
145  return dimensionedScalar
146  (
147  "sqrt(" + ds.name() + ')',
148  pow(ds.dimensions(), dimensionedScalar("0.5", dimless, 0.5)),
149  ::sqrt(ds.value())
150  );
151 }
152 
153 
155 {
156  return dimensionedScalar
157  (
158  "cbrt(" + ds.name() + ')',
159  pow(ds.dimensions(), dimensionedScalar("(1|3)", dimless, 1.0/3.0)),
160  ::cbrt(ds.value())
161  );
162 }
163 
164 
166 (
167  const dimensionedScalar& x,
168  const dimensionedScalar& y
169 )
170 {
171  return dimensionedScalar
172  (
173  "hypot(" + x.name() + ',' + y.name() + ')',
174  x.dimensions() + y.dimensions(),
175  ::hypot(x.value(), y.value())
176  );
177 }
178 
179 
181 {
182  return dimensionedScalar
183  (
184  "sign(" + ds.name() + ')',
185  sign(ds.dimensions()),
186  ::Foam::sign(ds.value())
187  );
188 }
189 
190 
192 {
193  return dimensionedScalar
194  (
195  "pos(" + ds.name() + ')',
196  pos(ds.dimensions()),
197  ::Foam::pos(ds.value())
198  );
199 }
200 
201 
203 {
204  return dimensionedScalar
205  (
206  "pos0(" + ds.name() + ')',
207  pos0(ds.dimensions()),
208  ::Foam::pos0(ds.value())
209  );
210 }
211 
212 
214 {
215  return dimensionedScalar
216  (
217  "neg(" + ds.name() + ')',
218  neg(ds.dimensions()),
219  ::Foam::neg(ds.value())
220  );
221 }
222 
223 
225 {
226  return dimensionedScalar
227  (
228  "neg0(" + ds.name() + ')',
229  neg0(ds.dimensions()),
230  ::Foam::neg0(ds.value())
231  );
232 }
233 
234 
236 {
237  return dimensionedScalar
238  (
239  "posPart(" + ds.name() + ')',
240  posPart(ds.dimensions()),
241  ::Foam::pos0(ds.value())
242  );
243 }
244 
245 
247 {
248  return dimensionedScalar
249  (
250  "negPart(" + ds.name() + ')',
251  negPart(ds.dimensions()),
252  ::Foam::neg(ds.value())
253  );
254 }
255 
256 
257 #define transFunc(func) \
258 dimensionedScalar func(const dimensionedScalar& ds) \
259 { \
260  if (!ds.dimensions().dimensionless()) \
261  { \
262  FatalErrorInFunction \
263  << "ds not dimensionless" \
264  << abort(FatalError); \
265  } \
266  \
267  return dimensionedScalar \
268  ( \
269  #func "(" + ds.name() + ')', \
270  dimless, \
271  ::func(ds.value()) \
272  ); \
273 }
274 
297 
298 #undef transFunc
299 
300 
301 #define transFunc(func) \
302 dimensionedScalar func(const int n, const dimensionedScalar& ds) \
303 { \
304  if (!ds.dimensions().dimensionless()) \
305  { \
306  FatalErrorInFunction \
307  << "ds not dimensionless" \
308  << abort(FatalError); \
309  } \
310  \
311  return dimensionedScalar \
312  ( \
313  #func "(" + name(n) + ',' + ds.name() + ')', \
314  dimless, \
315  ::func(n, ds.value()) \
316  ); \
317 }
318 
321 
322 #undef transFunc
323 
324 
326 (
327  const dimensionedScalar& x,
328  const dimensionedScalar& y
329 )
330 {
331  return dimensionedScalar
332  (
333  "atan2(" + x.name() + ',' + y.name() + ')',
334  atan2(x.dimensions(), y.dimensions()),
335  ::atan2(x.value(), y.value())
336  );
337 }
338 
339 
340 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
341 
342 } // End namespace Foam
343 
344 // ************************************************************************* //
Foam::dimless
const dimensionSet dimless(0, 0, 0, 0, 0, 0, 0)
Dimensionless.
Definition: dimensionSets.H:50
Foam::tan
dimensionedScalar tan(const dimensionedScalar &ds)
Definition: dimensionedScalar.C:280
Foam::cosh
dimensionedScalar cosh(const dimensionedScalar &ds)
Definition: dimensionedScalar.C:285
Foam::y1
dimensionedScalar y1(const dimensionedScalar &ds)
Definition: dimensionedScalar.C:296
Foam::sin
dimensionedScalar sin(const dimensionedScalar &ds)
Definition: dimensionedScalar.C:278
Foam::jn
dimensionedScalar jn(const int n, const dimensionedScalar &ds)
Definition: dimensionedScalar.C:319
Foam::posPart
dimensionedScalar posPart(const dimensionedScalar &ds)
Definition: dimensionedScalar.C:235
Foam::atan2
dimensionedScalar atan2(const dimensionedScalar &x, const dimensionedScalar &y)
Definition: dimensionedScalar.C:326
Foam::operator-
tmp< faMatrix< Type > > operator-(const faMatrix< Type > &)
Foam::neg0
dimensionedScalar neg0(const dimensionedScalar &ds)
Definition: dimensionedScalar.C:224
Foam::dimensioned::name
const word & name() const
Return const reference to name.
Definition: dimensionedType.C:336
Foam::dimensioned::value
const Type & value() const
Return const reference to value.
Definition: dimensionedType.C:364
Foam::pos0
dimensionedScalar pos0(const dimensionedScalar &ds)
Definition: dimensionedScalar.C:202
Foam::exp
dimensionedScalar exp(const dimensionedScalar &ds)
Definition: dimensionedScalar.C:275
Foam::sign
dimensionedScalar sign(const dimensionedScalar &ds)
Definition: dimensionedScalar.C:180
Foam::erf
dimensionedScalar erf(const dimensionedScalar &ds)
Definition: dimensionedScalar.C:290
Foam::pow025
dimensionedScalar pow025(const dimensionedScalar &ds)
Definition: dimensionedScalar.C:132
Foam::hypot
dimensionedScalar hypot(const dimensionedScalar &x, const dimensionedScalar &y)
Definition: dimensionedScalar.C:166
Foam::atanh
dimensionedScalar atanh(const dimensionedScalar &ds)
Definition: dimensionedScalar.C:289
Foam::lgamma
dimensionedScalar lgamma(const dimensionedScalar &ds)
Definition: dimensionedScalar.C:292
Foam::pow4
dimensionedScalar pow4(const dimensionedScalar &ds)
Definition: dimensionedScalar.C:99
Foam::pow6
dimensionedScalar pow6(const dimensionedScalar &ds)
Definition: dimensionedScalar.C:121
Foam::tanh
dimensionedScalar tanh(const dimensionedScalar &ds)
Definition: dimensionedScalar.C:286
Foam::pow3
dimensionedScalar pow3(const dimensionedScalar &ds)
Definition: dimensionedScalar.C:88
Foam::log10
dimensionedScalar log10(const dimensionedScalar &ds)
Definition: dimensionedScalar.C:277
Foam::y0
dimensionedScalar y0(const dimensionedScalar &ds)
Definition: dimensionedScalar.C:295
Foam::erfc
dimensionedScalar erfc(const dimensionedScalar &ds)
Definition: dimensionedScalar.C:291
Foam::dimensionedScalar
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
Definition: dimensionedScalarFwd.H:43
Foam::asinh
dimensionedScalar asinh(const dimensionedScalar &ds)
Definition: dimensionedScalar.C:287
Foam::pow
dimensionedScalar pow(const dimensionedScalar &ds, const dimensionedScalar &expt)
Definition: dimensionedScalar.C:74
Foam::pow5
dimensionedScalar pow5(const dimensionedScalar &ds)
Definition: dimensionedScalar.C:110
Foam::log
dimensionedScalar log(const dimensionedScalar &ds)
Definition: dimensionedScalar.C:276
Foam::dimensioned
Generic dimensioned Type class.
Definition: dimensionedScalarFwd.H:43
Foam
Namespace for OpenFOAM.
Definition: atmBoundaryLayer.C:32
Foam::operator/
dimensionedScalar operator/(const scalar s1, const dimensionedScalar &ds2)
Definition: dimensionedScalar.C:67
Foam::yn
dimensionedScalar yn(const int n, const dimensionedScalar &ds)
Definition: dimensionedScalar.C:320
dimensionedScalar.H
Foam::acosh
dimensionedScalar acosh(const dimensionedScalar &ds)
Definition: dimensionedScalar.C:288
Foam::negPart
dimensionedScalar negPart(const dimensionedScalar &ds)
Definition: dimensionedScalar.C:246
Foam::sqrt
dimensionedScalar sqrt(const dimensionedScalar &ds)
Definition: dimensionedScalar.C:143
Foam::j0
dimensionedScalar j0(const dimensionedScalar &ds)
Definition: dimensionedScalar.C:293
Foam::acos
dimensionedScalar acos(const dimensionedScalar &ds)
Definition: dimensionedScalar.C:282
Foam::operator+
tmp< faMatrix< Type > > operator+(const faMatrix< Type > &, const faMatrix< Type > &)
Foam::operator*
tmp< faMatrix< Type > > operator*(const areaScalarField &, const faMatrix< Type > &)
x
x
Definition: LISASMDCalcMethod2.H:52
Foam::transFunc
transFunc(sqrt) transFunc(cbrt) transFunc(exp) transFunc(log) transFunc(log10) transFunc(sin) transFunc(cos) transFunc(tan) transFunc(asin) transFunc(acos) transFunc(atan) transFunc(sinh) transFunc(cosh) transFunc(tanh) transFunc(asinh) transFunc(acosh) transFunc(atanh) transFunc(erf) transFunc(erfc) transFunc(lgamma) transFunc(tgamma) besselFunc(j0) besselFunc(j1) besselFunc(y0) besselFunc(y1) besselFunc2(jn) besselFunc2(yn) inline Scalar &setComponent(Scalar &s
Foam::atan
dimensionedScalar atan(const dimensionedScalar &ds)
Definition: dimensionedScalar.C:283
Foam::dimensioned::dimensions
const dimensionSet & dimensions() const
Return const reference to dimensions.
Definition: dimensionedType.C:350
Foam::cbrt
dimensionedScalar cbrt(const dimensionedScalar &ds)
Definition: dimensionedScalar.C:154
Foam::j1
dimensionedScalar j1(const dimensionedScalar &ds)
Definition: dimensionedScalar.C:294
Foam::neg
dimensionedScalar neg(const dimensionedScalar &ds)
Definition: dimensionedScalar.C:213
Foam::asin
dimensionedScalar asin(const dimensionedScalar &ds)
Definition: dimensionedScalar.C:281
y
scalar y
Definition: LISASMDCalcMethod1.H:14
Foam::cos
dimensionedScalar cos(const dimensionedScalar &ds)
Definition: dimensionedScalar.C:279
Foam::pos
dimensionedScalar pos(const dimensionedScalar &ds)
Definition: dimensionedScalar.C:191
Foam::sinh
dimensionedScalar sinh(const dimensionedScalar &ds)
Definition: dimensionedScalar.C:284