uniformInterpolationTable.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 
27 #include "Time.H"
28 
29 // * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * * //
30 
31 template<class Type>
33 {
34  if (size() < 2)
35  {
36  FatalErrorIn("uniformInterpolationTable<Type>::checkTable()")
37  << "Table " << name() << ": must have at least 2 values." << nl
38  << "Table size = " << size() << nl
39  << " min, interval width = " << x0_ << ", " << dx_ << nl
40  << exit(FatalError);
41  }
42 }
43 
44 
45 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
46 
47 template<class Type>
49 (
50  const IOobject& io,
51  bool readFields
52 )
53 :
54  IOobject(io),
55  List<scalar>(2, 0.0),
56  x0_(0.0),
57  dx_(1.0),
58  log10_(false),
59  bound_(false)
60 {
61  if (readFields)
62  {
63  IOdictionary dict(io);
64 
65  dict.lookup("data") >> *this;
66  dict.lookup("x0") >> x0_;
67  dict.lookup("dx") >> dx_;
68  dict.readIfPresent("log10", log10_);
69  dict.readIfPresent("bound", bound_);
70  }
71 
72  checkTable();
73 }
74 
75 
76 template<class Type>
78 (
79  const word& tableName,
80  const objectRegistry& db,
81  const dictionary& dict,
82  const bool initialiseOnly
83 )
84 :
85  IOobject
86  (
87  tableName,
88  db.time().constant(),
89  db,
90  IOobject::NO_READ,
91  IOobject::NO_WRITE,
92  false // if used in BCs, could be used by multiple patches
93  ),
94  List<scalar>(2, 0.0),
95  x0_(readScalar(dict.lookup("x0"))),
96  dx_(readScalar(dict.lookup("dx"))),
97  log10_(dict.lookupOrDefault<Switch>("log10", false)),
98  bound_(dict.lookupOrDefault<Switch>("bound", false))
99 {
100  if (initialiseOnly)
101  {
102  const scalar xMax = readScalar(dict.lookup("xMax"));
103  const label nIntervals = static_cast<label>(xMax - x0_)/dx_ + 1;
104  this->setSize(nIntervals);
105  }
106  else
107  {
108  dict.lookup("data") >> *this;
109  }
110 
111  checkTable();
112 }
113 
114 
115 template<class Type>
117 (
118  const uniformInterpolationTable& uit
119 )
120 :
121  IOobject(uit),
122  List<scalar>(uit),
123  x0_(uit.x0_),
124  dx_(uit.dx_),
125  log10_(uit.log10_),
126  bound_(uit.bound_)
127 {
128  checkTable();
129 }
130 
131 
132 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
133 
134 template<class Type>
136 {}
137 
138 
139 // * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * * //
140 
141 template<class Type>
143 {
144  if (bound_)
145  {
146  x = max(min(xMax() - SMALL*dx_, x), x0_);
147  }
148  else
149  {
150  if (x < x0_)
151  {
153  (
154  "uniformInterpolationTable<Type>::interpolate(scalar x)"
155  ) << "Supplied value is less than minimum table value:" << nl
156  << "xMin=" << x0_ << ", xMax=" << xMax() << ", x=" << x << nl
157  << exit(FatalError);
158  }
159 
160  if (x > xMax())
161  {
163  (
164  "uniformInterpolationTable<Type>::interpolate(scalar x)"
165  ) << "Supplied value is greater than maximum table value:" << nl
166  << "xMin=" << x0_ << ", xMax=" << xMax() << ", x=" << x << nl
167  << exit(FatalError);
168  }
169  }
170 
171  const label i = static_cast<label>((x - x0_)/dx_);
172 
173  const scalar xLo = x0_ + i*dx_;
174 
175  Type fx = (x - xLo)/dx_*(operator[](i+1) - operator[](i)) + operator[](i);
176 
177  if (debug)
178  {
179  Info<< "Table: " << name() << ", x=" << x
180  << ", x_lo=" << xLo << ", x_hi=" << xLo + dx_
181  << ", f(x_lo)=" << operator[](i) << ", f(x_hi)=" << operator[](i+1)
182  << ", f(x)=" << fx << endl;
183  }
184 
185  return fx;
186 }
187 
188 
189 template<class Type>
191 (
192  scalar x
193 ) const
194 {
195  if (log10_)
196  {
197  if (x > 0)
198  {
199  x = ::log10(x);
200  }
201  else if (bound_ && (x <= 0))
202  {
203  x = x0_;
204  }
205  else
206  {
208  (
209  "uniformInterpolationTable<Type>::interpolateLog10(scalar x)"
210  ) << "Table " << name() << nl
211  << "Supplied value must be greater than 0 when in log10 mode"
212  << nl << "x=" << x << nl << exit(FatalError);
213  }
214  }
215 
216  return interpolate(x);
217 }
218 
219 
220 template<class Type>
222 {
223  IOdictionary dict(*this);
224 
225  dict.add("data", static_cast<const List<scalar>&>(*this));
226  dict.add("x0", x0_);
227  dict.add("dx", dx_);
228  if (log10_)
229  {
230  dict.add("log10", log10_);
231  }
232  if (bound_)
233  {
234  dict.add("bound", bound_);
235  }
236 
237  dict.regIOobject::write();
238 }
239 
240 
241 // ************************************************************************* //
#define readScalar
Definition: doubleScalar.C:38
word name(const complex &)
Return a string representation of a complex.
Definition: complex.C:47
Type interpolateLog10(scalar x) const
Interpolate - takes log10 flag into account.
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
uniformInterpolationTable(const IOobject &, const bool readFields)
Construct from IOobject and readFields flag.
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:124
dimensionedScalar log10(const dimensionedScalar &ds)
A simple wrapper around bool so that it can be read as a word: true/false, on/off, yes/no, y/n, t/f, or none.
Definition: Switch.H:60
Type interpolate(scalar x) const
Interpolate.
messageStream Info
tmp< surfaceScalarField > interpolate(const RhoType &rho)
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:137
const word & constant() const
Return constant name.
Definition: TimePaths.H:124
dictionary dict
static const char nl
Definition: Ostream.H:260
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:251
IOobject defines the attributes of an object for which implicit objectRegistry management is supporte...
Definition: IOobject.H:91
dimensioned< Type > max(const dimensioned< Type > &, const dimensioned< Type > &)
IOdictionary is derived from dictionary and IOobject to give the dictionary automatic IO functionalit...
Definition: IOdictionary.H:53
Table with uniform interval in independant variable, with linear interpolation.
const label nIntervals(readLabel(pdfDictionary.lookup("nIntervals")))
bool add(entry *, bool mergeEntry=false)
Add a new entry.
Definition: dictionary.C:739
#define FatalErrorIn(functionName)
Report an error message using Foam::FatalError.
Definition: error.H:314
points setSize(newPointi)
ITstream & lookup(const word &, bool recursive=false, bool patternMatch=true) const
Find and return an entry data stream.
Definition: dictionary.C:452
error FatalError
T lookupOrDefault(const word &, const T &, bool recursive=false, bool patternMatch=true) const
Find and return a T,.
Registry of regIOobjects.
bool readIfPresent(const word &, T &, bool recursive=false, bool patternMatch=true) const
Find an entry if present, and assign to T.
const Time & time() const
Return time.
const scalar xMax
Definition: createFields.H:35
This function object reads fields from the time directories and adds them to the mesh database for fu...
Definition: readFields.H:102
dimensioned< Type > min(const dimensioned< Type > &, const dimensioned< Type > &)