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

Comparing branches/development/src/math/Polynomial.hpp (file contents):
Revision 1808 by gezelter, Mon Oct 22 20:42:10 2012 UTC vs.
Revision 1874 by gezelter, Wed May 15 15:09:35 2013 UTC

# Line 35 | Line 35
35   *                                                                      
36   * [1]  Meineke, et al., J. Comp. Chem. 26, 252-271 (2005).            
37   * [2]  Fennell & Gezelter, J. Chem. Phys. 124, 234104 (2006).          
38 < * [3]  Sun, Lin & Gezelter, J. Chem. Phys. 128, 24107 (2008).          
38 > * [3]  Sun, Lin & Gezelter, J. Chem. Phys. 128, 234107 (2008).          
39   * [4]  Kuang & Gezelter,  J. Chem. Phys. 133, 164101 (2010).
40   * [5]  Vardeman, Stocker & Gezelter, J. Chem. Theory Comput. 7, 834 (2011).
41   */
# Line 273 | Line 273 | namespace OpenMD {
273       * Returns the first derivative of this polynomial.
274       * @return the first derivative of this polynomial
275       */
276 <    PolynomialType & getDerivative() {
277 <      Polynomial<Real> p;
276 >    PolynomialType* getDerivative() {
277 >      Polynomial<Real>* p = new Polynomial<Real>();
278        
279        typename Polynomial<Real>::const_iterator i;
280        ExponentType exponent;
# Line 283 | Line 283 | namespace OpenMD {
283        for (i =  this->begin(); i  != this->end(); ++i) {
284          exponent = i->first;
285          coefficient = i->second;
286 <        p.setCoefficient(exponent-1, coefficient * exponent);
286 >        p->setCoefficient(exponent-1, coefficient * exponent);
287        }
288      
289        return p;
# Line 338 | Line 338 | namespace OpenMD {
338          roots.push_back( -fC0 / fC1);
339          return roots;
340        }
341 –        break;      
341        case 2: {
342          Real fC2 = getCoefficient(2);
343          Real fC1 = getCoefficient(1);
# Line 347 | Line 346 | namespace OpenMD {
346          if (abs(fDiscr) <= fEpsilon) {
347            fDiscr = (Real)0.0;
348          }
349 <      
349 >        
350          if (fDiscr < (Real)0.0) {  // complex roots only
351            return roots;
352          }
353 <      
353 >        
354          Real fTmp = ((Real)0.5)/fC2;
355 <      
355 >        
356          if (fDiscr > (Real)0.0) { // 2 real roots
357            fDiscr = sqrt(fDiscr);
358            roots.push_back(fTmp*(-fC1 - fDiscr));
# Line 362 | Line 361 | namespace OpenMD {
361            roots.push_back(-fTmp * fC1);  // 1 real root
362          }
363        }
364 <        return roots;
366 <        break;
367 <      
364 >        return roots;      
365        case 3: {
366          Real fC3 = getCoefficient(3);
367          Real fC2 = getCoefficient(2);
# Line 431 | Line 428 | namespace OpenMD {
428          }
429        }
430          return roots;
434 –        break;
431        case 4: {
432          Real fC4 = getCoefficient(4);
433          Real fC3 = getCoefficient(3);
# Line 517 | Line 513 | namespace OpenMD {
513          }
514        }
515          return roots;
520 –        break;
516        default: {
517          DynamicRectMatrix<Real> companion = CreateCompanion();
518          JAMA::Eigenvalue<Real> eig(companion);
# Line 531 | Line 526 | namespace OpenMD {
526          }      
527        }
528          return roots;
534 –        break;
529        }
536 –
537 –      return roots; // should be empty if you got here
530      }
531 <  
531 >    
532    private:
533          
534      PolynomialPairMap polyPairMap_;
# Line 630 | Line 622 | namespace OpenMD {
622     * @return the first derivative of this polynomial
623     */
624    template<typename Real>
625 <  Polynomial<Real> getDerivative(const Polynomial<Real>& p1) {
626 <    Polynomial<Real> p;
625 >  Polynomial<Real> * getDerivative(const Polynomial<Real>& p1) {
626 >    Polynomial<Real> * p = new Polynomial<Real>();
627      
628      typename Polynomial<Real>::const_iterator i;
629      int exponent;
# Line 640 | Line 632 | namespace OpenMD {
632      for (i =  p1.begin(); i  != p1.end(); ++i) {
633        exponent = i->first;
634        coefficient = i->second;
635 <      p.setCoefficient(exponent-1, coefficient * exponent);
635 >      p->setCoefficient(exponent-1, coefficient * exponent);
636      }
637      
638      return p;

Diff Legend

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