FreeFOAM The Cross-Platform CFD Toolkit
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) 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 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
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 dimensionedScalar operator+(const scalar s1, const dimensionedScalar& ds2)
41 {
42  return dimensionedScalar(s1) + ds2;
43 }
44 
45 dimensionedScalar operator-(const dimensionedScalar& ds1, const scalar s2)
46 {
47  return ds1 - dimensionedScalar(s2);
48 }
49 
50 dimensionedScalar operator-(const scalar s1, const dimensionedScalar& ds2)
51 {
52  return dimensionedScalar(s1) - ds2;
53 }
54 
55 dimensionedScalar operator*(const dimensionedScalar& ds1, const scalar s2)
56 {
57  return ds1 * dimensionedScalar(s2);
58 }
59 
60 dimensionedScalar operator/(const scalar s1, const dimensionedScalar& ds2)
61 {
62  return dimensionedScalar(s1)/ds2;
63 }
64 
65 
67 (
68  const dimensionedScalar& ds,
69  const dimensionedScalar& expt
70 )
71 {
72  return dimensionedScalar
73  (
74  "pow(" + ds.name() + ',' + expt.name() + ')',
75  pow(ds.dimensions(), expt),
76  ::pow(ds.value(), expt.value())
77  );
78 }
79 
81 {
82  return dimensionedScalar
83  (
84  "pow3(" + ds.name() + ')',
85  pow3(ds.dimensions()),
86  pow3(ds.value())
87  );
88 }
89 
91 {
92  return dimensionedScalar
93  (
94  "pow4(" + ds.name() + ')',
95  pow4(ds.dimensions()),
96  pow4(ds.value())
97  );
98 }
99 
101 {
102  return dimensionedScalar
103  (
104  "pow5(" + ds.name() + ')',
105  pow5(ds.dimensions()),
106  pow5(ds.value())
107  );
108 }
109 
111 {
112  return dimensionedScalar
113  (
114  "pow6(" + ds.name() + ')',
115  pow6(ds.dimensions()),
116  pow6(ds.value())
117  );
118 }
119 
121 {
122  return dimensionedScalar
123  (
124  "sqrt(" + ds.name() + ')',
125  pow(ds.dimensions(), dimensionedScalar("0.5", dimless, 0.5)),
126  ::sqrt(ds.value())
127  );
128 }
129 
131 {
132  return dimensionedScalar
133  (
134  "cbrt(" + ds.name() + ')',
135  pow(ds.dimensions(), dimensionedScalar("(1|3)", dimless, 1.0/3.0)),
136  ::cbrt(ds.value())
137  );
138 }
139 
141 (
142  const dimensionedScalar& x,
143  const dimensionedScalar& y
144 )
145 {
146  return dimensionedScalar
147  (
148  "hypot(" + x.name() + ',' + y.name() + ')',
149  x.dimensions() + y.dimensions(),
150  ::hypot(x.value(), y.value())
151  );
152 }
153 
155 {
156  return dimensionedScalar
157  (
158  "sign(" + ds.name() + ')',
159  sign(ds.dimensions()),
160  ::Foam::sign(ds.value())
161  );
162 }
163 
165 {
166  return dimensionedScalar
167  (
168  "pos(" + ds.name() + ')',
169  pos(ds.dimensions()),
170  ::Foam::pos(ds.value())
171  );
172 }
173 
175 {
176  return dimensionedScalar
177  (
178  "neg(" + ds.name() + ')',
179  neg(ds.dimensions()),
180  ::Foam::neg(ds.value())
181  );
182 }
183 
184 
185 #define transFunc(func) \
186 dimensionedScalar func(const dimensionedScalar& ds) \
187 { \
188  if (!ds.dimensions().dimensionless()) \
189  { \
190  FatalErrorIn(#func "(const dimensionedScalar& ds)") \
191  << "ds not dimensionless" \
192  << abort(FatalError); \
193  } \
194  \
195  return dimensionedScalar \
196  ( \
197  #func "(" + ds.name() + ')', \
198  dimless, \
199  ::func(ds.value()) \
200  ); \
201 }
202 
225 
226 #undef transFunc
227 
228 
229 #define transFunc(func) \
230 dimensionedScalar func(const int n, const dimensionedScalar& ds) \
231 { \
232  if (!ds.dimensions().dimensionless()) \
233  { \
234  FatalErrorIn(#func "(const int n, const dimensionedScalar& ds)") \
235  << "ds not dimensionless" \
236  << abort(FatalError); \
237  } \
238  \
239  return dimensionedScalar \
240  ( \
241  #func "(" + name(n) + ',' + ds.name() + ')', \
242  dimless, \
243  ::func(n, ds.value()) \
244  ); \
245 }
246 
249 
250 #undef transFunc
251 
252 
253 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
254 
255 } // End namespace Foam
256 
257 // ************************ vim: set sw=4 sts=4 et: ************************ //