polynomialFunction.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-2015 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 #include "polynomialFunction.H"
28 
29 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
30 
31 namespace Foam
32 {
33  defineTypeNameAndDebug(polynomialFunction, 0);
34 }
35 
36 
37 // * * * * * * * * * * * * * Static Member Functions * * * * * * * * * * * * //
38 
39 
40 Foam::polynomialFunction Foam::polynomialFunction::cloneIntegral
41 (
42  const polynomialFunction& poly,
43  const scalar intConstant
44 )
45 {
46  polynomialFunction newPoly(poly.size()+1);
47 
48  newPoly[0] = intConstant;
49  forAll(poly, i)
50  {
51  newPoly[i+1] = poly[i]/(i + 1);
52  }
53 
54  return newPoly;
55 }
56 
57 
58 Foam::polynomialFunction Foam::polynomialFunction::cloneIntegralMinus1
59 (
60  const polynomialFunction& poly,
61  const scalar intConstant
62 )
63 {
64  polynomialFunction newPoly(poly.size()+1);
65 
66  if (poly[0] > VSMALL)
67  {
68  newPoly.logActive_ = true;
69  newPoly.logCoeff_ = poly[0];
70  }
71 
72  newPoly[0] = intConstant;
73  for (label i=1; i < poly.size(); ++i)
74  {
75  newPoly[i] = poly[i]/i;
76  }
77 
78  return newPoly;
79 }
80 
81 
82 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
83 
85 :
86  scalarList(order, 0.0),
87  logActive_(false),
88  logCoeff_(0.0)
89 {
90  if (this->empty())
91  {
93  << "polynomialFunction coefficients are invalid (empty)"
94  << nl << exit(FatalError);
95  }
96 }
97 
98 
100 :
101  scalarList(poly),
102  logActive_(poly.logActive_),
103  logCoeff_(poly.logCoeff_)
104 {}
105 
106 
108 :
109  scalarList(coeffs),
110  logActive_(false),
111  logCoeff_(0.0)
112 {
113  if (this->empty())
114  {
116  << "polynomialFunction coefficients are invalid (empty)"
117  << nl << exit(FatalError);
118  }
119 }
120 
121 
123 :
124  scalarList(is),
125  logActive_(false),
126  logCoeff_(0.0)
127 {
128  if (this->empty())
129  {
131  << "polynomialFunction coefficients are invalid (empty)"
132  << nl << exit(FatalError);
133  }
134 }
135 
136 
137 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
138 
140 {}
141 
142 
143 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
144 
146 {
147  return logActive_;
148 }
149 
150 
152 {
153  return logCoeff_;
154 }
155 
156 
157 Foam::scalar Foam::polynomialFunction::value(const scalar x) const
158 {
159  const scalarList& coeffs = *this;
160  scalar val = coeffs[0];
161 
162  // avoid costly pow() in calculation
163  scalar powX = x;
164  for (label i=1; i<coeffs.size(); ++i)
165  {
166  val += coeffs[i]*powX;
167  powX *= x;
168  }
169 
170  if (logActive_)
171  {
172  val += this->logCoeff_*log(x);
173  }
174 
175  return val;
176 }
177 
178 
180 (
181  const scalar x1,
182  const scalar x2
183 ) const
184 {
185  const scalarList& coeffs = *this;
186 
187  if (logActive_)
188  {
190  << "Cannot integrate polynomial with logarithmic coefficients"
191  << nl << abort(FatalError);
192  }
193 
194  // avoid costly pow() in calculation
195  scalar powX1 = x1;
196  scalar powX2 = x2;
197 
198  scalar val = coeffs[0]*(powX2 - powX1);
199  for (label i=1; i<coeffs.size(); ++i)
200  {
201  val += coeffs[i]/(i + 1)*(powX2 - powX1);
202  powX1 *= x1;
203  powX2 *= x2;
204  }
205 
206  return val;
207 }
208 
209 
211 Foam::polynomialFunction::integral(const scalar intConstant) const
212 {
213  return cloneIntegral(*this, intConstant);
214 }
215 
216 
218 Foam::polynomialFunction::integralMinus1(const scalar intConstant) const
219 {
220  return cloneIntegralMinus1(*this, intConstant);
221 }
222 
223 
224 // * * * * * * * * * * * * * * Member Operators * * * * * * * * * * * * * * //
225 
228 {
229  scalarList& coeffs = *this;
230 
231  if (coeffs.size() > poly.size())
232  {
233  forAll(poly, i)
234  {
235  coeffs[i] += poly[i];
236  }
237  }
238  else
239  {
240  coeffs.setSize(poly.size(), 0.0);
241 
242  forAll(coeffs, i)
243  {
244  coeffs[i] += poly[i];
245  }
246  }
247 
248  return *this;
249 }
250 
251 
254 {
255  scalarList& coeffs = *this;
256 
257  if (coeffs.size() > poly.size())
258  {
259  forAll(poly, i)
260  {
261  coeffs[i] -= poly[i];
262  }
263  }
264  else
265  {
266  coeffs.setSize(poly.size(), 0.0);
267 
268  forAll(coeffs, i)
269  {
270  coeffs[i] -= poly[i];
271  }
272  }
273 
274  return *this;
275 }
276 
277 
280 {
281  scalarList& coeffs = *this;
282  forAll(coeffs, i)
283  {
284  coeffs[i] *= s;
285  }
286 
287  return *this;
288 }
289 
290 
293 {
294  scalarList& coeffs = *this;
295  forAll(coeffs, i)
296  {
297  coeffs[i] /= s;
298  }
299 
300  return *this;
301 }
302 
303 
304 // * * * * * * * * * * * * * * * IOstream Operators * * * * * * * * * * * * //
305 
307 {
308  // output like VectorSpace
309  os << token::BEGIN_LIST;
310 
311  if (!poly.empty())
312  {
313  for (int i=0; i<poly.size()-1; i++)
314  {
315  os << poly[i] << token::SPACE;
316  }
317  os << poly.last();
318  }
319  os << token::END_LIST;
320 
321 
322  // Check state of Ostream
323  os.check("operator<<(Ostream&, const polynomialFunction&)");
324 
325  return os;
326 }
327 
328 
329 // * * * * * * * * * * * * * * * Global Operators * * * * * * * * * * * * * //
330 
332 Foam::operator+
333 (
334  const polynomialFunction& p1,
335  const polynomialFunction& p2
336 )
337 {
338  polynomialFunction poly(p1);
339  return poly += p2;
340 }
341 
342 
344 Foam::operator-
345 (
346  const polynomialFunction& p1,
347  const polynomialFunction& p2
348 )
349 {
350  polynomialFunction poly(p1);
351  return poly -= p2;
352 }
353 
354 
356 Foam::operator*
357 (
358  const scalar s,
359  const polynomialFunction& p
360 )
361 {
362  polynomialFunction poly(p);
363  return poly *= s;
364 }
365 
366 
368 Foam::operator/
369 (
370  const scalar s,
371  const polynomialFunction& p
372 )
373 {
374  polynomialFunction poly(p);
375  return poly /= s;
376 }
377 
378 
380 Foam::operator*
381 (
382  const polynomialFunction& p,
383  const scalar s
384 )
385 {
386  polynomialFunction poly(p);
387  return poly *= s;
388 }
389 
390 
392 Foam::operator/
393 (
394  const polynomialFunction& p,
395  const scalar s
396 )
397 {
398  polynomialFunction poly(p);
399  return poly /= s;
400 }
401 
402 
403 // ************************************************************************* //
polynomialFunction & operator/=(const scalar)
polynomialFunction & operator+=(const polynomialFunction &)
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:428
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
dimensionedScalar log(const dimensionedScalar &ds)
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:124
error FatalError
scalar value(const scalar x) const
Return polynomial value.
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:319
bool empty() const
Return true if the UList is empty (ie, size() is zero)
Definition: UListI.H:313
An Istream is an abstract base class for all input systems (streams, files, token lists etc)...
Definition: Istream.H:57
void size(const label)
Override size to be inconsistent with allocated storage.
Definition: ListI.H:76
virtual bool check(const char *operation) const
Check IOstream status for given operation.
Definition: IOstream.C:92
Macros for easy insertion into run-time selection tables.
polynomialFunction(const label)
Construct a particular size, with all coefficients = 0.0.
scalar integrate(const scalar x1, const scalar x2) const
Integrate between two values.
gmvFile<< "tracers "<< particles.size()<< nl;forAllConstIter(Cloud< passiveParticle >, particles, iter){gmvFile<< iter().position().x()<< " ";}gmvFile<< nl;forAllConstIter(Cloud< passiveParticle >, particles, iter){gmvFile<< iter().position().y()<< " ";}gmvFile<< nl;forAllConstIter(Cloud< passiveParticle >, particles, iter){gmvFile<< iter().position().z()<< " ";}gmvFile<< nl;forAll(lagrangianScalarNames, i){word name=lagrangianScalarNames[i];IOField< scalar > s(IOobject(name, runTime.timeName(), cloud::prefix, mesh, IOobject::MUST_READ, IOobject::NO_WRITE))
scalar logCoeff() const
Return the log coefficient.
polynomialFunction & operator*=(const scalar)
Polynomial function representation.
polynomialFunction integral(const scalar intConstant=0.0) const
Return integral coefficients.
polynomialFunction integralMinus1(const scalar intConstant=0.0) const
Return integral coefficients when lowest order is -1.
errorManip< error > abort(error &err)
Definition: errorManip.H:131
A 1D vector of objects of type <T>, where the size of the vector is known and can be used for subscri...
Definition: HashTable.H:60
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
bool logActive() const
Return true if the log term is active.
defineTypeNameAndDebug(combustionModel, 0)
void setSize(const label)
Reset size of List.
Definition: List.C:295
polynomialFunction & operator-=(const polynomialFunction &)
Ostream & operator<<(Ostream &, const ensightPart &)
virtual ~polynomialFunction()
Destructor.
volScalarField & p
T & last()
Return the last element of the list.
Definition: UListI.H:128
Namespace for OpenFOAM.