--- branches/development/src/math/CubicSpline.hpp	2010/07/26 19:00:48	1479
+++ trunk/src/math/CubicSpline.hpp	2015/03/03 15:22:26	2057
@@ -35,8 +35,9 @@
  *                                                                      
  * [1]  Meineke, et al., J. Comp. Chem. 26, 252-271 (2005).             
  * [2]  Fennell & Gezelter, J. Chem. Phys. 124, 234104 (2006).          
- * [3]  Sun, Lin & Gezelter, J. Chem. Phys. 128, 24107 (2008).          
- * [4]  Vardeman & Gezelter, in progress (2009).                        
+ * [3]  Sun, Lin & Gezelter, J. Chem. Phys. 128, 234107 (2008).          
+ * [4]  Kuang & Gezelter,  J. Chem. Phys. 133, 164101 (2010).
+ * [5]  Vardeman, Stocker & Gezelter, J. Chem. Theory Comput. 7, 834 (2011).
  */
  
 #ifndef MATH_CUBICSPLINE_HPP
@@ -53,22 +54,38 @@ namespace OpenMD {
   public:    
     CubicSpline();
     virtual ~CubicSpline() {} 
-    void addPoint(RealType xp, RealType yp);
+    void addPoint(const RealType xp, const RealType yp);
     void addPoints(const vector<RealType>& xps, const vector<RealType>& yps);
-    RealType getValueAt(RealType t);
-    pair<RealType, RealType> getValueAndDerivativeAt(RealType t);
+    RealType getValueAt(const RealType& t);
+    pair<RealType, RealType> getLimits();
+    void getValueAt(const RealType& t, RealType& v);
+    void getValueAndDerivativeAt(const RealType& t, RealType& v, RealType& d);
     
   private:
     void generate();
+    std::vector<int> sort_permutation(std::vector<RealType>& v);
+    std::vector<RealType> apply_permutation(std::vector<RealType> const& v,
+                                            std::vector<int> const& p);
+    
     bool isUniform;
     bool generated;
-    RealType dx;
-    int n;
-    vector<pair<RealType, RealType> > data;
+    RealType dx, dt, yval, dydx;
+    int n, j;
+    vector<RealType> x_;
+    vector<RealType> y_;
     vector<RealType> b;
     vector<RealType> c;
     vector<RealType> d;    
   };
+
+  class Comparator{
+    std::vector<RealType>& _v;
+  public:
+    Comparator(std::vector<RealType>& v) : _v(v) {}
+    bool operator()(size_t i, size_t j){
+      return _v[i] < _v[j];
+    }
+  };
 }
 
 #endif