ViewVC Help
View File | Revision Log | Show Annotations | View Changeset | Root Listing
root/OpenMD/branches/development/src/math/Polynomial.hpp
(Generate patch)

Comparing trunk/src/math/Polynomial.hpp (file contents):
Revision 385 by tim, Tue Mar 1 20:10:14 2005 UTC vs.
Revision 1290 by cli2, Wed Sep 10 19:51:45 2008 UTC

# Line 1 | Line 1
1 < /*
1 > /*
2   * Copyright (c) 2005 The University of Notre Dame. All Rights Reserved.
3   *
4   * The University of Notre Dame grants you ("Licensee") a
# Line 53 | Line 53
53   #include <list>
54   #include <map>
55   #include <utility>
56 <
56 > #include "config.h"
57   namespace oopse {
58  
59 < template<typename ElemType> ElemType pow(ElemType x, int N) {
59 >  template<typename ElemType> ElemType pow(ElemType x, int N) {
60      ElemType result(1);
61  
62      for (int i = 0; i < N; ++i) {
63 <        result *= x;
63 >      result *= x;
64      }
65  
66      return result;
67 < }
67 >  }
68  
69 < /**
70 < * @class Polynomial Polynomial.hpp "math/Polynomial.hpp"
71 < * A generic Polynomial class
72 < */
73 < template<typename ElemType>
74 < class Polynomial {
69 >  /**
70 >   * @class Polynomial Polynomial.hpp "math/Polynomial.hpp"
71 >   * A generic Polynomial class
72 >   */
73 >  template<typename ElemType>
74 >  class Polynomial {
75  
76 <    public:
77 <        
78 <        typedef int ExponentType;
79 <        typedef ElemType CoefficientType;
80 <        typedef std::map<ExponentType, CoefficientType> PolynomialPairMap;
81 <        typedef typename PolynomialPairMap::iterator iterator;
82 <        typedef typename PolynomialPairMap::const_iterator const_iterator;
83 <        /**
84 <         * Calculates the value of this Polynomial evaluated at the given x value.
85 <         * @return The value of this Polynomial evaluates at the given x value
86 <         * @param x the value of the independent variable for this Polynomial function
87 <         */
88 <        ElemType evaluate(const ElemType& x) {
89 <            ElemType result = ElemType();
90 <            ExponentType exponent;
91 <            CoefficientType coefficient;
76 >  public:
77 >    typedef Polynomial<ElemType> PolynomialType;    
78 >    typedef int ExponentType;
79 >    typedef ElemType CoefficientType;
80 >    typedef std::map<ExponentType, CoefficientType> PolynomialPairMap;
81 >    typedef typename PolynomialPairMap::iterator iterator;
82 >    typedef typename PolynomialPairMap::const_iterator const_iterator;
83 >
84 >    Polynomial() {}
85 >    Polynomial(ElemType v) {setCoefficient(0, v);}
86 >    /**
87 >     * Calculates the value of this Polynomial evaluated at the given x value.
88 >     * @return The value of this Polynomial evaluates at the given x value
89 >     * @param x the value of the independent variable for this Polynomial function
90 >     */
91 >    ElemType evaluate(const ElemType& x) {
92 >      ElemType result = ElemType();
93 >      ExponentType exponent;
94 >      CoefficientType coefficient;
95              
96 <            for (iterator i = polyPairMap_.begin(); i != polyPairMap_.end(); ++i) {
97 <                exponent = i->first;
98 <                coefficient = i->second;
99 <                result  += pow(x, exponent) * coefficient;
100 <            }
96 >      for (iterator i = polyPairMap_.begin(); i != polyPairMap_.end(); ++i) {
97 >        exponent = i->first;
98 >        coefficient = i->second;
99 >        result  += pow(x, exponent) * coefficient;
100 >      }
101  
102 <            return result;
103 <        }
102 >      return result;
103 >    }
104  
105 <        /**
106 <         * Returns the first derivative of this polynomial.
107 <         * @return the first derivative of this polynomial
108 <         * @param x
109 <         */
110 <        ElemType evaluateDerivative(const ElemType& x) {
111 <            ElemType result = ElemType();
112 <            ExponentType exponent;
113 <            CoefficientType coefficient;
105 >    /**
106 >     * Returns the first derivative of this polynomial.
107 >     * @return the first derivative of this polynomial
108 >     * @param x
109 >     */
110 >    ElemType evaluateDerivative(const ElemType& x) {
111 >      ElemType result = ElemType();
112 >      ExponentType exponent;
113 >      CoefficientType coefficient;
114              
115 <            for (iterator i = polyPairMap_.begin(); i != polyPairMap_.end(); ++i) {
116 <                exponent = i->first;
117 <                coefficient = i->second;
118 <                result  += pow(x, exponent - 1) * coefficient * exponent;
119 <            }
115 >      for (iterator i = polyPairMap_.begin(); i != polyPairMap_.end(); ++i) {
116 >        exponent = i->first;
117 >        coefficient = i->second;
118 >        result  += pow(x, exponent - 1) * coefficient * exponent;
119 >      }
120  
121 <            return result;
122 <        }
121 >      return result;
122 >    }
123  
124 <        /**
125 <         * Set the coefficent of the specified exponent, if the coefficient is already there, it
126 <         * will be overwritten.
127 <         * @param exponent exponent of a term in this Polynomial
128 <         * @param coefficient multiplier of a term in this Polynomial
129 <         */
124 >    /**
125 >     * Set the coefficent of the specified exponent, if the coefficient is already there, it
126 >     * will be overwritten.
127 >     * @param exponent exponent of a term in this Polynomial
128 >     * @param coefficient multiplier of a term in this Polynomial
129 >     */
130          
131 <        void setCoefficient(int exponent, const ElemType& coefficient) {
132 <            polyPairMap_.insert(typename PolynomialPairMap::value_type(exponent, coefficient));
133 <        }
131 >    void setCoefficient(int exponent, const ElemType& coefficient) {
132 >      polyPairMap_[exponent] = coefficient;
133 >    }
134  
135 <        /**
136 <         * Set the coefficent of the specified exponent. If the coefficient is already there,  just add the
137 <         * new coefficient to the old one, otherwise,  just call setCoefficent
138 <         * @param exponent exponent of a term in this Polynomial
139 <         * @param coefficient multiplier of a term in this Polynomial
140 <         */
135 >    /**
136 >     * Set the coefficent of the specified exponent. If the coefficient is already there,  just add the
137 >     * new coefficient to the old one, otherwise,  just call setCoefficent
138 >     * @param exponent exponent of a term in this Polynomial
139 >     * @param coefficient multiplier of a term in this Polynomial
140 >     */
141          
142 <        void addCoefficient(int exponent, const ElemType& coefficient) {
143 <            iterator i = polyPairMap_.find(exponent);
142 >    void addCoefficient(int exponent, const ElemType& coefficient) {
143 >      iterator i = polyPairMap_.find(exponent);
144  
145 <            if (i != end()) {
146 <                i->second += coefficient;
147 <            } else {
148 <                setCoefficient(exponent, coefficient);
149 <            }
150 <        }
145 >      if (i != end()) {
146 >        i->second += coefficient;
147 >      } else {
148 >        setCoefficient(exponent, coefficient);
149 >      }
150 >    }
151  
152 +    /**
153 +     * Returns the coefficient associated with the given power for this Polynomial.
154 +     * @return the coefficient associated with the given power for this Polynomial
155 +     * @exponent exponent of any term in this Polynomial
156 +     */
157 +    ElemType getCoefficient(ExponentType exponent) {
158 +      iterator i = polyPairMap_.find(exponent);
159  
160 <        /**
161 <         * Returns the coefficient associated with the given power for this Polynomial.
162 <         * @return the coefficient associated with the given power for this Polynomial
163 <         * @exponent exponent of any term in this Polynomial
164 <         */
165 <        ElemType getCoefficient(ExponentType exponent) {
156 <            iterator i = polyPairMap_.find(exponent);
160 >      if (i != end()) {
161 >        return i->second;
162 >      } else {
163 >        return ElemType(0);
164 >      }
165 >    }
166  
167 <            if (i != end()) {
168 <                return i->second;
169 <            } else {
161 <                return ElemType(0);
162 <            }
163 <        }
167 >    iterator begin() {
168 >      return polyPairMap_.begin();
169 >    }
170  
171 <        iterator begin() {
172 <            return polyPairMap_.begin();
173 <        }
168 <
169 <        const_iterator begin() const{
170 <            return polyPairMap_.begin();
171 <        }
171 >    const_iterator begin() const{
172 >      return polyPairMap_.begin();
173 >    }
174          
175 <        iterator end() {
176 <            return polyPairMap_.end();
177 <        }
175 >    iterator end() {
176 >      return polyPairMap_.end();
177 >    }
178  
179 <        const_iterator end() const{
180 <            return polyPairMap_.end();
181 <        }
179 >    const_iterator end() const{
180 >      return polyPairMap_.end();
181 >    }
182  
183 <        iterator find(ExponentType exponent) {
184 <            return polyPairMap_.find(exponent);
183 >    iterator find(ExponentType exponent) {
184 >      return polyPairMap_.find(exponent);
185 >    }
186 >
187 >    size_t size() {
188 >      return polyPairMap_.size();
189 >    }
190 >
191 >    PolynomialType& operator = (const PolynomialType& p) {
192 >
193 >      if (this != &p)  // protect against invalid self-assignment
194 >      {
195 >        typename Polynomial<ElemType>::const_iterator i;
196 >
197 >        polyPairMap_.clear();  // clear out the old map
198 >      
199 >        for (i =  p.begin(); i != p.end(); ++i) {
200 >          this->setCoefficient(i->first, i->second);
201          }
202 +      }
203 +      // by convention, always return *this
204 +      return *this;
205 +    }
206  
207 <        size_t size() {
208 <            return polyPairMap_.size();
207 >    PolynomialType& operator += (const PolynomialType& p) {
208 >        typename Polynomial<ElemType>::const_iterator i;
209 >
210 >        for (i =  p.begin(); i  != p.end(); ++i) {
211 >          this->addCoefficient(i->first, i->second);
212          }
213 +
214 +        return *this;        
215 +    }
216 +
217 +    PolynomialType& operator -= (const PolynomialType& p) {
218 +        typename Polynomial<ElemType>::const_iterator i;
219 +        for (i =  p.begin(); i  != p.end(); ++i) {
220 +          this->addCoefficient(i->first, -i->second);
221 +        }        
222 +        return *this;
223 +    }
224 +
225 +    PolynomialType& operator *= (const PolynomialType& p) {
226 +    typename Polynomial<ElemType>::const_iterator i;
227 +    typename Polynomial<ElemType>::const_iterator j;
228 +    Polynomial<ElemType> p2(*this);
229 +  
230 +    polyPairMap_.clear();  // clear out old map
231 +    for (i = p2.begin(); i !=p2.end(); ++i) {
232 +      for (j = p.begin(); j !=p.end(); ++j) {
233 +        this->addCoefficient( i->first + j->first, i->second * j->second);
234 +      }
235 +    }
236 +    return *this;
237 +    }
238 +
239 +    //PolynomialType& operator *= (const ElemType v)
240 +    PolynomialType& operator *= (const ElemType v) {
241 +    typename Polynomial<ElemType>::const_iterator i;
242 +    //Polynomial<ElemType> result;
243 +
244 +    for (i = this->begin(); i != this->end(); ++i) {
245 +        this->setCoefficient( i->first, i->second*v);
246 +    }
247 +
248 +    return *this;
249 +    }
250 +
251 +    PolynomialType& operator += (const ElemType v) {    
252 +    this->addCoefficient( 0, v);
253 +    return *this;
254 +    }
255 +  
256 +  private:
257          
258 <    private:
259 <        
191 <        PolynomialPairMap polyPairMap_;
192 < };
258 >    PolynomialPairMap polyPairMap_;
259 >  };
260  
261  
262 < /**
263 < * Generates and returns the product of two given Polynomials.
264 < * @return A Polynomial containing the product of the two given Polynomial parameters
265 < */
266 < template<typename ElemType>
267 < Polynomial<ElemType> operator *(const Polynomial<ElemType>& p1, const Polynomial<ElemType>& p2) {
262 >  /**
263 >   * Generates and returns the product of two given Polynomials.
264 >   * @return A Polynomial containing the product of the two given Polynomial parameters
265 >   */
266 >  template<typename ElemType>
267 >  Polynomial<ElemType> operator *(const Polynomial<ElemType>& p1, const Polynomial<ElemType>& p2) {
268      typename Polynomial<ElemType>::const_iterator i;
269      typename Polynomial<ElemType>::const_iterator j;
270      Polynomial<ElemType> p;
271      
272      for (i = p1.begin(); i !=p1.end(); ++i) {
273 <        for (j = p2.begin(); j !=p2.end(); ++j) {
274 <            p.addCoefficient( i->first + j->first, i->second * j->second);
275 <        }
273 >      for (j = p2.begin(); j !=p2.end(); ++j) {
274 >        p.addCoefficient( i->first + j->first, i->second * j->second);
275 >      }
276      }
277  
278      return p;
279 < }
279 >  }
280  
281 < /**
282 < * Generates and returns the sum of two given Polynomials.
283 < * @param p1 the first polynomial
284 < * @param p2 the second polynomial
285 < */
286 < template<typename ElemType>
287 < Polynomial<ElemType> operator +(const Polynomial<ElemType>& p1, const Polynomial<ElemType>& p2) {
281 >  template<typename ElemType>
282 >  Polynomial<ElemType> operator *(const Polynomial<ElemType>& p, const ElemType v) {
283 >    typename Polynomial<ElemType>::const_iterator i;
284 >    Polynomial<ElemType> result;
285 >    
286 >    for (i = p.begin(); i !=p.end(); ++i) {
287 >        result.setCoefficient( i->first , i->second * v);
288 >    }
289 >
290 >    return result;
291 >  }
292 >
293 >  template<typename ElemType>
294 >  Polynomial<ElemType> operator *( const ElemType v, const Polynomial<ElemType>& p) {
295 >    typename Polynomial<ElemType>::const_iterator i;
296 >    Polynomial<ElemType> result;
297 >    
298 >    for (i = p.begin(); i !=p.end(); ++i) {
299 >        result.setCoefficient( i->first , i->second * v);
300 >    }
301 >
302 >    return result;
303 >  }
304 >  
305 >  /**
306 >   * Generates and returns the sum of two given Polynomials.
307 >   * @param p1 the first polynomial
308 >   * @param p2 the second polynomial
309 >   */
310 >  template<typename ElemType>
311 >  Polynomial<ElemType> operator +(const Polynomial<ElemType>& p1, const Polynomial<ElemType>& p2) {
312      Polynomial<ElemType> p(p1);
313  
314      typename Polynomial<ElemType>::const_iterator i;
315  
316      for (i =  p2.begin(); i  != p2.end(); ++i) {
317 <        p.addCoefficient(i->first, i->second);
317 >      p.addCoefficient(i->first, i->second);
318      }
319  
320      return p;
321  
322 < }
322 >  }
323  
324 < /**
325 < * Generates and returns the difference of two given Polynomials.
326 < * @return
327 < * @param p1 the first polynomial
328 < * @param p2 the second polynomial
329 < */
330 < template<typename ElemType>
331 < Polynomial<ElemType> operator -(const Polynomial<ElemType>& p1, const Polynomial<ElemType>& p2) {
324 >  /**
325 >   * Generates and returns the difference of two given Polynomials.
326 >   * @return
327 >   * @param p1 the first polynomial
328 >   * @param p2 the second polynomial
329 >   */
330 >  template<typename ElemType>
331 >  Polynomial<ElemType> operator -(const Polynomial<ElemType>& p1, const Polynomial<ElemType>& p2) {
332      Polynomial<ElemType> p(p1);
333  
334      typename Polynomial<ElemType>::const_iterator i;
335  
336      for (i =  p2.begin(); i  != p2.end(); ++i) {
337 <        p.addCoefficient(i->first, -i->second);
337 >      p.addCoefficient(i->first, -i->second);
338      }
339  
340      return p;
341  
342 < }
342 >  }
343  
344 < /**
345 < * Tests if two polynomial have the same exponents
346 < * @return true if these all of the exponents in these Polynomial are identical
347 < * @param p1 the first polynomial
348 < * @param p2 the second polynomial
349 < * @note this function does not compare the coefficient
350 < */
351 < template<typename ElemType>
352 < bool equal(const Polynomial<ElemType>& p1, const Polynomial<ElemType>& p2) {
344 >  /**
345 >   * Tests if two polynomial have the same exponents
346 >   * @return true if all of the exponents in these Polynomial are identical
347 >   * @param p1 the first polynomial
348 >   * @param p2 the second polynomial
349 >   * @note this function does not compare the coefficient
350 >   */
351 >  template<typename ElemType>
352 >  bool equal(const Polynomial<ElemType>& p1, const Polynomial<ElemType>& p2) {
353  
354      typename Polynomial<ElemType>::const_iterator i;
355      typename Polynomial<ElemType>::const_iterator j;
356  
357      if (p1.size() != p2.size() ) {
358 <        return false;
358 >      return false;
359      }
360      
361      for (i =  p1.begin(), j = p2.begin(); i  != p1.end() && j != p2.end(); ++i, ++j) {
362 <        if (i->first != j->first) {
363 <            return false;
364 <        }
362 >      if (i->first != j->first) {
363 >        return false;
364 >      }
365      }
366  
367      return true;
368 < }
368 >  }
369  
370 < typedef Polynomial<double> DoublePolynomial;
370 >  typedef Polynomial<RealType> DoublePolynomial;
371  
372   } //end namespace oopse
373   #endif //MATH_POLYNOMIAL_HPP

Diff Legend

– Removed lines
+ Added lines
< Changed lines
> Changed lines