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-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 
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  << "table is empty" << nl << exit(FatalError);
222  }
223 }
224 
225 
226 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
227 
228 template<class Type>
230 :
231  List<scalarField>(),
232  fileName_("fileNameIsUndefined")
233 {}
234 
235 
236 template<class Type>
238 (
239  const fileName& fn,
240  const word& instance,
241  const objectRegistry& obr
242 )
243 :
245  fileName_(fn),
246  dim_(0),
247  min_(0),
248  delta_(0.0),
249  max_(0.0),
250  entries_(0),
251  output_(0),
252  entryIndices_(0),
253  outputIndices_(0),
254  interpOutput_(0)
255 {
256  readTable(instance, obr);
257 }
258 
259 
260 template<class Type>
262 (
263  const interpolationLookUpTable& interpTable
264 )
265 :
266  List<scalarField>(interpTable),
267  fileName_(interpTable.fileName_),
268  entryIndices_(interpTable.entryIndices_),
269  outputIndices_(interpTable.outputIndices_),
270  dim_(interpTable.dim_),
271  min_(interpTable.min_),
272  delta_(interpTable.delta_),
273  max_(interpTable.max_),
274  entries_(0),
275  output_(0),
276  interpOutput_(interpTable.interpOutput_)
277 {}
278 
279 
280 template<class Type>
282 (
283  const dictionary& dict
284 )
285 :
287  fileName_(fileName(dict.lookup("file")).expand()),
288  dim_(0),
289  min_(0.0),
290  delta_(0.0),
291  max_(0.0),
292  entries_(dict.lookup("fields")),
293  output_(dict.lookup("output")),
294  entryIndices_(0),
295  outputIndices_(0),
296  fieldIndices_(0),
297  interpOutput_(0)
298 {
299  dimensionTable();
300 }
301 
302 
303 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
304 
305 template<class Type>
307 {
308  // check order in the first dimension.
309  scalar prevValue = List<scalarField>::operator[](0).operator[](0);
310  label dim = 1;
311  for (int j = 1; j < dim_.size(); j++)
312  {
313  dim *= dim_[j] + 1;
314  }
315 
316  for (label i = 1; i < dim_[0]; i++)
317  {
318  label index = i*dim;
319  const scalar currValue =
320  List<scalarField>::operator[](0).operator[](index);
321 
322  // avoid duplicate values (divide-by-zero error)
323  if (currValue <= prevValue)
324  {
326  << "out-of-order value: " << currValue
327  << " at index " << index << nl << exit(FatalError);
328  }
329  prevValue = currValue;
330  }
331 }
332 
333 
334 template<class Type>
336 (
337  Ostream& os,
338  const fileName& fn,
339  const word& instance,
340  const objectRegistry& obr
341 ) const
342 {
343  IOdictionary control
344  (
345  IOobject
346  (
347  fn,
348  instance,
349  obr,
350  IOobject::NO_READ,
351  IOobject::NO_WRITE
352  )
353  );
354 
355  control.writeHeader(os);
356 
357  os.writeKeyword("fields")
358  << entries_ << token::END_STATEMENT << nl;
359 
360  os.writeKeyword("output")
361  << output_ << token::END_STATEMENT << nl;
362 
363  if (this->size() == 0)
364  {
366  << "table is empty" << nl << exit(FatalError);
367  }
368  os.writeKeyword("values")
369  << *this << token::END_STATEMENT << nl;
370 }
371 
372 
373 // * * * * * * * * * * * * * * * Member Operators * * * * * * * * * * * * * //
374 
375 template<class Type>
378 {
379  const label n = this->size();
380 
381  if (n <= 1)
382  {
384  << "table has (" << n << ") columns" << nl << exit(FatalError);
385  }
386  else if (i < 0)
387  {
389  << "index (" << i << ") underflow" << nl << exit(FatalError);
390  }
391  else if (i >= n)
392  {
394  << "index (" << i << ") overflow" << nl << exit(FatalError);
395  }
396 
398 }
399 
400 
401 template<class Type>
402 const Foam::scalarField&
404 {
405  const label n = this->size();
406 
407  if (n <= 1)
408  {
410  << "table has (" << n << ") columns" << nl << exit(FatalError);
411  }
412  else if (i < 0)
413  {
415  << "index (" << i << ") underflow" << nl << exit(FatalError);
416  }
417  else if (i >= n)
418  {
420  << "index (" << i << ") overflow" << nl
421  << exit(FatalError);
422  }
423 
425 }
426 
427 
428 template<class Type>
430 {
431  return fieldIndices_.found(fieldName);
432 }
433 
434 
435 template<class Type>
436 const Foam::scalarList&
438 {
439  const label lo = index(retvals);
440  findHi(lo, retvals);
441  return interpOutput_;
442 }
443 
444 
445 template<class Type>
447 (
448  const label lo,
449  const scalar retvals
450 )
451 {
452  forAll(outputIndices_,j)
453  {
454  scalar tmp = 0;
455  label ofield = outputIndices_[j];
456  scalar baseValue = List<scalarField>::operator[](ofield).operator[](lo);
457 
458  forAll(entryIndices_,i)
459  {
460  if (checkRange(retvals, entryIndices_[i]))
461  {
462  label dim = 1;
463 
464  label hi = Foam::min(lo + dim, (*this)[0].size() - 1);
465 
466  tmp += interpolate(lo, hi, retvals, ofield, entryIndices_[i])
467  - baseValue;
468  }
469  interpOutput_[entryIndices_[i]] = retvals;
470  }
471 
472  tmp += baseValue;
473  interpOutput_[outputIndices_[j]] = tmp;
474  }
475 }
476 
477 
478 // ************************************************************************* //
#define readScalar
Definition: doubleScalar.C:38
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:428
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:75
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
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:137
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:319
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
IOdictionary is derived from dictionary and IOobject to give the dictionary automatic IO functionalit...
Definition: IOdictionary.H:52
void write(Ostream &, const fileName &, const word &instance, const objectRegistry &) const
Write lookup table to filename.
stressControl lookup("compactNormalStress") >> compactNormalStress
A class for handling words, derived from string.
Definition: word.H:59
A list of lists. Interpolates based on the first dimension. The values must be positive and monotonic...
const scalarField & operator[](const label) const
Return an element of constant List<scalar, Type>
const List< scalar > & lookUp(const scalar)
Return the output list given a single input scalar.
label readLabel(Istream &is)
Definition: label.H:64
bool found(const word &fieldName) const
Return true if the field exists in the table.
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
dimensioned< Type > min(const dimensioned< Type > &, const dimensioned< Type > &)
static tmp< GeometricField< Type, fvsPatchField, surfaceMesh > > interpolate(const GeometricField< Type, fvPatchField, volMesh > &tvf, const surfaceScalarField &faceFlux, Istream &schemeData)
Interpolate field onto faces using scheme given by Istream.
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.
IOobject defines the attributes of an object for which implicit objectRegistry management is supporte...
Definition: IOobject.H:92
ITstream & lookup(const word &, bool recursive=false, bool patternMatch=true) const
Find and return an entry data stream.
Definition: dictionary.C:576