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