dimensionSet.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 "dimensionSet.H"
27 #include "dimensionedScalar.H"
28 #include "OStringStream.H"
29 
30 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
31 
32 namespace Foam
33 {
34  defineTypeNameAndDebug(dimensionSet, 1);
35  const scalar dimensionSet::smallExponent = small;
36 }
37 
38 
39 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
40 
42 (
43  const scalar mass,
44  const scalar length,
45  const scalar time,
46  const scalar temperature,
47  const scalar moles,
48  const scalar current,
49  const scalar luminousIntensity
50 )
51 {
52  exponents_[MASS] = mass;
53  exponents_[LENGTH] = length;
54  exponents_[TIME] = time;
55  exponents_[TEMPERATURE] = temperature;
56  exponents_[MOLES] = moles;
57  exponents_[CURRENT] = current;
58  exponents_[LUMINOUS_INTENSITY] = luminousIntensity;
59 }
60 
61 
63 (
64  const scalar mass,
65  const scalar length,
66  const scalar time,
67  const scalar temperature,
68  const scalar moles
69 )
70 {
71  exponents_[MASS] = mass;
72  exponents_[LENGTH] = length;
73  exponents_[TIME] = time;
74  exponents_[TEMPERATURE] = temperature;
75  exponents_[MOLES] = moles;
76  exponents_[CURRENT] = 0;
77  exponents_[LUMINOUS_INTENSITY] = 0;
78 }
79 
80 
82 {
83  reset(ds);
84 }
85 
86 
87 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
88 
90 {
91  for (int Dimension=0; Dimension<nDimensions; ++Dimension)
92  {
93  // ie, mag(exponents_[Dimension]) > smallExponent
94  if
95  (
96  exponents_[Dimension] > smallExponent
97  || exponents_[Dimension] < -smallExponent
98  )
99  {
100  return false;
101  }
102  }
103 
104  return true;
105 }
106 
107 
109 {
110  for (int Dimension=0; Dimension<nDimensions; ++Dimension)
111  {
112  exponents_[Dimension] = ds.exponents_[Dimension];
113  }
114 }
115 
116 
117 // * * * * * * * * * * * * * * * Member Operators * * * * * * * * * * * * * //
118 
120 {
121  return exponents_[type];
122 }
123 
124 
126 {
127  return exponents_[type];
128 }
129 
130 
131 Foam::scalar Foam::dimensionSet::operator[](const label type) const
132 {
133  return exponents_[type];
134 }
135 
136 
138 {
139  return exponents_[type];
140 }
141 
142 
144 {
145  for (int Dimension=0; Dimension < nDimensions; ++Dimension)
146  {
147  if
148  (
149  mag(exponents_[Dimension] - ds.exponents_[Dimension])
150  > smallExponent
151  )
152  {
153  return false;
154  }
155  }
156 
157  return true;
158 }
159 
160 
162 {
163  return !(operator==(ds));
164 }
165 
166 
168 {
169  if (dimensionSet::debug && *this != ds)
170  {
172  << "Different dimensions for =" << endl
173  << " dimensions : " << *this << " = " << ds << endl
174  << abort(FatalError);
175  }
176 
177  return true;
178 }
179 
180 
182 {
183  if (dimensionSet::debug && *this != ds)
184  {
186  << "Different dimensions for +=" << endl
187  << " dimensions : " << *this << " = " << ds << endl
188  << abort(FatalError);
189  }
190 
191  return true;
192 }
193 
194 
196 {
197  if (dimensionSet::debug && *this != ds)
198  {
200  << "Different dimensions for -=" << endl
201  << " dimensions : " << *this << " = " << ds << endl
202  << abort(FatalError);
203  }
204 
205  return true;
206 }
207 
208 
210 {
211  reset((*this)*ds);
212 
213  return true;
214 }
215 
216 
218 {
219  reset((*this)/ds);
220 
221  return true;
222 }
223 
224 
225 // * * * * * * * * * * * * * * * Friend functions * * * * * * * * * * * * * * //
226 
228 {
229  if (dimensionSet::debug && ds1 != ds2)
230  {
232  << "Arguments of max have different dimensions" << endl
233  << " dimensions : " << ds1 << " and " << ds2 << endl
234  << abort(FatalError);
235  }
236 
237  return ds1;
238 }
239 
240 
242 {
243  if (dimensionSet::debug && ds1 != ds2)
244  {
246  << "Arguments of min have different dimensions" << endl
247  << " dimensions : " << ds1 << " and " << ds2 << endl
248  << abort(FatalError);
249  }
250 
251  return ds1;
252 }
253 
254 
256 (
257  const dimensionSet& ds1,
258  const dimensionSet& ds2
259 )
260 {
261  return ds1*ds2;
262 }
263 
264 
266 (
267  const dimensionSet& ds1,
268  const dimensionSet& ds2
269 )
270 {
271  return ds1/ds2;
272 }
273 
274 
275 Foam::dimensionSet Foam::pow(const dimensionSet& ds, const scalar p)
276 {
277  dimensionSet dimPow
278  (
279  ds[dimensionSet::MASS]*p,
280  ds[dimensionSet::LENGTH]*p,
281  ds[dimensionSet::TIME]*p,
283  ds[dimensionSet::MOLES]*p,
284  ds[dimensionSet::CURRENT]*p,
286  );
287 
288  return dimPow;
289 }
290 
291 
293 (
294  const dimensionSet& ds,
295  const dimensionedScalar& dS
296 )
297 {
298  if (dimensionSet::debug && !dS.dimensions().dimensionless())
299  {
301  << "Exponent of pow is not dimensionless"
302  << abort(FatalError);
303  }
304 
305  dimensionSet dimPow
306  (
307  ds[dimensionSet::MASS]*dS.value(),
308  ds[dimensionSet::LENGTH]*dS.value(),
309  ds[dimensionSet::TIME]*dS.value(),
311  ds[dimensionSet::MOLES]*dS.value(),
312  ds[dimensionSet::CURRENT]*dS.value(),
314  );
315 
316  return dimPow;
317 }
318 
319 
321 (
322  const dimensionedScalar& dS,
323  const dimensionSet& ds
324 )
325 {
326  if
327  (
328  dimensionSet::debug
329  && !dS.dimensions().dimensionless()
330  && !ds.dimensionless()
331  )
332  {
334  << "Argument or exponent of pow not dimensionless" << endl
335  << abort(FatalError);
336  }
337 
338  return ds;
339 }
340 
341 
343 {
344  return pow(ds, 2);
345 }
346 
347 
349 {
350  return pow(ds, 3);
351 }
352 
353 
355 {
356  return pow(ds, 4);
357 }
358 
359 
361 {
362  return pow(ds, 5);
363 }
364 
365 
367 {
368  return pow(ds, 6);
369 }
370 
371 
373 {
374  return sqrt(sqrt(ds));
375 }
376 
377 
379 {
380  return pow(ds, 0.5);
381 }
382 
383 
385 {
386  return pow(ds, 1.0/3.0);
387 }
388 
389 
391 {
392  return pow(ds, 2);
393 }
394 
395 
397 {
398  return ds;
399 }
400 
401 
403 {
404  return dimless;
405 }
406 
407 
409 {
410  return dimless;
411 }
412 
413 
415 {
416  return dimless;
417 }
418 
419 
421 {
422  return dimless;
423 }
424 
425 
427 {
428  return dimless;
429 }
430 
431 
433 {
434  return ds;
435 }
436 
437 
439 {
440  return ds;
441 }
442 
443 
445 {
446  return dimless/ds;
447 }
448 
449 
451 {
452  if (dimensionSet::debug && !ds.dimensionless())
453  {
455  << "Argument of trancendental function not dimensionless"
456  << abort(FatalError);
457  }
458 
459  return ds;
460 }
461 
462 
464 {
465  if (dimensionSet::debug && ds1 != ds2)
466  {
468  << "Arguments of atan2 have different dimensions" << endl
469  << " dimensions : " << ds1 << " and " << ds2 << endl
470  << abort(FatalError);
471  }
472 
473  return dimless;
474 }
475 
476 
478 {
479  return ds;
480 }
481 
482 
483 // * * * * * * * * * * * * * * * Friend Operators * * * * * * * * * * * * * //
484 
486 {
487  return ds;
488 }
489 
490 
491 Foam::dimensionSet Foam::operator+
492 (
493  const dimensionSet& ds1,
494  const dimensionSet& ds2
495 )
496 {
497  dimensionSet dimSum(ds1);
498 
499  if (dimensionSet::debug && ds1 != ds2)
500  {
502  << "LHS and RHS of + have different dimensions" << endl
503  << " dimensions : " << ds1 << " + " << ds2 << endl
504  << abort(FatalError);
505  }
506 
507  return dimSum;
508 }
509 
510 
511 Foam::dimensionSet Foam::operator-
512 (
513  const dimensionSet& ds1,
514  const dimensionSet& ds2
515 )
516 {
517  dimensionSet dimDifference(ds1);
518 
519  if (dimensionSet::debug && ds1 != ds2)
520  {
522  << "LHS and RHS of - have different dimensions" << endl
523  << " dimensions : " << ds1 << " - " << ds2 << endl
524  << abort(FatalError);
525  }
526 
527  return dimDifference;
528 }
529 
530 
531 Foam::dimensionSet Foam::operator*
532 (
533  const dimensionSet& ds1,
534  const dimensionSet& ds2
535 )
536 {
537  dimensionSet dimProduct(ds1);
538 
539  for (int Dimension=0; Dimension<dimensionSet::nDimensions; Dimension++)
540  {
541  dimProduct.exponents_[Dimension] += ds2.exponents_[Dimension];
542  }
543 
544  return dimProduct;
545 }
546 
547 
548 Foam::dimensionSet Foam::operator/
549 (
550  const dimensionSet& ds1,
551  const dimensionSet& ds2
552 )
553 {
554  dimensionSet dimQuotient(ds1);
555 
556  for (int Dimension=0; Dimension<dimensionSet::nDimensions; Dimension++)
557  {
558  dimQuotient.exponents_[Dimension] -= ds2.exponents_[Dimension];
559  }
560 
561  return dimQuotient;
562 }
563 
564 
565 Foam::dimensionSet Foam::operator&
566 (
567  const dimensionSet& ds1,
568  const dimensionSet& ds2
569 )
570 {
571  return ds1*ds2;
572 }
573 
574 
575 Foam::dimensionSet Foam::operator^
576 (
577  const dimensionSet& ds1,
578  const dimensionSet& ds2
579 )
580 {
581  return ds1*ds2;
582 }
583 
584 
585 Foam::dimensionSet Foam::operator&&
586 (
587  const dimensionSet& ds1,
588  const dimensionSet& ds2
589 )
590 {
591  return ds1*ds2;
592 }
593 
594 
595 // ************************************************************************* //
dimensionedScalar sign(const dimensionedScalar &ds)
dimensionSet trans(const dimensionSet &)
Definition: dimensionSet.C:450
intWM_LABEL_SIZE_t label
A label is an int32_t or int64_t as specified by the pre-processor macro WM_LABEL_SIZE.
Definition: label.H:59
bool operator/=(const dimensionSet &)
Definition: dimensionSet.C:217
error FatalError
dimensioned< Type > max(const dimensioned< Type > &, const dimensioned< Type > &)
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:323
dimensionedSymmTensor sqr(const dimensionedVector &dv)
dimensionedSphericalTensor inv(const dimensionedSphericalTensor &dt)
bool operator-=(const dimensionSet &) const
Definition: dimensionSet.C:195
dimensionedScalar sqrt(const dimensionedScalar &ds)
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:251
dimensionedScalar pow025(const dimensionedScalar &ds)
bool operator=(const dimensionSet &) const
Definition: dimensionSet.C:167
dimensionSet(const scalar mass, const scalar length, const scalar time, const scalar temperature, const scalar moles, const scalar current, const scalar luminousIntensity)
Construct given individual dimension exponents for all.
Definition: dimensionSet.C:42
const dimensionSet dimless
dimensionedScalar posPart(const dimensionedScalar &ds)
dimensionedScalar neg(const dimensionedScalar &ds)
dimensionedScalar pow5(const dimensionedScalar &ds)
dimensioned< Type > cmptDivide(const dimensioned< Type > &, const dimensioned< Type > &)
bool operator*=(const dimensionSet &)
Definition: dimensionSet.C:209
dimensionedScalar pos(const dimensionedScalar &ds)
Dimension set for the base types.
Definition: dimensionSet.H:120
tmp< fvMatrix< Type > > operator==(const fvMatrix< Type > &, const fvMatrix< Type > &)
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 > &)
bool operator+=(const dimensionSet &) const
Definition: dimensionSet.C:181
errorManip< error > abort(error &err)
Definition: errorManip.H:131
bool operator!=(const dimensionSet &) const
Definition: dimensionSet.C:161
dimensioned< Type > cmptMultiply(const dimensioned< Type > &, const dimensioned< Type > &)
dimensioned< scalar > magSqr(const dimensioned< Type > &)
dimensionedScalar pos0(const dimensionedScalar &ds)
dimensionedScalar atan2(const dimensionedScalar &x, const dimensionedScalar &y)
defineTypeNameAndDebug(combustionModel, 0)
dimensioned< Type > min(const dimensioned< Type > &, const dimensioned< Type > &)
dimensionedScalar pow(const dimensionedScalar &ds, const dimensionedScalar &expt)
void reset(const dimensionSet &)
Definition: dimensionSet.C:108
dimensionedScalar pow3(const dimensionedScalar &ds)
scalar operator[](const dimensionType) const
Definition: dimensionSet.C:119
fileType type(const fileName &, const bool checkVariants=true, const bool followLink=true)
Return the file type: directory or file.
Definition: POSIX.C:488
dimensionedScalar pow4(const dimensionedScalar &ds)
const dimensionSet & dimensions() const
Return const reference to dimensions.
dimensionedScalar pow6(const dimensionedScalar &ds)
dimensionType
Define an enumeration for the names of the dimension exponents.
Definition: dimensionSet.H:133
dimensioned< scalar > mag(const dimensioned< Type > &)
volScalarField & p
bool operator==(const dimensionSet &) const
Definition: dimensionSet.C:143
static const scalar smallExponent
Definition: dimensionSet.H:147
bool dimensionless() const
Return true if it is dimensionless.
Definition: dimensionSet.C:89
Namespace for OpenFOAM.
dimensionedScalar negPart(const dimensionedScalar &ds)
dimensionSet transform(const dimensionSet &)
Definition: dimensionSet.C:477