From ffd20770ddf3a9e236d2217ad62b41c6430edc34 Mon Sep 17 00:00:00 2001
From: Mark Olesen <Mark.Olesen@Germany>
Date: Fri, 4 Feb 2011 18:59:45 +0100
Subject: [PATCH] ENH: add (non-templated) polynomialFunction
STYLE: integrateLimits() -> integrate() for consistency with
DataEntry/polynomial. Rename old 'integrate', 'integrateMinus1'
methods to 'integral', 'integralMinus1' (ie, substantive instead of
verb for clarity).
---
.../test/Polynomial/Test-Polynomial.C | 65 ++-
src/OpenFOAM/Make/files | 1 +
.../functions/Polynomial/Polynomial.C | 109 +++--
.../functions/Polynomial/Polynomial.H | 45 +-
.../functions/Polynomial/polynomialFunction.C | 415 ++++++++++++++++++
.../functions/Polynomial/polynomialFunction.H | 246 +++++++++++
.../thermo/hPolynomial/hPolynomialThermo.C | 8 +-
7 files changed, 808 insertions(+), 81 deletions(-)
create mode 100644 src/OpenFOAM/primitives/functions/Polynomial/polynomialFunction.C
create mode 100644 src/OpenFOAM/primitives/functions/Polynomial/polynomialFunction.H
diff --git a/applications/test/Polynomial/Test-Polynomial.C b/applications/test/Polynomial/Test-Polynomial.C
index fbfbcc20006..c251bade747 100644
--- a/applications/test/Polynomial/Test-Polynomial.C
+++ b/applications/test/Polynomial/Test-Polynomial.C
@@ -2,7 +2,7 @@
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
- \\ / A nd | Copyright (C) 2009-2010 OpenCFD Ltd.
+ \\ / A nd | Copyright (C) 2009-2011 OpenCFD Ltd.
\\/ M anipulation |
-------------------------------------------------------------------------------
License
@@ -31,6 +31,7 @@ Description
#include "IStringStream.H"
#include "Polynomial.H"
+#include "polynomialFunction.H"
#include "Random.H"
#include "cpuTime.H"
@@ -38,7 +39,7 @@ using namespace Foam;
const int PolySize = 8;
const scalar coeff[] = { 0.11, 0.45, -0.94, 1.58, -2.58, 0.08, 3.15, -4.78 };
-const char* polyDef = "testPoly (0.11 0.45 -0.94 1.58 -2.58 0.08 3.15 -4.78)";
+const char* polyDef = "(0.11 0.45 -0.94 1.58 -2.58 0.08 3.15 -4.78)";
scalar polyValue(const scalar x)
@@ -58,7 +59,7 @@ scalar polyValue(const scalar x)
scalar intPolyValue(const scalar x)
{
- // Hard-coded integrated form of above polynomial
+ // Hard-coded integral form of above polynomial
return
coeff[0]*x
+ coeff[1]/2.0*sqr(x)
@@ -109,27 +110,44 @@ int main(int argc, char *argv[])
const label nIters = 1000;
scalar sum = 0.0;
- Info<< "null poly = " << (Polynomial<8>()) << nl << endl;
+ Info<< "null poly = " << (Polynomial<8>()) << nl
+ << "null poly = " << (polynomialFunction(8)) << nl
+ << endl;
- // Polynomial<8> poly("testPoly", IStringStream(polyDef)());
Polynomial<8> poly(coeff);
- Polynomial<9> intPoly(poly.integrate(0.0));
+ Polynomial<9> intPoly(poly.integral(0.0));
- Info<< "poly = " << poly << endl;
- Info<< "intPoly = " << intPoly << nl << endl;
+ IStringStream is(polyDef);
+ polynomialFunction polyfunc(is);
- Info<< "2*poly = " << 2*poly << endl;
- Info<< "poly+poly = " << poly + poly << nl << endl;
+ Info<< "poly = " << poly << nl
+ << "intPoly = " << intPoly << nl
+ << endl;
- Info<< "3*poly = " << 3*poly << endl;
- Info<< "poly+poly+poly = " << poly + poly + poly << nl << endl;
+ Info<< "2*poly = " << 2*poly << nl
+ << "poly+poly = " << poly + poly << nl
+ << "3*poly = " << 3*poly << nl
+ << "poly+poly+poly = " << poly + poly + poly << nl
+ << "3*poly - 2*poly = " << 3*poly - 2*poly << nl
+ << endl;
+
+ Info<< nl << "--- as polynomialFunction" << nl << endl;
+ Info<< "polyf = " << polyfunc << nl
+ << "intPoly = " << poly.integral(0.0) << nl
+ << endl;
+
+ Info<< "2*polyf = " << 2*polyfunc << nl
+ << "polyf+polyf = " << polyfunc + polyfunc << nl
+ << "3*polyf = " << 3*polyfunc << nl
+ << "polyf+polyf+polyf = " << polyfunc + polyfunc + polyfunc << nl
+ << "3*polyf - 2*polyf = " << 3*polyfunc - 2*polyfunc << nl
+ << endl;
- Info<< "3*poly - 2*poly = " << 3*poly - 2*poly << nl << endl;
Polynomial<8> polyCopy = poly;
Info<< "poly, polyCopy = " << poly << ", " << polyCopy << nl << endl;
polyCopy = 2.5*poly;
- Info<< "2.5*polyCopy = " << polyCopy << nl << endl;
+ Info<< "2.5*poly = " << polyCopy << nl << endl;
Random rnd(123456);
for (int i=0; i<10; i++)
@@ -139,8 +157,8 @@ int main(int argc, char *argv[])
scalar px = polyValue(x);
scalar ipx = intPolyValue(x);
- scalar pxTest = poly.evaluate(x);
- scalar ipxTest = intPoly.evaluate(x);
+ scalar pxTest = poly.value(x);
+ scalar ipxTest = intPoly.value(x);
Info<<"\nx = " << x << endl;
Info<< " px, pxTest = " << px << ", " << pxTest << endl;
@@ -173,10 +191,21 @@ int main(int argc, char *argv[])
sum = 0.0;
for (label iter = 0; iter < nIters; ++iter)
{
- sum += poly.evaluate(loop+iter);
+ sum += poly.value(loop+iter);
+ }
+ }
+ Info<< "value: " << sum
+ << " in " << timer.cpuTimeIncrement() << " s\n";
+
+ for (int loop = 0; loop < n; ++loop)
+ {
+ sum = 0.0;
+ for (label iter = 0; iter < nIters; ++iter)
+ {
+ sum += polyfunc.value(loop+iter);
}
}
- Info<< "evaluate: " << sum
+ Info<< "via function: " << sum
<< " in " << timer.cpuTimeIncrement() << " s\n";
diff --git a/src/OpenFOAM/Make/files b/src/OpenFOAM/Make/files
index 7e293fc7eb1..ba3c157b106 100644
--- a/src/OpenFOAM/Make/files
+++ b/src/OpenFOAM/Make/files
@@ -66,6 +66,7 @@ primitives/functions/DataEntry/makeDataEntries.C
primitives/functions/DataEntry/polynomial/polynomial.C
primitives/functions/DataEntry/polynomial/polynomialIO.C
+primitives/functions/Polynomial/polynomialFunction.C
strings = primitives/strings
$(strings)/string/string.C
diff --git a/src/OpenFOAM/primitives/functions/Polynomial/Polynomial.C b/src/OpenFOAM/primitives/functions/Polynomial/Polynomial.C
index 7f9ae731f6b..3c642de2d3a 100644
--- a/src/OpenFOAM/primitives/functions/Polynomial/Polynomial.C
+++ b/src/OpenFOAM/primitives/functions/Polynomial/Polynomial.C
@@ -41,6 +41,18 @@ Foam::Polynomial<PolySize>::Polynomial()
}
+template<int PolySize>
+Foam::Polynomial<PolySize>::Polynomial
+(
+ const Polynomial<PolySize>& poly
+)
+:
+ VectorSpace<Polynomial<PolySize>, scalar, PolySize>(poly),
+ logActive_(poly.logActive_),
+ logCoeff_(poly.logCoeff_)
+{}
+
+
template<int PolySize>
Foam::Polynomial<PolySize>::Polynomial(const scalar coeffs[PolySize])
:
@@ -68,7 +80,7 @@ Foam::Polynomial<PolySize>::Polynomial(const UList<scalar>& coeffs)
(
"Polynomial<PolySize>::Polynomial(const UList<scalar>&)"
) << "Size mismatch: Needed " << PolySize
- << " but got " << coeffs.size()
+ << " but given " << coeffs.size()
<< nl << exit(FatalError);
}
@@ -79,6 +91,39 @@ Foam::Polynomial<PolySize>::Polynomial(const UList<scalar>& coeffs)
}
+// template<int PolySize>
+// Foam::Polynomial<PolySize>::Polynomial(const polynomialFunction& poly)
+// :
+// VectorSpace<Polynomial<PolySize>, scalar, PolySize>(),
+// logActive_(poly.logActive()),
+// logCoeff_(poly.logCoeff())
+// {
+// if (poly.size() != PolySize)
+// {
+// FatalErrorIn
+// (
+// "Polynomial<PolySize>::Polynomial(const polynomialFunction&)"
+// ) << "Size mismatch: Needed " << PolySize
+// << " but given " << poly.size()
+// << nl << exit(FatalError);
+// }
+//
+// for (int i = 0; i < PolySize; ++i)
+// {
+// this->v_[i] = poly[i];
+// }
+// }
+
+
+template<int PolySize>
+Foam::Polynomial<PolySize>::Polynomial(Istream& is)
+:
+ VectorSpace<Polynomial<PolySize>, scalar, PolySize>(is),
+ logActive_(false),
+ logCoeff_(0.0)
+{}
+
+
template<int PolySize>
Foam::Polynomial<PolySize>::Polynomial(const word& name, Istream& is)
:
@@ -111,38 +156,17 @@ Foam::Polynomial<PolySize>::Polynomial(const word& name, Istream& is)
}
-template<int PolySize>
-Foam::Polynomial<PolySize>::Polynomial(Istream& is)
-:
- VectorSpace<Polynomial<PolySize>, scalar, PolySize>(is),
- logActive_(false),
- logCoeff_(0.0)
-{}
-
-
-template<int PolySize>
-Foam::Polynomial<PolySize>::Polynomial
-(
- const Polynomial<PolySize>& poly
-)
-:
- VectorSpace<Polynomial<PolySize>, scalar, PolySize>(poly),
- logActive_(poly.logActive_),
- logCoeff_(poly.logCoeff_)
-{}
-
-
// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
template<int PolySize>
-bool& Foam::Polynomial<PolySize>::logActive()
+bool Foam::Polynomial<PolySize>::logActive() const
{
return logActive_;
}
template<int PolySize>
-Foam::scalar& Foam::Polynomial<PolySize>::logCoeff()
+Foam::scalar Foam::Polynomial<PolySize>::logCoeff() const
{
return logCoeff_;
}
@@ -151,27 +175,27 @@ Foam::scalar& Foam::Polynomial<PolySize>::logCoeff()
template<int PolySize>
Foam::scalar Foam::Polynomial<PolySize>::value(const scalar x) const
{
- scalar y = this->v_[0];
+ scalar val = this->v_[0];
// avoid costly pow() in calculation
scalar powX = x;
for (label i=1; i<PolySize; ++i)
{
- y += this->v_[i]*powX;
+ val += this->v_[i]*powX;
powX *= x;
}
if (logActive_)
{
- y += logCoeff_*log(x);
+ val += logCoeff_*log(x);
}
- return y;
+ return val;
}
template<int PolySize>
-Foam::scalar Foam::Polynomial<PolySize>::integrateLimits
+Foam::scalar Foam::Polynomial<PolySize>::integrate
(
const scalar x1,
const scalar x2
@@ -181,7 +205,7 @@ Foam::scalar Foam::Polynomial<PolySize>::integrateLimits
{
FatalErrorIn
(
- "scalar Polynomial<PolySize>::integrateLimits"
+ "scalar Polynomial<PolySize>::integrate"
"("
"const scalar, "
"const scalar"
@@ -190,22 +214,33 @@ Foam::scalar Foam::Polynomial<PolySize>::integrateLimits
<< nl << abort(FatalError);
}
- intPolyType poly = this->integrate();
- return poly.value(x2) - poly.value(x1);
+ // avoid costly pow() in calculation
+ scalar powX1 = x1;
+ scalar powX2 = x2;
+
+ scalar val = this->v_[0]*(powX2 - powX1);
+ for (label i=1; i<PolySize; ++i)
+ {
+ val += this->v_[i]/(i + 1) * (powX2 - powX1);
+ powX1 *= x1;
+ powX2 *= x2;
+ }
+
+ return val;
}
template<int PolySize>
typename Foam::Polynomial<PolySize>::intPolyType
-Foam::Polynomial<PolySize>::integrate(const scalar intConstant)
+Foam::Polynomial<PolySize>::integral(const scalar intConstant) const
{
intPolyType newCoeffs;
newCoeffs[0] = intConstant;
forAll(*this, i)
{
- newCoeffs[i + 1] = this->v_[i]/(i + 1);
+ newCoeffs[i+1] = this->v_[i]/(i + 1);
}
return newCoeffs;
@@ -214,14 +249,14 @@ Foam::Polynomial<PolySize>::integrate(const scalar intConstant)
template<int PolySize>
typename Foam::Polynomial<PolySize>::polyType
-Foam::Polynomial<PolySize>::integrateMinus1(const scalar intConstant)
+Foam::Polynomial<PolySize>::integralMinus1(const scalar intConstant) const
{
polyType newCoeffs;
if (this->v_[0] > VSMALL)
{
- newCoeffs.logActive() = true;
- newCoeffs.logCoeff() = this->v_[0];
+ newCoeffs.logActive_ = true;
+ newCoeffs.logCoeff_ = this->v_[0];
}
newCoeffs[0] = intConstant;
diff --git a/src/OpenFOAM/primitives/functions/Polynomial/Polynomial.H b/src/OpenFOAM/primitives/functions/Polynomial/Polynomial.H
index 06511b0003e..5dfa271c8de 100644
--- a/src/OpenFOAM/primitives/functions/Polynomial/Polynomial.H
+++ b/src/OpenFOAM/primitives/functions/Polynomial/Polynomial.H
@@ -29,14 +29,14 @@ Description
poly = logCoeff*log(x) + sum(coeff_[i]*x^i)
- where 0 \<= i \<= n
+ where 0 \<= i \<= N
- integer powers, starting at zero
- value(x) to evaluate the poly for a given value
- integrate(x1, x2) between two scalar values
- - integrate() to return a new, intergated coeff polynomial
+ - integral() to return a new, integral coeff polynomial
- increases the size (order)
- - integrateMinus1() to return a new, integrated coeff polynomial where
+ - integralMinus1() to return a new, integral coeff polynomial where
the base poly starts at order -1
SourceFiles
@@ -85,10 +85,10 @@ class Polynomial
// Private data
- //- Include the log term? - only activated using integrateMinus1()
+ //- Include the log term? - only activated using integralMinus1()
bool logActive_;
- //- Log coefficient - only activated using integrateMinus1()
+ //- Log coefficient - only activated using integralMinus1()
scalar logCoeff_;
@@ -104,6 +104,9 @@ public:
//- Construct null, with all coefficients = 0.0
Polynomial();
+ //- Copy constructor
+ Polynomial(const Polynomial&);
+
//- Construct from C-array of coefficients
explicit Polynomial(const scalar coeffs[PolySize]);
@@ -111,24 +114,21 @@ public:
explicit Polynomial(const UList<scalar>& coeffs);
//- Construct from Istream
- Polynomial(Istream& is);
+ Polynomial(Istream&);
//- Construct from name and Istream
- Polynomial(const word& name, Istream& is);
-
- //- Copy constructor
- Polynomial(const Polynomial& poly);
+ Polynomial(const word& name, Istream&);
// Member Functions
// Access
- //- Return access to the log term active flag
- bool& logActive();
+ //- Return true if the log term is active
+ bool logActive() const;
- //- Return access to the log coefficient
- scalar& logCoeff();
+ //- Return the log coefficient
+ scalar logCoeff() const;
// Evaluation
@@ -136,16 +136,17 @@ public:
//- Return polynomial value
scalar value(const scalar x) const;
- //- Return integrated polynomial coefficients
- // argument becomes zeroth element (constant of integration)
- intPolyType integrate(const scalar intConstant = 0.0);
+ //- Integrate between two values
+ scalar integrate(const scalar x1, const scalar x2) const;
- //- Return integrated polynomial coefficients when lowest order
- // is -1. Argument added to zeroth element
- polyType integrateMinus1(const scalar intConstant = 0.0);
- //- Integrate between two values
- scalar integrateLimits(const scalar x1, const scalar x2) const;
+ //- Return integral coefficients.
+ // Argument becomes zeroth element (constant of integration)
+ intPolyType integral(const scalar intConstant = 0.0) const;
+
+ //- Return integral coefficients when lowest order is -1.
+ // Argument becomes zeroth element (constant of integration)
+ polyType integralMinus1(const scalar intConstant = 0.0) const;
//- Ostream Operator
diff --git a/src/OpenFOAM/primitives/functions/Polynomial/polynomialFunction.C b/src/OpenFOAM/primitives/functions/Polynomial/polynomialFunction.C
new file mode 100644
index 00000000000..44e82ac26d6
--- /dev/null
+++ b/src/OpenFOAM/primitives/functions/Polynomial/polynomialFunction.C
@@ -0,0 +1,415 @@
+/*---------------------------------------------------------------------------*\
+ ========= |
+ \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
+ \\ / O peration |
+ \\ / A nd | Copyright (C) 2011-2011 OpenCFD Ltd.
+ \\/ M anipulation |
+-------------------------------------------------------------------------------
+License
+ This file is part of OpenFOAM.
+
+ OpenFOAM is free software: you can redistribute it and/or modify it
+ under the terms of the GNU General Public License as published by
+ the Free Software Foundation, either version 3 of the License, or
+ (at your option) any later version.
+
+ OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
+ ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
+ FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
+ for more details.
+
+ You should have received a copy of the GNU General Public License
+ along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
+
+\*---------------------------------------------------------------------------*/
+
+#include "polynomialFunction.H"
+#include "addToRunTimeSelectionTable.H"
+
+// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
+
+namespace Foam
+{
+ defineTypeNameAndDebug(polynomialFunction, 0);
+}
+
+
+// * * * * * * * * * * * * * Static Member Functions * * * * * * * * * * * * //
+
+
+Foam::polynomialFunction Foam::polynomialFunction::cloneIntegral
+(
+ const polynomialFunction& poly,
+ const scalar intConstant
+)
+{
+ polynomialFunction newPoly(poly.size()+1);
+
+ newPoly[0] = intConstant;
+ forAll(poly, i)
+ {
+ newPoly[i+1] = poly[i]/(i + 1);
+ }
+
+ return newPoly;
+}
+
+
+Foam::polynomialFunction Foam::polynomialFunction::cloneIntegralMinus1
+(
+ const polynomialFunction& poly,
+ const scalar intConstant
+)
+{
+ polynomialFunction newPoly(poly.size()+1);
+
+ if (poly[0] > VSMALL)
+ {
+ newPoly.logActive_ = true;
+ newPoly.logCoeff_ = poly[0];
+ }
+
+ newPoly[0] = intConstant;
+ for (label i=1; i < poly.size(); ++i)
+ {
+ newPoly[i] = poly[i]/i;
+ }
+
+ return newPoly;
+}
+
+
+// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
+
+Foam::polynomialFunction::polynomialFunction(const label order)
+:
+ scalarList(order, 0.0),
+ logActive_(false),
+ logCoeff_(0.0)
+{
+ if (this->empty())
+ {
+ FatalErrorIn
+ (
+ "polynomialFunction::polynomialFunction(const label order)"
+ ) << "polynomialFunction coefficients are invalid (empty)"
+ << nl << exit(FatalError);
+ }
+}
+
+
+Foam::polynomialFunction::polynomialFunction(const polynomialFunction& poly)
+:
+ scalarList(poly),
+ logActive_(poly.logActive_),
+ logCoeff_(poly.logCoeff_)
+{}
+
+
+Foam::polynomialFunction::polynomialFunction(const UList<scalar>& coeffs)
+:
+ scalarList(coeffs),
+ logActive_(false),
+ logCoeff_(0.0)
+{
+ if (this->empty())
+ {
+ FatalErrorIn
+ (
+ "polynomialFunction::polynomialFunction(const UList<scalar>&)"
+ ) << "polynomialFunction coefficients are invalid (empty)"
+ << nl << exit(FatalError);
+ }
+}
+
+
+Foam::polynomialFunction::polynomialFunction(Istream& is)
+:
+ scalarList(is),
+ logActive_(false),
+ logCoeff_(0.0)
+{
+ if (this->empty())
+ {
+ FatalErrorIn
+ (
+ "polynomialFunction::polynomialFunction(Istream&)"
+ ) << "polynomialFunction coefficients are invalid (empty)"
+ << nl << exit(FatalError);
+ }
+}
+
+
+// * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
+
+Foam::polynomialFunction::~polynomialFunction()
+{}
+
+
+// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
+
+bool Foam::polynomialFunction::logActive() const
+{
+ return logActive_;
+}
+
+
+Foam::scalar Foam::polynomialFunction::logCoeff() const
+{
+ return logCoeff_;
+}
+
+
+Foam::scalar Foam::polynomialFunction::value(const scalar x) const
+{
+ const scalarList& coeffs = *this;
+ scalar val = coeffs[0];
+
+ // avoid costly pow() in calculation
+ scalar powX = x;
+ for (label i=1; i<coeffs.size(); ++i)
+ {
+ val += coeffs[i]*powX;
+ powX *= x;
+ }
+
+ if (logActive_)
+ {
+ val += this->logCoeff_*log(x);
+ }
+
+ return val;
+}
+
+
+Foam::scalar Foam::polynomialFunction::integrate
+(
+ const scalar x1,
+ const scalar x2
+) const
+{
+ const scalarList& coeffs = *this;
+
+ if (logActive_)
+ {
+ FatalErrorIn
+ (
+ "scalar polynomialFunction::integrate"
+ "("
+ "const scalar, "
+ "const scalar"
+ ") const"
+ ) << "Cannot integrate polynomial with logarithmic coefficients"
+ << nl << abort(FatalError);
+ }
+
+ // avoid costly pow() in calculation
+ scalar powX1 = x1;
+ scalar powX2 = x2;
+
+ scalar val = coeffs[0]*(powX2 - powX1);
+ for (label i=1; i<coeffs.size(); ++i)
+ {
+ val += coeffs[i]/(i + 1)*(powX2 - powX1);
+ powX1 *= x1;
+ powX2 *= x2;
+ }
+
+ return val;
+}
+
+
+Foam::polynomialFunction
+Foam::polynomialFunction::integral(const scalar intConstant) const
+{
+ return cloneIntegral(*this, intConstant);
+}
+
+
+Foam::polynomialFunction
+Foam::polynomialFunction::integralMinus1(const scalar intConstant) const
+{
+ return cloneIntegralMinus1(*this, intConstant);
+}
+
+
+// * * * * * * * * * * * * * * Member Operators * * * * * * * * * * * * * * //
+
+Foam::polynomialFunction&
+Foam::polynomialFunction::operator+=(const polynomialFunction& poly)
+{
+ scalarList& coeffs = *this;
+
+ if (coeffs.size() > poly.size())
+ {
+ forAll(poly, i)
+ {
+ coeffs[i] += poly[i];
+ }
+ }
+ else
+ {
+ coeffs.setSize(poly.size(), 0.0);
+
+ forAll(coeffs, i)
+ {
+ coeffs[i] += poly[i];
+ }
+ }
+
+ return *this;
+}
+
+
+Foam::polynomialFunction&
+Foam::polynomialFunction::operator-=(const polynomialFunction& poly)
+{
+ scalarList& coeffs = *this;
+
+ if (coeffs.size() > poly.size())
+ {
+ forAll(poly, i)
+ {
+ coeffs[i] -= poly[i];
+ }
+ }
+ else
+ {
+ coeffs.setSize(poly.size(), 0.0);
+
+ forAll(coeffs, i)
+ {
+ coeffs[i] -= poly[i];
+ }
+ }
+
+ return *this;
+}
+
+
+Foam::polynomialFunction&
+Foam::polynomialFunction::operator*=(const scalar s)
+{
+ scalarList& coeffs = *this;
+ forAll(coeffs, i)
+ {
+ coeffs[i] *= s;
+ }
+
+ return *this;
+}
+
+
+Foam::polynomialFunction&
+Foam::polynomialFunction::operator/=(const scalar s)
+{
+ scalarList& coeffs = *this;
+ forAll(coeffs, i)
+ {
+ coeffs[i] /= s;
+ }
+
+ return *this;
+}
+
+
+// * * * * * * * * * * * * * * * IOstream Operators * * * * * * * * * * * * //
+
+Foam::Ostream& Foam::operator<<(Ostream& os, const polynomialFunction& poly)
+{
+ // output like VectorSpace
+ os << token::BEGIN_LIST;
+
+ if (!poly.empty())
+ {
+ for (int i=0; i<poly.size()-1; i++)
+ {
+ os << poly[i] << token::SPACE;
+ }
+ os << poly.last();
+ }
+ os << token::END_LIST;
+
+
+ // Check state of Ostream
+ os.check("operator<<(Ostream&, const polynomialFunction&)");
+
+ return os;
+}
+
+
+// * * * * * * * * * * * * * * * Global Operators * * * * * * * * * * * * * //
+
+Foam::polynomialFunction
+Foam::operator+
+(
+ const polynomialFunction& p1,
+ const polynomialFunction& p2
+)
+{
+ polynomialFunction poly(p1);
+ return poly += p2;
+}
+
+
+Foam::polynomialFunction
+Foam::operator-
+(
+ const polynomialFunction& p1,
+ const polynomialFunction& p2
+)
+{
+ polynomialFunction poly(p1);
+ return poly -= p2;
+}
+
+
+Foam::polynomialFunction
+Foam::operator*
+(
+ const scalar s,
+ const polynomialFunction& p
+)
+{
+ polynomialFunction poly(p);
+ return poly *= s;
+}
+
+
+Foam::polynomialFunction
+Foam::operator/
+(
+ const scalar s,
+ const polynomialFunction& p
+)
+{
+ polynomialFunction poly(p);
+ return poly /= s;
+}
+
+
+Foam::polynomialFunction
+Foam::operator*
+(
+ const polynomialFunction& p,
+ const scalar s
+)
+{
+ polynomialFunction poly(p);
+ return poly *= s;
+}
+
+
+Foam::polynomialFunction
+Foam::operator/
+(
+ const polynomialFunction& p,
+ const scalar s
+)
+{
+ polynomialFunction poly(p);
+ return poly /= s;
+}
+
+
+// ************************************************************************* //
diff --git a/src/OpenFOAM/primitives/functions/Polynomial/polynomialFunction.H b/src/OpenFOAM/primitives/functions/Polynomial/polynomialFunction.H
new file mode 100644
index 00000000000..b9ce2c7d090
--- /dev/null
+++ b/src/OpenFOAM/primitives/functions/Polynomial/polynomialFunction.H
@@ -0,0 +1,246 @@
+/*---------------------------------------------------------------------------*\
+ ========= |
+ \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
+ \\ / O peration |
+ \\ / A nd | Copyright (C) 2011-2011 OpenCFD Ltd.
+ \\/ M anipulation |
+-------------------------------------------------------------------------------
+License
+ This file is part of OpenFOAM.
+
+ OpenFOAM is free software: you can redistribute it and/or modify it
+ under the terms of the GNU General Public License as published by
+ the Free Software Foundation, either version 3 of the License, or
+ (at your option) any later version.
+
+ OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
+ ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
+ FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
+ for more details.
+
+ You should have received a copy of the GNU General Public License
+ along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
+
+Class
+ Foam::polynomialFunction
+
+Description
+ Polynomial function representation
+
+ poly = logCoeff*log(x) + sum(coeff_[i]*x^i)
+
+ where 0 \<= i \<= N
+
+ - integer powers, starting at zero
+ - value(x) to evaluate the poly for a given value
+ - integrate(x1, x2) between two scalar values
+ - integral() to return a new, integral coeff polynomial
+ - increases the size (order)
+ - integralMinus1() to return a new, integral coeff polynomial where
+ the base poly starts at order -1
+
+SeeAlso
+ Foam::Polynomial for a templated implementation
+
+SourceFiles
+ polynomialFunction.C
+
+\*---------------------------------------------------------------------------*/
+
+#ifndef polynomialFunction_H
+#define polynomialFunction_H
+
+#include "scalarList.H"
+#include "Ostream.H"
+#include "runTimeSelectionTables.H"
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+namespace Foam
+{
+
+// Forward declaration of classes
+class polynomialFunction;
+
+// Forward declaration of friend functions
+Ostream& operator<<(Ostream&, const polynomialFunction&);
+
+
+/*---------------------------------------------------------------------------*\
+ Class polynomialFunction Declaration
+\*---------------------------------------------------------------------------*/
+
+class polynomialFunction
+:
+ private scalarList
+{
+ // Private data
+
+ //- Include the log term? - only activated using integralMinus1()
+ bool logActive_;
+
+ //- Log coefficient - only activated using integralMinus1()
+ scalar logCoeff_;
+
+
+ // Private Member Functions
+
+ //- Return integral coefficients.
+ // Argument becomes zeroth element (constant of integration)
+ static polynomialFunction cloneIntegral
+ (
+ const polynomialFunction&,
+ const scalar intConstant = 0.0
+ );
+
+ //- Return integral coefficients when lowest order is -1.
+ // Argument becomes zeroth element (constant of integration)
+ static polynomialFunction cloneIntegralMinus1
+ (
+ const polynomialFunction&,
+ const scalar intConstant = 0.0
+ );
+
+
+ //- Disallow default bitwise assignment
+ void operator=(const polynomialFunction&);
+
+
+
+public:
+
+ //- Runtime type information
+ TypeName("polynomialFunction");
+
+
+ // Constructors
+
+ //- Construct a particular size, with all coefficients = 0.0
+ explicit polynomialFunction(const label);
+
+ //- Copy constructor
+ polynomialFunction(const polynomialFunction&);
+
+ //- Construct from a list of coefficients
+ explicit polynomialFunction(const UList<scalar>& coeffs);
+
+ //- Construct from Istream
+ polynomialFunction(Istream&);
+
+
+ //- Destructor
+ virtual ~polynomialFunction();
+
+
+ // Member Functions
+
+ //- Return the number of coefficients
+ using scalarList::size;
+
+ //- Return coefficient
+ using scalarList::operator[];
+
+
+ // Access
+
+
+ //- Return true if the log term is active
+ bool logActive() const;
+
+ //- Return the log coefficient
+ scalar logCoeff() const;
+
+
+ // Evaluation
+
+ //- Return polynomial value
+ scalar value(const scalar x) const;
+
+ //- Integrate between two values
+ scalar integrate(const scalar x1, const scalar x2) const;
+
+
+ //- Return integral coefficients.
+ // Argument becomes zeroth element (constant of integration)
+ polynomialFunction integral
+ (
+ const scalar intConstant = 0.0
+ ) const;
+
+ //- Return integral coefficients when lowest order is -1.
+ // Argument becomes zeroth element (constant of integration)
+ polynomialFunction integralMinus1
+ (
+ const scalar intConstant = 0.0
+ ) const;
+
+
+ // Member Operators
+
+ polynomialFunction& operator+=(const polynomialFunction&);
+ polynomialFunction& operator-=(const polynomialFunction&);
+
+ polynomialFunction& operator*=(const scalar);
+ polynomialFunction& operator/=(const scalar);
+
+
+ //- Ostream Operator
+ friend Ostream& operator<<(Ostream&, const polynomialFunction&);
+};
+
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+// * * * * * * * * * * * * * * * Global Operators * * * * * * * * * * * * * //
+
+polynomialFunction operator+
+(
+ const polynomialFunction&,
+ const polynomialFunction&
+);
+
+
+polynomialFunction operator-
+(
+ const polynomialFunction&,
+ const polynomialFunction&
+);
+
+
+polynomialFunction operator*
+(
+ const scalar,
+ const polynomialFunction&
+);
+
+
+polynomialFunction operator/
+(
+ const scalar,
+ const polynomialFunction&
+);
+
+
+polynomialFunction operator*
+(
+ const polynomialFunction&,
+ const scalar
+);
+
+
+polynomialFunction operator/
+(
+ const polynomialFunction&,
+ const scalar
+);
+
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+} // End namespace Foam
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+#endif
+
+// ************************************************************************* //
diff --git a/src/thermophysicalModels/specie/thermo/hPolynomial/hPolynomialThermo.C b/src/thermophysicalModels/specie/thermo/hPolynomial/hPolynomialThermo.C
index 15269f6ce2c..c86db56ae7f 100644
--- a/src/thermophysicalModels/specie/thermo/hPolynomial/hPolynomialThermo.C
+++ b/src/thermophysicalModels/specie/thermo/hPolynomial/hPolynomialThermo.C
@@ -45,8 +45,8 @@ Foam::hPolynomialThermo<EquationOfState, PolySize>::hPolynomialThermo
Sf_ *= this->W();
CpCoeffs_ *= this->W();
- hCoeffs_ = CpCoeffs_.integrate();
- sCoeffs_ = CpCoeffs_.integrateMinus1();
+ hCoeffs_ = CpCoeffs_.integral();
+ sCoeffs_ = CpCoeffs_.integralMinus1();
// Offset h poly so that it is relative to the enthalpy at Tstd
hCoeffs_[0] += Hf_ - hCoeffs_.value(specie::Tstd);
@@ -73,8 +73,8 @@ Foam::hPolynomialThermo<EquationOfState, PolySize>::hPolynomialThermo
Sf_ *= this->W();
CpCoeffs_ *= this->W();
- hCoeffs_ = CpCoeffs_.integrate();
- sCoeffs_ = CpCoeffs_.integrateMinus1();
+ hCoeffs_ = CpCoeffs_.integral();
+ sCoeffs_ = CpCoeffs_.integralMinus1();
// Offset h poly so that it is relative to the enthalpy at Tstd
hCoeffs_[0] += Hf_ - hCoeffs_.value(specie::Tstd);
--
GitLab