uniformInterpolationTable.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-2018 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 #include "IOdictionary.H"
29 
30 // * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * * //
31 
32 template<class Type>
34 {
35  if (size() < 2)
36  {
38  << "Table " << name() << ": must have at least 2 values." << nl
39  << "Table size = " << size() << nl
40  << " min, interval width = " << x0_ << ", " << dx_ << nl
41  << exit(FatalError);
42  }
43 }
44 
45 
46 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
47 
48 template<class Type>
50 (
51  const IOobject& io,
52  bool readFields
53 )
54 :
55  IOobject(io),
56  List<scalar>(2, 0.0),
57  x0_(0.0),
58  dx_(1.0),
59  log10_(false),
60  bound_(false)
61 {
62  if (readFields)
63  {
64  IOdictionary dict(io);
65 
66  dict.lookup("data") >> *this;
67  dict.lookup("x0") >> x0_;
68  dict.lookup("dx") >> dx_;
69  dict.readIfPresent("log10", log10_);
70  dict.readIfPresent("bound", bound_);
71  }
72 
73  checkTable();
74 }
75 
76 
77 template<class Type>
79 (
80  const word& tableName,
81  const objectRegistry& db,
82  const dictionary& dict,
83  const bool initialiseOnly
84 )
85 :
86  IOobject
87  (
88  tableName,
89  db.time().constant(),
90  db,
91  IOobject::NO_READ,
92  IOobject::NO_WRITE,
93  false // if used in BCs, could be used by multiple patches
94  ),
95  List<scalar>(2, 0.0),
96  x0_(readScalar(dict.lookup("x0"))),
97  dx_(readScalar(dict.lookup("dx"))),
98  log10_(dict.lookupOrDefault<Switch>("log10", false)),
99  bound_(dict.lookupOrDefault<Switch>("bound", false))
100 {
101  if (initialiseOnly)
102  {
103  const scalar xMax = readScalar(dict.lookup("xMax"));
104  const label nIntervals = static_cast<label>(xMax - x0_)/dx_ + 1;
105  this->setSize(nIntervals);
106  }
107  else
108  {
109  dict.lookup("data") >> *this;
110  }
111 
112  checkTable();
113 }
114 
115 
116 template<class Type>
118 (
119  const uniformInterpolationTable& uit
120 )
121 :
122  IOobject(uit),
123  List<scalar>(uit),
124  x0_(uit.x0_),
125  dx_(uit.dx_),
126  log10_(uit.log10_),
127  bound_(uit.bound_)
128 {
129  checkTable();
130 }
131 
132 
133 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
134 
135 template<class Type>
137 {}
138 
139 
140 // * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * * //
141 
142 template<class Type>
144 {
145  if (bound_)
146  {
147  x = max(min(xMax() - small*dx_, x), x0_);
148  }
149  else
150  {
151  if (x < x0_)
152  {
154  << "Supplied value is less than minimum table value:" << nl
155  << "xMin=" << x0_ << ", xMax=" << xMax() << ", x=" << x << nl
156  << exit(FatalError);
157  }
158 
159  if (x > xMax())
160  {
162  << "Supplied value is greater than maximum table value:" << nl
163  << "xMin=" << x0_ << ", xMax=" << xMax() << ", x=" << x << nl
164  << exit(FatalError);
165  }
166  }
167 
168  const label i = static_cast<label>((x - x0_)/dx_);
169 
170  const scalar xLo = x0_ + i*dx_;
171 
172  Type fx = (x - xLo)/dx_*(operator[](i+1) - operator[](i)) + operator[](i);
173 
174  if (debug)
175  {
176  Info<< "Table: " << name() << ", x=" << x
177  << ", x_lo=" << xLo << ", x_hi=" << xLo + dx_
178  << ", f(x_lo)=" << operator[](i) << ", f(x_hi)=" << operator[](i+1)
179  << ", f(x)=" << fx << endl;
180  }
181 
182  return fx;
183 }
184 
185 
186 template<class Type>
188 (
189  scalar x
190 ) const
191 {
192  if (log10_)
193  {
194  if (x > 0)
195  {
196  x = ::log10(x);
197  }
198  else if (bound_ && (x <= 0))
199  {
200  x = x0_;
201  }
202  else
203  {
205  << "Table " << name() << nl
206  << "Supplied value must be greater than 0 when in log10 mode"
207  << nl << "x=" << x << nl << exit(FatalError);
208  }
209  }
210 
211  return interpolate(x);
212 }
213 
214 
215 template<class Type>
217 {
218  IOdictionary dict(*this);
219 
220  dict.add("data", static_cast<const List<scalar>&>(*this));
221  dict.add("x0", x0_);
222  dict.add("dx", dx_);
223  if (log10_)
224  {
225  dict.add("log10", log10_);
226  }
227  if (bound_)
228  {
229  dict.add("bound", bound_);
230  }
231 
232  dict.regIOobject::write();
233 }
234 
235 
236 // ************************************************************************* //
#define readScalar
Definition: doubleScalar.C:38
const label nIntervals(readLabel(pdfDictionary.lookup("nIntervals")))
dictionary dict
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
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:124
void readFields(const typename GeoFieldType::Mesh &mesh, const IOobjectList &objects, const HashSet< word > &selectedFields, LIFOStack< regIOobject *> &storedObjects)
Read the selected GeometricFields of the specified type.
Definition: ReadFields.C:244
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 > &)
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:319
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:256
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
Table with uniform interval in independent variable, with linear interpolation.
uniformInterpolationTable(const IOobject &, const bool readFields)
Construct from IOobject and readFields flag.
bool add(entry *, bool mergeEntry=false)
Add a new entry.
Definition: dictionary.C:814
IOdictionary is derived from dictionary and IOobject to give the dictionary automatic IO functionalit...
Definition: IOdictionary.H:52
points setSize(newPointi)
Type interpolateLog10(scalar x) const
Interpolate - takes log10 flag into account.
A class for handling words, derived from string.
Definition: word.H:59
const word & constant() const
Return constant name.
Definition: TimePaths.H:124
bool readIfPresent(const word &, T &, bool recursive=false, bool patternMatch=true) const
Find an entry if present, and assign to T.
static const char nl
Definition: Ostream.H:265
const Time & time() const
Return time.
dimensioned< Type > min(const dimensioned< Type > &, const dimensioned< Type > &)
word name(const complex &)
Return a string representation of a complex.
Definition: complex.C:47
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.
const scalar xMax
Definition: createFields.H:35
T lookupOrDefault(const word &, const T &, bool recursive=false, bool patternMatch=true) const
Find and return a T,.
messageStream Info
Registry of regIOobjects.
IOobject defines the attributes of an object for which implicit objectRegistry management is supporte...
Definition: IOobject.H:92
dimensionedScalar log10(const dimensionedScalar &ds)
ITstream & lookup(const word &, bool recursive=false, bool patternMatch=true) const
Find and return an entry data stream.
Definition: dictionary.C:576
Type interpolate(scalar x) const
Interpolate.