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