dimensionedScalar.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 "dimensionedScalar.H"
27 
28 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
29 
30 namespace Foam
31 {
32 
33 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
34 
35 dimensionedScalar operator+(const dimensionedScalar& ds1, const scalar s2)
36 {
37  return ds1 + dimensionedScalar(s2);
38 }
39 
40 
41 dimensionedScalar operator+(const scalar s1, const dimensionedScalar& ds2)
42 {
43  return dimensionedScalar(s1) + ds2;
44 }
45 
46 
47 dimensionedScalar operator-(const dimensionedScalar& ds1, const scalar s2)
48 {
49  return ds1 - dimensionedScalar(s2);
50 }
51 
52 
53 dimensionedScalar operator-(const scalar s1, const dimensionedScalar& ds2)
54 {
55  return dimensionedScalar(s1) - ds2;
56 }
57 
58 
59 dimensionedScalar operator*(const dimensionedScalar& ds1, const scalar s2)
60 {
61  return ds1 * dimensionedScalar(s2);
62 }
63 
64 
65 dimensionedScalar operator/(const scalar s1, const dimensionedScalar& ds2)
66 {
67  return dimensionedScalar(s1)/ds2;
68 }
69 
70 
71 
73 (
74  const dimensionedScalar& ds,
75  const dimensionedScalar& expt
76 )
77 {
78  return dimensionedScalar
79  (
80  "pow(" + ds.name() + ',' + expt.name() + ')',
81  pow(ds.dimensions(), expt),
82  ::pow(ds.value(), expt.value())
83  );
84 }
85 
86 
88 {
89  return dimensionedScalar
90  (
91  "pow3(" + ds.name() + ')',
92  pow3(ds.dimensions()),
93  pow3(ds.value())
94  );
95 }
96 
97 
99 {
100  return dimensionedScalar
101  (
102  "pow4(" + ds.name() + ')',
103  pow4(ds.dimensions()),
104  pow4(ds.value())
105  );
106 }
107 
108 
110 {
111  return dimensionedScalar
112  (
113  "pow5(" + ds.name() + ')',
114  pow5(ds.dimensions()),
115  pow5(ds.value())
116  );
117 }
118 
119 
121 {
122  return dimensionedScalar
123  (
124  "pow6(" + ds.name() + ')',
125  pow6(ds.dimensions()),
126  pow6(ds.value())
127  );
128 }
129 
130 
132 {
133  return dimensionedScalar
134  (
135  "pow025(" + ds.name() + ')',
136  pow025(ds.dimensions()),
137  pow025(ds.value())
138  );
139 }
140 
141 
143 {
144  return dimensionedScalar
145  (
146  "sqrt(" + ds.name() + ')',
147  pow(ds.dimensions(), dimensionedScalar("0.5", dimless, 0.5)),
148  ::sqrt(ds.value())
149  );
150 }
151 
152 
154 {
155  return dimensionedScalar
156  (
157  "cbrt(" + ds.name() + ')',
158  pow(ds.dimensions(), dimensionedScalar("(1|3)", dimless, 1.0/3.0)),
159  ::cbrt(ds.value())
160  );
161 }
162 
163 
165 (
166  const dimensionedScalar& x,
167  const dimensionedScalar& y
168 )
169 {
170  return dimensionedScalar
171  (
172  "hypot(" + x.name() + ',' + y.name() + ')',
173  x.dimensions() + y.dimensions(),
174  ::hypot(x.value(), y.value())
175  );
176 }
177 
178 
180 {
181  return dimensionedScalar
182  (
183  "sign(" + ds.name() + ')',
184  sign(ds.dimensions()),
185  ::Foam::sign(ds.value())
186  );
187 }
188 
189 
191 {
192  return dimensionedScalar
193  (
194  "pos(" + ds.name() + ')',
195  pos(ds.dimensions()),
196  ::Foam::pos(ds.value())
197  );
198 }
199 
200 
202 {
203  return dimensionedScalar
204  (
205  "pos0(" + ds.name() + ')',
206  pos0(ds.dimensions()),
207  ::Foam::pos0(ds.value())
208  );
209 }
210 
211 
213 {
214  return dimensionedScalar
215  (
216  "neg(" + ds.name() + ')',
217  neg(ds.dimensions()),
218  ::Foam::neg(ds.value())
219  );
220 }
221 
222 
224 {
225  return dimensionedScalar
226  (
227  "neg0(" + ds.name() + ')',
228  neg0(ds.dimensions()),
229  ::Foam::neg0(ds.value())
230  );
231 }
232 
233 
235 {
236  return dimensionedScalar
237  (
238  "posPart(" + ds.name() + ')',
239  posPart(ds.dimensions()),
240  ::Foam::pos0(ds.value())
241  );
242 }
243 
244 
246 {
247  return dimensionedScalar
248  (
249  "negPart(" + ds.name() + ')',
250  negPart(ds.dimensions()),
251  ::Foam::neg(ds.value())
252  );
253 }
254 
255 
256 #define transFunc(func) \
257 dimensionedScalar func(const dimensionedScalar& ds) \
258 { \
259  if (!ds.dimensions().dimensionless()) \
260  { \
261  FatalErrorInFunction \
262  << "ds not dimensionless" \
263  << abort(FatalError); \
264  } \
265  \
266  return dimensionedScalar \
267  ( \
268  #func "(" + ds.name() + ')', \
269  dimless, \
270  ::func(ds.value()) \
271  ); \
272 }
273 
296 
297 #undef transFunc
298 
299 
300 #define transFunc(func) \
301 dimensionedScalar func(const int n, const dimensionedScalar& ds) \
302 { \
303  if (!ds.dimensions().dimensionless()) \
304  { \
305  FatalErrorInFunction \
306  << "ds not dimensionless" \
307  << abort(FatalError); \
308  } \
309  \
310  return dimensionedScalar \
311  ( \
312  #func "(" + name(n) + ',' + ds.name() + ')', \
313  dimless, \
314  ::func(n, ds.value()) \
315  ); \
316 }
317 
320 
321 #undef transFunc
322 
323 
325 (
326  const dimensionedScalar& x,
327  const dimensionedScalar& y
328 )
329 {
330  return dimensionedScalar
331  (
332  "atan2(" + x.name() + ',' + y.name() + ')',
333  atan2(x.dimensions(), y.dimensions()),
334  ::atan2(x.value(), y.value())
335  );
336 }
337 
338 
339 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
340 
341 } // End namespace Foam
342 
343 // ************************************************************************* //
dimensionedScalar sign(const dimensionedScalar &ds)
dimensionedScalar tanh(const dimensionedScalar &ds)
dimensionedScalar acos(const dimensionedScalar &ds)
dimensionedScalar log(const dimensionedScalar &ds)
dimensionedScalar yn(const int n, const dimensionedScalar &ds)
dimensionedScalar sqrt(const dimensionedScalar &ds)
dimensionedScalar pow025(const dimensionedScalar &ds)
dimensionedScalar operator/(const scalar s1, const dimensionedScalar &ds2)
dimensionedScalar y0(const dimensionedScalar &ds)
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)
dimensionedScalar pos(const dimensionedScalar &ds)
dimensionedScalar acosh(const dimensionedScalar &ds)
tmp< fvMatrix< Type > > operator*(const volScalarField::Internal &, const fvMatrix< Type > &)
dimensionedScalar cos(const dimensionedScalar &ds)
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.
tmp< fvMatrix< Type > > operator-(const fvMatrix< Type > &)
dimensionedScalar atanh(const dimensionedScalar &ds)
tmp< fvMatrix< Type > > operator+(const fvMatrix< Type > &, const fvMatrix< Type > &)
dimensionedScalar y1(const dimensionedScalar &ds)
dimensionedScalar pos0(const dimensionedScalar &ds)
dimensionedScalar sin(const dimensionedScalar &ds)
dimensionedScalar erf(const dimensionedScalar &ds)
dimensionedScalar atan2(const dimensionedScalar &x, const dimensionedScalar &y)
const word & name() const
Return const reference to name.
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) transFunc(j0) transFunc(j1) transFunc(y0) transFunc(y1) inline Scalar &setComponent(Scalar &s
dimensionedScalar lgamma(const dimensionedScalar &ds)
dimensionedScalar pow(const dimensionedScalar &ds, const dimensionedScalar &expt)
dimensionedScalar pow3(const dimensionedScalar &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)
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)
dimensionedScalar tan(const dimensionedScalar &ds)
dimensionedScalar j0(const dimensionedScalar &ds)
dimensionedScalar hypot(const dimensionedScalar &x, const dimensionedScalar &y)
dimensionedScalar log10(const dimensionedScalar &ds)
Namespace for OpenFOAM.
dimensionedScalar negPart(const dimensionedScalar &ds)