1.1 --- /dev/null Thu Jan 01 00:00:00 1970 +0000
1.2 +++ b/lemon/polynomial.h Tue May 16 16:59:57 2006 +0000
1.3 @@ -0,0 +1,349 @@
1.4 +/* -*- C++ -*-
1.5 + *
1.6 + * This file is a part of LEMON, a generic C++ optimization library
1.7 + *
1.8 + * Copyright (C) 2003-2006
1.9 + * Egervary Jeno Kombinatorikus Optimalizalasi Kutatocsoport
1.10 + * (Egervary Research Group on Combinatorial Optimization, EGRES).
1.11 + *
1.12 + * Permission to use, modify and distribute this software is granted
1.13 + * provided that this copyright notice appears in all copies. For
1.14 + * precise terms see the accompanying LICENSE file.
1.15 + *
1.16 + * This software is provided "AS IS" with no warranty of any kind,
1.17 + * express or implied, and with no claim as to its suitability for any
1.18 + * purpose.
1.19 + *
1.20 + */
1.21 +
1.22 +#ifndef LEMON_BEZIER_H
1.23 +#define LEMON_BEZIER_H
1.24 +
1.25 +///\ingroup misc
1.26 +///\file
1.27 +///\brief A simple class implementing polynomials.
1.28 +///
1.29 +///\author Alpar Juttner
1.30 +
1.31 +#include<vector>
1.32 +
1.33 +namespace lemon {
1.34 +
1.35 + /// \addtogroup misc
1.36 + /// @{
1.37 +
1.38 + ///Simple polinomial class
1.39 +
1.40 + ///This class implements a polynomial where the coefficients are of
1.41 + ///type \c T.
1.42 + ///
1.43 + ///The coefficients are stored in an std::vector.
1.44 + template<class T>
1.45 + class Polynomial
1.46 + {
1.47 + std::vector<T> _coeff;
1.48 + public:
1.49 + ///Construct a polynomial of degree \c d.
1.50 + explicit Polynomial(int d=0) : _coeff(d+1) {}
1.51 + ///\e
1.52 + template<class U> Polynomial(const U &u) : _coeff(1,u) {}
1.53 + ///\e
1.54 + template<class U> Polynomial(const Polynomial<U> &u) : _coeff(u.deg()+1)
1.55 + {
1.56 + for(int i=0;i<(int)_coeff.size();i++) _coeff[i]=u[i];
1.57 + }
1.58 + ///Query the degree of the polynomial.
1.59 +
1.60 + ///Query the degree of the polynomial.
1.61 + ///\warning This number differs from real degree of the polinomial if
1.62 + ///the coefficient of highest degree is 0.
1.63 + int deg() const { return _coeff.size()-1; }
1.64 + ///Set the degree of the polynomial.
1.65 +
1.66 + ///Set the degree of the polynomial. In fact it resizes the
1.67 + ///coefficient vector.
1.68 + void deg(int d) { _coeff.resize(d+1);}
1.69 +
1.70 + ///Returns (as a reference) the coefficient of degree \c d.
1.71 + typename std::vector<T>::reference operator[](int d) { return _coeff[d]; }
1.72 + ///Returns (as a const reference) the coefficient of degree \c d.
1.73 + typename std::vector<T>::const_reference
1.74 + operator[](int d) const {return _coeff[d];}
1.75 +
1.76 + ///Substitute the value u into the polinomial.
1.77 +
1.78 + ///Substitute the value u into the polinomial.
1.79 + ///The calculation will be done using type \c R.
1.80 + ///The following examples shows the usage of the template parameter \c R.
1.81 + ///\code
1.82 + /// Polynomial<xy<double> > line(1);
1.83 + /// line[0]=xy<double>(12,25);
1.84 + /// line[1]=xy<double>(2,7);
1.85 + /// ...
1.86 + /// xy<double> d = line.subst<xy<double> >(23.2);
1.87 + ///\endcode
1.88 + ///
1.89 + ///\code
1.90 + /// Polynomial<double> p;
1.91 + /// Polynomial<double> q;
1.92 + /// ...
1.93 + /// Polynomial<double> s = p.subst<Polynomial<double> >(q);
1.94 + ///\endcode
1.95 + template<class R,class U>
1.96 + R subst(const U &u) const
1.97 + {
1.98 + typename std::vector<T>::const_reverse_iterator i=_coeff.rbegin();
1.99 + R v=*i;
1.100 + for(++i;i!=_coeff.rend();++i) v=v*u+*i;
1.101 + return v;
1.102 + }
1.103 + ///Substitute the value u into the polinomial.
1.104 +
1.105 + ///Substitute the value u into the polinomial.
1.106 + ///The calculation will be done using type \c T
1.107 + ///(i.e. using the type of the coefficients.)
1.108 + template<class U>
1.109 + T operator()(const U &u) const
1.110 + {
1.111 + return subst<T>(u);
1.112 + }
1.113 +
1.114 + ///Derivate the polynomial (in place)
1.115 + Polynomial &derivateMyself()
1.116 + {
1.117 + for(int i=1;i<(int)_coeff.size();i++) _coeff[i-1]=i*_coeff[i];
1.118 + _coeff.pop_back();
1.119 + return *this;
1.120 + }
1.121 +
1.122 + ///Return the derivate of the polynomial
1.123 + Polynomial derivate() const
1.124 + {
1.125 + Polynomial tmp(deg()-1);
1.126 + for(int i=1;i<(int)_coeff.size();i++) tmp[i-1]=i*_coeff[i];
1.127 + return tmp;
1.128 + }
1.129 +
1.130 + ///Integrate the polynomial (in place)
1.131 + Polynomial &integrateMyself()
1.132 + {
1.133 + _coeff.push_back(T());
1.134 + for(int i=_coeff.size()-1;i>=0;i--) _coeff[i]=_coeff[i-1]/i;
1.135 + _coeff[0]=0;
1.136 + return *this;
1.137 + }
1.138 +
1.139 + ///Return the integrate of the polynomial
1.140 + Polynomial integrate() const
1.141 + {
1.142 + Polynomial tmp(deg()+1);
1.143 + tmp[0]=0;
1.144 + for(int i=0;i<(int)_coeff.size();i++) tmp[i+1]=_coeff[i]/(i+1);
1.145 + return tmp;
1.146 + }
1.147 +
1.148 + ///\e
1.149 + template<class U>
1.150 + Polynomial &operator+=(const Polynomial<U> &p)
1.151 + {
1.152 + if(p.deg()>deg()) _coeff.resize(p.deg()+1);
1.153 + for(int i=0;i<=(int)std::min(deg(),p.deg());i++)
1.154 + _coeff[i]+=p[i];
1.155 + return *this;
1.156 + }
1.157 + ///\e
1.158 + template<class U>
1.159 + Polynomial &operator-=(const Polynomial<U> &p)
1.160 + {
1.161 + if(p.deg()>deg()) _coeff.resize(p.deg()+1);
1.162 + for(int i=0;i<=std::min(deg(),p.deg());i++) _coeff[i]-=p[i];
1.163 + return *this;
1.164 + }
1.165 + ///\e
1.166 + template<class U>
1.167 + Polynomial &operator+=(const U &u)
1.168 + {
1.169 + _coeff[0]+=u;
1.170 + return *this;
1.171 + }
1.172 + ///\e
1.173 + template<class U>
1.174 + Polynomial &operator-=(const U &u)
1.175 + {
1.176 + _coeff[0]+=u;
1.177 + return *this;
1.178 + }
1.179 + ///\e
1.180 + template<class U>
1.181 + Polynomial &operator*=(const U &u)
1.182 + {
1.183 + for(typename std::vector<T>::iterator i=_coeff.begin();i!=_coeff.end();++i)
1.184 + *i*=u;
1.185 + return *this;
1.186 + }
1.187 + ///\e
1.188 + template<class U>
1.189 + Polynomial &operator/=(const U &u)
1.190 + {
1.191 + for(typename std::vector<T>::iterator i=_coeff.begin();i!=_coeff.end();++i)
1.192 + *i/=u;
1.193 + return *this;
1.194 + }
1.195 +
1.196 + };
1.197 +
1.198 + ///Equality comparison
1.199 +
1.200 + ///\relates Polynomial
1.201 + ///\warning Two polynomials are defined to be unequal if their degrees differ,
1.202 + ///even if the non-zero coefficients are the same.
1.203 + template<class U,class V>
1.204 + bool operator==(const Polynomial<U> &u,const Polynomial<V> &v)
1.205 + {
1.206 + if(u.deg()!=v.deg()) return false;
1.207 + for(int i=0;i<=u.deg();i++) if(u[i]!=v[i]) return false;
1.208 + return true;
1.209 + }
1.210 +
1.211 + ///Non-equality comparison
1.212 +
1.213 + ///\relates Polynomial
1.214 + ///\warning Two polynomials are defined to be unequal if their degrees differ,
1.215 + ///even if the non-zero coefficients are the same.
1.216 + template<class U,class V>
1.217 + bool operator!=(const Polynomial<U> &u,const Polynomial<V> &v)
1.218 + {
1.219 + return !(u==v);
1.220 + }
1.221 +
1.222 + ///\e
1.223 +
1.224 + ///\relates Polynomial
1.225 + ///
1.226 + template<class U,class V>
1.227 + Polynomial<U> operator+(const Polynomial<U> &u,const Polynomial<V> &v)
1.228 + {
1.229 + Polynomial<U> tmp=u;
1.230 + tmp+=v;
1.231 + return tmp;
1.232 + }
1.233 +
1.234 + ///\e
1.235 +
1.236 + ///\relates Polynomial
1.237 + ///
1.238 + template<class U,class V>
1.239 + Polynomial<U> operator-(const Polynomial<U> &u,const Polynomial<V> &v)
1.240 + {
1.241 + Polynomial<U> tmp=u;
1.242 + tmp-=v;
1.243 + return tmp;
1.244 + }
1.245 +
1.246 + ///\e
1.247 +
1.248 + ///\relates Polynomial
1.249 + ///
1.250 + template<class U,class V>
1.251 + Polynomial<U> operator*(const Polynomial<U> &u,const Polynomial<V> &v)
1.252 + {
1.253 + Polynomial<U> tmp(u.deg()+v.deg());
1.254 + for(int i=0;i<=v.deg();i++)
1.255 + for(int j=0;j<=u.deg();j++)
1.256 + tmp[i+j]+=v[i]*u[j];
1.257 + return tmp;
1.258 + }
1.259 + ///\e
1.260 +
1.261 + ///\relates Polynomial
1.262 + ///
1.263 + template<class U,class V>
1.264 + Polynomial<U> operator+(const Polynomial<U> &u,const V &v)
1.265 + {
1.266 + Polynomial<U> tmp=u;
1.267 + tmp+=v;
1.268 + return tmp;
1.269 + }
1.270 + ///\e
1.271 +
1.272 + ///\relates Polynomial
1.273 + ///
1.274 + template<class U,class V>
1.275 + Polynomial<U> operator+(const V &v,const Polynomial<U> &u)
1.276 + {
1.277 + Polynomial<U> tmp=u;
1.278 + tmp+=v;
1.279 + return tmp;
1.280 + }
1.281 + ///\e
1.282 +
1.283 + ///\relates Polynomial
1.284 + ///
1.285 + template<class U,class V>
1.286 + Polynomial<U> operator-(const Polynomial<U> &u,const V &v)
1.287 + {
1.288 + Polynomial<U> tmp=u;
1.289 + tmp-=v;
1.290 + return tmp;
1.291 + }
1.292 + ///\e
1.293 +
1.294 + ///\relates Polynomial
1.295 + ///
1.296 + template<class U>
1.297 + Polynomial<U> operator-(const Polynomial<U> &u)
1.298 + {
1.299 + Polynomial<U> tmp(u.deg());
1.300 + for(int i=0;i<=u.deg();i++) tmp[i]=-u[i];
1.301 + return tmp;
1.302 + }
1.303 + ///\e
1.304 +
1.305 + ///\relates Polynomial
1.306 + ///
1.307 + template<class U,class V>
1.308 + Polynomial<U> operator-(const V &v,const Polynomial<U> &u)
1.309 + {
1.310 + Polynomial<U> tmp=-u;
1.311 + tmp+=v;
1.312 + return tmp;
1.313 + }
1.314 + ///\e
1.315 +
1.316 + ///\relates Polynomial
1.317 + ///
1.318 + template<class U,class V>
1.319 + Polynomial<U> operator*(const Polynomial<U> &u,const V &v)
1.320 + {
1.321 + Polynomial<U> tmp=u;
1.322 + tmp*=v;
1.323 + return tmp;
1.324 + }
1.325 + ///\e
1.326 +
1.327 + ///\relates Polynomial
1.328 + ///
1.329 + template<class U,class V>
1.330 + Polynomial<U> operator*(const V &v,const Polynomial<U> &u)
1.331 + {
1.332 + Polynomial<U> tmp=u;
1.333 + tmp*=v;
1.334 + return tmp;
1.335 + }
1.336 + ///\e
1.337 +
1.338 + ///\relates Polynomial
1.339 + ///
1.340 + template<class U,class V>
1.341 + Polynomial<U> operator/(const Polynomial<U> &u,const V &v)
1.342 + {
1.343 + Polynomial<U> tmp=u;
1.344 + tmp/=v;
1.345 + return tmp;
1.346 + }
1.347 +
1.348 + /// @}
1.349 +
1.350 +} //END OF NAMESPACE LEMON
1.351 +
1.352 +#endif // LEMON_POLYNOMIAL_H