FreeFOAM The Cross-Platform CFD Toolkit
dimensionSet.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 
26 #include "dimensionSet.H"
28 
29 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
30 
32 const Foam::scalar Foam::dimensionSet::smallExponent = SMALL;
33 
34 
35 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
36 
38 (
39  const scalar mass,
40  const scalar length,
41  const scalar time,
42  const scalar temperature,
43  const scalar moles,
44  const scalar current,
45  const scalar luminousIntensity
46 )
47 {
48  exponents_[MASS] = mass;
49  exponents_[LENGTH] = length;
50  exponents_[TIME] = time;
51  exponents_[TEMPERATURE] = temperature;
52  exponents_[MOLES] = moles;
53  exponents_[CURRENT] = current;
54  exponents_[LUMINOUS_INTENSITY] = luminousIntensity;
55 }
56 
57 
59 (
60  const scalar mass,
61  const scalar length,
62  const scalar time,
63  const scalar temperature,
64  const scalar moles
65 )
66 {
67  exponents_[MASS] = mass;
68  exponents_[LENGTH] = length;
69  exponents_[TIME] = time;
70  exponents_[TEMPERATURE] = temperature;
71  exponents_[MOLES] = moles;
72  exponents_[CURRENT] = 0;
73  exponents_[LUMINOUS_INTENSITY] = 0;
74 }
75 
76 
77 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
78 
80 {
81  bool Dimensionless = true;
82 
83  for (int Dimension=0; Dimension<nDimensions; Dimension++)
84  {
85  Dimensionless = Dimensionless &&
86  (
87  exponents_[Dimension] < smallExponent
88  && exponents_[Dimension] > -smallExponent
89  );
90  }
91 
92  return Dimensionless;
93 }
94 
95 
97 {
98  for (int Dimension=0; Dimension<nDimensions; Dimension++)
99  {
100  exponents_[Dimension] = ds.exponents_[Dimension];
101  }
102 }
103 
104 
105 // * * * * * * * * * * * * * * * Member Operators * * * * * * * * * * * * * //
106 
108 {
109  return exponents_[type];
110 }
111 
113 {
114  return exponents_[type];
115 }
116 
117 
119 {
120  bool equall = true;
121 
122  for (int Dimension=0; Dimension<nDimensions; Dimension++)
123  {
124  equall = equall &&
125  (mag(exponents_[Dimension] - ds.exponents_[Dimension])
126  < smallExponent);
127  }
128 
129  return equall;
130 }
131 
133 {
134  return !(operator==(ds));
135 }
136 
137 
139 {
140  if (dimensionSet::debug && *this != ds)
141  {
142  FatalErrorIn("dimensionSet::operator=(const dimensionSet& ds) const")
143  << "Different dimensions for =" << endl
144  << " dimensions : " << *this << " = " << ds << endl
145  << abort(FatalError);
146  }
147 
148  return true;
149 }
150 
151 
153 {
154  if (dimensionSet::debug && *this != ds)
155  {
156  FatalErrorIn("dimensionSet::operator+=(const dimensionSet& ds) const")
157  << "Different dimensions for +=" << endl
158  << " dimensions : " << *this << " = " << ds << endl
159  << abort(FatalError);
160  }
161 
162  return true;
163 }
164 
166 {
167  if (dimensionSet::debug && *this != ds)
168  {
169  FatalErrorIn("dimensionSet::operator-=(const dimensionSet& ds) const")
170  << "Different dimensions for -=" << endl
171  << " dimensions : " << *this << " = " << ds << endl
172  << abort(FatalError);
173  }
174 
175  return true;
176 }
177 
179 {
180  reset((*this)*ds);
181 
182  return true;
183 }
184 
186 {
187  reset((*this)/ds);
188 
189  return true;
190 }
191 
192 
193 // * * * * * * * * * * * * * * * Friend functions * * * * * * * * * * * * * * //
194 
196 {
197  if (dimensionSet::debug && ds1 != ds2)
198  {
199  FatalErrorIn("max(const dimensionSet& ds1, const dimensionSet& ds2)")
200  << "Arguments of max have different dimensions" << endl
201  << " dimensions : " << ds1 << " and " << ds2 << endl
202  << abort(FatalError);
203  }
204 
205  return ds1;
206 }
207 
209 {
210  if (dimensionSet::debug && ds1 != ds2)
211  {
212  FatalErrorIn("min(const dimensionSet& ds1, const dimensionSet& ds2)")
213  << "Arguments of min have different dimensions" << endl
214  << " dimensions : " << ds1 << " and " << ds2 << endl
215  << abort(FatalError);
216  }
217 
218  return ds1;
219 }
220 
221 
223 (
224  const dimensionSet& ds1,
225  const dimensionSet& ds2
226 )
227 {
228  return ds1*ds2;
229 }
230 
231 
233 (
234  const dimensionSet& ds1,
235  const dimensionSet& ds2
236 )
237 {
238  return ds1/ds2;
239 }
240 
241 
242 Foam::dimensionSet Foam::pow(const dimensionSet& ds, const scalar p)
243 {
244  dimensionSet dimPow
245  (
246  ds[dimensionSet::MASS]*p,
247  ds[dimensionSet::LENGTH]*p,
248  ds[dimensionSet::TIME]*p,
250  ds[dimensionSet::MOLES]*p,
251  ds[dimensionSet::CURRENT]*p,
253  );
254 
255  return dimPow;
256 }
257 
259 (
260  const dimensionSet& ds,
261  const dimensionedScalar& dS
262 )
263 {
264  if (dimensionSet::debug && !dS.dimensions().dimensionless())
265  {
266  FatalErrorIn("pow(const dimensionSet& ds, const dimensionedScalar& dS)")
267  << "Exponent of pow are not dimensionless"
268  << abort(FatalError);
269  }
270 
271  dimensionSet dimPow
272  (
273  ds[dimensionSet::MASS]*dS.value(),
274  ds[dimensionSet::LENGTH]*dS.value(),
275  ds[dimensionSet::TIME]*dS.value(),
277  ds[dimensionSet::MOLES]*dS.value(),
278  ds[dimensionSet::CURRENT]*dS.value(),
280  );
281 
282  return dimPow;
283 }
284 
286 (
287  const dimensionedScalar& dS,
288  const dimensionSet& ds
289 )
290 {
291  if
292  (
293  dimensionSet::debug
294  && !dS.dimensions().dimensionless()
295  && !ds.dimensionless())
296  {
297  FatalErrorIn("pow(const dimensionedScalar& dS, const dimensionSet& ds)")
298  << "Argument or exponent of pow not dimensionless" << endl
299  << abort(FatalError);
300  }
301 
302  return ds;
303 }
304 
305 
307 {
308  return pow(ds, 2);
309 }
310 
312 {
313  return pow(ds, 3);
314 }
315 
317 {
318  return pow(ds, 4);
319 }
320 
322 {
323  return pow(ds, 5);
324 }
325 
327 {
328  return pow(ds, 6);
329 }
330 
332 {
333  return pow(ds, 0.5);
334 }
335 
337 {
338  return pow(ds, 2);
339 }
340 
342 {
343  return ds;
344 }
345 
347 {
348  return dimless;
349 }
350 
352 {
353  return dimless;
354 }
355 
357 {
358  return dimless;
359 }
360 
362 {
363  return dimless/ds;
364 }
365 
367 {
368  if (dimensionSet::debug && !ds.dimensionless())
369  {
370  FatalErrorIn("trans(const dimensionSet& ds)")
371  << "Argument of trancendental function not dimensionless"
372  << abort(FatalError);
373  }
374 
375  return ds;
376 }
377 
379 {
380  return ds;
381 }
382 
383 
384 // * * * * * * * * * * * * * * * Friend Operators * * * * * * * * * * * * * //
385 
387 {
388  return ds;
389 }
390 
391 Foam::dimensionSet Foam::operator+
392 (
393  const dimensionSet& ds1,
394  const dimensionSet& ds2
395 )
396 {
397  dimensionSet dimSum(ds1);
398 
399  if (dimensionSet::debug && ds1 != ds2)
400  {
402  ("operator+(const dimensionSet& ds1, const dimensionSet& ds2)")
403  << "LHS and RHS of + have different dimensions" << endl
404  << " dimensions : " << ds1 << " + " << ds2 << endl
405  << abort(FatalError);
406  }
407 
408  return dimSum;
409 }
410 
411 Foam::dimensionSet Foam::operator-
412 (
413  const dimensionSet& ds1,
414  const dimensionSet& ds2
415 )
416 {
417  dimensionSet dimDifference(ds1);
418 
419  if (dimensionSet::debug && ds1 != ds2)
420  {
422  ("operator-(const dimensionSet& ds1, const dimensionSet& ds2)")
423  << "LHS and RHS of - have different dimensions" << endl
424  << " dimensions : " << ds1 << " - " << ds2 << endl
425  << abort(FatalError);
426  }
427 
428  return dimDifference;
429 }
430 
431 Foam::dimensionSet Foam::operator*
432 (
433  const dimensionSet& ds1,
434  const dimensionSet& ds2
435 )
436 {
437  dimensionSet dimProduct(ds1);
438 
439  for (int Dimension=0; Dimension<dimensionSet::nDimensions; Dimension++)
440  {
441  dimProduct.exponents_[Dimension] += ds2.exponents_[Dimension];
442  }
443 
444  return dimProduct;
445 }
446 
447 Foam::dimensionSet Foam::operator/
448 (
449  const dimensionSet& ds1,
450  const dimensionSet& ds2
451 )
452 {
453  dimensionSet dimQuotient(ds1);
454 
455  for (int Dimension=0; Dimension<dimensionSet::nDimensions; Dimension++)
456  {
457  dimQuotient.exponents_[Dimension] -= ds2.exponents_[Dimension];
458  }
459 
460  return dimQuotient;
461 }
462 
463 
464 Foam::dimensionSet Foam::operator&
465 (
466  const dimensionSet& ds1,
467  const dimensionSet& ds2
468 )
469 {
470  return ds1*ds2;
471 }
472 
473 Foam::dimensionSet Foam::operator^
474 (
475  const dimensionSet& ds1,
476  const dimensionSet& ds2
477 )
478 {
479  return ds1*ds2;
480 }
481 
482 Foam::dimensionSet Foam::operator&&
483 (
484  const dimensionSet& ds1,
485  const dimensionSet& ds2
486 )
487 {
488  return ds1*ds2;
489 }
490 
491 
492 // ************************ vim: set sw=4 sts=4 et: ************************ //