Converted CUDA and OpenCL platforms to process tabulated functions with the TabulatedFunction API

platforms/cuda/include/CudaExpressionUtilities.h

@@ -9,7 +9,7 @@
  * Biological Structures at Stanford, funded under the NIH Roadmap for        *
  * Medical Research, grant U54 GM072970. See https://simtk.org.               *
  *                                                                            *
- * Portions copyright (c) 2009-2012 Stanford University and the Authors.      *
+ * Portions copyright (c) 2009-2014 Stanford University and the Authors.      *
  * Authors: Peter Eastman                                                     *
  * Contributors:                                                              *
  *                                                                            *
@@ -28,6 +28,7 @@
  * -------------------------------------------------------------------------- */
 
 #include "CudaContext.h"
+#include "openmm/TabulatedFunction.h"
 #include "lepton/CustomFunction.h"
 #include "lepton/ExpressionTreeNode.h"
 #include "lepton/ParsedExpression.h"
@@ -76,12 +77,17 @@ public:
     /**
      * Calculate the spline coefficients for a tabulated function that appears in expressions.
      *
-     * @param values         the tabulated values of the function
-     * @param min            the value of the independent variable corresponding to the first element of values
-     * @param max            the value of the independent variable corresponding to the last element of values
+     * @param function   the function for which to compute coefficients
      * @return the spline coefficients
      */
-    std::vector<float4> computeFunctionCoefficients(const std::vector<double>& values, double min, double max);
+    std::vector<float> computeFunctionCoefficients(const TabulatedFunction& function);
+    /**
+     * Given the list of TabulatedFunctions used by a Force, create the parameter array describing them.
+     *
+     * @param functions   the list of functions to include in the array
+     * @return the parameter array
+     */
+    std::vector<float4> computeFunctionParameters(const std::vector<const TabulatedFunction*>& functions);
     class FunctionPlaceholder;
 private:
     void processExpression(std::stringstream& out, const Lepton::ExpressionTreeNode& node,

platforms/cuda/src/CudaExpressionUtilities.cpp

@@ -6,7 +6,7 @@
  * Biological Structures at Stanford, funded under the NIH Roadmap for        *
  * Medical Research, grant U54 GM072970. See https://simtk.org.               *
  *                                                                            *
- * Portions copyright (c) 2009-2012 Stanford University and the Authors.      *
+ * Portions copyright (c) 2009-2014 Stanford University and the Authors.      *
  * Authors: Peter Eastman                                                     *
  * Contributors:                                                              *
  *                                                                            *
@@ -341,16 +341,43 @@ void CudaExpressionUtilities::findRelatedPowers(const ExpressionTreeNode& node,
             findRelatedPowers(node, searchNode.getChildren()[i], powers);
 }
 
-vector<float4> CudaExpressionUtilities::computeFunctionCoefficients(const vector<double>& values, double min, double max) {
+vector<float> CudaExpressionUtilities::computeFunctionCoefficients(const TabulatedFunction& function) {
     // Compute the spline coefficients.
 
+    if (dynamic_cast<const Continuous1DFunction*>(&function) != NULL) {
+        const Continuous1DFunction& fn = dynamic_cast<const Continuous1DFunction&>(function);
+        vector<double> values;
+        double min, max;
+        fn.getFunctionParameters(values, min, max);
         int numValues = values.size();
         vector<double> x(numValues), derivs;
         for (int i = 0; i < numValues; i++)
             x[i] = min+i*(max-min)/(numValues-1);
         SplineFitter::createNaturalSpline(x, values, derivs);
-    vector<float4> f(numValues-1);
-    for (int i = 0; i < (int) values.size()-1; i++)
-        f[i] = make_float4((float) values[i], (float) values[i+1], (float) (derivs[i]/6.0), (float) (derivs[i+1]/6.0));
+        vector<float> f(4*(numValues-1));
+        for (int i = 0; i < (int) values.size()-1; i++) {
+            f[4*i] = (float) values[i];
+            f[4*i+1] = (float) values[i+1];
+            f[4*i+2] = (float) (derivs[i]/6.0);
+            f[4*i+3] = (float) (derivs[i+1]/6.0);
+        }
         return f;
+    }
+    throw OpenMMException("computeFunctionCoefficients: Unknown function type");
+}
+
+vector<float4> CudaExpressionUtilities::computeFunctionParameters(const vector<const TabulatedFunction*>& functions) {
+    vector<float4> params(functions.size());
+    for (int i = 0; i < (int) functions.size(); i++) {
+        if (dynamic_cast<const Continuous1DFunction*>(functions[i]) != NULL) {
+            const Continuous1DFunction& fn = dynamic_cast<const Continuous1DFunction&>(*functions[i]);
+            vector<double> values;
+            double min, max;
+            fn.getFunctionParameters(values, min, max);
+            params[i] = make_float4((float) min, (float) max, (float) ((values.size()-1)/(max-min)), (float) values.size()-2);
+        }
+        else
+            throw OpenMMException("computeFunctionParameters: Unknown function type");
+    }
+    return params;
 }

platforms/cuda/src/CudaKernels.cpp

@@ -6,7 +6,7 @@
  * Biological Structures at Stanford, funded under the NIH Roadmap for        *
  * Medical Research, grant U54 GM072970. See https://simtk.org.               *
  *                                                                            *
- * Portions copyright (c) 2008-2013 Stanford University and the Authors.      *
+ * Portions copyright (c) 2008-2014 Stanford University and the Authors.      *
  * Authors: Peter Eastman                                                     *
  * Contributors:                                                              *
  *                                                                            *
@@ -1958,21 +1958,19 @@ void CudaCalcCustomNonbondedForceKernel::initialize(const System& system, const
     CudaExpressionUtilities::FunctionPlaceholder fp;
     map<string, Lepton::CustomFunction*> functions;
     vector<pair<string, string> > functionDefinitions;
-    vector<float4> tabulatedFunctionParamsVec(force.getNumFunctions());
+    vector<const TabulatedFunction*> functionList;
     for (int i = 0; i < force.getNumFunctions(); i++) {
-        string name;
-        vector<double> values;
-        double min, max;
-        force.getFunctionParameters(i, name, values, min, max);
+        functionList.push_back(&force.getFunction(i));
+        string name = force.getFunctionName(i);
         string arrayName = prefix+"table"+cu.intToString(i);
         functionDefinitions.push_back(make_pair(name, arrayName));
         functions[name] = &fp;
-        tabulatedFunctionParamsVec[i] = make_float4((float) min, (float) max, (float) ((values.size()-1)/(max-min)), (float) values.size()-2);
-        vector<float4> f = cu.getExpressionUtilities().computeFunctionCoefficients(values, min, max);
-        tabulatedFunctions.push_back(CudaArray::create<float4>(cu, values.size()-1, "TabulatedFunction"));
+        vector<float> f = cu.getExpressionUtilities().computeFunctionCoefficients(force.getFunction(i));
+        tabulatedFunctions.push_back(CudaArray::create<float>(cu, f.size(), "TabulatedFunction"));
         tabulatedFunctions[tabulatedFunctions.size()-1]->upload(f);
         cu.getNonbondedUtilities().addArgument(CudaNonbondedUtilities::ParameterInfo(arrayName, "float", 4, sizeof(float4), tabulatedFunctions[tabulatedFunctions.size()-1]->getDevicePointer()));
     }
+    vector<float4> tabulatedFunctionParamsVec = cu.getExpressionUtilities().computeFunctionParameters(functionList);
     if (force.getNumFunctions() > 0) {
         tabulatedFunctionParams = CudaArray::create<float4>(cu, tabulatedFunctionParamsVec.size(), "tabulatedFunctionParameters");
         tabulatedFunctionParams->upload(tabulatedFunctionParamsVec);
@@ -2672,23 +2670,21 @@ void CudaCalcCustomGBForceKernel::initialize(const System& system, const CustomG
     CudaExpressionUtilities::FunctionPlaceholder fp;
     map<string, Lepton::CustomFunction*> functions;
     vector<pair<string, string> > functionDefinitions;
-    vector<float4> tabulatedFunctionParamsVec(force.getNumFunctions());
+    vector<const TabulatedFunction*> functionList;
     stringstream tableArgs;
     for (int i = 0; i < force.getNumFunctions(); i++) {
-        string name;
-        vector<double> values;
-        double min, max;
-        force.getFunctionParameters(i, name, values, min, max);
+        functionList.push_back(&force.getFunction(i));
+        string name = force.getFunctionName(i);
         string arrayName = prefix+"table"+cu.intToString(i);
         functionDefinitions.push_back(make_pair(name, arrayName));
         functions[name] = &fp;
-        tabulatedFunctionParamsVec[i] = make_float4((float) min, (float) max, (float) ((values.size()-1)/(max-min)), (float) values.size()-2);
-        vector<float4> f = cu.getExpressionUtilities().computeFunctionCoefficients(values, min, max);
-        tabulatedFunctions.push_back(CudaArray::create<float4>(cu, values.size()-1, "TabulatedFunction"));
+        vector<float> f = cu.getExpressionUtilities().computeFunctionCoefficients(force.getFunction(i));
+        tabulatedFunctions.push_back(CudaArray::create<float>(cu, f.size(), "TabulatedFunction"));
         tabulatedFunctions[tabulatedFunctions.size()-1]->upload(f);
         cu.getNonbondedUtilities().addArgument(CudaNonbondedUtilities::ParameterInfo(arrayName, "float", 4, sizeof(float4), tabulatedFunctions[tabulatedFunctions.size()-1]->getDevicePointer()));
         tableArgs << ", const float4* __restrict__ " << arrayName;
     }
+    vector<float4> tabulatedFunctionParamsVec = cu.getExpressionUtilities().computeFunctionParameters(functionList);
     if (force.getNumFunctions() > 0) {
         tabulatedFunctionParams = CudaArray::create<float4>(cu, tabulatedFunctionParamsVec.size(), "tabulatedFunctionParameters");
         tabulatedFunctionParams->upload(tabulatedFunctionParamsVec);
@@ -3787,22 +3783,20 @@ void CudaCalcCustomHbondForceKernel::initialize(const System& system, const Cust
     CudaExpressionUtilities::FunctionPlaceholder fp;
     map<string, Lepton::CustomFunction*> functions;
     vector<pair<string, string> > functionDefinitions;
-    vector<float4> tabulatedFunctionParamsVec(force.getNumFunctions());
+    vector<const TabulatedFunction*> functionList;
     stringstream tableArgs;
     for (int i = 0; i < force.getNumFunctions(); i++) {
-        string name;
-        vector<double> values;
-        double min, max;
-        force.getFunctionParameters(i, name, values, min, max);
+        functionList.push_back(&force.getFunction(i));
+        string name = force.getFunctionName(i);
         string arrayName = "table"+cu.intToString(i);
         functionDefinitions.push_back(make_pair(name, arrayName));
         functions[name] = &fp;
-        tabulatedFunctionParamsVec[i] = make_float4((float) min, (float) max, (float) ((values.size()-1)/(max-min)), (float) values.size()-2);
-        vector<float4> f = cu.getExpressionUtilities().computeFunctionCoefficients(values, min, max);
-        tabulatedFunctions.push_back(CudaArray::create<float4>(cu, values.size()-1, "TabulatedFunction"));
+        vector<float> f = cu.getExpressionUtilities().computeFunctionCoefficients(force.getFunction(i));
+        tabulatedFunctions.push_back(CudaArray::create<float>(cu, f.size(), "TabulatedFunction"));
         tabulatedFunctions[tabulatedFunctions.size()-1]->upload(f);
         tableArgs << ", const float4* __restrict__ " << arrayName;
     }
+    vector<float4> tabulatedFunctionParamsVec = cu.getExpressionUtilities().computeFunctionParameters(functionList);
     if (force.getNumFunctions() > 0) {
         tabulatedFunctionParams = CudaArray::create<float4>(cu, tabulatedFunctionParamsVec.size(), "tabulatedFunctionParameters");
         tabulatedFunctionParams->upload(tabulatedFunctionParamsVec);
@@ -4181,22 +4175,20 @@ void CudaCalcCustomCompoundBondForceKernel::initialize(const System& system, con
     CudaExpressionUtilities::FunctionPlaceholder fp;
     map<string, Lepton::CustomFunction*> functions;
     vector<pair<string, string> > functionDefinitions;
-    vector<float4> tabulatedFunctionParamsVec(force.getNumFunctions());
+    vector<const TabulatedFunction*> functionList;
     stringstream tableArgs;
     for (int i = 0; i < force.getNumFunctions(); i++) {
-        string name;
-        vector<double> values;
-        double min, max;
-        force.getFunctionParameters(i, name, values, min, max);
+        functionList.push_back(&force.getFunction(i));
+        string name = force.getFunctionName(i);
         functions[name] = &fp;
-        tabulatedFunctionParamsVec[i] = make_float4((float) min, (float) max, (float) ((values.size()-1)/(max-min)), (float) values.size()-2);
-        vector<float4> f = cu.getExpressionUtilities().computeFunctionCoefficients(values, min, max);
-        CudaArray* array = CudaArray::create<float4>(cu, values.size()-1, "TabulatedFunction");
+        vector<float> f = cu.getExpressionUtilities().computeFunctionCoefficients(force.getFunction(i));
+        CudaArray* array = CudaArray::create<float>(cu, f.size(), "TabulatedFunction");
         tabulatedFunctions.push_back(array);
         array->upload(f);
         string arrayName = cu.getBondedUtilities().addArgument(array->getDevicePointer(), "float4");
         functionDefinitions.push_back(make_pair(name, arrayName));
     }
+    vector<float4> tabulatedFunctionParamsVec = cu.getExpressionUtilities().computeFunctionParameters(functionList);
     string functionParamsName;
     if (force.getNumFunctions() > 0) {
         tabulatedFunctionParams = CudaArray::create<float4>(cu, tabulatedFunctionParamsVec.size(), "tabulatedFunctionParameters");

platforms/opencl/include/OpenCLExpressionUtilities.h

@@ -9,7 +9,7 @@
  * Biological Structures at Stanford, funded under the NIH Roadmap for        *
  * Medical Research, grant U54 GM072970. See https://simtk.org.               *
  *                                                                            *
- * Portions copyright (c) 2009-2011 Stanford University and the Authors.      *
+ * Portions copyright (c) 2009-2014 Stanford University and the Authors.      *
  * Authors: Peter Eastman                                                     *
  * Contributors:                                                              *
  *                                                                            *
@@ -28,6 +28,7 @@
  * -------------------------------------------------------------------------- */
 
 #include "OpenCLContext.h"
+#include "openmm/TabulatedFunction.h"
 #include "lepton/CustomFunction.h"
 #include "lepton/ExpressionTreeNode.h"
 #include "lepton/ParsedExpression.h"
@@ -76,12 +77,17 @@ public:
     /**
      * Calculate the spline coefficients for a tabulated function that appears in expressions.
      *
-     * @param values         the tabulated values of the function
-     * @param min            the value of the independent variable corresponding to the first element of values
-     * @param max            the value of the independent variable corresponding to the last element of values
+     * @param function   the function for which to compute coefficients
      * @return the spline coefficients
      */
-    std::vector<mm_float4> computeFunctionCoefficients(const std::vector<double>& values, double min, double max);
+    std::vector<float> computeFunctionCoefficients(const TabulatedFunction& function);
+    /**
+     * Given the list of TabulatedFunctions used by a Force, create the parameter array describing them.
+     *
+     * @param functions   the list of functions to include in the array
+     * @return the parameter array
+     */
+    std::vector<mm_float4> computeFunctionParameters(const std::vector<const TabulatedFunction*>& functions);
     class FunctionPlaceholder;
 private:
     void processExpression(std::stringstream& out, const Lepton::ExpressionTreeNode& node,

platforms/opencl/src/OpenCLExpressionUtilities.cpp

@@ -6,7 +6,7 @@
  * Biological Structures at Stanford, funded under the NIH Roadmap for        *
  * Medical Research, grant U54 GM072970. See https://simtk.org.               *
  *                                                                            *
- * Portions copyright (c) 2009-2011 Stanford University and the Authors.      *
+ * Portions copyright (c) 2009-2014 Stanford University and the Authors.      *
  * Authors: Peter Eastman                                                     *
  * Contributors:                                                              *
  *                                                                            *
@@ -341,16 +341,43 @@ void OpenCLExpressionUtilities::findRelatedPowers(const ExpressionTreeNode& node
             findRelatedPowers(node, searchNode.getChildren()[i], powers);
 }
 
-vector<mm_float4> OpenCLExpressionUtilities::computeFunctionCoefficients(const vector<double>& values, double min, double max) {
+vector<float> OpenCLExpressionUtilities::computeFunctionCoefficients(const TabulatedFunction& function) {
     // Compute the spline coefficients.
 
+    if (dynamic_cast<const Continuous1DFunction*>(&function) != NULL) {
+        const Continuous1DFunction& fn = dynamic_cast<const Continuous1DFunction&>(function);
+        vector<double> values;
+        double min, max;
+        fn.getFunctionParameters(values, min, max);
         int numValues = values.size();
         vector<double> x(numValues), derivs;
         for (int i = 0; i < numValues; i++)
             x[i] = min+i*(max-min)/(numValues-1);
         SplineFitter::createNaturalSpline(x, values, derivs);
-    vector<mm_float4> f(numValues-1);
-    for (int i = 0; i < (int) values.size()-1; i++)
-        f[i] = mm_float4((cl_float) values[i], (cl_float) values[i+1], (cl_float) (derivs[i]/6.0), (cl_float) (derivs[i+1]/6.0));
+        vector<float> f(4*(numValues-1));
+        for (int i = 0; i < (int) values.size()-1; i++) {
+            f[4*i] = (float) values[i];
+            f[4*i+1] = (float) values[i+1];
+            f[4*i+2] = (float) (derivs[i]/6.0);
+            f[4*i+3] = (float) (derivs[i+1]/6.0);
+        }
         return f;
+    }
+    throw OpenMMException("computeFunctionCoefficients: Unknown function type");
+}
+
+vector<mm_float4> OpenCLExpressionUtilities::computeFunctionParameters(const vector<const TabulatedFunction*>& functions) {
+    vector<mm_float4> params(functions.size());
+    for (int i = 0; i < (int) functions.size(); i++) {
+        if (dynamic_cast<const Continuous1DFunction*>(functions[i]) != NULL) {
+            const Continuous1DFunction& fn = dynamic_cast<const Continuous1DFunction&>(*functions[i]);
+            vector<double> values;
+            double min, max;
+            fn.getFunctionParameters(values, min, max);
+            params[i] = mm_float4((float) min, (float) max, (float) ((values.size()-1)/(max-min)), (float) values.size()-2);
+        }
+        else
+            throw OpenMMException("computeFunctionParameters: Unknown function type");
+    }
+    return params;
 }

platforms/opencl/src/OpenCLKernels.cpp

@@ -6,7 +6,7 @@
  * Biological Structures at Stanford, funded under the NIH Roadmap for        *
  * Medical Research, grant U54 GM072970. See https://simtk.org.               *
  *                                                                            *
- * Portions copyright (c) 2008-2013 Stanford University and the Authors.      *
+ * Portions copyright (c) 2008-2014 Stanford University and the Authors.      *
  * Authors: Peter Eastman                                                     *
  * Contributors:                                                              *
  *                                                                            *
@@ -1968,21 +1968,19 @@ void OpenCLCalcCustomNonbondedForceKernel::initialize(const System& system, cons
     OpenCLExpressionUtilities::FunctionPlaceholder fp;
     map<string, Lepton::CustomFunction*> functions;
     vector<pair<string, string> > functionDefinitions;
-    vector<mm_float4> tabulatedFunctionParamsVec(force.getNumFunctions());
+    vector<const TabulatedFunction*> functionList;
     for (int i = 0; i < force.getNumFunctions(); i++) {
-        string name;
-        vector<double> values;
-        double min, max;
-        force.getFunctionParameters(i, name, values, min, max);
+        functionList.push_back(&force.getFunction(i));
+        string name = force.getFunctionName(i);
         string arrayName = prefix+"table"+cl.intToString(i);
         functionDefinitions.push_back(make_pair(name, arrayName));
         functions[name] = &fp;
-        tabulatedFunctionParamsVec[i] = mm_float4((float) min, (float) max, (float) ((values.size()-1)/(max-min)), (float) values.size()-2);
-        vector<mm_float4> f = cl.getExpressionUtilities().computeFunctionCoefficients(values, min, max);
-        tabulatedFunctions.push_back(OpenCLArray::create<mm_float4>(cl, values.size()-1, "TabulatedFunction"));
+        vector<float> f = cl.getExpressionUtilities().computeFunctionCoefficients(force.getFunction(i));
+        tabulatedFunctions.push_back(OpenCLArray::create<float>(cl, f.size(), "TabulatedFunction"));
         tabulatedFunctions[tabulatedFunctions.size()-1]->upload(f);
         cl.getNonbondedUtilities().addArgument(OpenCLNonbondedUtilities::ParameterInfo(arrayName, "float", 4, sizeof(cl_float4), tabulatedFunctions[tabulatedFunctions.size()-1]->getDeviceBuffer()));
     }
+    vector<mm_float4> tabulatedFunctionParamsVec = cl.getExpressionUtilities().computeFunctionParameters(functionList);
     if (force.getNumFunctions() > 0) {
         tabulatedFunctionParams = OpenCLArray::create<mm_float4>(cl, tabulatedFunctionParamsVec.size(), "tabulatedFunctionParameters", CL_MEM_READ_ONLY);
         tabulatedFunctionParams->upload(tabulatedFunctionParamsVec);
@@ -2725,23 +2723,21 @@ void OpenCLCalcCustomGBForceKernel::initialize(const System& system, const Custo
     OpenCLExpressionUtilities::FunctionPlaceholder fp;
     map<string, Lepton::CustomFunction*> functions;
     vector<pair<string, string> > functionDefinitions;
-    vector<mm_float4> tabulatedFunctionParamsVec(force.getNumFunctions());
+    vector<const TabulatedFunction*> functionList;
     stringstream tableArgs;
     for (int i = 0; i < force.getNumFunctions(); i++) {
-        string name;
-        vector<double> values;
-        double min, max;
-        force.getFunctionParameters(i, name, values, min, max);
+        functionList.push_back(&force.getFunction(i));
+        string name = force.getFunctionName(i);
         string arrayName = prefix+"table"+cl.intToString(i);
         functionDefinitions.push_back(make_pair(name, arrayName));
         functions[name] = &fp;
-        tabulatedFunctionParamsVec[i] = mm_float4((float) min, (float) max, (float) ((values.size()-1)/(max-min)), (float) values.size()-2);
-        vector<mm_float4> f = cl.getExpressionUtilities().computeFunctionCoefficients(values, min, max);
-        tabulatedFunctions.push_back(OpenCLArray::create<mm_float4>(cl, values.size()-1, "TabulatedFunction"));
+        vector<float> f = cl.getExpressionUtilities().computeFunctionCoefficients(force.getFunction(i));
+        tabulatedFunctions.push_back(OpenCLArray::create<float>(cl, f.size(), "TabulatedFunction"));
         tabulatedFunctions[tabulatedFunctions.size()-1]->upload(f);
         cl.getNonbondedUtilities().addArgument(OpenCLNonbondedUtilities::ParameterInfo(arrayName, "float", 4, sizeof(cl_float4), tabulatedFunctions[tabulatedFunctions.size()-1]->getDeviceBuffer()));
         tableArgs << ", __global const float4* restrict " << arrayName;
     }
+    vector<mm_float4> tabulatedFunctionParamsVec = cl.getExpressionUtilities().computeFunctionParameters(functionList);
     if (force.getNumFunctions() > 0) {
         tabulatedFunctionParams = OpenCLArray::create<mm_float4>(cl, tabulatedFunctionParamsVec.size(), "tabulatedFunctionParameters", CL_MEM_READ_ONLY);
         tabulatedFunctionParams->upload(tabulatedFunctionParamsVec);
@@ -3950,22 +3946,20 @@ void OpenCLCalcCustomHbondForceKernel::initialize(const System& system, const Cu
     OpenCLExpressionUtilities::FunctionPlaceholder fp;
     map<string, Lepton::CustomFunction*> functions;
     vector<pair<string, string> > functionDefinitions;
-    vector<mm_float4> tabulatedFunctionParamsVec(force.getNumFunctions());
+    vector<const TabulatedFunction*> functionList;
     stringstream tableArgs;
     for (int i = 0; i < force.getNumFunctions(); i++) {
-        string name;
-        vector<double> values;
-        double min, max;
-        force.getFunctionParameters(i, name, values, min, max);
+        functionList.push_back(&force.getFunction(i));
+        string name = force.getFunctionName(i);
         string arrayName = "table"+cl.intToString(i);
         functionDefinitions.push_back(make_pair(name, arrayName));
         functions[name] = &fp;
-        tabulatedFunctionParamsVec[i] = mm_float4((float) min, (float) max, (float) ((values.size()-1)/(max-min)), (float) values.size()-2);
-        vector<mm_float4> f = cl.getExpressionUtilities().computeFunctionCoefficients(values, min, max);
-        tabulatedFunctions.push_back(OpenCLArray::create<mm_float4>(cl, values.size()-1, "TabulatedFunction"));
+        vector<float> f = cl.getExpressionUtilities().computeFunctionCoefficients(force.getFunction(i));
+        tabulatedFunctions.push_back(OpenCLArray::create<float>(cl, f.size(), "TabulatedFunction"));
         tabulatedFunctions[tabulatedFunctions.size()-1]->upload(f);
         tableArgs << ", __global const float4* restrict " << arrayName;
     }
+    vector<mm_float4> tabulatedFunctionParamsVec = cl.getExpressionUtilities().computeFunctionParameters(functionList);
     if (force.getNumFunctions() > 0) {
         tabulatedFunctionParams = OpenCLArray::create<mm_float4>(cl, tabulatedFunctionParamsVec.size(), "tabulatedFunctionParameters", CL_MEM_READ_ONLY);
         tabulatedFunctionParams->upload(tabulatedFunctionParamsVec);
@@ -4346,22 +4340,20 @@ void OpenCLCalcCustomCompoundBondForceKernel::initialize(const System& system, c
     OpenCLExpressionUtilities::FunctionPlaceholder fp;
     map<string, Lepton::CustomFunction*> functions;
     vector<pair<string, string> > functionDefinitions;
-    vector<mm_float4> tabulatedFunctionParamsVec(force.getNumFunctions());
+    vector<const TabulatedFunction*> functionList;
     stringstream tableArgs;
     for (int i = 0; i < force.getNumFunctions(); i++) {
-        string name;
-        vector<double> values;
-        double min, max;
-        force.getFunctionParameters(i, name, values, min, max);
+        functionList.push_back(&force.getFunction(i));
+        string name = force.getFunctionName(i);
         functions[name] = &fp;
-        tabulatedFunctionParamsVec[i] = mm_float4((float) min, (float) max, (float) ((values.size()-1)/(max-min)), (float) values.size()-2);
-        vector<mm_float4> f = cl.getExpressionUtilities().computeFunctionCoefficients(values, min, max);
-        OpenCLArray* array = OpenCLArray::create<mm_float4>(cl, values.size()-1, "TabulatedFunction");
+        vector<float> f = cl.getExpressionUtilities().computeFunctionCoefficients(force.getFunction(i));
+        OpenCLArray* array = OpenCLArray::create<float>(cl, f.size(), "TabulatedFunction");
         tabulatedFunctions.push_back(array);
         array->upload(f);
         string arrayName = cl.getBondedUtilities().addArgument(array->getDeviceBuffer(), "float4");
         functionDefinitions.push_back(make_pair(name, arrayName));
     }
+    vector<mm_float4> tabulatedFunctionParamsVec = cl.getExpressionUtilities().computeFunctionParameters(functionList);
     string functionParamsName;
     if (force.getNumFunctions() > 0) {
         tabulatedFunctionParams = OpenCLArray::create<mm_float4>(cl, tabulatedFunctionParamsVec.size(), "tabulatedFunctionParameters", CL_MEM_READ_ONLY);