Minor cleanup, make chain propagation private

openmmapi/include/openmm/NoseHooverChain.h

@@ -253,19 +253,6 @@ private:
     std::vector<int> thermostatedAtoms, parentAtoms;
 };
 
-/**
- * Check if two Nose-Hoover chains are identical (checks all member variables but the chain id).
- */
-//inline bool operator==(const NoseHooverChain& lhs, const NoseHooverChain& rhs){
-//    if (lhs.getDefaultTemperature() != rhs.getDefaultTemperature()) return false;
-//    if (lhs.getDefaultCollisionFrequency() != rhs.getDefaultCollisionFrequency()) return false;
-//    if (lhs.getDefaultNumDegreesOfFreedom() != rhs.getDefaultNumDegreesOfFreedom()) return false;
-//    if (lhs.getDefaultNumMultiTimeSteps() != rhs.getDefaultNumMultiTimeSteps()) return false;
-//    if (lhs.getDefaultNumYoshidaSuzukiTimeSteps() != rhs.getDefaultNumYoshidaSuzukiTimeSteps()) return false;
-//    if (lhs.getThermostatedAtoms() != rhs.getThermostatedAtoms()) return false;
-//    if (lhs.getParentAtoms() != rhs.getParentAtoms()) return false;
-//}
-
 } // namespace OpenMM
 
 #endif /*OPENMM_NOSEHOOVERCHAIN_H_*/

openmmapi/include/openmm/VelocityVerletIntegrator.h

@@ -63,13 +63,6 @@ public:
      * @param steps   the number of time steps to take
      */
     void step(int steps);
-   /**
-     * Advance a simulation through time by taking a series of time steps.
-     * 
-     * @param kineticEnergy  the kinetic energy of the system that the chain is thermostating
-     * @param chainID        id of the Nose-Hoover-Chain
-     */
-    double propagateChain(double kineticEnergy, int chainID=0);
     /**
      * Add a Nose-Hoover Chain thermostat to control the temperature of the system
      *
@@ -157,6 +150,15 @@ public:
     }
 
 protected:
+   /**
+     * Advance any Nose-Hoover chains associated with this integrator and determine
+     * scale factor for the velocities.
+     * 
+     * @param kineticEnergy  the kinetic energy of the system that the chain is thermostating
+     * @param chainID        id of the Nose-Hoover-Chain
+     * @return the scale factor to be applied to the velocities of the particles thermostated by the chain.
+     */
+    double propagateChain(double kineticEnergy, int chainID=0);
     /**
      * This will be called by the Context when it is created.  It informs the Integrator
      * of what context it will be integrating, and gives it a chance to do any necessary initialization.

openmmapi/src/ContextImpl.cpp

@@ -52,11 +52,6 @@ using namespace OpenMM;
 using namespace std;
 const static char CHECKPOINT_MAGIC_BYTES[] = "OpenMM Binary Checkpoint\n";
 
-void printvars(std::map<std::string, double> & map) {
-    for(auto & v: map) {
-        std::cout << v.first << "  " << v.second << std::endl;
-    }
-}
 
 ContextImpl::ContextImpl(Context& owner, const System& system, Integrator& integrator, Platform* platform, const map<string, string>& properties, ContextImpl* originalContext) :
         owner(owner), system(system), integrator(integrator), hasInitializedForces(false), hasSetPositions(false), integratorIsDeleted(false),

platforms/cuda/src/kernels/noseHooverChain.cu

 
 #include <initializer_list>
 
-/**
- * Propagate the Nose-Hoover chain with one yoshida-suzuki term
- */
 extern "C" __global__ void propagateNoseHooverChain(mixed2* __restrict__ chainData, const mixed * __restrict__ energySum, mixed* __restrict__ scaleFactor,
                                                     mixed* __restrict__ chainMasses, mixed* __restrict__ chainForces, 
                                                     int chainLength, int numMTS, int numDOFs, float timeStep,

platforms/reference/src/SimTKReference/ReferenceVerletDynamics.cpp

@@ -24,7 +24,6 @@
 
 #include <cstring>
 #include <sstream>
-#include <iostream>
 
 #include "SimTKOpenMMUtilities.h"
 #include "ReferenceVerletDynamics.h"

platforms/reference/tests/TestReferenceNoseHooverThermostat.cpp

@@ -33,6 +33,4 @@
 #include "TestNoseHooverThermostat.h"
 
 void runPlatformTests() {
-    testNHCPropagation();
-    testPropagateChainConsistentWithPythonReference();
 }

plugins/drude/openmmapi/include/OpenMMDrude.h

@@ -35,5 +35,6 @@
 #include "openmm/DrudeForce.h"
 #include "openmm/DrudeLangevinIntegrator.h"
 #include "openmm/DrudeSCFIntegrator.h"
+#include "openmm/DrudeVelocityVerletIntegrator.h"
 
 #endif /*OPENMM_DRUDE_H_*/

tests/TestNoseHooverThermostat.h

@@ -210,60 +210,6 @@ void testDimerBox(bool constrain=true) {
     ASSERT_EQUAL_TOL(relative_std, 0, 1e-4);
 }
 
-void testNHCPropagation() {
-    // test if the velocity scale factor goes to one for a single particle
-    // with no forces in the system
-    for (int numMTS = 1; numMTS < 5; numMTS++){
-        for (int numYS = 1; numYS <=5; numYS+=2){
-            for (int chainLength = 2; chainLength < 6; chainLength += 2){
-            double temperature = 300; // kelvin
-            double collisionFrequency = 10; // 1/ps
-            VelocityVerletIntegrator integrator(0.001);
-            // make system
-            System system;
-            double mass = 1;
-            system.addParticle(mass);
-            int chainID = integrator.addNoseHooverChainThermostat(system, temperature, collisionFrequency,
-                                                                  chainLength, numMTS, numYS);
-            Context context(system, integrator, platform);
-            // propagate the chain
-            double velocity = 1;
-            double temp, scale, kineticEnergy;
-            double mean_temp=0, mean_scale=0;
-            for(int i = 0; i < 10000; ++i) {
-                kineticEnergy = 3 * 0.5 * mass * velocity * velocity;
-                scale = integrator.propagateChain(kineticEnergy, chainID);
-                velocity *= scale;
-                temp = 2* kineticEnergy / BOLTZ / 3;
-                mean_temp = (i*mean_temp + temp)/(i+1);
-                mean_scale = (i*mean_scale + scale)/(i+1);
-            }
-            //std::cout << mean_scale <<  " " << mean_temp << std::endl;
-            ASSERT_EQUAL_TOL(1, mean_scale,  1e-2);
-            ASSERT_EQUAL_TOL(temperature, mean_temp,  0.25);
-            }
-        }
-    }
-}
-
-void testPropagateChainConsistentWithPythonReference() {
-    VelocityVerletIntegrator integrator(0.001);
-    System system;
-    double mass = 1;
-    system.addParticle(mass);
-    double kineticEnergy = 1e6;
-    double temperature=300, collisionFrequency=1, chainLength=3, numMTS=3, numYS=3;
-    int chainID = integrator.addNoseHooverChainThermostat(system, temperature, collisionFrequency,
-                                                                  chainLength, numMTS, numYS);
-    Context context(system, integrator, platform);
-    double scale = integrator.propagateChain(kineticEnergy, chainID);
-#if DEBUG
-    std::cout << std::setw(12) << std::setprecision(10) << scale << std::endl;
-#endif
-    ASSERT_EQUAL_TOL(0.9674732261005896, scale, 1e-5)
-}
-
-
 void testCheckpoints() {
     double timeStep = 0.001;
     VelocityVerletIntegrator integrator(timeStep), newIntegrator(timeStep);