interpolationTable.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-2017 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 "interpolationTable.H"
27 #include "openFoamTableReader.H"
28 
29 // * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * * //
30 
31 template<class Type>
33 {
34  // preserve the original (unexpanded) fileName to avoid absolute paths
35  // appearing subsequently in the write() method
36  fileName fName(fileName_);
37 
38  fName.expand();
39 
40  // Read data from file
41  reader_()(fName, *this);
42 
43  if (this->empty())
44  {
46  << "table read from " << fName << " is empty" << nl
47  << exit(FatalError);
48  }
49 
50  // Check that the data are okay
51  check();
52 }
53 
54 
55 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
56 
57 template<class Type>
59 :
60  List<Tuple2<scalar, Type>>(),
61  boundsHandling_(interpolationTable::WARN),
62  fileName_("fileNameIsUndefined"),
63  reader_(nullptr)
64 {}
65 
66 
67 template<class Type>
69 (
70  const List<Tuple2<scalar, Type>>& values,
71  const boundsHandling bounds,
72  const fileName& fName
73 )
74 :
76  boundsHandling_(bounds),
77  fileName_(fName),
78  reader_(nullptr)
79 {}
80 
81 
82 template<class Type>
84 :
85  List<Tuple2<scalar, Type>>(),
86  boundsHandling_(interpolationTable::WARN),
87  fileName_(fName),
88  reader_(new openFoamTableReader<Type>(dictionary()))
89 {
90  readTable();
91 }
92 
93 
94 template<class Type>
96 :
97  List<Tuple2<scalar, Type>>(),
98  boundsHandling_(wordToBoundsHandling(dict.lookup("outOfBounds"))),
99  fileName_(dict.lookup("file")),
100  reader_(tableReader<Type>::New(dict))
101 {
102  readTable();
103 }
104 
105 
106 template<class Type>
108 (
109  const interpolationTable& interpTable
110 )
111 :
112  List<Tuple2<scalar, Type>>(interpTable),
113  boundsHandling_(interpTable.boundsHandling_),
114  fileName_(interpTable.fileName_),
115  reader_(interpTable.reader_) // note: steals reader. Used in write().
116 {}
117 
118 
119 
120 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
121 
122 template<class Type>
124 (
125  const boundsHandling& bound
126 ) const
127 {
128  word enumName("warn");
129 
130  switch (bound)
131  {
132  case interpolationTable::ERROR:
133  {
134  enumName = "error";
135  break;
136  }
137  case interpolationTable::WARN:
138  {
139  enumName = "warn";
140  break;
141  }
142  case interpolationTable::CLAMP:
143  {
144  enumName = "clamp";
145  break;
146  }
147  case interpolationTable::REPEAT:
148  {
149  enumName = "repeat";
150  break;
151  }
152  }
153 
154  return enumName;
155 }
156 
157 
158 template<class Type>
161 (
162  const word& bound
163 ) const
164 {
165  if (bound == "error")
166  {
167  return interpolationTable::ERROR;
168  }
169  else if (bound == "warn")
170  {
171  return interpolationTable::WARN;
172  }
173  else if (bound == "clamp")
174  {
175  return interpolationTable::CLAMP;
176  }
177  else if (bound == "repeat")
178  {
179  return interpolationTable::REPEAT;
180  }
181  else
182  {
184  << "bad outOfBounds specifier " << bound << " using 'warn'" << endl;
185 
186  return interpolationTable::WARN;
187  }
188 }
189 
190 
191 template<class Type>
194 (
195  const boundsHandling& bound
196 )
197 {
198  boundsHandling prev = boundsHandling_;
199  boundsHandling_ = bound;
200  return prev;
201 }
202 
203 
204 template<class Type>
206 {
207  label n = this->size();
208  scalar prevValue = List<Tuple2<scalar, Type>>::operator[](0).first();
209 
210  for (label i=1; i<n; ++i)
211  {
212  const scalar currValue =
213  List<Tuple2<scalar, Type>>::operator[](i).first();
214 
215  // avoid duplicate values (divide-by-zero error)
216  if (currValue <= prevValue)
217  {
219  << "out-of-order value: "
220  << currValue << " at index " << i << nl
221  << exit(FatalError);
222  }
223  prevValue = currValue;
224  }
225 }
226 
227 
228 template<class Type>
230 {
231  os.writeKeyword("file")
232  << fileName_ << token::END_STATEMENT << nl;
233  os.writeKeyword("outOfBounds")
234  << boundsHandlingToWord(boundsHandling_) << token::END_STATEMENT << nl;
235  if (reader_.valid())
236  {
237  reader_->write(os);
238  }
239 }
240 
241 
242 template<class Type>
243 Type Foam::interpolationTable<Type>::rateOfChange(const scalar value) const
244 {
245  label n = this->size();
246 
247  if (n <= 1)
248  {
249  // There are not enough entries to provide a rate of change
250  return 0;
251  }
252 
253  scalar minLimit = List<Tuple2<scalar, Type>>::operator[](0).first();
254  scalar maxLimit = List<Tuple2<scalar, Type>>::operator[](n-1).first();
255  scalar lookupValue = value;
256 
257  if (lookupValue < minLimit)
258  {
259  switch (boundsHandling_)
260  {
261  case interpolationTable::ERROR:
262  {
264  << "value (" << lookupValue << ") underflow" << nl
265  << exit(FatalError);
266  break;
267  }
268  case interpolationTable::WARN:
269  {
271  << "value (" << lookupValue << ") underflow" << nl
272  << " Zero rate of change."
273  << endl;
274  // fall-through to 'CLAMP'
275  }
276  case interpolationTable::CLAMP:
277  {
278  return 0;
279  break;
280  }
281  case interpolationTable::REPEAT:
282  {
283  // adjust lookupValue to >= minLimit
284  scalar span = maxLimit-minLimit;
285  lookupValue = fmod(lookupValue-minLimit, span) + minLimit;
286  break;
287  }
288  }
289  }
290  else if (lookupValue >= maxLimit)
291  {
292  switch (boundsHandling_)
293  {
294  case interpolationTable::ERROR:
295  {
297  << "value (" << lookupValue << ") overflow" << nl
298  << exit(FatalError);
299  break;
300  }
301  case interpolationTable::WARN:
302  {
304  << "value (" << lookupValue << ") overflow" << nl
305  << " Zero rate of change."
306  << endl;
307  // fall-through to 'CLAMP'
308  }
309  case interpolationTable::CLAMP:
310  {
311  return 0;
312  break;
313  }
314  case interpolationTable::REPEAT:
315  {
316  // adjust lookupValue <= maxLimit
317  scalar span = maxLimit-minLimit;
318  lookupValue = fmod(lookupValue-minLimit, span) + minLimit;
319  break;
320  }
321  }
322  }
323 
324  label lo = 0;
325  label hi = 0;
326 
327  // look for the correct range
328  for (label i = 0; i < n; ++i)
329  {
330  if (lookupValue >= List<Tuple2<scalar, Type>>::operator[](i).first())
331  {
332  lo = hi = i;
333  }
334  else
335  {
336  hi = i;
337  break;
338  }
339  }
340 
341  if (lo == hi)
342  {
343  // we are at the end of the table - or there is only a single entry
344  return 0;
345  }
346  else if (hi == 0)
347  {
348  // this treatment should should only occur under these conditions:
349  // -> the 'REPEAT' treatment
350  // -> (0 <= value <= minLimit)
351  // -> minLimit > 0
352  // Use the value at maxLimit as the value for value=0
353  lo = n - 1;
354 
355  return
356  (
357  (
358  List<Tuple2<scalar, Type>>::operator[](hi).second()
359  - List<Tuple2<scalar, Type>>::operator[](lo).second()
360  )
361  /(
362  List<Tuple2<scalar, Type>>::operator[](hi).first()
363  + minLimit
364  - List<Tuple2<scalar, Type>>::operator[](lo).first()
365  )
366  );
367  }
368  else
369  {
370  // normal rate of change
371  return
372  (
373  (
374  List<Tuple2<scalar, Type>>::operator[](hi).second()
375  - List<Tuple2<scalar, Type>>::operator[](lo).second()
376  )
377  /(
378  List<Tuple2<scalar, Type>>::operator[](hi).first()
379  - List<Tuple2<scalar, Type>>::operator[](lo).first()
380  )
381  );
382  }
383 }
384 
385 
386 // * * * * * * * * * * * * * * * Member Operators * * * * * * * * * * * * * //
387 
388 template<class Type>
391 {
392  label ii = i;
393  label n = this->size();
394 
395  if (n <= 1)
396  {
397  ii = 0;
398  }
399  else if (ii < 0)
400  {
401  switch (boundsHandling_)
402  {
403  case interpolationTable::ERROR:
404  {
406  << "index (" << ii << ") underflow" << nl
407  << exit(FatalError);
408  break;
409  }
410  case interpolationTable::WARN:
411  {
413  << "index (" << ii << ") underflow" << nl
414  << " Continuing with the first entry"
415  << endl;
416  // fall-through to 'CLAMP'
417  }
418  case interpolationTable::CLAMP:
419  {
420  ii = 0;
421  break;
422  }
423  case interpolationTable::REPEAT:
424  {
425  while (ii < 0)
426  {
427  ii += n;
428  }
429  break;
430  }
431  }
432  }
433  else if (ii >= n)
434  {
435  switch (boundsHandling_)
436  {
437  case interpolationTable::ERROR:
438  {
440  << "index (" << ii << ") overflow" << nl
441  << exit(FatalError);
442  break;
443  }
444  case interpolationTable::WARN:
445  {
447  << "index (" << ii << ") overflow" << nl
448  << " Continuing with the last entry"
449  << endl;
450  // fall-through to 'CLAMP'
451  }
452  case interpolationTable::CLAMP:
453  {
454  ii = n - 1;
455  break;
456  }
457  case interpolationTable::REPEAT:
458  {
459  while (ii >= n)
460  {
461  ii -= n;
462  }
463  break;
464  }
465  }
466  }
467 
468  return List<Tuple2<scalar, Type>>::operator[](ii);
469 }
470 
471 
472 template<class Type>
473 Type Foam::interpolationTable<Type>::operator()(const scalar value) const
474 {
475  label n = this->size();
476 
477  if (n <= 1)
478  {
479  return List<Tuple2<scalar, Type>>::operator[](0).second();
480  }
481 
482  scalar minLimit = List<Tuple2<scalar, Type>>::operator[](0).first();
483  scalar maxLimit = List<Tuple2<scalar, Type>>::operator[](n-1).first();
484  scalar lookupValue = value;
485 
486  if (lookupValue < minLimit)
487  {
488  switch (boundsHandling_)
489  {
490  case interpolationTable::ERROR:
491  {
493  << "value (" << lookupValue << ") underflow" << nl
494  << exit(FatalError);
495  break;
496  }
497  case interpolationTable::WARN:
498  {
500  << "value (" << lookupValue << ") underflow" << nl
501  << " Continuing with the first entry"
502  << endl;
503  // fall-through to 'CLAMP'
504  }
505  case interpolationTable::CLAMP:
506  {
507  return List<Tuple2<scalar, Type>>::operator[](0).second();
508  break;
509  }
510  case interpolationTable::REPEAT:
511  {
512  // adjust lookupValue to >= minLimit
513  scalar span = maxLimit-minLimit;
514  lookupValue = fmod(lookupValue-minLimit, span) + minLimit;
515  break;
516  }
517  }
518  }
519  else if (lookupValue >= maxLimit)
520  {
521  switch (boundsHandling_)
522  {
523  case interpolationTable::ERROR:
524  {
526  << "value (" << lookupValue << ") overflow" << nl
527  << exit(FatalError);
528  break;
529  }
530  case interpolationTable::WARN:
531  {
533  << "value (" << lookupValue << ") overflow" << nl
534  << " Continuing with the last entry"
535  << endl;
536  // fall-through to 'CLAMP'
537  }
538  case interpolationTable::CLAMP:
539  {
540  return List<Tuple2<scalar, Type>>::operator[](n-1).second();
541  break;
542  }
543  case interpolationTable::REPEAT:
544  {
545  // adjust lookupValue <= maxLimit
546  scalar span = maxLimit-minLimit;
547  lookupValue = fmod(lookupValue-minLimit, span) + minLimit;
548  break;
549  }
550  }
551  }
552 
553  label lo = 0;
554  label hi = 0;
555 
556  // look for the correct range
557  for (label i = 0; i < n; ++i)
558  {
559  if (lookupValue >= List<Tuple2<scalar, Type>>::operator[](i).first())
560  {
561  lo = hi = i;
562  }
563  else
564  {
565  hi = i;
566  break;
567  }
568  }
569 
570  if (lo == hi)
571  {
572  // we are at the end of the table - or there is only a single entry
573  return List<Tuple2<scalar, Type>>::operator[](hi).second();
574  }
575  else if (hi == 0)
576  {
577  // this treatment should should only occur under these conditions:
578  // -> the 'REPEAT' treatment
579  // -> (0 <= value <= minLimit)
580  // -> minLimit > 0
581  // Use the value at maxLimit as the value for value=0
582  lo = n - 1;
583 
584  return
585  (
586  List<Tuple2<scalar, Type>>::operator[](lo).second()
587  + (
588  List<Tuple2<scalar, Type>>::operator[](hi).second()
589  - List<Tuple2<scalar, Type>>::operator[](lo).second()
590  )
591  *(lookupValue / minLimit)
592  );
593  }
594  else
595  {
596  // normal interpolation
597  return
598  (
599  List<Tuple2<scalar, Type>>::operator[](lo).second()
600  + (
601  List<Tuple2<scalar, Type>>::operator[](hi).second()
602  - List<Tuple2<scalar, Type>>::operator[](lo).second()
603  )
604  *(
605  lookupValue
606  - List<Tuple2<scalar, Type>>::operator[](lo).first()
607  )
608  /(
609  List<Tuple2<scalar, Type>>::operator[](hi).first()
610  - List<Tuple2<scalar, Type>>::operator[](lo).first()
611  )
612  );
613  }
614 }
615 
616 
617 // ************************************************************************* //
dictionary dict
word boundsHandlingToWord(const boundsHandling &bound) const
Return the out-of-bounds handling as a word.
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
A class for handling file names.
Definition: fileName.H:69
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:124
const Tuple2< scalar, Type > & operator[](const label) const
Return an element of constant Tuple2<scalar, Type>
interpolationTable()
Construct null.
error FatalError
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:137
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:319
A 2-tuple for storing two objects of different types.
Definition: HashTable.H:66
A 1D array of objects of type <T>, where the size of the vector is known and used for subscript bound...
Definition: HashTable.H:60
Type rateOfChange(const scalar) const
Return the rate of change at the interpolation location.
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:253
Reads an interpolation table from a file - OpenFOAM-format.
Base class to read table data for the interpolationTable.
Definition: tableReader.H:57
boundsHandling outOfBounds(const boundsHandling &bound)
Set the out-of-bounds handling from enum, return previous setting.
tmp< DimensionedField< TypeR, GeoMesh > > New(const tmp< DimensionedField< TypeR, GeoMesh >> &tdf1, const word &name, const dimensionSet &dimensions)
stressControl lookup("compactNormalStress") >> compactNormalStress
A class for handling words, derived from string.
Definition: word.H:59
void check() const
Check that list is monotonically increasing.
An interpolation/look-up table of scalar vs <Type> values. The reference scalar values must be monoto...
Type operator()(const scalar) const
Return an interpolated value.
An Ostream is an abstract base class for all output systems (streams, files, token lists...
Definition: Ostream.H:53
static const char nl
Definition: Ostream.H:262
Ostream & writeKeyword(const keyType &)
Write the keyword followed by an appropriate indentation.
Definition: Ostream.C:54
void write(Ostream &os) const
Write.
boundsHandling
Enumeration for handling out-of-bound values.
volScalarField & bound(volScalarField &, const dimensionedScalar &lowerBound)
Bound the given scalar field if it has gone unbounded.
Definition: bound.C:33
#define WarningInFunction
Report a warning using Foam::Warning.
string & expand(const bool allowEmpty=false)
Expand initial tildes and all occurences of environment variables.
Definition: string.C:95
label n
boundsHandling wordToBoundsHandling(const word &bound) const
Return the out-of-bounds handling as an enumeration.