TableBase.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-2019 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 "TableBase.H"
27 #include "Time.H"
28 #include "interpolationWeights.H"
29 
30 // * * * * * * * * * * * * Protected Member Functions * * * * * * * * * * * //
31 
32 template<class Type>
35 {
36  if (interpolatorPtr_.empty())
37  {
38  // Re-work table into linear list
39  tableSamplesPtr_.reset(new scalarField(table_.size()));
40  scalarField& tableSamples = tableSamplesPtr_();
41  forAll(table_, i)
42  {
43  tableSamples[i] = table_[i].first();
44  }
45  interpolatorPtr_ = interpolationWeights::New
46  (
47  interpolationScheme_,
48  tableSamples
49  );
50  }
51 
52  return interpolatorPtr_();
53 }
54 
55 
56 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
57 
58 template<class Type>
60 (
61  const word& name,
62  const dictionary& dict
63 )
64 :
66  name_(name),
67  boundsHandling_
68  (
69  wordToBoundsHandling
70  (
71  dict.lookupOrDefault<word>("outOfBounds", "clamp")
72  )
73  ),
74  interpolationScheme_
75  (
76  dict.lookupOrDefault<word>("interpolationScheme", "linear")
77  ),
78  table_()
79 {}
80 
81 
82 template<class Type>
84 :
85  Function1<Type>(tbl),
86  name_(tbl.name_),
87  boundsHandling_(tbl.boundsHandling_),
88  interpolationScheme_(tbl.interpolationScheme_),
89  table_(tbl.table_),
90  tableSamplesPtr_(tbl.tableSamplesPtr_),
91  interpolatorPtr_(tbl.interpolatorPtr_)
92 {}
93 
94 
95 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
96 
97 template<class Type>
99 {}
100 
101 
102 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
103 
104 template<class Type>
106 (
107  const boundsHandling& bound
108 ) const
109 {
110  word enumName("warn");
111 
112  switch (bound)
113  {
114  case ERROR:
115  {
116  enumName = "error";
117  break;
118  }
119  case WARN:
120  {
121  enumName = "warn";
122  break;
123  }
124  case CLAMP:
125  {
126  enumName = "clamp";
127  break;
128  }
129  case REPEAT:
130  {
131  enumName = "repeat";
132  break;
133  }
134  }
135 
136  return enumName;
137 }
138 
139 
140 template<class Type>
143 (
144  const word& bound
145 ) const
146 {
147  if (bound == "error")
148  {
149  return ERROR;
150  }
151  else if (bound == "warn")
152  {
153  return WARN;
154  }
155  else if (bound == "clamp")
156  {
157  return CLAMP;
158  }
159  else if (bound == "repeat")
160  {
161  return REPEAT;
162  }
163  else
164  {
166  << "bad outOfBounds specifier " << bound << " using 'warn'"
167  << endl;
168 
169  return WARN;
170  }
171 }
172 
173 
174 template<class Type>
177 (
178  const boundsHandling& bound
179 )
180 {
183 
184  return prev;
185 }
186 
187 
188 template<class Type>
190 {
191  if (!table_.size())
192  {
194  << "Table for entry " << this->name_ << " is invalid (empty)"
195  << nl << exit(FatalError);
196  }
197 
198  label n = table_.size();
199  scalar prevValue = table_[0].first();
200 
201  for (label i = 1; i < n; ++i)
202  {
203  const scalar currValue = table_[i].first();
204 
205  // avoid duplicate values (divide-by-zero error)
206  if (currValue <= prevValue)
207  {
209  << "out-of-order value: " << currValue << " at index " << i
210  << exit(FatalError);
211  }
212  prevValue = currValue;
213  }
214 }
215 
216 
217 template<class Type>
219 (
220  const scalar x,
221  scalar& xDash
222 ) const
223 {
224  if (x < table_[0].first())
225  {
226  switch (boundsHandling_)
227  {
228  case ERROR:
229  {
231  << "value (" << x << ") underflow"
232  << exit(FatalError);
233  break;
234  }
235  case WARN:
236  {
238  << "value (" << x << ") underflow" << nl
239  << endl;
240 
241  // fall-through to 'CLAMP'
242  [[fallthrough]];
243  }
244  case CLAMP:
245  {
246  xDash = table_[0].first();
247  return true;
248  break;
249  }
250  case REPEAT:
251  {
252  // adjust x to >= minX
253  scalar span = table_.last().first() - table_[0].first();
254  xDash = fmod(x - table_[0].first(), span) + table_[0].first();
255  break;
256  }
257  }
258  }
259  else
260  {
261  xDash = x;
262  }
263 
264  return false;
265 }
266 
267 
268 template<class Type>
270 (
271  const scalar x,
272  scalar& xDash
273 ) const
274 {
275  if (x > table_.last().first())
276  {
277  switch (boundsHandling_)
278  {
279  case ERROR:
280  {
282  << "value (" << x << ") overflow"
283  << exit(FatalError);
284  break;
285  }
286  case WARN:
287  {
289  << "value (" << x << ") overflow" << nl
290  << endl;
291 
292  // fall-through to 'CLAMP'
293  [[fallthrough]];
294  }
295  case CLAMP:
296  {
297  xDash = table_.last().first();
298  return true;
299  break;
300  }
301  case REPEAT:
302  {
303  // adjust x to >= minX
304  scalar span = table_.last().first() - table_[0].first();
305  xDash = fmod(x - table_[0].first(), span) + table_[0].first();
306  break;
307  }
308  }
309  }
310  else
311  {
312  xDash = x;
313  }
314 
315  return false;
316 }
317 
318 
319 template<class Type>
321 {
322  forAll(table_, i)
323  {
324  scalar value = table_[i].first();
325  table_[i].first() = t.userTimeToTime(value);
326  }
327 
328  tableSamplesPtr_.clear();
329  interpolatorPtr_.clear();
330 }
331 
332 
333 template<class Type>
335 {
336  scalar xDash = x;
337 
338  if (checkMinBounds(x, xDash))
339  {
340  return table_[0].second();
341  }
342 
343  if (checkMaxBounds(xDash, xDash))
344  {
345  return table_.last().second();
346  }
347 
348  // Use interpolator
350 
351  Type t = currentWeights_[0]*table_[currentIndices_[0]].second();
352  for (label i = 1; i < currentIndices_.size(); i++)
353  {
354  t += currentWeights_[i]*table_[currentIndices_[i]].second();
355  }
356 
357  return t;
358 }
359 
360 
361 template<class Type>
363 (
364  const scalar x1,
365  const scalar x2
366 ) const
367 {
368  // Use interpolator
370 
371  Type sum = currentWeights_[0]*table_[currentIndices_[0]].second();
372  for (label i = 1; i < currentIndices_.size(); i++)
373  {
374  sum += currentWeights_[i]*table_[currentIndices_[i]].second();
375  }
376 
377  return sum;
378 }
379 
380 
381 template<class Type>
383 {
384  tmp<scalarField> tfld(new scalarField(table_.size(), 0.0));
385  scalarField& fld = tfld.ref();
386 
387  forAll(table_, i)
388  {
389  fld[i] = table_[i].first();
390  }
391 
392  return tfld;
393 }
394 
395 
396 template<class Type>
398 {
399  tmp<Field<Type>> tfld(new Field<Type>(table_.size(), Zero));
400  Field<Type>& fld = tfld.ref();
401 
402  forAll(table_, i)
403  {
404  fld[i] = table_[i].second();
405  }
406 
407  return tfld;
408 }
409 
410 
411 template<class Type>
413 {
414  if (boundsHandling_ != CLAMP)
415  {
416  writeEntry(os, "outOfBounds", boundsHandlingToWord(boundsHandling_));
417  }
418  if (interpolationScheme_ != "linear")
419  {
420  writeEntry(os, "interpolationScheme", interpolationScheme_);
421  }
422 }
423 
424 
425 template<class Type>
427 {
429  os << nl << indent << table_ << token::END_STATEMENT << nl;
430  writeEntries(os);
431 }
432 
433 
434 // ************************************************************************* //
Issue warning and clamp value (default)
Definition: TableBase.H:67
Top level data entry class for use in dictionaries. Provides a mechanism to specify a variable as a c...
Definition: Function1.H:53
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:434
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
boundsHandling outOfBounds(const boundsHandling &bound)
Set the out-of-bounds handling from enum, return previous setting.
Definition: TableBase.C:177
Ostream & indent(Ostream &os)
Indent stream.
Definition: Ostream.H:226
TableBase(const word &name, const dictionary &dict)
Construct from dictionary - note table is not populated.
Definition: TableBase.C:60
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:124
List< Tuple2< scalar, Type > > table_
Table data.
Definition: TableBase.H:87
error FatalError
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:158
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:319
const interpolationWeights & interpolator() const
Return (demand driven) interpolator.
Definition: TableBase.C:34
T & ref() const
Return non-const reference or generate a fatal error.
Definition: tmpI.H:174
Exit with a FatalError.
Definition: TableBase.H:66
void size(const label)
Override size to be inconsistent with allocated storage.
Definition: ListI.H:163
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:256
autoPtr< interpolationWeights > interpolatorPtr_
Interpolator method.
Definition: TableBase.H:93
T & first()
Return the first element of the list.
Definition: UListI.H:114
Abstract base class for interpolating in 1D.
Class to control time during OpenFOAM simulations that is also the top-level objectRegistry.
Definition: Time.H:68
virtual ~TableBase()
Destructor.
Definition: TableBase.C:98
dimensioned< Type > sum(const DimensionedField< Type, GeoMesh > &df)
tmp< DimensionedField< TypeR, GeoMesh > > New(const tmp< DimensionedField< TypeR, GeoMesh >> &tdf1, const word &name, const dimensionSet &dimensions)
virtual void writeData(Ostream &os) const
Write in dictionary format.
Definition: Function1.C:150
boundsHandling wordToBoundsHandling(const word &bound) const
Return the out-of-bounds handling as an enumeration.
Definition: TableBase.C:143
Treat as a repeating list.
Definition: TableBase.H:69
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){ const word &name=lagrangianScalarNames[i];IOField< scalar > fld(IOobject(name, runTime.timeName(), cloud::prefix, mesh, IOobject::MUST_READ, IOobject::NO_WRITE))
Base class for table with bounds handling, interpolation and integration.
Definition: TableBase.H:55
virtual scalar userTimeToTime(const scalar theta) const
Convert the user-time (e.g. CA deg) to real-time (s).
Definition: TimeState.C:52
const word name_
Table name.
Definition: TableBase.H:78
A class for handling words, derived from string.
Definition: word.H:59
Field< scalar > scalarField
Specialisation of Field<T> for scalar.
virtual Type value(const scalar x) const
Return Table value.
Definition: TableBase.C:334
bool checkMinBounds(const scalar x, scalar &xDash) const
Check minimum table bounds.
Definition: TableBase.C:219
void check() const
Check the table for size and consistency.
Definition: TableBase.C:189
virtual tmp< scalarField > x() const
Return the reference values.
Definition: TableBase.C:382
static const zero Zero
Definition: zero.H:97
word boundsHandlingToWord(const boundsHandling &bound) const
Return the out-of-bounds handling as a word.
Definition: TableBase.C:106
virtual Type integrate(const scalar x1, const scalar x2) const
Integrate between two (scalar) values.
Definition: TableBase.C:363
virtual bool integrationWeights(const scalar t1, const scalar t2, labelList &indices, scalarField &weights) const =0
Calculate weights and indices to calculate integrand of t1..t2.
An Ostream is an abstract base class for all output systems (streams, files, token lists...
Definition: Ostream.H:53
volScalarField scalarField(fieldObject, mesh)
virtual tmp< Field< Type > > y() const
Return the dependent values.
Definition: TableBase.C:397
static const char nl
Definition: Ostream.H:265
void writeEntry(Ostream &os, const HashTable< T, Key, Hash > &ht)
Definition: HashTableIO.C:96
word name(const complex &)
Return a string representation of a complex.
Definition: complex.C:47
const word interpolationScheme_
Interpolation type.
Definition: TableBase.H:84
virtual void convertTimeBase(const Time &t)
Convert time.
Definition: TableBase.C:320
autoPtr< scalarField > tableSamplesPtr_
Extracted values.
Definition: TableBase.H:90
const boundsHandling boundsHandling_
Enumeration for handling out-of-bound values.
Definition: TableBase.H:81
labelList currentIndices_
Cached indices and weights.
Definition: TableBase.H:96
volScalarField & bound(volScalarField &, const dimensionedScalar &lowerBound)
Bound the given scalar field if it has gone unbounded.
Definition: bound.C:33
T lookupOrDefault(const word &, const T &, bool recursive=false, bool patternMatch=true) const
Find and return a T,.
#define WarningInFunction
Report a warning using Foam::Warning.
bool checkMaxBounds(const scalar x, scalar &xDash) const
Check maximum table bounds.
Definition: TableBase.C:270
label n
A class for managing temporary objects.
Definition: PtrList.H:53
T & last()
Return the last element of the list.
Definition: UListI.H:128
virtual bool valueWeights(const scalar t, labelList &indices, scalarField &weights) const =0
Calculate weights and indices to calculate t from samples.
virtual void writeData(Ostream &os) const
Write all table data in dictionary format.
Definition: TableBase.C:426
Clamp value to the start/end value.
Definition: TableBase.H:68
boundsHandling
Enumeration for handling out-of-bound values.
Definition: TableBase.H:64
virtual void writeEntries(Ostream &os) const
Write keywords only in dictionary format. Used for non-inline.
Definition: TableBase.C:412