dimensionSetIO.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) 2011-2014 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 "IOstreams.H"
28 #include "dimensionedScalar.H"
29 #include <limits>
30 
31 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
32 
34 {
35  is >> *this;
36 }
37 
38 
39 Foam::dimensionSet::tokeniser::tokeniser(Istream& is)
40 :
41  is_(is),
42  tokens_(100),
43  start_(0),
44  size_(0)
45 {}
46 
47 
48 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
49 
50 void Foam::dimensionSet::tokeniser::push(const token& t)
51 {
52  label end = (start_+size_)%tokens_.size();
53  tokens_[end] = t;
54  if (size_ == tokens_.size())
55  {
56  start_ = tokens_.fcIndex(start_);
57  }
58  else
59  {
60  size_++;
61  }
62 }
63 
64 
65 Foam::token Foam::dimensionSet::tokeniser::pop()
66 {
67  token t = tokens_[start_];
68  start_ = tokens_.fcIndex(start_);
69  --size_;
70  return t;
71 }
72 
73 
74 void Foam::dimensionSet::tokeniser::unpop(const token& t)
75 {
76  ++size_;
77  start_ = tokens_.rcIndex(start_);
78  tokens_[start_] = t;
79 }
80 
81 
82 bool Foam::dimensionSet::tokeniser::hasToken() const
83 {
84  return size_ || is_.good();
85 }
86 
87 
88 bool Foam::dimensionSet::tokeniser::valid(char c)
89 {
90  return
91  (
92  !isspace(c)
93  && c != '"' // string quote
94  && c != '\'' // string quote
95  && c != '/' // div
96  && c != ';' // end statement
97  && c != '{' // beg subdict
98  && c != '}' // end subdict
99  && c != '(' // beg expr
100  && c != ')' // end expr
101  && c != '[' // beg dim
102  && c != ']' // end dim
103  && c != '^' // power
104  && c != '*' // mult
105  );
106 }
107 
108 
109 Foam::label Foam::dimensionSet::tokeniser::priority(const token& t)
110 {
111  if (!t.isPunctuation())
112  {
113  return 0;
114  }
115  else if
116  (
117  t.pToken() == token::MULTIPLY
118  || t.pToken() == token::DIVIDE
119  )
120  {
121  return 2;
122  }
123  else if (t.pToken() == '^')
124  {
125  return 3;
126  }
127  else
128  {
129  return 0;
130  }
131 }
132 
133 
134 void Foam::dimensionSet::tokeniser::splitWord(const word& w)
135 {
136  size_t start = 0;
137  for (size_t i=0; i<w.size(); ++i)
138  {
139  if (!valid(w[i]))
140  {
141  if (i > start)
142  {
143  word subWord = w(start, i-start);
144  if (isdigit(subWord[0]) || subWord[0] == token::SUBTRACT)
145  {
146  push(token(readScalar(IStringStream(subWord)())));
147  }
148  else
149  {
150  push(token(subWord));
151  }
152  }
153  if (w[i] != token::SPACE)
154  {
155  if (isdigit(w[i]))
156  {
157  push(token(readScalar(IStringStream(w[i])())));
158  }
159  else
160  {
161  push(token::punctuationToken(w[i]));
162  }
163  }
164  start = i+1;
165  }
166  }
167  if (start < w.size())
168  {
169  word subWord = w(start, w.size()-start);
170  if (isdigit(subWord[0]) || subWord[0] == token::SUBTRACT)
171  {
172  push(token(readScalar(IStringStream(subWord)())));
173  }
174  else
175  {
176  push(token(subWord));
177  }
178  }
179 }
180 
181 
182 Foam::token Foam::dimensionSet::tokeniser::nextToken()
183 {
184  if (size_ == 0)
185  {
186  token t(is_);
187  if (t.isWord())
188  {
189  splitWord(t.wordToken());
190  return pop();
191  }
192  else
193  {
194  return t;
195  }
196  }
197  else
198  {
199  return pop();
200  }
201 }
202 
203 
204 void Foam::dimensionSet::tokeniser::putBack(const token& t)
205 {
206  if (size_ == 0)
207  {
208  push(t);
209  }
210  else
211  {
212  unpop(t);
213  }
214 }
215 
216 
217 void Foam::dimensionSet::round(const scalar tol)
218 {
219  for (int i=0; i < dimensionSet::nDimensions; ++i)
220  {
221  scalar integralPart;
222  scalar fractionalPart = std::modf(exponents_[i], &integralPart);
223 
224  if (mag(fractionalPart-1.0) <= tol)
225  {
226  exponents_[i] = 1.0+integralPart;
227  }
228  else if (mag(fractionalPart+1.0) <= tol)
229  {
230  exponents_[i] = -1.0+integralPart;
231  }
232  else if (mag(fractionalPart) <= tol)
233  {
234  exponents_[i] = integralPart;
235  }
236  }
237 }
238 
239 
240 Foam::dimensionedScalar Foam::dimensionSet::parse
241 (
242  const label lastPrior,
243  tokeniser& tis,
244  const HashTable<dimensionedScalar>& readSet
245 ) const
246 {
247  dimensionedScalar ds("", dimless, 1.0);
248 
249  // Get initial token
250  token nextToken(tis.nextToken());
251 
252  // Store type of last token read. Used to detect two consecutive
253  // symbols and assume multiplication
254  bool haveReadSymbol = false;
255 
256 
257  while (true)
258  {
259  if (nextToken.isWord())
260  {
261  const word& unitName = nextToken.wordToken();
262  const dimensionedScalar& unitDim = readSet[unitName];
263  ds.dimensions() *= unitDim.dimensions();
264  ds.value() *= unitDim.value();
265  haveReadSymbol = true;
266  }
267  else if (nextToken.isNumber())
268  {
269  // no dimensions, just value
270  ds.value() *= nextToken.number();
271  haveReadSymbol = true;
272  }
273  else if (nextToken.isPunctuation())
274  {
275  label nextPrior = tokeniser::priority(nextToken);
276 
277  if (nextToken.pToken() == token::BEGIN_SQR)
278  {
279  // No idea when this will happen
280  tis.putBack(nextToken);
281  return ds;
282  }
283  else if (nextToken.pToken() == token::END_SQR)
284  {
285  tis.putBack(nextToken);
286  return ds;
287  }
288  else if (nextToken.pToken() == token::BEGIN_LIST)
289  {
290  dimensionedScalar sub(parse(nextPrior, tis, readSet));
291 
292  token t = tis.nextToken();
293  if (!t.isPunctuation() || t.pToken() != token::END_LIST)
294  {
296  (
297  "dimensionSet::parse"
298  "(const label, tokeniser&"
299  ", const HashTable<dimensionedScalar>&)",
300  tis.stream()
301  ) << "Illegal token " << t << exit(FatalIOError);
302  }
303 
304  ds.dimensions() *= sub.dimensions();
305  ds.value() *= sub.value();
306 
307  haveReadSymbol = true;
308  }
309  else if (nextToken.pToken() == token::END_LIST)
310  {
311  tis.putBack(nextToken);
312  return ds;
313  }
314  else if (nextToken.pToken() == token::MULTIPLY)
315  {
316  if (nextPrior > lastPrior)
317  {
318  dimensionedScalar sub(parse(nextPrior, tis, readSet));
319 
320  ds.dimensions() *= sub.dimensions();
321  ds.value() *= sub.value();
322  }
323  else
324  {
325  // Restore token
326  tis.putBack(nextToken);
327  return ds;
328  }
329  haveReadSymbol = false;
330  }
331  else if (nextToken.pToken() == token::DIVIDE)
332  {
333  if (nextPrior > lastPrior)
334  {
335  dimensionedScalar sub(parse(nextPrior, tis, readSet));
336 
337  ds.dimensions() /= sub.dimensions();
338  ds.value() /= sub.value();
339  }
340  else
341  {
342  tis.putBack(nextToken);
343  return ds;
344  }
345  haveReadSymbol = false;
346  }
347  else if (nextToken.pToken() == '^')
348  {
349  if (nextPrior > lastPrior)
350  {
351  dimensionedScalar exp(parse(nextPrior, tis, readSet));
352 
353  ds.dimensions().reset(pow(ds.dimensions(), exp.value()));
354  // Round to nearest integer if close to it
355  ds.dimensions().round(10*smallExponent);
356  ds.value() = Foam::pow(ds.value(), exp.value());
357  }
358  else
359  {
360  tis.putBack(nextToken);
361  return ds;
362  }
363  haveReadSymbol = false;
364  }
365  else
366  {
368  (
369  "dimensionSet::parse"
370  "(const label, tokeniser&"
371  ", const HashTable<dimensionedScalar>&)",
372  tis.stream()
373  ) << "Illegal token " << nextToken << exit(FatalIOError);
374  }
375  }
376  else
377  {
379  (
380  "dimensionSet::parse"
381  "(const label, tokeniser&"
382  ", const HashTable<dimensionedScalar>&)",
383  tis.stream()
384  ) << "Illegal token " << nextToken << exit(FatalIOError);
385  }
386 
387 
388  if (!tis.hasToken())
389  {
390  break;
391  }
392 
393  nextToken = tis.nextToken();
394  if (nextToken.error())
395  {
396  break;
397  }
398 
399  if (haveReadSymbol && (nextToken.isWord() || nextToken.isNumber()))
400  {
401  // Two consecutive symbols. Assume multiplication
402  tis.putBack(nextToken);
403  nextToken = token(token::MULTIPLY);
404  }
405  }
406 
407  return ds;
408 }
409 
410 
412 (
413  Istream& is,
414  scalar& multiplier,
415  const HashTable<dimensionedScalar>& readSet
416 )
417 {
418  multiplier = 1.0;
419 
420  // Read begining of dimensionSet
421  token startToken(is);
422 
423  if (startToken != token::BEGIN_SQR)
424  {
426  (
427  "dimensionSet::read"
428  "(Istream&, scalar&, const HashTable<dimensionedScalar>&)",
429  is
430  ) << "expected a " << token::BEGIN_SQR << " in dimensionSet"
431  << endl << "in stream " << is.info()
432  << exit(FatalIOError);
433  }
434 
435  // Read next token
436  token nextToken(is);
437 
438  if (!nextToken.isNumber())
439  {
440  is.putBack(nextToken);
441 
442  tokeniser tis(is);
443 
444  dimensionedScalar ds(parse(0, tis, readSet));
445 
446  multiplier = ds.value();
447  for (int i=0; i < dimensionSet::nDimensions; ++i)
448  {
449  exponents_[i] = ds.dimensions()[i];
450  }
451  }
452  else
453  {
454  // Read first five dimensions
455  exponents_[dimensionSet::MASS] = nextToken.number();
456  for (int Dimension=1; Dimension<dimensionSet::CURRENT; Dimension++)
457  {
458  is >> exponents_[Dimension];
459  }
460 
461  // Read next token
462  token nextToken(is);
463 
464  // If next token is another number
465  // read last two dimensions
466  // and then read another token for the end of the dimensionSet
467  if (nextToken.isNumber())
468  {
469  exponents_[dimensionSet::CURRENT] = nextToken.number();
470  is >> nextToken;
471  exponents_[dimensionSet::LUMINOUS_INTENSITY] = nextToken.number();
472  is >> nextToken;
473  }
474  else
475  {
476  exponents_[dimensionSet::CURRENT] = 0;
477  exponents_[dimensionSet::LUMINOUS_INTENSITY] = 0;
478  }
479 
480  // Check end of dimensionSet
481  if (nextToken != token::END_SQR)
482  {
484  (
485  "dimensionSet::read"
486  "(Istream&, scalar&, const HashTable<dimensionedScalar>&)",
487  is
488  ) << "expected a " << token::END_SQR << " in dimensionSet "
489  << endl << "in stream " << is.info()
490  << exit(FatalIOError);
491  }
492  }
493  // Check state of Istream
494  is.check("Istream& operator>>(Istream&, dimensionSet&)");
495 
496  return is;
497 }
498 
499 
501 (
502  Istream& is,
503  scalar& multiplier
504 )
505 {
506  return read(is, multiplier, unitSet());
507 }
508 
509 
511 (
512  Istream& is,
513  scalar& multiplier,
514  const dictionary& readSet
515 )
516 {
517  multiplier = 1.0;
518 
519  // Read begining of dimensionSet
520  token startToken(is);
521 
522  if (startToken != token::BEGIN_SQR)
523  {
525  (
526  "dimensionSet::read"
527  "(Istream&, scalar&, const dictionary&)",
528  is
529  ) << "expected a " << token::BEGIN_SQR << " in dimensionSet"
530  << endl << "in stream " << is.info()
531  << exit(FatalIOError);
532  }
533 
534  // Read next token
535  token nextToken(is);
536 
537  if (nextToken.isWord())
538  {
539  bool continueParsing = true;
540  do
541  {
542  word symbolPow = nextToken.wordToken();
543  if (symbolPow[symbolPow.size()-1] == token::END_SQR)
544  {
545  symbolPow = symbolPow(0, symbolPow.size()-1);
546  continueParsing = false;
547  }
548 
549 
550  // Parse unit
551  dimensionSet symbolSet(dimless);
552 
553  size_t index = symbolPow.find('^');
554  if (index != string::npos)
555  {
556  word symbol = symbolPow(0, index);
557  word exp = symbolPow(index+1, symbolPow.size()-index+1);
558  scalar exponent = readScalar(IStringStream(exp)());
559 
561  s.read(readSet[symbol], readSet);
562 
563  symbolSet.reset(pow(s.dimensions(), exponent));
564  // Round to nearest integer if close to it
565  symbolSet.round(10*smallExponent);
566  multiplier *= Foam::pow(s.value(), exponent);
567  }
568  else
569  {
571  s.read(readSet[symbolPow], readSet);
572 
573  symbolSet.reset(s.dimensions());
574  multiplier *= s.value();
575  }
576 
577  // Add dimensions without checking
578  for (int i=0; i < dimensionSet::nDimensions; ++i)
579  {
580  exponents_[i] += symbolSet[i];
581  }
582 
583  if (continueParsing)
584  {
585  nextToken = token(is);
586 
587  if (!nextToken.isWord() || nextToken == token::END_SQR)
588  {
589  continueParsing = false;
590  }
591  }
592  }
593  while (continueParsing);
594  }
595  else
596  {
597  // Read first five dimensions
598  exponents_[dimensionSet::MASS] = nextToken.number();
599  for (int Dimension=1; Dimension<dimensionSet::CURRENT; Dimension++)
600  {
601  is >> exponents_[Dimension];
602  }
603 
604  // Read next token
605  token nextToken(is);
606 
607  // If next token is another number
608  // read last two dimensions
609  // and then read another token for the end of the dimensionSet
610  if (nextToken.isNumber())
611  {
612  exponents_[dimensionSet::CURRENT] = nextToken.number();
613  is >> nextToken;
614  exponents_[dimensionSet::LUMINOUS_INTENSITY] = nextToken.number();
615  is >> nextToken;
616  }
617  else
618  {
619  exponents_[dimensionSet::CURRENT] = 0;
620  exponents_[dimensionSet::LUMINOUS_INTENSITY] = 0;
621  }
622 
623  // Check end of dimensionSet
624  if (nextToken != token::END_SQR)
625  {
627  (
628  "dimensionSet::read"
629  "(Istream&, scalar&, const dictionary&)",
630  is
631  ) << "expected a " << token::END_SQR << " in dimensionSet "
632  << endl << "in stream " << is.info()
633  << exit(FatalIOError);
634  }
635  }
636 
637  // Check state of Istream
638  is.check("Istream& operator>>(Istream&, dimensionSet&)");
639 
640  return is;
641 }
642 
643 
645 (
646  Ostream& os,
647  scalar& multiplier,
648  const dimensionSets& writeUnits
649 ) const
650 {
651  multiplier = 1.0;
652 
653  os << token::BEGIN_SQR;
654 
655  if (writeUnits.valid() && os.format() == IOstream::ASCII)
656  {
658  for (int d=0; d<dimensionSet::nDimensions; d++)
659  {
660  exponents[d] = exponents_[d];
661  }
662  writeUnits.coefficients(exponents);
663 
664  bool hasPrinted = false;
665 
666  // Set precision to lots
667  std::streamsize oldPrecision = os.precision
668  (
669  std::numeric_limits<scalar>::digits10
670  );
671 
672  forAll(exponents, i)
673  {
674  if (mag(exponents[i]) > smallExponent)
675  {
676  const dimensionedScalar& ds = writeUnits.units()[i];
677 
678  if (hasPrinted)
679  {
680  os << token::SPACE;
681  }
682  hasPrinted = true;
683  os << ds.name();
684  if (mag(exponents[i]-1) > smallExponent)
685  {
686  os << '^' << exponents[i];
687 
688  multiplier *= Foam::pow(ds.value(), exponents[i]);
689  }
690  else
691  {
692  multiplier *= ds.value();
693  }
694  }
695  }
696 
697  // Reset precision
698  os.precision(oldPrecision);
699  }
700  else
701  {
702  for (int d=0; d<dimensionSet::nDimensions-1; d++)
703  {
704  os << exponents_[d] << token::SPACE;
705  }
706  os << exponents_[dimensionSet::nDimensions-1];
707  }
708 
709  os << token::END_SQR;
710 
711  // Check state of Ostream
712  os.check("Ostream& operator<<(Ostream&, const dimensionSet&)");
713 
714  return os;
715 }
716 
717 
719 (
720  Ostream& os,
721  scalar& multiplier
722 ) const
723 {
724  return write(os, multiplier, writeUnitSet());
725 }
726 
727 
728 // * * * * * * * * * * * * * * * IOstream Operators * * * * * * * * * * * * //
729 
731 {
732  scalar multiplier;
733  dset.read(is, multiplier);
734 
735  if (mag(multiplier-1.0) > dimensionSet::smallExponent)
736  {
737  FatalIOErrorIn("Foam::operator>>(Istream&, dimensionSet&)", is)
738  << "Cannot use scaled units in dimensionSet"
739  << exit(FatalIOError);
740  }
741 
742  // Check state of Istream
743  is.check("Istream& operator>>(Istream&, dimensionSet&)");
744 
745  return is;
746 }
747 
748 
750 {
751  scalar multiplier;
752  dset.write(os, multiplier);
753 
754  // Check state of Ostream
755  os.check("Ostream& operator<<(Ostream&, const dimensionSet&)");
756 
757  return os;
758 }
759 
760 
761 // ************************************************************************* //
Input from memory buffer stream.
Definition: IStringStream.H:49
InfoProxy< IOstream > info() const
Return info proxy.
Definition: IOstream.H:531
Ostream & write(Ostream &os, scalar &multiplier, const dimensionSets &) const
Write using provided units.
scalar number() const
Definition: tokenI.H:345
bool isWord() const
Definition: tokenI.H:221
dimensionedScalar pow(const dimensionedScalar &ds, const dimensionedScalar &expt)
gmvFile<< "tracers "<< particles.size()<< nl;forAllConstIter(Cloud< passiveParticle >, particles, iter){gmvFile<< iter().position().x()<< " ";}gmvFile<< nl;forAllConstIter(Cloud< passiveParticle >, particles, iter){gmvFile<< iter().position().y()<< " ";}gmvFile<< nl;forAllConstIter(Cloud< passiveParticle >, particles, iter){gmvFile<< iter().position().z()<< " ";}gmvFile<< nl;forAll(lagrangianScalarNames, i){word name=lagrangianScalarNames[i];IOField< scalar > s(IOobject( name, runTime.timeName(), cloud::prefix, mesh, IOobject::MUST_READ, IOobject::NO_WRITE ))
An STL-conforming hash table.
Definition: HashTable.H:61
An Istream is an abstract base class for all input systems (streams, files, token lists etc)...
Definition: Istream.H:57
friend dimensionSet pow(const dimensionSet &, const scalar)
A class for handling words, derived from string.
Definition: word.H:59
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 isPunctuation() const
Definition: tokenI.H:203
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:124
dimensionedScalar exp(const dimensionedScalar &ds)
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:137
const word & wordToken() const
Definition: tokenI.H:226
Ostream & operator<<(Ostream &, const edgeMesh &)
Definition: edgeMeshIO.C:133
const HashTable< dimensionedScalar > & unitSet()
Set of all dimensions.
Definition: dimensionSets.C:97
const dimensionSet & dimensions() const
Return const reference to dimensions.
bool isNumber() const
Definition: tokenI.H:340
virtual bool check(const char *operation) const
Check IOstream status for given operation.
Definition: IOstream.C:92
punctuationToken
Standard punctuation tokens.
Definition: token.H:92
Useful combination of include files which define Sin, Sout and Serr and the use of IO streams general...
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:251
IOerror FatalIOError
void putBack(const token &)
Put back token.
Definition: Istream.C:30
const word & name() const
Return const reference to name.
friend dimensionSet mag(const dimensionSet &)
bool valid() const
Is there a valid inverse of the selected unit.
Dimension set for the base types.
Definition: dimensionSet.H:116
#define forAll(list, i)
Definition: UList.H:421
bool readScalar(const char *buf, doubleScalar &s)
Read whole of buf as a scalar. Return true if succesful.
Definition: doubleScalar.H:63
static const scalar smallExponent
Definition: dimensionSet.H:143
streamFormat format() const
Return current stream format.
Definition: IOstream.H:377
Istream & read(Istream &is, scalar &multiplier, const dictionary &)
Read using provided units. Used only in initial parsing.
A token holds items read from Istream.
Definition: token.H:67
const dimensionSets & writeUnitSet()
Set of units.
virtual int precision() const =0
Get precision of output field.
void coefficients(scalarField &) const
(if valid) obtain set of coefficients of unitNames
Istream & operator>>(Istream &, edgeMesh &)
Definition: edgeMeshIO.C:144
const PtrList< dimensionedScalar > & units() const
Return the units.
An Ostream is an abstract base class for all output systems (streams, files, token lists...
Definition: Ostream.H:53
const dimensionedScalar c
Speed of light in a vacuum.
const dimensionSet dimless(0, 0, 0, 0, 0, 0, 0)
Definition: dimensionSets.H:47
bool good() const
Definition: tokenI.H:188
void reset(const dimensionSet &)
Definition: dimensionSet.C:108
punctuationToken pToken() const
Definition: tokenI.H:208
bool isspace(char c)
Definition: char.H:53
#define FatalIOErrorIn(functionName, ios)
Report an error message using Foam::FatalIOError.
Definition: error.H:325
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
void read(const dictionary &)
Update the value of dimensioned<Type>
const Type & value() const
Return const reference to value.