FreeFOAM The Cross-Platform CFD Toolkit
Hosted by SourceForge:
Get FreeFOAM at SourceForge.net.
            Fast, secure and Free Open Source software downloads

Polynomial.C

Go to the documentation of this file.
00001 /*---------------------------------------------------------------------------*\
00002   =========                 |
00003   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
00004    \\    /   O peration     |
00005     \\  /    A nd           | Copyright (C) 2008-2010 OpenCFD Ltd.
00006      \\/     M anipulation  |
00007 -------------------------------------------------------------------------------
00008 License
00009     This file is part of OpenFOAM.
00010 
00011     OpenFOAM is free software: you can redistribute it and/or modify it
00012     under the terms of the GNU General Public License as published by
00013     the Free Software Foundation, either version 3 of the License, or
00014     (at your option) any later version.
00015 
00016     OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
00017     ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
00018     FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
00019     for more details.
00020 
00021     You should have received a copy of the GNU General Public License
00022     along with OpenFOAM.  If not, see <http://www.gnu.org/licenses/>.
00023 
00024 \*---------------------------------------------------------------------------*/
00025 
00026 #include "Polynomial.H"
00027 
00028 // * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
00029 
00030 template<int PolySize>
00031 Foam::Polynomial<PolySize>::Polynomial()
00032 :
00033     VectorSpace<Polynomial<PolySize>, scalar, PolySize>(),
00034     logActive_(false),
00035     logCoeff_(0.0)
00036 {}
00037 
00038 
00039 template<int PolySize>
00040 Foam::Polynomial<PolySize>::Polynomial(const word& name, Istream& is)
00041 :
00042     VectorSpace<Polynomial<PolySize>, scalar, PolySize>(),
00043     logActive_(false),
00044     logCoeff_(0.0)
00045 {
00046     word isName(is);
00047 
00048     if (isName != name)
00049     {
00050         FatalErrorIn
00051         (
00052             "Polynomial<PolySize>::Polynomial(const word&, Istream&)"
00053         )   << "Expected polynomial name " << name << " but read " << isName
00054             << nl << exit(FatalError);
00055     }
00056 
00057     VectorSpace<Polynomial<PolySize>, scalar, PolySize>::
00058         operator=(VectorSpace<Polynomial<PolySize>, scalar, PolySize>(is));
00059 
00060     if (this->size() == 0)
00061     {
00062         FatalErrorIn
00063         (
00064             "Polynomial<PolySize>::Polynomial(const word&, Istream&)"
00065         )   << "Polynomial coefficients for entry " << isName
00066             << " are invalid (empty)" << nl << exit(FatalError);
00067     }
00068 }
00069 
00070 
00071 template<int PolySize>
00072 Foam::Polynomial<PolySize>::Polynomial
00073 (
00074     const Polynomial<PolySize>& poly
00075 )
00076 :
00077     VectorSpace<Polynomial<PolySize>, scalar, PolySize>(poly),
00078     logActive_(poly.logActive_),
00079     logCoeff_(poly.logCoeff_)
00080 {}
00081 
00082 
00083 // * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
00084 
00085 template<int PolySize>
00086 bool& Foam::Polynomial<PolySize>::logActive()
00087 {
00088     return logActive_;
00089 }
00090 
00091 
00092 template<int PolySize>
00093 Foam::scalar& Foam::Polynomial<PolySize>::logCoeff()
00094 {
00095     return logCoeff_;
00096 }
00097 
00098 
00099 template<int PolySize>
00100 Foam::scalar Foam::Polynomial<PolySize>::evaluate(const scalar x) const
00101 {
00102     scalar y = this->v_[0];
00103 
00104     for (label i=1; i<PolySize; i++)
00105     {
00106         y += this->v_[i]*pow(x, i);
00107     }
00108 
00109     if (logActive_)
00110     {
00111         y += logCoeff_*log(x);
00112     }
00113 
00114     return y;
00115 }
00116 
00117 
00118 template<int PolySize>
00119 Foam::scalar Foam::Polynomial<PolySize>::integrateLimits
00120 (
00121     const scalar x1,
00122     const scalar x2
00123 ) const
00124 {
00125     if (logActive_)
00126     {
00127         FatalErrorIn
00128         (
00129             "scalar Polynomial<PolySize>::integrateLimits"
00130             "("
00131                 "const scalar, "
00132                 "const scalar"
00133             ") const"
00134         )   << "Cannot integrate polynomial with logarithmic coefficients"
00135             << nl << abort(FatalError);
00136     }
00137 
00138     intPolyType poly = this->integrate();
00139 
00140     return poly.evaluate(x2) - poly.evaluate(x1);
00141 }
00142 
00143 
00144 template<int PolySize>
00145 typename Foam::Polynomial<PolySize>::intPolyType
00146 Foam::Polynomial<PolySize>::integrate(const scalar intConstant)
00147 {
00148     intPolyType newCoeffs;
00149 
00150     newCoeffs[0] = intConstant;
00151     forAll(*this, i)
00152     {
00153         newCoeffs[i + 1] = this->v_[i]/(i + 1);
00154     }
00155 
00156     return newCoeffs;
00157 }
00158 
00159 
00160 template<int PolySize>
00161 typename Foam::Polynomial<PolySize>::polyType
00162 Foam::Polynomial<PolySize>::integrateMinus1(const scalar intConstant)
00163 {
00164     polyType newCoeffs;
00165 
00166     if (this->v_[0] > VSMALL)
00167     {
00168         newCoeffs.logActive() = true;
00169         newCoeffs.logCoeff() = this->v_[0];
00170     }
00171 
00172     newCoeffs[0] = intConstant;
00173 
00174     if (PolySize > 0)
00175     {
00176         for (label i=1; i<PolySize; i++)
00177         {
00178             newCoeffs[i] = this->v_[i]/i;
00179         }
00180     }
00181 
00182     return newCoeffs;
00183 }
00184 
00185 
00186 // ************************ vim: set sw=4 sts=4 et: ************************ //
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Defines