Added updateParametersInContext() to CMAPTorsionForce, and implemented in reference platform

olla/include/openmm/kernels.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) 2008-2012 Stanford University and the Authors.      *
+ * Portions copyright (c) 2008-2015 Stanford University and the Authors.      *
  * Authors: Peter Eastman                                                     *
  * Contributors:                                                              *
  *                                                                            *
@@ -494,6 +494,13 @@ public:
      * @return the potential energy due to the force
      */
     virtual double execute(ContextImpl& context, bool includeForces, bool includeEnergy) = 0;
+    /**
+     * Copy changed parameters over to a context.
+     *
+     * @param context    the context to copy parameters to
+     * @param force      the CMAPTorsionForce to copy the parameters from
+     */
+    virtual void copyParametersToContext(ContextImpl& context, const CMAPTorsionForce& force) = 0;
 };
 
 /**

openmmapi/include/openmm/CMAPTorsionForce.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) 2010 Stanford University and the Authors.           *
+ * Portions copyright (c) 2010-2015 Stanford University and the Authors.      *
  * Authors: Peter Eastman                                                     *
  * Contributors:                                                              *
  *                                                                            *
@@ -148,6 +148,17 @@ public:
      * @param b4    the index of the fourth particle forming the second torsion
      */
     void setTorsionParameters(int index, int map, int a1, int a2, int a3, int a4, int b1, int b2, int b3, int b4);
+    /**
+     * Update the map and torsion parameters in a Context to match those stored in this Force object.  This method provides
+     * an efficient method to update certain parameters in an existing Context without needing to reinitialize it.
+     * Simply call setMapParameters() and setTorsionParameters() to modify this object's parameters, then call updateParametersInContext()
+     * to copy them over to the Context.
+     *
+     * The only information that can be updated with this method is the energy values for a map, and the map index
+     * for a torsion.  The size of a map and the set of particles involved in a torsion cannot be changed.  Also,
+     * new bonds and torsions cannot be added.
+     */
+    void updateParametersInContext(Context& context);
     /**
      * Returns whether or not this force makes use of periodic boundary
      * conditions.

openmmapi/include/openmm/internal/CMAPTorsionForceImpl.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) 2008-2010 Stanford University and the Authors.      *
+ * Portions copyright (c) 2008-2015 Stanford University and the Authors.      *
  * Authors: Peter Eastman                                                     *
  * Contributors:                                                              *
  *                                                                            *
@@ -60,6 +60,7 @@ public:
         return std::map<std::string, double>(); // This force field doesn't define any parameters.
     }
     std::vector<std::string> getKernelNames();
+    void updateParametersInContext(ContextImpl& context);
     /**
      * Given the energy values for a map, compute the spline coefficients at each point of the map.
      */

openmmapi/src/CMAPTorsionForce.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) 2010 Stanford University and the Authors.           *
+ * Portions copyright (c) 2010-2015 Stanford University and the Authors.      *
  * Authors: Peter Eastman                                                     *
  * Contributors:                                                              *
  *                                                                            *
@@ -96,3 +96,6 @@ ForceImpl* CMAPTorsionForce::createImpl() const {
     return new CMAPTorsionForceImpl(*this);
 }
 
+void CMAPTorsionForce::updateParametersInContext(Context& context) {
+    dynamic_cast<CMAPTorsionForceImpl&>(getImplInContext(context)).updateParametersInContext(getContextImpl(context));
+}

openmmapi/src/CMAPTorsionForceImpl.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) 2010 Stanford University and the Authors.           *
+ * Portions copyright (c) 2010-2015 Stanford University and the Authors.      *
  * Authors: Peter Eastman                                                     *
  * Contributors:                                                              *
  *                                                                            *
@@ -155,3 +155,7 @@ void CMAPTorsionForceImpl::calcMapDerivatives(int size, const vector<double>& en
         }
     }
 }
+
+void CMAPTorsionForceImpl::updateParametersInContext(ContextImpl& context) {
+    kernel.getAs<CalcCMAPTorsionForceKernel>().copyParametersToContext(context, owner);
+}

platforms/cuda/include/CudaKernels.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) 2008-2013 Stanford University and the Authors.      *
+ * Portions copyright (c) 2008-2015 Stanford University and the Authors.      *
  * Authors: Peter Eastman                                                     *
  * Contributors:                                                              *
  *                                                                            *
@@ -520,6 +520,13 @@ public:
      * @return the potential energy due to the force
      */
     double execute(ContextImpl& context, bool includeForces, bool includeEnergy);
+    /**
+     * Copy changed parameters over to a context.
+     *
+     * @param context    the context to copy parameters to
+     * @param force      the CMAPTorsionForce to copy the parameters from
+     */
+    void copyParametersToContext(ContextImpl& context, const CMAPTorsionForce& force);
 private:
     int numTorsions;
     bool hasInitializedKernel;

platforms/cuda/include/CudaParallelKernels.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) 2011-2013 Stanford University and the Authors.      *
+ * Portions copyright (c) 2011-2015 Stanford University and the Authors.      *
  * Authors: Peter Eastman                                                     *
  * Contributors:                                                              *
  *                                                                            *
@@ -345,6 +345,13 @@ public:
      * @return the potential energy due to the force
      */
     double execute(ContextImpl& context, bool includeForces, bool includeEnergy);
+    /**
+     * Copy changed parameters over to a context.
+     *
+     * @param context    the context to copy parameters to
+     * @param force      the CMAPTorsionForce to copy the parameters from
+     */
+    void copyParametersToContext(ContextImpl& context, const CMAPTorsionForce& force);
 private:
     class Task;
     CudaPlatform::PlatformData& data;

platforms/cuda/src/CudaKernels.cpp

@@ -1166,6 +1166,9 @@ double CudaCalcCMAPTorsionForceKernel::execute(ContextImpl& context, bool includ
     return 0.0;
 }
 
+void CudaCalcCMAPTorsionForceKernel::copyParametersToContext(ContextImpl& context, const CMAPTorsionForce& force) {
+}
+
 class CudaCustomTorsionForceInfo : public CudaForceInfo {
 public:
     CudaCustomTorsionForceInfo(const CustomTorsionForce& force) : force(force) {

platforms/cuda/src/CudaParallelKernels.cpp

@@ -536,6 +536,11 @@ double CudaParallelCalcCMAPTorsionForceKernel::execute(ContextImpl& context, boo
     return 0.0;
 }
 
+void CudaParallelCalcCMAPTorsionForceKernel::copyParametersToContext(ContextImpl& context, const CMAPTorsionForce& force) {
+    for (int i = 0; i < (int) kernels.size(); i++)
+        getKernel(i).copyParametersToContext(context, force);
+}
+
 class CudaParallelCalcCustomTorsionForceKernel::Task : public CudaContext::WorkTask {
 public:
     Task(ContextImpl& context, CudaCalcCustomTorsionForceKernel& kernel, bool includeForce,

platforms/opencl/include/OpenCLKernels.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) 2008-2013 Stanford University and the Authors.      *
+ * Portions copyright (c) 2008-2015 Stanford University and the Authors.      *
  * Authors: Peter Eastman                                                     *
  * Contributors:                                                              *
  *                                                                            *
@@ -498,6 +498,13 @@ public:
      * @return the potential energy due to the force
      */
     double execute(ContextImpl& context, bool includeForces, bool includeEnergy);
+    /**
+     * Copy changed parameters over to a context.
+     *
+     * @param context    the context to copy parameters to
+     * @param force      the CMAPTorsionForce to copy the parameters from
+     */
+    void copyParametersToContext(ContextImpl& context, const CMAPTorsionForce& force);
 private:
     int numTorsions;
     bool hasInitializedKernel;

platforms/opencl/include/OpenCLParallelKernels.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) 2011-2013 Stanford University and the Authors.      *
+ * Portions copyright (c) 2011-2015 Stanford University and the Authors.      *
  * Authors: Peter Eastman                                                     *
  * Contributors:                                                              *
  *                                                                            *
@@ -344,6 +344,13 @@ public:
      * @return the potential energy due to the force
      */
     double execute(ContextImpl& context, bool includeForces, bool includeEnergy);
+    /**
+     * Copy changed parameters over to a context.
+     *
+     * @param context    the context to copy parameters to
+     * @param force      the CMAPTorsionForce to copy the parameters from
+     */
+    void copyParametersToContext(ContextImpl& context, const CMAPTorsionForce& force);
 private:
     class Task;
     OpenCLPlatform::PlatformData& data;

platforms/opencl/src/OpenCLKernels.cpp

@@ -1164,6 +1164,9 @@ double OpenCLCalcCMAPTorsionForceKernel::execute(ContextImpl& context, bool incl
     return 0.0;
 }
 
+void OpenCLCalcCMAPTorsionForceKernel::copyParametersToContext(ContextImpl& context, const CMAPTorsionForce& force) {
+}
+
 class OpenCLCustomTorsionForceInfo : public OpenCLForceInfo {
 public:
     OpenCLCustomTorsionForceInfo(const CustomTorsionForce& force) : OpenCLForceInfo(0), force(force) {

platforms/opencl/src/OpenCLParallelKernels.cpp

@@ -492,6 +492,11 @@ double OpenCLParallelCalcCMAPTorsionForceKernel::execute(ContextImpl& context, b
     return 0.0;
 }
 
+void OpenCLParallelCalcCMAPTorsionForceKernel::copyParametersToContext(ContextImpl& context, const CMAPTorsionForce& force) {
+    for (int i = 0; i < (int) kernels.size(); i++)
+        getKernel(i).copyParametersToContext(context, force);
+}
+
 class OpenCLParallelCalcCustomTorsionForceKernel::Task : public OpenCLContext::WorkTask {
 public:
     Task(ContextImpl& context, OpenCLCalcCustomTorsionForceKernel& kernel, bool includeForce,

platforms/reference/include/ReferenceKernels.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) 2008-2013 Stanford University and the Authors.      *
+ * Portions copyright (c) 2008-2015 Stanford University and the Authors.      *
  * Authors: Peter Eastman                                                     *
  * Contributors:                                                              *
  *                                                                            *
@@ -493,6 +493,13 @@ public:
      * @return the potential energy due to the force
      */
     double execute(ContextImpl& context, bool includeForces, bool includeEnergy);
+    /**
+     * Copy changed parameters over to a context.
+     *
+     * @param context    the context to copy parameters to
+     * @param force      the CMAPTorsionForce to copy the parameters from
+     */
+    void copyParametersToContext(ContextImpl& context, const CMAPTorsionForce& force);
 private:
     std::vector<std::vector<std::vector<RealOpenMM> > > coeff;
     std::vector<int> torsionMaps;

platforms/reference/src/ReferenceKernels.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-2015 Stanford University and the Authors.      *
  * Authors: Peter Eastman                                                     *
  * Contributors:                                                              *
  *                                                                            *
@@ -710,6 +710,40 @@ double ReferenceCalcCMAPTorsionForceKernel::execute(ContextImpl& context, bool i
     return totalEnergy;
 }
 
+void ReferenceCalcCMAPTorsionForceKernel::copyParametersToContext(ContextImpl& context, const CMAPTorsionForce& force) {
+    int numMaps = force.getNumMaps();
+    int numTorsions = force.getNumTorsions();
+    if (coeff.size() != numMaps)
+        throw OpenMMException("updateParametersInContext: The number of maps has changed");
+    if (torsionMaps.size() != numTorsions)
+        throw OpenMMException("updateParametersInContext: The number of CMAP torsions has changed");
+
+    // Update the maps.
+
+    vector<double> energy;
+    vector<vector<double> > c;
+    for (int i = 0; i < numMaps; i++) {
+        int size;
+        force.getMapParameters(i, size, energy);
+        if (coeff[i].size() != size*size)
+            throw OpenMMException("updateParametersInContext: The size of a map has changed");
+        CMAPTorsionForceImpl::calcMapDerivatives(size, energy, c);
+        for (int j = 0; j < size*size; j++)
+            for (int k = 0; k < 16; k++)
+                coeff[i][j][k] = c[j][k];
+    }
+
+    // Update the indices.
+
+    for (int i = 0; i < numTorsions; i++) {
+        int index[8];
+        force.getTorsionParameters(i, torsionMaps[i], index[0], index[1], index[2], index[3], index[4], index[5], index[6], index[7]);
+        for (int j = 0; j < 8; j++)
+            if (index[j] != torsionIndices[i][j])
+                throw OpenMMException("updateParametersInContext: The set of particles in a CMAP torsion has changed");
+    }
+}
+
 ReferenceCalcCustomTorsionForceKernel::~ReferenceCalcCustomTorsionForceKernel() {
     disposeIntArray(torsionIndexArray, numTorsions);
     disposeRealArray(torsionParamArray, numTorsions);

platforms/reference/tests/TestReferenceCMAPTorsionForce.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) 2010 Stanford University and the Authors.           *
+ * Portions copyright (c) 2010-2015 Stanford University and the Authors.      *
  * Authors: Peter Eastman                                                     *
  * Contributors:                                                              *
  *                                                                            *
@@ -48,6 +48,8 @@
 using namespace OpenMM;
 using namespace std;
 
+ReferencePlatform platform;
+
 const double TOL = 1e-5;
 
 void testCMAPTorsions() {
@@ -56,7 +58,6 @@ void testCMAPTorsions() {
     // Create two systems: one with a pair of periodic torsions, and one with a CMAP torsion
     // that approximates the same force.
 
-    ReferencePlatform platform;
     System system1;
     for (int i = 0; i < 5; i++)
         system1.addParticle(1.0);
@@ -108,9 +109,66 @@ void testCMAPTorsions() {
     }
 }
 
+void testChangingParameters() {
+    // Create a system with two maps and one torsion.
+
+    const int mapSize = 8;
+    System system;
+    for (int i = 0; i < 5; i++)
+        system.addParticle(1.0);
+    CMAPTorsionForce* cmap = new CMAPTorsionForce();
+    vector<double> mapEnergy1(mapSize*mapSize);
+    vector<double> mapEnergy2(mapSize*mapSize);
+    for (int i = 0; i < mapSize; i++) {
+        double angle1 = i*2*M_PI/mapSize;
+        double energy1 = cos(angle1);
+        for (int j = 0; j < mapSize; j++) {
+            double angle2 = j*2*M_PI/mapSize;
+            double energy2 = 10*sin(angle2);
+            mapEnergy1[i+j*mapSize] = energy1+energy2;
+            mapEnergy2[i+j*mapSize] = energy1-energy2;
+        }
+    }
+    cmap->addMap(mapSize, mapEnergy1);
+    cmap->addMap(mapSize, mapEnergy2);
+    cmap->addTorsion(0, 0, 1, 2, 3, 1, 2, 3, 4);
+    system.addForce(cmap);
+
+    // Set particle positions so angle1=0 and angle2=PI/4.
+
+    vector<Vec3> positions(5);
+    positions[0] = Vec3(0, 0, 1);
+    positions[1] = Vec3(0, 0, 0);
+    positions[2] = Vec3(1, 0, 0);
+    positions[3] = Vec3(1, 0, 1);
+    positions[4] = Vec3(0.5, -0.5, 1);
+    VerletIntegrator integrator(0.01);
+    Context context(system, integrator, platform);
+    context.setPositions(positions);
+
+    // Check that the energy is correct.
+
+    double energy = context.getState(State::Energy).getPotentialEnergy();
+    ASSERT_EQUAL_TOL(1+10*sin(M_PI/4), energy, 1e-5);
+
+    // Modify the parameters.
+
+    cmap->setTorsionParameters(0, 1, 0, 1, 2, 3, 1, 2, 3, 4);
+    for (int i = 0; i < mapSize*mapSize; i++)
+        mapEnergy2[i] *= 2.0;
+    cmap->setMapParameters(1, mapSize, mapEnergy2);
+    cmap->updateParametersInContext(context);
+
+    // See if the results are correct.
+
+    energy = context.getState(State::Energy).getPotentialEnergy();
+    ASSERT_EQUAL_TOL(2-20*sin(M_PI/4), energy, 1e-5);
+}
+
 int main() {
     try {
         testCMAPTorsions();
+        testChangingParameters();
     }
     catch(const exception& e) {
         cout << "exception: " << e.what() << endl;