Implementation of ContinuousPeriodic1DFunction: Reference platform

openmmapi/include/openmm/TabulatedFunction.h

@@ -111,6 +111,55 @@ private:
     double min, max;
 };
 
+/**
+ * This is a TabulatedFunction that computes a continuous, periodic one dimensional function.
+ */
+class OPENMM_EXPORT ContinuousPeriodic1DFunction : public TabulatedFunction {
+public:
+    /**
+     * Create a ContinuousPeriodic1DFunction f(x) based on a set of tabulated values.
+     *
+     * @param values         the tabulated values of the function f(x) at uniformly spaced values of x between min
+     *                       and max.  A periodic cubic spline is used to interpolate between the tabulated values.
+     *                       The function is assumed to be periodic with period L=max-min.  The first and last
+     *                       elements must have the same value.
+     * @param min            the value of x corresponding to the first element of values
+     * @param max            the value of x corresponding to the last element of values
+     */
+    ContinuousPeriodic1DFunction(const std::vector<double>& values, double min, double max);
+    /**
+     * Get the parameters for the tabulated function.
+     *
+     * @param values         the tabulated values of the function f(x) at uniformly spaced values of x between min
+     *                       and max.  A periodic cubic spline is used to interpolate between the tabulated values.
+     *                       The function is assumed to be periodic with period L=max-min.  The first and last
+     *                       elements must have the same value.
+     * @param[out] min       the value of x corresponding to the first element of values
+     * @param[out] max       the value of x corresponding to the last element of values
+     */
+    void getFunctionParameters(std::vector<double>& values, double& min, double& max) const;
+    /**
+     * Set the parameters for the tabulated function.
+     *
+     * @param values         the tabulated values of the function f(x) at uniformly spaced values of x between min
+     *                       and max.  A periodic cubic spline is used to interpolate between the tabulated values.
+     *                       The function is assumed to be periodic with period L=max-min.  The first and last
+     *                       elements must have the same value.
+     * @param min            the value of x corresponding to the first element of values
+     * @param max            the value of x corresponding to the last element of values
+     */
+    void setFunctionParameters(const std::vector<double>& values, double min, double max);
+    /**
+     * Create a deep copy of the tabulated function.
+     *
+     * @deprecated This will be removed in a future release.
+     */
+    ContinuousPeriodic1DFunction* Copy() const;
+private:
+    std::vector<double> values;
+    double min, max;
+};
+
 /**
  * This is a TabulatedFunction that computes a continuous two dimensional function.
  */

openmmapi/src/TabulatedFunction.cpp

@@ -68,6 +68,45 @@ Continuous1DFunction* Continuous1DFunction::Copy() const {
     return new Continuous1DFunction(new_vec, min, max);
 }
 
+ContinuousPeriodic1DFunction::ContinuousPeriodic1DFunction(const vector<double>& values, double min, double max) {
+    if (max <= min)
+        throw OpenMMException("ContinuousPeriodic1DFunction: max <= min for a tabulated function.");
+    int n = values.size();
+    if (n < 3)
+        throw OpenMMException("ContinuousPeriodic1DFunction: a periodic tabulated function must have at least three points");
+    if (values[0] != values[n-1])
+        throw OpenMMException("ContinuousPeriodic1DFunction: the first and last points must have the same value");
+    this->values = values;
+    this->min = min;
+    this->max = max;
+}
+
+void ContinuousPeriodic1DFunction::getFunctionParameters(vector<double>& values, double& min, double& max) const {
+    values = this->values;
+    min = this->min;
+    max = this->max;
+}
+
+void ContinuousPeriodic1DFunction::setFunctionParameters(const vector<double>& values, double min, double max) {
+    if (max <= min)
+        throw OpenMMException("ContinuousPeriodic1DFunction: max <= min for a tabulated function.");
+    int n = values.size();
+    if (n < 3)
+        throw OpenMMException("ContinuousPeriodic1DFunction: a periodic tabulated function must have at least three points");
+    if (values[0] != values[n-1])
+        throw OpenMMException("ContinuousPeriodic1DFunction: the first and last points must have the same value");
+    this->values = values;
+    this->min = min;
+    this->max = max;
+}
+
+ContinuousPeriodic1DFunction* ContinuousPeriodic1DFunction::Copy() const {
+    vector<double> new_vec(values.size());
+    for (size_t i = 0; i < values.size(); i++)
+        new_vec[i] = values[i];
+    return new ContinuousPeriodic1DFunction(new_vec, min, max);
+}
+
 Continuous2DFunction::Continuous2DFunction(int xsize, int ysize, const vector<double>& values, double xmin, double xmax, double ymin, double ymax) {
     if (xsize < 2 || ysize < 2)
         throw OpenMMException("Continuous2DFunction: must have at least two points along each axis");

platforms/reference/include/ReferenceTabulatedFunction.h

@@ -62,6 +62,23 @@ private:
     std::vector<double> x, values, derivs;
 };
 
+/**
+ * This class adapts a ContinuousPeriodic1DFunction into a Lepton::CustomFunction.
+ */
+class OPENMM_EXPORT ReferenceContinuousPeriodic1DFunction : public Lepton::CustomFunction {
+public:
+    ReferenceContinuousPeriodic1DFunction(const ContinuousPeriodic1DFunction& function);
+    int getNumArguments() const;
+    double evaluate(const double* arguments) const;
+    double evaluateDerivative(const double* arguments, const int* derivOrder) const;
+    CustomFunction* clone() const;
+private:
+    ReferenceContinuousPeriodic1DFunction(const ReferenceContinuousPeriodic1DFunction& other);
+    const ContinuousPeriodic1DFunction& function;
+    double min, max;
+    std::vector<double> x, values, derivs;
+};
+
 /**
  * This class adapts a Continuous2DFunction into a Lepton::CustomFunction.
  */

platforms/reference/src/ReferenceTabulatedFunction.cpp

@@ -59,6 +59,8 @@ extern "C" OPENMM_EXPORT CustomFunction* createReferenceTabulatedFunction(const
     CustomFunction* fn;
     if (dynamic_cast<const Continuous1DFunction*>(&function) != NULL)
         fn = new ReferenceContinuous1DFunction(dynamic_cast<const Continuous1DFunction&>(function));
+    else if (dynamic_cast<const ContinuousPeriodic1DFunction*>(&function) != NULL)
+        fn = new ReferenceContinuousPeriodic1DFunction(dynamic_cast<const ContinuousPeriodic1DFunction&>(function));
     else if (dynamic_cast<const Continuous2DFunction*>(&function) != NULL)
         fn = new ReferenceContinuous2DFunction(dynamic_cast<const Continuous2DFunction&>(function));
     else if (dynamic_cast<const Continuous3DFunction*>(&function) != NULL)
@@ -112,6 +114,44 @@ CustomFunction* ReferenceContinuous1DFunction::clone() const {
     return new ReferenceContinuous1DFunction(*this);
 }
 
+ReferenceContinuousPeriodic1DFunction::ReferenceContinuousPeriodic1DFunction(const ContinuousPeriodic1DFunction& function) : function(function) {
+    function.getFunctionParameters(values, min, max);
+    int numValues = values.size();
+    x.resize(numValues);
+    for (int i = 0; i < numValues; i++)
+        x[i] = min+i*(max-min)/(numValues-1);
+    SplineFitter::createPeriodicSpline(x, values, derivs);
+}
+
+ReferenceContinuousPeriodic1DFunction::ReferenceContinuousPeriodic1DFunction(const ReferenceContinuousPeriodic1DFunction& other) : function(other.function) {
+    function.getFunctionParameters(values, min, max);
+    x = other.x;
+    values = other.values;
+    derivs = other.derivs;
+}
+
+int ReferenceContinuousPeriodic1DFunction::getNumArguments() const {
+    return 1;
+}
+
+double ReferenceContinuousPeriodic1DFunction::evaluate(const double* arguments) const {
+    double t = arguments[0];
+    if (t < min || t > max)
+        return 0.0;
+    return SplineFitter::evaluateSpline(x, values, derivs, t);
+}
+
+double ReferenceContinuousPeriodic1DFunction::evaluateDerivative(const double* arguments, const int* derivOrder) const {
+    double t = arguments[0];
+    if (t < min || t > max)
+        return 0.0;
+    return SplineFitter::evaluateSplineDerivative(x, values, derivs, t);
+}
+
+CustomFunction* ReferenceContinuousPeriodic1DFunction::clone() const {
+    return new ReferenceContinuousPeriodic1DFunction(*this);
+}
+
 ReferenceContinuous2DFunction::ReferenceContinuous2DFunction(const Continuous2DFunction& function) : function(function) {
     function.getFunctionParameters(xsize, ysize, values, xmin, xmax, ymin, ymax);
     x.resize(xsize);