interpolationLookUpTable.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-2013 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 
27 // * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * * //
28 
29 template<class Type>
31 (
32  const List<scalar>& indices,
33  const bool lastDim
34 ) const
35 {
36  label totalIndex = 0;
37 
38  forAll(dim_, i)
39  {
40  label dim = 1;
41  for (int j = i + 1; j < dim_.size(); j++)
42  {
43  dim *= dim_[j] + 1;
44  }
45 
46  totalIndex +=
47  dim
48  *min
49  (
50  max(label((indices[i] - min_[i])/delta_[i]), 0),
51  dim_[i]
52  );
53  }
54 
55  if (lastDim)
56  {
57  label iLastdim = dim_.size() - 1;
58  totalIndex += Foam::min
59  (
60  max
61  (
62  label((indices[iLastdim] - min_[iLastdim])/delta_[iLastdim]),
63  0
64  ),
65  dim_[iLastdim]
66  );
67  }
68 
69  return totalIndex;
70 }
71 
72 
73 template<class Type>
75 (
76  const scalar indice
77 ) const
78 {
79  label i = 0;
80  label totalIndex =
81  Foam::min
82  (
83  Foam::max
84  (
85  label((indice - min_[i])/delta_[i]),
86  0
87  ),
88  dim_[i]
89  );
90 
91  return totalIndex;
92 }
93 
94 
95 template<class Type>
97 (
98  const scalar lookUpValue,
99  const label interfield
100 ) const
101 {
102  return lookUpValue >= min_[interfield] && lookUpValue <= max_[interfield];
103 }
104 
105 
106 template<class Type>
108 (
109  const label lo,
110  const label hi,
111  const scalar lookUpValue,
112  const label ofield,
113  const label interfield
114 ) const
115 {
116  if
117  (
118  List<scalarField>::operator[](interfield).operator[](hi)
119  != List<scalarField>::operator[](interfield).operator[](lo)
120  )
121  {
122  scalar output
123  (
124  List<scalarField>::operator[](ofield).operator[](lo)
125  + (
126  List<scalarField>::operator[](ofield).operator[](hi)
127  - List<scalarField>::operator[](ofield).operator[](lo)
128  )
129  *(
130  lookUpValue
131  - List<scalarField>::operator[](interfield).operator[](lo)
132  )
133  /(
134  List<scalarField>::operator[](interfield).operator[](hi)
135  - List<scalarField>::operator[](interfield).operator[](lo)
136  )
137  );
138  return output;
139  }
140  else
141  {
142  return List<scalarField>::operator[](ofield).operator[](lo);
143  }
144 }
145 
146 
147 template<class Type>
149 {
150  min_.setSize(entries_.size());
151  dim_.setSize(entries_.size());
152  delta_.setSize(entries_.size());
153  max_.setSize(entries_.size());
154  entryIndices_.setSize(entries_.size());
155  outputIndices_.setSize(output_.size());
156  label index = 0;
157  label tableDim = 1;
158 
159  forAll(entries_,i)
160  {
161  dim_[i] = readLabel(entries_[i].lookup("N"));
162  max_[i] = readScalar(entries_[i].lookup("max"));
163  min_[i] = readScalar(entries_[i].lookup("min"));
164  delta_[i] = (max_[i] - min_[i])/dim_[i];
165  tableDim *= dim_[i] + 1;
166  fieldIndices_.insert(entries_[i].lookup("name"), index);
167  entryIndices_[i] = index;
168  index++;
169  }
170 
171  forAll(output_,i)
172  {
173  fieldIndices_.insert(output_[i].lookup("name"), index);
174  outputIndices_[i] = index;
175  index++;
176  }
177 
178  List<scalarField>& internal = *this;
179 
180  internal.setSize(entries_.size() + output_.size());
181 
182  interpOutput_.setSize(entries_.size() + output_.size());
183 
184  forAll(internal, i)
185  {
186  internal[i].setSize(tableDim);
187  }
188 }
189 
190 
191 template<class Type>
193 (
194  const word& instance,
195  const objectRegistry& obr
196 )
197 {
198  IOdictionary control
199  (
200  IOobject
201  (
202  fileName_,
203  instance,
204  obr,
205  IOobject::MUST_READ_IF_MODIFIED,
206  IOobject::NO_WRITE
207  )
208  );
209 
210  control.lookup("fields") >> entries_;
211  control.lookup("output") >> output_;
212  control.lookup("values") >> *this;
213 
214  dimensionTable();
215 
216  check();
217 
218  if (this->size() == 0)
219  {
221  (
222  "Foam::interpolationLookUpTable<Type>::readTable()"
223  ) << "table is empty" << nl << exit(FatalError);
224  }
225 }
226 
227 
228 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
229 
230 template<class Type>
232 :
233  List<scalarField>(),
234  fileName_("fileNameIsUndefined")
235 {}
236 
237 
238 template<class Type>
240 (
241  const fileName& fn,
242  const word& instance,
243  const objectRegistry& obr
244 )
245 :
247  fileName_(fn),
248  dim_(0),
249  min_(0),
250  delta_(0.0),
251  max_(0.0),
252  entries_(0),
253  output_(0),
254  entryIndices_(0),
255  outputIndices_(0),
256  interpOutput_(0)
257 {
258  readTable(instance, obr);
259 }
260 
261 
262 template<class Type>
264 (
265  const interpolationLookUpTable& interpTable
266 )
267 :
268  List<scalarField>(interpTable),
269  fileName_(interpTable.fileName_),
270  entryIndices_(interpTable.entryIndices_),
271  outputIndices_(interpTable.outputIndices_),
272  dim_(interpTable.dim_),
273  min_(interpTable.min_),
274  delta_(interpTable.delta_),
275  max_(interpTable.max_),
276  entries_(0),
277  output_(0),
278  interpOutput_(interpTable.interpOutput_)
279 {}
280 
281 
282 template<class Type>
284 (
285  const dictionary& dict
286 )
287 :
289  fileName_(fileName(dict.lookup("fileName")).expand()),
290  dim_(0),
291  min_(0.0),
292  delta_(0.0),
293  max_(0.0),
294  entries_(dict.lookup("fields")),
295  output_(dict.lookup("output")),
296  entryIndices_(0),
297  outputIndices_(0),
298  fieldIndices_(0),
299  interpOutput_(0)
300 {
301  dimensionTable();
302 }
303 
304 
305 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
306 
307 template<class Type>
309 {
310  // check order in the first dimension.
311  scalar prevValue = List<scalarField>::operator[](0).operator[](0);
312  label dim = 1;
313  for (int j = 1; j < dim_.size(); j++)
314  {
315  dim *= dim_[j] + 1;
316  }
317 
318  for (label i = 1; i < dim_[0]; i++)
319  {
320  label index = i*dim;
321  const scalar currValue =
322  List<scalarField>::operator[](0).operator[](index);
323 
324  // avoid duplicate values (divide-by-zero error)
325  if (currValue <= prevValue)
326  {
328  (
329  "Foam::interpolationLookUpTable<Type>::checkOrder() const"
330  ) << "out-of-order value: " << currValue
331  << " at index " << index << nl << exit(FatalError);
332  }
333  prevValue = currValue;
334  }
335 }
336 
337 
338 template<class Type>
340 (
341  Ostream& os,
342  const fileName& fn,
343  const word& instance,
344  const objectRegistry& obr
345 ) const
346 {
347  IOdictionary control
348  (
349  IOobject
350  (
351  fn,
352  instance,
353  obr,
356  )
357  );
358 
359  control.writeHeader(os);
360 
361  os.writeKeyword("fields")
362  << entries_ << token::END_STATEMENT << nl;
363 
364  os.writeKeyword("output")
365  << output_ << token::END_STATEMENT << nl;
366 
367  if (this->size() == 0)
368  {
370  (
371  "Foam::interpolationTable<Type>::write()"
372  ) << "table is empty" << nl << exit(FatalError);
373  }
374  os.writeKeyword("values")
375  << *this << token::END_STATEMENT << nl;
376 }
377 
378 
379 // * * * * * * * * * * * * * * * Member Operators * * * * * * * * * * * * * //
380 
381 template<class Type>
384 {
385  const label n = this->size();
386 
387  if (n <= 1)
388  {
390  (
391  "Foam::interpolationLookUpTable<Type>::operator[](const label)"
392  ) << "table has (" << n << ") columns" << nl << exit(FatalError);
393  }
394  else if (i < 0)
395  {
397  (
398  "Foam::interpolationLookUpTable<Type>::operator[](const label)"
399  ) << "index (" << i << ") underflow" << nl << exit(FatalError);
400  }
401  else if (i >= n)
402  {
404  (
405  "Foam::interpolationLookUpTable<Type>::operator[](const label)"
406  ) << "index (" << i << ") overflow" << nl << exit(FatalError);
407  }
408 
410 }
411 
412 
413 template<class Type>
414 const Foam::scalarField&
416 {
417  const label n = this->size();
418 
419  if (n <= 1)
420  {
422  (
423  "Foam::interpolationLookUpTable<Type>::operator[]"
424  "(const label) const"
425  ) << "table has (" << n << ") columns" << nl << exit(FatalError);
426  }
427  else if (i < 0)
428  {
430  (
431  "Foam::interpolationLookUpTable<Type>::operator[]"
432  "(const label) const"
433  ) << "index (" << i << ") underflow" << nl << exit(FatalError);
434  }
435  else if (i >= n)
436  {
438  (
439  "Foam::interpolationLookUpTable<Type>::operator[]"
440  "(const label) const"
441  ) << "index (" << i << ") overflow" << nl
442  << exit(FatalError);
443  }
444 
446 }
447 
448 
449 template<class Type>
451 {
452  return fieldIndices_.found(fieldName);
453 }
454 
455 
456 template<class Type>
457 const Foam::scalarList&
459 {
460  const label lo = index(retvals);
461  findHi(lo, retvals);
462  return interpOutput_;
463 }
464 
465 
466 template<class Type>
468 (
469  const label lo,
470  const scalar retvals
471 )
472 {
473  forAll(outputIndices_,j)
474  {
475  scalar tmp = 0;
476  label ofield = outputIndices_[j];
477  scalar baseValue = List<scalarField>::operator[](ofield).operator[](lo);
478 
479  forAll(entryIndices_,i)
480  {
481  if (checkRange(retvals, entryIndices_[i]))
482  {
483  label dim = 1;
484 
485  label hi = Foam::min(lo + dim, (*this)[0].size() - 1);
486 
487  tmp += interpolate(lo, hi, retvals, ofield, entryIndices_[i])
488  - baseValue;
489  }
490  interpOutput_[entryIndices_[i]] = retvals;
491  }
492 
493  tmp += baseValue;
494  interpOutput_[outputIndices_[j]] = tmp;
495  }
496 }
497 
498 
499 // ************************************************************************* //
#define readScalar
Definition: doubleScalar.C:38
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
friend Ostream & operator(Ostream &, const UList< T > &)
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:137
label readLabel(Istream &is)
Definition: label.H:64
string expand(const string &, const HashTable< string, word, string::hash > &mapping, const char sigil= '$')
Expand occurences of variables according to the mapping.
Definition: stringOps.C:74
label n
A 1D array of objects of type <T>, where the size of the vector is known and used for subscript bound...
Definition: HashTable.H:59
T & operator[](const label)
Return element of UList.
Definition: UListI.H:163
static const char nl
Definition: Ostream.H:260
void setSize(const label)
Reset size of List.
Definition: List.C:318
stressControl lookup("compactNormalStress") >> compactNormalStress
IOobject defines the attributes of an object for which implicit objectRegistry management is supporte...
Definition: IOobject.H:91
bool found(const word &fieldName) const
Return true if the field exists in the table.
dimensioned< Type > max(const dimensioned< Type > &, const dimensioned< Type > &)
#define forAll(list, i)
Definition: UList.H:421
IOdictionary is derived from dictionary and IOobject to give the dictionary automatic IO functionalit...
Definition: IOdictionary.H:53
const List< label > & dim() const
Return const access to the list of dimensions.
Ostream & writeKeyword(const keyType &)
Write the keyword followed by an appropriate indentation.
Definition: Ostream.C:59
#define FatalErrorIn(functionName)
Report an error message using Foam::FatalError.
Definition: error.H:314
const scalarField & operator[](const label) const
Return an element of constant List<scalar, Type>
ITstream & lookup(const word &, bool recursive=false, bool patternMatch=true) const
Find and return an entry data stream.
Definition: dictionary.C:452
error FatalError
A list of lists. Interpolates based on the first dimension. The values must be positive and monotonic...
Registry of regIOobjects.
A class for handling file names.
Definition: fileName.H:69
label size() const
Return the number of elements in the UList.
An Ostream is an abstract base class for all output systems (streams, files, token lists...
Definition: Ostream.H:53
bool writeHeader(Ostream &) const
Write header.
dimensioned< Type > min(const dimensioned< Type > &, const dimensioned< Type > &)
void write(Ostream &, const fileName &, const word &instance, const objectRegistry &) const
Write lookup table to filename.
bool found(const Key &) const
Return true if hashedEntry is found in table.
Definition: HashTable.C:109
A class for managing temporary objects.
Definition: PtrList.H:118
const List< scalar > & lookUp(const scalar)
Return the output list given a single input scalar.