Trying to use an Integrator for multiple Contexts produces an exception

openmmapi/include/openmm/BrownianIntegrator.h

@@ -124,6 +124,7 @@ private:
     double temperature, friction;
     int randomNumberSeed;
     ContextImpl* context;
+    Context* owner;
     Kernel kernel;
 };
 

openmmapi/include/openmm/LangevinIntegrator.h

@@ -124,6 +124,7 @@ private:
     double temperature, friction;
     int randomNumberSeed;
     ContextImpl* context;
+    Context* owner;
     Kernel kernel;
 };
 

openmmapi/include/openmm/VariableLangevinIntegrator.h

@@ -155,6 +155,7 @@ private:
     double temperature, friction, errorTol;
     int randomNumberSeed;
     ContextImpl* context;
+    Context* owner;
     Kernel kernel;
 };
 

openmmapi/include/openmm/VariableVerletIntegrator.h

@@ -107,6 +107,7 @@ protected:
 private:
     double errorTol;
     ContextImpl* context;
+    Context* owner;
     Kernel kernel;
 };
 

openmmapi/include/openmm/VerletIntegrator.h

@@ -69,6 +69,7 @@ protected:
     std::vector<std::string> getKernelNames();
 private:
     ContextImpl* context;
+    Context* owner;
     Kernel kernel;
 };
 

openmmapi/src/BrownianIntegrator.cpp

@@ -31,6 +31,7 @@
 
 #include "openmm/BrownianIntegrator.h"
 #include "openmm/Context.h"
+#include "openmm/OpenMMException.h"
 #include "openmm/internal/ContextImpl.h"
 #include "openmm/kernels.h"
 #include <ctime>
@@ -40,7 +41,7 @@ using namespace OpenMM;
 using std::string;
 using std::vector;
 
-BrownianIntegrator::BrownianIntegrator(double temperature, double frictionCoeff, double stepSize) {
+BrownianIntegrator::BrownianIntegrator(double temperature, double frictionCoeff, double stepSize) : owner(NULL) {
     setTemperature(temperature);
     setFriction(frictionCoeff);
     setStepSize(stepSize);
@@ -49,7 +50,10 @@ BrownianIntegrator::BrownianIntegrator(double temperature, double frictionCoeff,
 }
 
 void BrownianIntegrator::initialize(ContextImpl& contextRef) {
+    if (owner != NULL && &contextRef.getOwner() != owner)
+        throw OpenMMException("This Integrator is already bound to a context");
     context = &contextRef;
+    owner = &contextRef.getOwner();
     kernel = context->getPlatform().createKernel(IntegrateBrownianStepKernel::Name(), contextRef);
     dynamic_cast<IntegrateBrownianStepKernel&>(kernel.getImpl()).initialize(contextRef.getSystem(), *this);
 }

openmmapi/src/LangevinIntegrator.cpp

@@ -31,6 +31,7 @@
 
 #include "openmm/LangevinIntegrator.h"
 #include "openmm/Context.h"
+#include "openmm/OpenMMException.h"
 #include "openmm/internal/ContextImpl.h"
 #include "openmm/kernels.h"
 #include <ctime>
@@ -40,7 +41,7 @@ using namespace OpenMM;
 using std::string;
 using std::vector;
 
-LangevinIntegrator::LangevinIntegrator(double temperature, double frictionCoeff, double stepSize) {
+LangevinIntegrator::LangevinIntegrator(double temperature, double frictionCoeff, double stepSize) : owner(NULL) {
     setTemperature(temperature);
     setFriction(frictionCoeff);
     setStepSize(stepSize);
@@ -49,7 +50,10 @@ LangevinIntegrator::LangevinIntegrator(double temperature, double frictionCoeff,
 }
 
 void LangevinIntegrator::initialize(ContextImpl& contextRef) {
+    if (owner != NULL && &contextRef.getOwner() != owner)
+        throw OpenMMException("This Integrator is already bound to a context");
     context = &contextRef;
+    owner = &contextRef.getOwner();
     kernel = context->getPlatform().createKernel(IntegrateLangevinStepKernel::Name(), contextRef);
     dynamic_cast<IntegrateLangevinStepKernel&>(kernel.getImpl()).initialize(contextRef.getSystem(), *this);
 }

openmmapi/src/VariableLangevinIntegrator.cpp

@@ -31,6 +31,7 @@
 
 #include "openmm/VariableLangevinIntegrator.h"
 #include "openmm/Context.h"
+#include "openmm/OpenMMException.h"
 #include "openmm/internal/ContextImpl.h"
 #include "openmm/kernels.h"
 #include <limits>
@@ -41,7 +42,7 @@ using namespace OpenMM;
 using std::string;
 using std::vector;
 
-VariableLangevinIntegrator::VariableLangevinIntegrator(double temperature, double frictionCoeff, double errorTol) {
+VariableLangevinIntegrator::VariableLangevinIntegrator(double temperature, double frictionCoeff, double errorTol) : owner(NULL) {
     setTemperature(temperature);
     setFriction(frictionCoeff);
     setErrorTolerance(errorTol);
@@ -50,7 +51,10 @@ VariableLangevinIntegrator::VariableLangevinIntegrator(double temperature, doubl
 }
 
 void VariableLangevinIntegrator::initialize(ContextImpl& contextRef) {
+    if (owner != NULL && &contextRef.getOwner() != owner)
+        throw OpenMMException("This Integrator is already bound to a context");
     context = &contextRef;
+    owner = &contextRef.getOwner();
     kernel = context->getPlatform().createKernel(IntegrateVariableLangevinStepKernel::Name(), contextRef);
     dynamic_cast<IntegrateVariableLangevinStepKernel&>(kernel.getImpl()).initialize(contextRef.getSystem(), *this);
 }

openmmapi/src/VariableVerletIntegrator.cpp

@@ -31,6 +31,7 @@
 
 #include "openmm/VariableVerletIntegrator.h"
 #include "openmm/Context.h"
+#include "openmm/OpenMMException.h"
 #include "openmm/internal/ContextImpl.h"
 #include "openmm/kernels.h"
 #include <limits>
@@ -40,12 +41,15 @@ using namespace OpenMM;
 using std::string;
 using std::vector;
 
-VariableVerletIntegrator::VariableVerletIntegrator(double errorTol) : errorTol(errorTol) {
+VariableVerletIntegrator::VariableVerletIntegrator(double errorTol) : errorTol(errorTol), owner(NULL) {
     setConstraintTolerance(1e-4);
 }
 
 void VariableVerletIntegrator::initialize(ContextImpl& contextRef) {
+    if (owner != NULL && &contextRef.getOwner() != owner)
+        throw OpenMMException("This Integrator is already bound to a context");
     context = &contextRef;
+    owner = &contextRef.getOwner();
     kernel = context->getPlatform().createKernel(IntegrateVariableVerletStepKernel::Name(), contextRef);
     dynamic_cast<IntegrateVariableVerletStepKernel&>(kernel.getImpl()).initialize(contextRef.getSystem(), *this);
 }

openmmapi/src/VerletIntegrator.cpp

@@ -31,6 +31,7 @@
 
 #include "openmm/VerletIntegrator.h"
 #include "openmm/Context.h"
+#include "openmm/OpenMMException.h"
 #include "openmm/internal/ContextImpl.h"
 #include "openmm/kernels.h"
 #include <string>
@@ -39,13 +40,16 @@ using namespace OpenMM;
 using std::string;
 using std::vector;
 
-VerletIntegrator::VerletIntegrator(double stepSize) {
+VerletIntegrator::VerletIntegrator(double stepSize) : owner(NULL) {
     setStepSize(stepSize);
     setConstraintTolerance(1e-4);
 }
 
 void VerletIntegrator::initialize(ContextImpl& contextRef) {
+    if (owner != NULL && &contextRef.getOwner() != owner)
+        throw OpenMMException("This Integrator is already bound to a context");
     context = &contextRef;
+    owner = &contextRef.getOwner();
     kernel = context->getPlatform().createKernel(IntegrateVerletStepKernel::Name(), contextRef);
     dynamic_cast<IntegrateVerletStepKernel&>(kernel.getImpl()).initialize(contextRef.getSystem(), *this);
 }

platforms/cuda/tests/TestCudaCustomAngleForce.cpp

@@ -79,7 +79,6 @@ void testAngles() {
     harmonicSystem.addParticle(1.0);
     harmonicSystem.addParticle(1.0);
     harmonicSystem.addParticle(1.0);
-    VerletIntegrator integrator(0.01);
     HarmonicAngleForce* harmonic = new HarmonicAngleForce();
     harmonic->addAngle(0, 1, 2, 1.5, 0.8);
     harmonic->addAngle(1, 2, 3, 2.0, 0.5);
@@ -91,22 +90,22 @@ void testAngles() {
     init_gen_rand(0, sfmt);
 
     vector<Vec3> positions(4);
-    VerletIntegrator integrator1(0.01);
-    VerletIntegrator integrator2(0.01);
     for (int i = 0; i < 10; i++) {
         for (int j = 0; j < (int) positions.size(); j++)
             positions[j] = Vec3(5.0*genrand_real2(sfmt), 5.0*genrand_real2(sfmt), 5.0*genrand_real2(sfmt));
         double energy1, energy2;
         vector<Vec3> forces1, forces2;
         {
-            Context c(customSystem, integrator1, platform);
+            VerletIntegrator integrator(0.01);
+            Context c(customSystem, integrator, platform);
             c.setPositions(positions);
             State s = c.getState(State::Forces | State::Energy);
             energy1 = s.getPotentialEnergy();
             forces1 = s.getForces();
         }
         {
-            Context c(harmonicSystem, integrator1, platform);
+            VerletIntegrator integrator(0.01);
+            Context c(harmonicSystem, integrator, platform);
             c.setPositions(positions);
             State s = c.getState(State::Forces | State::Energy);
             energy2 = s.getPotentialEnergy();

platforms/cuda/tests/TestCudaCustomTorsionForce.cpp

@@ -96,22 +96,22 @@ void testTorsions() {
     init_gen_rand(0, sfmt);
 
     vector<Vec3> positions(5);
-    VerletIntegrator integrator1(0.01);
-    VerletIntegrator integrator2(0.01);
     for (int i = 0; i < 10; i++) {
         for (int j = 0; j < (int) positions.size(); j++)
             positions[j] = Vec3(5.0*genrand_real2(sfmt), 5.0*genrand_real2(sfmt), 5.0*genrand_real2(sfmt));
         double energy1, energy2;
         vector<Vec3> forces1, forces2;
         {
-            Context c(customSystem, integrator1, platform);
+            VerletIntegrator integrator(0.01);
+            Context c(customSystem, integrator, platform);
             c.setPositions(positions);
             State s = c.getState(State::Forces | State::Energy);
             energy1 = s.getPotentialEnergy();
             forces1 = s.getForces();
         }
         {
-            Context c(harmonicSystem, integrator1, platform);
+            VerletIntegrator integrator(0.01);
+            Context c(harmonicSystem, integrator, platform);
             c.setPositions(positions);
             State s = c.getState(State::Forces | State::Energy);
             energy2 = s.getPotentialEnergy();

platforms/cuda/tests/TestCudaEwald.cpp

@@ -65,7 +65,6 @@ void testEwaldPME(bool includeExceptions) {
     CudaPlatform cuda;
     ReferencePlatform reference;
     System system;
-    VerletIntegrator integrator(0.01);
     NonbondedForce* nonbonded = new NonbondedForce();
     nonbonded->setNonbondedMethod(NonbondedForce::Ewald);
     nonbonded->setCutoffDistance(cutoff);
@@ -96,8 +95,10 @@ void testEwaldPME(bool includeExceptions) {
 
 //    (1)  Check whether the Reference and Cuda platforms agree when using Ewald Method 
  
-    Context cudaContext(system, integrator, cuda);
-    Context referenceContext(system, integrator, reference);
+    VerletIntegrator integrator1(0.01);
+    VerletIntegrator integrator2(0.01);
+    Context cudaContext(system, integrator1, cuda);
+    Context referenceContext(system, integrator2, reference);
     cudaContext.setPositions(positions);
     referenceContext.setPositions(positions);
     State cudaState = cudaContext.getState(State::Forces | State::Energy);

platforms/cuda/tests/TestCudaGBSAOBCForce.cpp

@@ -136,7 +136,6 @@ void testForce(int numParticles, NonbondedForce::NonbondedMethod method, GBSAOBC
     CudaPlatform cuda;
     ReferencePlatform reference;
     System system;
-    LangevinIntegrator integrator(0, 0.1, 0.01);
     GBSAOBCForce* gbsa = new GBSAOBCForce();
     NonbondedForce* nonbonded = new NonbondedForce();
     for (int i = 0; i < numParticles; ++i) {
@@ -156,8 +155,10 @@ void testForce(int numParticles, NonbondedForce::NonbondedMethod method, GBSAOBC
     }
     system.addForce(gbsa);
     system.addForce(nonbonded);
-    Context context(system, integrator, cuda);
-    Context refContext(system, integrator, reference);
+    LangevinIntegrator integrator1(0, 0.1, 0.01);
+    LangevinIntegrator integrator2(0, 0.1, 0.01);
+    Context context(system, integrator1, cuda);
+    Context refContext(system, integrator2, reference);
     
     // Set random (but uniformly distributed) positions for all the particles.
     

platforms/cuda/tests/TestCudaNonbondedForce.cpp

@@ -105,7 +105,6 @@ void testLJ() {
 void testExclusionsAnd14() {
     CudaPlatform platform;
     System system;
-    LangevinIntegrator integrator(0.0, 0.1, 0.01);
     NonbondedForce* nonbonded = new NonbondedForce();
     for (int i = 0; i < 5; ++i) {
         system.addParticle(1.0);
@@ -146,6 +145,7 @@ void testExclusionsAnd14() {
         // The following is in its own block, because CUDA can't deal with multiple Contexts
         // existing on the same thread at the same time.
         {
+            LangevinIntegrator integrator(0.0, 0.1, 0.01);
             Context context(system, integrator, platform);
             context.setPositions(positions);
             State state = context.getState(State::Forces | State::Energy);
@@ -174,6 +174,7 @@ void testExclusionsAnd14() {
             nonbonded->setParticleParameters(i, 2, 1.5, 0);
             nonbonded->setExceptionParameters(first14, 0, 3, i == 3 ? 4/1.2 : 0, 1.5, 0);
             nonbonded->setExceptionParameters(second14, 1, 4, 0, 1.5, 0);
+            LangevinIntegrator integrator(0.0, 0.1, 0.01);
             Context context(system, integrator, platform);
             context.setPositions(positions);
             State state = context.getState(State::Forces | State::Energy);
@@ -375,7 +376,6 @@ void testLargeSystem() {
     System system;
     for (int i = 0; i < numParticles; i++)
         system.addParticle(1.0);
-    VerletIntegrator integrator(0.01);
     NonbondedForce* nonbonded = new NonbondedForce();
     HarmonicBondForce* bonds = new HarmonicBondForce();
     vector<Vec3> positions(numParticles);
@@ -407,8 +407,10 @@ void testLargeSystem() {
     nonbonded->setCutoffDistance(cutoff);
     system.addForce(nonbonded);
     system.addForce(bonds);
-    Context cudaContext(system, integrator, cuda);
-    Context referenceContext(system, integrator, reference);
+    VerletIntegrator integrator1(0.01);
+    VerletIntegrator integrator2(0.01);
+    Context cudaContext(system, integrator1, cuda);
+    Context referenceContext(system, integrator2, reference);
     cudaContext.setPositions(positions);
     cudaContext.setVelocities(velocities);
     referenceContext.setPositions(positions);

platforms/opencl/tests/TestOpenCLCustomAngleForce.cpp

@@ -79,7 +79,6 @@ void testAngles() {
     harmonicSystem.addParticle(1.0);
     harmonicSystem.addParticle(1.0);
     harmonicSystem.addParticle(1.0);
-    VerletIntegrator integrator(0.01);
     HarmonicAngleForce* harmonic = new HarmonicAngleForce();
     harmonic->addAngle(0, 1, 2, 1.5, 0.8);
     harmonic->addAngle(1, 2, 3, 2.0, 0.5);
@@ -91,9 +90,9 @@ void testAngles() {
     init_gen_rand(0, sfmt);
 
     vector<Vec3> positions(4);
+    for (int i = 0; i < 10; i++) {
         VerletIntegrator integrator1(0.01);
         VerletIntegrator integrator2(0.01);
-    for (int i = 0; i < 10; i++) {
         Context c1(customSystem, integrator1, platform);
         Context c2(harmonicSystem, integrator2, platform);
         for (int j = 0; j < (int) positions.size(); j++)

platforms/opencl/tests/TestOpenCLCustomHbondForce.cpp

@@ -109,9 +109,9 @@ void testHbond() {
     init_gen_rand(0, sfmt);
 
     vector<Vec3> positions(5);
+    for (int i = 0; i < 10; i++) {
         VerletIntegrator integrator1(0.01);
         VerletIntegrator integrator2(0.01);
-    for (int i = 0; i < 10; i++) {
         Context c1(customSystem, integrator1, platform);
         Context c2(standardSystem, integrator2, platform);
         for (int j = 0; j < (int) positions.size(); j++)

platforms/opencl/tests/TestOpenCLEwald.cpp

@@ -64,7 +64,6 @@ void testEwaldPME(bool includeExceptions) {
     OpenCLPlatform cl;
     ReferencePlatform reference;
     System system;
-    VerletIntegrator integrator(0.01);
     NonbondedForce* nonbonded = new NonbondedForce();
     nonbonded->setNonbondedMethod(NonbondedForce::Ewald);
     nonbonded->setCutoffDistance(cutoff);
@@ -95,8 +94,10 @@ void testEwaldPME(bool includeExceptions) {
 
 //    (1)  Check whether the Reference and OpenCL platforms agree when using Ewald Method
 
-    Context clContext(system, integrator, cl);
-    Context referenceContext(system, integrator, reference);
+    VerletIntegrator integrator1(0.01);
+    VerletIntegrator integrator2(0.01);
+    Context clContext(system, integrator1, cl);
+    Context referenceContext(system, integrator2, reference);
     clContext.setPositions(positions);
     referenceContext.setPositions(positions);
     State clState = clContext.getState(State::Forces | State::Energy);
@@ -124,7 +125,8 @@ void testEwaldPME(bool includeExceptions) {
         Vec3 f = clState.getForces()[i];
         positions[i] = Vec3(p[0]-f[0]*step, p[1]-f[1]*step, p[2]-f[2]*step);
     }
-    Context clContext2(system, integrator, cl);
+    VerletIntegrator integrator3(0.01);
+    Context clContext2(system, integrator3, cl);
     clContext2.setPositions(positions);
 
     tol = 1e-2;
@@ -162,7 +164,8 @@ void testEwaldPME(bool includeExceptions) {
         Vec3 f = clState.getForces()[i];
         positions[i] = Vec3(p[0]-f[0]*step, p[1]-f[1]*step, p[2]-f[2]*step);
     }
-    Context clContext3(system, integrator, cl);
+    VerletIntegrator integrator4(0.01);
+    Context clContext3(system, integrator4, cl);
     clContext3.setPositions(positions);
 
     tol = 1e-2;

platforms/opencl/tests/TestOpenCLGBSAOBCForce.cpp

@@ -136,7 +136,6 @@ void testForce(int numParticles, NonbondedForce::NonbondedMethod method, GBSAOBC
     OpenCLPlatform cl;
     ReferencePlatform reference;
     System system;
-    LangevinIntegrator integrator(0, 0.1, 0.01);
     GBSAOBCForce* gbsa = new GBSAOBCForce();
     NonbondedForce* nonbonded = new NonbondedForce();
     for (int i = 0; i < numParticles; ++i) {
@@ -156,8 +155,10 @@ void testForce(int numParticles, NonbondedForce::NonbondedMethod method, GBSAOBC
     }
     system.addForce(gbsa);
     system.addForce(nonbonded);
-    Context context(system, integrator, cl);
-    Context refContext(system, integrator, reference);
+    LangevinIntegrator integrator1(0, 0.1, 0.01);
+    LangevinIntegrator integrator2(0, 0.1, 0.01);
+    Context context(system, integrator1, cl);
+    Context refContext(system, integrator2, reference);
 
     // Set random (but uniformly distributed) positions for all the particles.
 

platforms/opencl/tests/TestOpenCLNonbondedForce.cpp

@@ -106,7 +106,6 @@ void testLJ() {
 void testExclusionsAnd14() {
     OpenCLPlatform platform;
     System system;
-    LangevinIntegrator integrator(0.0, 0.1, 0.01);
     NonbondedForce* nonbonded = new NonbondedForce();
     for (int i = 0; i < 5; ++i) {
         system.addParticle(1.0);
@@ -144,6 +143,7 @@ void testExclusionsAnd14() {
         nonbonded->setExceptionParameters(first14, 0, 3, 0, 1.5, i == 3 ? 0.5 : 0.0);
         nonbonded->setExceptionParameters(second14, 1, 4, 0, 1.5, 0.0);
         positions[i] = Vec3(r, 0, 0);
+        LangevinIntegrator integrator(0.0, 0.1, 0.01);
         Context context(system, integrator, platform);
         context.setPositions(positions);
         State state = context.getState(State::Forces | State::Energy);
@@ -170,7 +170,8 @@ void testExclusionsAnd14() {
         nonbonded->setParticleParameters(i, 2, 1.5, 0);
         nonbonded->setExceptionParameters(first14, 0, 3, i == 3 ? 4/1.2 : 0, 1.5, 0);
         nonbonded->setExceptionParameters(second14, 1, 4, 0, 1.5, 0);
-        Context context2(system, integrator, platform);
+        LangevinIntegrator integrator2(0.0, 0.1, 0.01);
+        Context context2(system, integrator2, platform);
         context2.setPositions(positions);
         state = context2.getState(State::Forces | State::Energy);
         const vector<Vec3>& forces2 = state.getForces();
@@ -370,7 +371,6 @@ void testLargeSystem() {
     System system;
     for (int i = 0; i < numParticles; i++)
         system.addParticle(1.0);
-    VerletIntegrator integrator(0.01);
     NonbondedForce* nonbonded = new NonbondedForce();
     HarmonicBondForce* bonds = new HarmonicBondForce();
     vector<Vec3> positions(numParticles);
@@ -402,8 +402,10 @@ void testLargeSystem() {
     nonbonded->setCutoffDistance(cutoff);
     system.addForce(nonbonded);
     system.addForce(bonds);
-    Context clContext(system, integrator, cl);
-    Context referenceContext(system, integrator, reference);
+    VerletIntegrator integrator1(0.01);
+    VerletIntegrator integrator2(0.01);
+    Context clContext(system, integrator1, cl);
+    Context referenceContext(system, integrator2, reference);
     clContext.setPositions(positions);
     clContext.setVelocities(velocities);
     referenceContext.setPositions(positions);