Teuchos - Trilinos Tools Package Version of the Day
Teuchos_Polynomial.hpp
00001 // @HEADER
00002 // ***********************************************************************
00003 //
00004 //                    Teuchos: Common Tools Package
00005 //                 Copyright (2004) Sandia Corporation
00006 //
00007 // Under terms of Contract DE-AC04-94AL85000, there is a non-exclusive
00008 // license for use of this work by or on behalf of the U.S. Government.
00009 //
00010 // Redistribution and use in source and binary forms, with or without
00011 // modification, are permitted provided that the following conditions are
00012 // met:
00013 //
00014 // 1. Redistributions of source code must retain the above copyright
00015 // notice, this list of conditions and the following disclaimer.
00016 //
00017 // 2. Redistributions in binary form must reproduce the above copyright
00018 // notice, this list of conditions and the following disclaimer in the
00019 // documentation and/or other materials provided with the distribution.
00020 //
00021 // 3. Neither the name of the Corporation nor the names of the
00022 // contributors may be used to endorse or promote products derived from
00023 // this software without specific prior written permission.
00024 //
00025 // THIS SOFTWARE IS PROVIDED BY SANDIA CORPORATION "AS IS" AND ANY
00026 // EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
00027 // IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR
00028 // PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL SANDIA CORPORATION OR THE
00029 // CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL,
00030 // EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO,
00031 // PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR
00032 // PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF
00033 // LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING
00034 // NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
00035 // SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
00036 //
00037 // Questions? Contact Michael A. Heroux (maherou@sandia.gov)
00038 //
00039 // ***********************************************************************
00040 // @HEADER
00041 
00042 #ifndef TEUCHOS_POLYNOMIAL_HPP
00043 #define TEUCHOS_POLYNOMIAL_HPP
00044 
00045 #include "Teuchos_PolynomialDecl.hpp"
00046 #include "Teuchos_ScalarTraits.hpp"
00047 
00048 template <typename CoeffT>
00049 Teuchos::Polynomial<CoeffT>::Polynomial(unsigned int deg,
00050           const CoeffT& cloneCoeff, 
00051           unsigned int reserve) :
00052   d(deg)
00053 {
00054   if (reserve > d)
00055     sz = reserve+1;
00056   else
00057     sz = d+1;
00058 
00059   coeff.resize(sz);
00060   for (unsigned int i=0; i<sz; i++)
00061     coeff[i] = PolynomialTraits<CoeffT>::clone(cloneCoeff);
00062 }
00063 
00064 template <typename CoeffT>
00065 Teuchos::Polynomial<CoeffT>::Polynomial(unsigned int deg, 
00066           unsigned int reserve) :
00067   d(deg)
00068 {
00069   if (reserve > d)
00070     sz = reserve+1;
00071   else
00072     sz = d+1;
00073 
00074   coeff.resize(sz);
00075 }
00076 
00077 template <typename CoeffT>
00078 Teuchos::Polynomial<CoeffT>::~Polynomial()
00079 {
00080 }
00081 
00082 template <typename CoeffT>
00083 void 
00084 Teuchos::Polynomial<CoeffT>::setDegree(unsigned int deg) 
00085 {
00086   d = deg;
00087   if (d+1 > sz) {
00088     coeff.resize(d+1);
00089     if (coeff[0] != Teuchos::null) {
00090       for (unsigned int i=sz; i<d+1; i++)
00091   coeff[i] = PolynomialTraits<CoeffT>::clone(*coeff[0]);
00092     }
00093     sz = d+1;
00094   }
00095 }
00096 
00097 template <typename CoeffT>
00098 Teuchos::RCP<CoeffT>
00099 Teuchos::Polynomial<CoeffT>::getCoefficient(unsigned int i) 
00100 {
00101 #ifdef TEUCHOS_DEBUG
00102   TEUCHOS_TEST_FOR_EXCEPTION(i > d, 
00103          std::out_of_range,
00104          "Polynomial<CoeffT>::getCoefficient(i): " << 
00105          "Error, coefficient i = " << i << 
00106          " is not in range, degree = " << d << "." );
00107 #endif
00108   return coeff[i];
00109 }
00110 
00111 template <typename CoeffT>
00112 Teuchos::RCP<const CoeffT>
00113 Teuchos::Polynomial<CoeffT>::getCoefficient(unsigned int i) const
00114 {
00115 #ifdef TEUCHOS_DEBUG
00116   TEUCHOS_TEST_FOR_EXCEPTION(i > d, 
00117          std::out_of_range,
00118          "Polynomial<CoeffT>::getCoefficient(i): " << 
00119          "Error, coefficient i = " << i << 
00120          " is not in range, degree = " << d << "." );
00121 #endif
00122   return coeff[i];
00123 }
00124 
00125 template <typename CoeffT>
00126 void
00127 Teuchos::Polynomial<CoeffT>::setCoefficient(unsigned int i, const CoeffT& v) 
00128 {
00129 #ifdef TEUCHOS_DEBUG
00130   TEUCHOS_TEST_FOR_EXCEPTION(i > d, 
00131          std::out_of_range,
00132          "Polynomial<CoeffT>::setCoefficient(i,v): " << 
00133          "Error, coefficient i = " << i << 
00134          " is not in range, degree = " << d << "." );
00135   TEUCHOS_TEST_FOR_EXCEPTION(coeff[i] == Teuchos::null, 
00136          std::runtime_error,
00137          "Polynomial<CoeffT>::setCoefficient(i,v): " << 
00138          "Error, coefficient i = " << i << " is null!");
00139 #endif
00140   PolynomialTraits<CoeffT>::copy(v, coeff[i].get());
00141 }
00142 
00143 template <typename CoeffT>
00144 void
00145 Teuchos::Polynomial<CoeffT>::setCoefficientPtr(
00146                             unsigned int i,
00147                             const Teuchos::RCP<CoeffT>& v)
00148 {
00149 #ifdef TEUCHOS_DEBUG
00150   TEUCHOS_TEST_FOR_EXCEPTION(i > d, 
00151          std::out_of_range,
00152          "Polynomial<CoeffT>::setCoefficientPtr(i,v): " << 
00153          "Error, coefficient i = " << i << 
00154          " is not in range, degree = " << d << "." );
00155 #endif
00156   coeff[i] = v;
00157 }
00158 
00159 template <typename CoeffT>
00160 void
00161 Teuchos::Polynomial<CoeffT>::evaluate(
00162         typename Teuchos::Polynomial<CoeffT>::scalar_type t, 
00163         CoeffT* x, CoeffT* xdot) const
00164 {
00165   bool evaluate_xdot = (xdot != NULL);
00166 
00167 #ifdef TEUCHOS_DEBUG
00168   for (unsigned int i=0; i<=d; i++)
00169     TEUCHOS_TEST_FOR_EXCEPTION(coeff[i] == Teuchos::null, 
00170            std::runtime_error,
00171            "Polynomial<CoeffT>::evaluate(): " << 
00172            "Error, coefficient i = " << i << " is null!");
00173 #endif
00174 
00175   // Initialize x, xdot with coeff[d]
00176   PolynomialTraits<CoeffT>::copy(*coeff[d], x);
00177   if (evaluate_xdot) {
00178     if (d > 0)
00179       PolynomialTraits<CoeffT>::copy(*coeff[d], xdot);
00180     else
00181       PolynomialTraits<CoeffT>::assign(
00182          xdot, 
00183          Teuchos::ScalarTraits<scalar_type>::zero());
00184   }
00185 
00186   // If this is a degree 0 polynomial, we're done
00187   if (d == 0)
00188     return;
00189 
00190   for (int k=d-1; k>=0; --k) {
00191     // compute x = coeff[k] + t*x
00192     PolynomialTraits<CoeffT>::update(x, *coeff[k], t);
00193 
00194     // compute xdot = x + t*xdot
00195     if (evaluate_xdot && k > 0)
00196       PolynomialTraits<CoeffT>::update(xdot, *x, t);
00197   }
00198 }
00199 
00200 #endif  // TEUCHOS_VECTOR_POLYNOMIAL_HPP
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Defines