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 | www.openfoam.com
6 \\/ M anipulation |
7-------------------------------------------------------------------------------
8 Copyright (C) 2011-2017 OpenFOAM Foundation
9 Copyright (C) 2019-2022 OpenCFD Ltd.
10-------------------------------------------------------------------------------
11License
12 This file is part of OpenFOAM.
13
14 OpenFOAM is free software: you can redistribute it and/or modify it
15 under the terms of the GNU General Public License as published by
16 the Free Software Foundation, either version 3 of the License, or
17 (at your option) any later version.
18
19 OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
20 ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
21 FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
22 for more details.
23
24 You should have received a copy of the GNU General Public License
25 along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
26
27\*---------------------------------------------------------------------------*/
28
29#include "dimensionedScalar.H"
30
31// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
32
33namespace Foam
34{
35
36// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
37
38dimensionedScalar operator+(const dimensionedScalar& ds1, const scalar s2)
39{
40 return ds1 + dimensionedScalar(s2);
41}
42
43
44dimensionedScalar operator+(const scalar s1, const dimensionedScalar& ds2)
45{
46 return dimensionedScalar(s1) + ds2;
47}
48
49
50dimensionedScalar operator-(const dimensionedScalar& ds1, const scalar s2)
51{
52 return ds1 - dimensionedScalar(s2);
53}
54
55
56dimensionedScalar operator-(const scalar s1, const dimensionedScalar& ds2)
57{
58 return dimensionedScalar(s1) - ds2;
59}
60
61
62dimensionedScalar operator*(const dimensionedScalar& ds1, const scalar s2)
63{
64 return ds1 * dimensionedScalar(s2);
65}
66
67
68dimensionedScalar operator/(const scalar s1, const dimensionedScalar& ds2)
69{
70 return dimensionedScalar(s1)/ds2;
71}
72
73
75(
76 const dimensionedScalar& ds,
77 const dimensionedScalar& expt
78)
79{
81 (
82 "pow(" + ds.name() + ',' + expt.name() + ')',
83 pow(ds.dimensions(), expt),
84 ::pow(ds.value(), expt.value())
85 );
86}
87
88
90{
92 (
93 "pow3(" + ds.name() + ')',
94 pow3(ds.dimensions()),
95 pow3(ds.value())
96 );
97}
98
99
101{
102 return dimensionedScalar
103 (
104 "pow4(" + ds.name() + ')',
105 pow4(ds.dimensions()),
106 pow4(ds.value())
107 );
108}
109
110
112{
113 return dimensionedScalar
114 (
115 "pow5(" + ds.name() + ')',
116 pow5(ds.dimensions()),
117 pow5(ds.value())
118 );
119}
120
121
123{
124 return dimensionedScalar
125 (
126 "pow6(" + ds.name() + ')',
127 pow6(ds.dimensions()),
128 pow6(ds.value())
129 );
130}
131
132
134{
135 return dimensionedScalar
136 (
137 "pow025(" + ds.name() + ')',
138 pow025(ds.dimensions()),
139 pow025(ds.value())
140 );
141}
142
143
145{
146 return dimensionedScalar
147 (
148 "sqrt(" + ds.name() + ')',
149 pow(ds.dimensions(), dimensionedScalar("0.5", dimless, 0.5)),
150 ::sqrt(ds.value())
151 );
152}
153
154
156{
157 return dimensionedScalar
158 (
159 "cbrt(" + ds.name() + ')',
160 pow(ds.dimensions(), dimensionedScalar("(1|3)", dimless, 1.0/3.0)),
161 ::cbrt(ds.value())
162 );
163}
164
165
167{
168 return dimensionedScalar
169 (
170 "sign(" + ds.name() + ')',
171 sign(ds.dimensions()),
172 ::Foam::sign(ds.value())
173 );
174}
175
176
178{
179 return dimensionedScalar
180 (
181 "pos(" + ds.name() + ')',
182 pos(ds.dimensions()),
183 ::Foam::pos(ds.value())
184 );
185}
186
187
189{
190 return dimensionedScalar
191 (
192 "pos0(" + ds.name() + ')',
193 pos0(ds.dimensions()),
194 ::Foam::pos0(ds.value())
195 );
196}
197
198
200{
201 return dimensionedScalar
202 (
203 "neg(" + ds.name() + ')',
204 neg(ds.dimensions()),
205 ::Foam::neg(ds.value())
206 );
207}
208
209
211{
212 return dimensionedScalar
213 (
214 "neg0(" + ds.name() + ')',
215 neg0(ds.dimensions()),
216 ::Foam::neg0(ds.value())
217 );
218}
219
220
222{
223 return dimensionedScalar
224 (
225 "posPart(" + ds.name() + ')',
226 posPart(ds.dimensions()),
227 ::Foam::pos0(ds.value())
228 );
229}
230
231
233{
234 return dimensionedScalar
235 (
236 "negPart(" + ds.name() + ')',
237 negPart(ds.dimensions()),
238 ::Foam::neg(ds.value())
239 );
240}
241
242
243#define transFunc(func) \
244dimensionedScalar func(const dimensionedScalar& ds) \
245{ \
246 if (!ds.dimensions().dimensionless()) \
247 { \
248 FatalErrorInFunction \
249 << "ds not dimensionless" \
250 << abort(FatalError); \
251 } \
252 \
253 return dimensionedScalar \
254 ( \
255 #func "(" + ds.name() + ')', \
256 dimless, \
257 ::func(ds.value()) \
258 ); \
259}
260
283
284#undef transFunc
285
286
287#define transFunc(func) \
288dimensionedScalar func(const int n, const dimensionedScalar& ds) \
289{ \
290 if (!ds.dimensions().dimensionless()) \
291 { \
292 FatalErrorInFunction \
293 << "ds not dimensionless" \
294 << abort(FatalError); \
295 } \
296 \
297 return dimensionedScalar \
298 ( \
299 #func "(" + name(n) + ',' + ds.name() + ')', \
300 dimless, \
301 ::func(n, ds.value()) \
302 ); \
303}
304
307
308#undef transFunc
309
310
312(
313 const dimensionedScalar& x,
314 const dimensionedScalar& y
315)
316{
317 return dimensionedScalar
318 (
319 "atan2(" + x.name() + ',' + y.name() + ')',
320 atan2(x.dimensions(), y.dimensions()),
321 ::atan2(x.value(), y.value())
322 );
323}
324
325
327(
328 const dimensionedScalar& x,
329 const dimensionedScalar& y
330)
331{
332 return dimensionedScalar
333 (
334 "hypot(" + x.name() + ',' + y.name() + ')',
335 hypot(x.dimensions(), y.dimensions()),
336 ::hypot(x.value(), y.value())
337 );
338}
339
340
342(
343 const dimensionedScalar& x,
344 const dimensionedScalar& y
345)
346{
347 return dimensionedScalar
348 (
349 "stabilise(" + x.name() + ',' + y.name() + ')',
350 stabilise(x.dimensions(), y.dimensions()),
351 stabilise(x.value(), y.value())
352 );
353}
354
355
356// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
357
358} // End namespace Foam
359
360// ************************************************************************* //
scalar y
const dimensionSet & dimensions() const
Return const reference to dimensions.
const Type & value() const
Return const reference to value.
const word & name() const
Return const reference to name.
#define transFunc(func)
Namespace for OpenFOAM.
tmp< faMatrix< Type > > operator-(const faMatrix< Type > &)
Unary negation.
dimensionedScalar pow6(const dimensionedScalar &ds)
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
dimensionedScalar pos(const dimensionedScalar &ds)
dimensionedScalar pow5(const dimensionedScalar &ds)
dimensionedScalar erfc(const dimensionedScalar &ds)
dimensionedScalar asin(const dimensionedScalar &ds)
dimensionedScalar exp(const dimensionedScalar &ds)
const dimensionSet dimless
Dimensionless.
dimensionedScalar tan(const dimensionedScalar &ds)
dimensionedScalar pos0(const dimensionedScalar &ds)
dimensionedScalar sign(const dimensionedScalar &ds)
dimensionedScalar lgamma(const dimensionedScalar &ds)
dimensionedScalar j1(const dimensionedScalar &ds)
tmp< faMatrix< Type > > operator+(const faMatrix< Type > &, const faMatrix< Type > &)
dimensionedScalar pow3(const dimensionedScalar &ds)
dimensionedScalar y0(const dimensionedScalar &ds)
dimensionedScalar cosh(const dimensionedScalar &ds)
dimensionedScalar sin(const dimensionedScalar &ds)
dimensionedScalar tanh(const dimensionedScalar &ds)
dimensionedScalar erf(const dimensionedScalar &ds)
dimensionedScalar sinh(const dimensionedScalar &ds)
dimensionedScalar log10(const dimensionedScalar &ds)
tmp< faMatrix< Type > > operator*(const areaScalarField::Internal &, const faMatrix< Type > &)
dimensionedScalar yn(const int n, const dimensionedScalar &ds)
dimensionedScalar log(const dimensionedScalar &ds)
dimensionedScalar pow(const dimensionedScalar &ds, const dimensionedScalar &expt)
dimensionedScalar operator/(const scalar s1, const dimensionedScalar &ds2)
dimensionedScalar atan2(const dimensionedScalar &x, const dimensionedScalar &y)
dimensionedScalar y1(const dimensionedScalar &ds)
dimensionedScalar negPart(const dimensionedScalar &ds)
dimensionedScalar acosh(const dimensionedScalar &ds)
dimensionedScalar sqrt(const dimensionedScalar &ds)
dimensionedScalar hypot(const dimensionedScalar &x, const dimensionedScalar &y)
dimensionedScalar pow4(const dimensionedScalar &ds)
dimensionedScalar jn(const int n, const dimensionedScalar &ds)
dimensionedScalar neg(const dimensionedScalar &ds)
dimensionedScalar atanh(const dimensionedScalar &ds)
dimensionedScalar stabilise(const dimensionedScalar &x, const dimensionedScalar &y)
dimensionedScalar neg0(const dimensionedScalar &ds)
dimensionedScalar cbrt(const dimensionedScalar &ds)
dimensionedScalar atan(const dimensionedScalar &ds)
dimensionedScalar cos(const dimensionedScalar &ds)
dimensionedScalar posPart(const dimensionedScalar &ds)
dimensionedScalar acos(const dimensionedScalar &ds)
dimensionedScalar j0(const dimensionedScalar &ds)
dimensionedScalar pow025(const dimensionedScalar &ds)
dimensionedScalar asinh(const dimensionedScalar &ds)