Hack particular gccxml location for one particular nightly build machine.

Incorporate Normal Mode Langevin openmm code from Chris Sweet as a plugin.  (compiles on windows, but no actual testing done)

CMakeLists.txt

@@ -268,6 +268,10 @@ set(cmv "${CMAKE_MAJOR_VERSION}.${CMAKE_MINOR_VERSION}")
 # Java and gccxml are required for API wrappers
 FIND_PACKAGE(Java)
 MARK_AS_ADVANCED(CLEAR JAVA_RUNTIME)
+# One particular location to try first for nightly builds
+FIND_PROGRAM(GCCXML_PATH gccxml PATH
+    "C:/Program Files/gccxml_sherm/bin")
+# In case first attempt fails
 FIND_PROGRAM(GCCXML_PATH gccxml PATH
     /usr/local/bin
     "C:/Program Files/gccxml 0.9/bin"
@@ -356,6 +360,12 @@ IF(OPENMM_BUILD_FREE_ENERGY_PLUGIN)
    ADD_SUBDIRECTORY(plugins/freeEnergy)
 ENDIF(OPENMM_BUILD_FREE_ENERGY_PLUGIN)
 
+# Normal mode langevin plugin
+set(OPENMM_BUILD_NML_PLUGIN ON CACHE BOOL "Build normal-mode langevin plugin")
+if(OPENMM_BUILD_NML_PLUGIN)
+    add_subdirectory(plugins/normalModeLangevin)
+endif(OPENMM_BUILD_NML_PLUGIN)
+
 INSTALL_TARGETS(/lib RUNTIME_DIRECTORY /lib ${SHARED_TARGET})
 INSTALL_TARGETS(/lib RUNTIME_DIRECTORY /lib ${STATIC_TARGET})
 FILE(GLOB CORE_HEADERS     include/*.h          */include/*.h)

platforms/reference/src/SimTKReference/ReferenceCCMAAlgorithm.h

@@ -32,7 +32,7 @@
 
 // ---------------------------------------------------------------------------------------
 
-class ReferenceCCMAAlgorithm : public ReferenceConstraintAlgorithm {
+class OPENMM_EXPORT ReferenceCCMAAlgorithm : public ReferenceConstraintAlgorithm {
 
    protected:
 

platforms/reference/src/SimTKReference/ReferenceDynamics.h

@@ -41,7 +41,7 @@ Main method (virtual) is update()
 
 --------------------------------------------------------------------------------------- */
 
-class ReferenceDynamics {
+class OPENMM_EXPORT ReferenceDynamics {
 
    public:
 

platforms/reference/src/SimTKUtilities/SimTKOpenMMLog.h

@@ -28,6 +28,7 @@
 #include <stdio.h>
 #include <sstream>
 #include "SimTKOpenMMCommon.h"
+#include "openmm/internal/windowsExport.h"
 
 /** ---------------------------------------------------------------------------------------
 
@@ -35,7 +36,7 @@
 
    --------------------------------------------------------------------------------------- */
 
-class SimTKOpenMMLog {
+class OPENMM_EXPORT SimTKOpenMMLog {
 
    public:
 

plugins/normalModeLangevin/CMakeLists.txt

+#/* -------------------------------------------------------------------------- *
+# *                                   OpenMM                                   *
+# * -------------------------------------------------------------------------- *
+# * This is part of the OpenMM molecular simulation toolkit originating from   *
+# * Simbios, the NIH National Center for Physics-Based Simulation of           *
+# * Biological Structures at Stanford, funded under the NIH Roadmap for        *
+# * Medical Research, grant U54 GM072970. See https://simtk.org.               *
+# *                                                                            *
+# * Portions copyright (c) 2009 Stanford University and the Authors.           *
+# * Authors: Chris Sweet, Christopher Bruns                                    *
+# * Contributors:                                                              *
+# *                                                                            *
+# * Permission is hereby granted, free of charge, to any person obtaining a    *
+# * copy of this software and associated documentation files (the "Software"), *
+# * to deal in the Software without restriction, including without limitation  *
+# * the rights to use, copy, modify, merge, publish, distribute, sublicense,   *
+# * and/or sell copies of the Software, and to permit persons to whom the      *
+# * Software is furnished to do so, subject to the following conditions:       *
+# *                                                                            *
+# * The above copyright notice and this permission notice shall be included in *
+# * all copies or substantial portions of the Software.                        *
+# *                                                                            *
+# * THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR *
+# * IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,   *
+# * FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL    *
+# * THE AUTHORS, CONTRIBUTORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM,    *
+# * DAMAGES OR OTHER LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR      *
+# * OTHERWISE, ARISING FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE  *
+# * USE OR OTHER DEALINGS IN THE SOFTWARE.                                     *
+# * -------------------------------------------------------------------------- */
+
+# Create Normal Mode Langevin (NML) plugin for OpenMM
+
+include_directories(${CMAKE_CURRENT_SOURCE_DIR}/include)
+file(GLOB nml_headers "include/*.h")
+file(GLOB nml_sources "src/*.cpp")
+
+# Reference platform implementation
+include_directories(${CMAKE_CURRENT_SOURCE_DIR}/platforms/reference/include)
+include_directories(${CMAKE_SOURCE_DIR}/platforms/reference/src)
+file(GLOB headers "${CMAKE_CURRENT_SOURCE_DIR}/platforms/reference/include/*.h")
+file(GLOB sources "${CMAKE_CURRENT_SOURCE_DIR}/platforms/reference/src/*.cpp")
+set(nml_headers ${nml_headers} ${headers})
+set(nml_sources ${nml_sources} ${sources})
+
+add_library(NormalModeLangevin SHARED ${nml_sources} ${nml_headers})
+set_target_properties(NormalModeLangevin PROPERTIES COMPILE_FLAGS "-DOPENMM_BUILDING_SHARED_LIBRARY")
+# On Unix or Cygwin we have to add the debug suffix manually
+if(UNIX)
+    set_target_properties(NormalModeLangevin PROPERTIES DEBUG_OUTPUT_NAME NormalModeLangevin_d)
+endif(UNIX)
+target_link_libraries( NormalModeLangevin OpenMM )
+
+if(BUILD_TESTING)
+    add_subdirectory(test)
+endif(BUILD_TESTING)
+
+install(TARGETS NormalModeLangevin RUNTIME DESTINATION lib/plugins)
+

plugins/normalModeLangevin/include/EigenvectorInfo.h

+/*  -*- c++ -*-  */
+#ifndef EIGENVECTORINFO_H
+#define EIGENVECTORINFO_H
+
+/* -------------------------------------------------------------------------- *
+ *                                   OpenMM                                   *
+ * -------------------------------------------------------------------------- *
+ * This is part of the OpenMM molecular simulation toolkit originating from   *
+ * Simbios, the NIH National Center for Physics-Based Simulation of           *
+ * Biological Structures at Stanford, funded under the NIH Roadmap for        *
+ * Medical Research, grant U54 GM072970. See https://simtk.org.               *
+ *                                                                            *
+ * Portions copyright (c) 2009 Stanford University and the Authors.           *
+ * Authors: Chris Sweet                                                       *
+ * Contributors: Christopher Bruns                                            *
+ *                                                                            *
+ * Permission is hereby granted, free of charge, to any person obtaining a    *
+ * copy of this software and associated documentation files (the "Software"), *
+ * to deal in the Software without restriction, including without limitation  *
+ * the rights to use, copy, modify, merge, publish, distribute, sublicense,   *
+ * and/or sell copies of the Software, and to permit persons to whom the      *
+ * Software is furnished to do so, subject to the following conditions:       *
+ *                                                                            *
+ * The above copyright notice and this permission notice shall be included in *
+ * all copies or substantial portions of the Software.                        *
+ *                                                                            *
+ * THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR *
+ * IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,   *
+ * FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL    *
+ * THE AUTHORS, CONTRIBUTORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM,    *
+ * DAMAGES OR OTHER LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR      *
+ * OTHERWISE, ARISING FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE  *
+ * USE OR OTHER DEALINGS IN THE SOFTWARE.                                     *
+ * -------------------------------------------------------------------------- */
+
+#include <string>
+#include <iostream>
+#include <vector>
+
+using std::cout;
+using std::endl;
+namespace ProtoMol
+{
+  //_________________________________________________________________XYZ
+  /**
+   * Container holding coordinates/Vector3D and names
+   */
+  struct EigenvectorInfo {
+    //~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
+    // Constructors, destructors, assignment
+    //~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
+    EigenvectorInfo() : myEigenvectorLength ( 0 ), myNumEigenvectors ( 0 ), myNumUsedEigenvectors( 0 ),
+        myEigenvectors ( 0 ), myOrigCEigval( 0.0 ), myNewCEigval( 0.0 ), myOrigTimestep( 0.0 ), 
+        reDiagonalize( false ),
+        mySingleEigs( 0 ), myEigVecChanged( true ), myMinimumLimit( 0.5 ), currentMode( -1 ) {};
+
+    EigenvectorInfo( unsigned int n, unsigned int m ) : myEigenvectorLength( n ), myNumEigenvectors( m ), myNumUsedEigenvectors( m ),
+        myMaxEigenvalue( 0.0 ), myEigenvectors ( new double[n * m * 3] ), 
+        myOrigCEigval( 0.0 ), myNewCEigval( 0.0 ), myOrigTimestep( 0.0 ), reDiagonalize( false ),
+        mySingleEigs( 0 ),
+        myEigVecChanged( true ), myMinimumLimit( 0.5 ), currentMode( -1 ) {}
+
+    ~EigenvectorInfo() {
+      if ( myEigenvectors ) {
+        delete [] myEigenvectors;
+        myEigenvectors = 0;
+      }
+      if ( mySingleEigs ) {
+        delete [] mySingleEigs;
+        mySingleEigs = 0;
+      }
+    };
+
+    //~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
+    // New methods of class EigenvectorInfo
+    //~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
+    bool initializeEigenvectors() {
+      try {
+
+        myEigenvectors = new double[myEigenvectorLength * myNumEigenvectors * 3];
+
+      } catch ( std::bad_alloc& ) {
+
+        return false;
+
+      }
+
+      return true;
+
+    }
+
+    float* getFloatEigPointer() {
+
+      //no eigenvectors assigned?
+      if(!myEigenvectors) return (float*)0;
+
+      const unsigned int arrayLen = myEigenvectorLength * myNumEigenvectors * 3;
+
+      //assign storage if required
+      if(!mySingleEigs) {  
+
+        try {
+          mySingleEigs = new float[arrayLen];
+        } catch ( std::bad_alloc& ) {
+          return (float*)0;
+        }
+
+      }
+
+      //update values if double array updated
+      if(myEigVecChanged) {
+
+        for( unsigned int i=0; i<arrayLen; i++)
+          mySingleEigs[i] = (float)myEigenvectors[i];
+
+        myEigVecChanged = false;
+
+      }
+
+      return mySingleEigs;
+
+    }
+
+    //~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
+    // My data members
+    //~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
+    // Eigenvector information
+    unsigned int myEigenvectorLength;
+    unsigned int myNumEigenvectors;
+    unsigned int myNumUsedEigenvectors;
+    double myMaxEigenvalue;
+    double *myEigenvectors;
+    
+    //Current and original max subspace eigenvalue, 
+    //for adaptive timestep
+    double myOrigCEigval, myNewCEigval;
+    double myOrigTimestep;
+
+    //re-diagonalization flag
+    bool reDiagonalize;
+    
+    //OpenMM single precision interface
+    float *mySingleEigs;
+    bool myEigVecChanged;
+
+    double myMinimumLimit;
+
+    //Analytic integrator
+    std::vector< double > myEigenvalues;
+    int currentMode;
+
+
+  };
+}
+#endif /* EIGENVECTORINFO_H */

plugins/normalModeLangevin/include/IntegrateNMLStepKernel.h

+#ifndef INTEGRATE_NML_STEP_KERNEL_H
+#define INTEGRATE_NML_STEP_KERNEL_H
+
+/* -------------------------------------------------------------------------- *
+ *                                   OpenMM                                   *
+ * -------------------------------------------------------------------------- *
+ * This is part of the OpenMM molecular simulation toolkit originating from   *
+ * Simbios, the NIH National Center for Physics-Based Simulation of           *
+ * Biological Structures at Stanford, funded under the NIH Roadmap for        *
+ * Medical Research, grant U54 GM072970. See https://simtk.org.               *
+ *                                                                            *
+ * Portions copyright (c) 2009 Stanford University and the Authors.           *
+ * Authors: Chris Sweet                                                       *
+ * Contributors: Christopher Bruns                                            *
+ *                                                                            *
+ * Permission is hereby granted, free of charge, to any person obtaining a    *
+ * copy of this software and associated documentation files (the "Software"), *
+ * to deal in the Software without restriction, including without limitation  *
+ * the rights to use, copy, modify, merge, publish, distribute, sublicense,   *
+ * and/or sell copies of the Software, and to permit persons to whom the      *
+ * Software is furnished to do so, subject to the following conditions:       *
+ *                                                                            *
+ * The above copyright notice and this permission notice shall be included in *
+ * all copies or substantial portions of the Software.                        *
+ *                                                                            *
+ * THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR *
+ * IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,   *
+ * FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL    *
+ * THE AUTHORS, CONTRIBUTORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM,    *
+ * DAMAGES OR OTHER LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR      *
+ * OTHERWISE, ARISING FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE  *
+ * USE OR OTHER DEALINGS IN THE SOFTWARE.                                     *
+ * -------------------------------------------------------------------------- */
+
+
+#include "openmm/KernelImpl.h"
+#include "openmm/System.h"
+#include "NMLIntegrator.h"
+
+namespace OpenMM {
+
+/**
+ * This kernel is invoked by NMLIntegrator to take one time step.
+ */
+class IntegrateNMLStepKernel : public KernelImpl {
+    public:
+        static std::string Name() {
+            return "IntegrateNMLStep";
+        }
+        IntegrateNMLStepKernel(std::string name, const Platform& platform) : KernelImpl(name, platform) {
+        }
+        /**
+         * Initialize the kernel.
+         *
+         * @param system     the System this kernel will be applied to
+         * @param integrator the NMLIntegrator this kernel will be used for
+         */
+        virtual void initialize(const System& system, const NMLIntegrator& integrator) = 0;
+        /**
+         * Execute the kernel.
+         *
+         * @param context    the context in which to execute this kernel
+         * @param integrator the NMLIntegrator this kernel is being used for
+         */
+        virtual void execute(ContextImpl& context, const NMLIntegrator& integrator, const double currentPE, const int stepType) = 0;
+};
+
+} /* namespace OpenMM */
+
+#endif /* INTEGRATE_NML_STEP_KERNEL_H */
+

plugins/normalModeLangevin/include/NMLIntegrator.h

+#ifndef OPENMM_NMLIntegrator_H_
+#define OPENMM_NMLIntegrator_H_
+
+/* -------------------------------------------------------------------------- *
+ *                                   OpenMM                                   *
+ * -------------------------------------------------------------------------- *
+ * This is part of the OpenMM molecular simulation toolkit originating from   *
+ * Simbios, the NIH National Center for Physics-Based Simulation of           *
+ * Biological Structures at Stanford, funded under the NIH Roadmap for        *
+ * Medical Research, grant U54 GM072970. See https://simtk.org.               *
+ *                                                                            *
+ * Portions copyright (c) 2009 Stanford University and the Authors.           *
+ * Authors: Chris Sweet                                                       *
+ * Contributors: Christopher Bruns                                            *
+ *                                                                            *
+ * Permission is hereby granted, free of charge, to any person obtaining a    *
+ * copy of this software and associated documentation files (the "Software"), *
+ * to deal in the Software without restriction, including without limitation  *
+ * the rights to use, copy, modify, merge, publish, distribute, sublicense,   *
+ * and/or sell copies of the Software, and to permit persons to whom the      *
+ * Software is furnished to do so, subject to the following conditions:       *
+ *                                                                            *
+ * The above copyright notice and this permission notice shall be included in *
+ * all copies or substantial portions of the Software.                        *
+ *                                                                            *
+ * THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR *
+ * IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,   *
+ * FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL    *
+ * THE AUTHORS, CONTRIBUTORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM,    *
+ * DAMAGES OR OTHER LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR      *
+ * OTHERWISE, ARISING FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE  *
+ * USE OR OTHER DEALINGS IN THE SOFTWARE.                                     *
+ * -------------------------------------------------------------------------- */
+
+#include "openmm/Integrator.h"
+#include "openmm/Kernel.h"
+#include "openmm/internal/windowsExport.h"
+#include "EigenvectorInfo.h"
+
+namespace OpenMM {
+
+/**
+ * This is an Integrator which simulates a System using Langevin dynamics.
+ */
+
+class OPENMM_EXPORT NMLIntegrator : public Integrator {
+public:
+    /**
+     * Create a NMLIntegrator.
+     *
+     * @param temperature    the temperature of the heat bath (in Kelvin)
+     * @param frictionCoeff  the friction coefficient which couples the system to the heat bath
+     * @param stepSize       the step size with which to integrator the system (in picoseconds)
+     */
+    NMLIntegrator(double temperature, double frictionCoeff, double stepSize, ProtoMol::EigenvectorInfo* projectionVectorInfo );
+    /**
+     * Get the temperature of the heat bath (in Kelvin).
+     */
+    double getTemperature() const {
+        return temperature;
+    }
+    /**
+     * Set the temperature of the heat bath (in Kelvin).
+     */
+    void setTemperature(double temp) {
+        temperature = temp;
+    }
+    /**
+     * Get the friction coefficient which determines how strongly the system is coupled to
+     * the heat bath.
+     */
+    double getFriction() const {
+        return friction;
+    }
+    /**
+     * Set the friction coefficient which determines how strongly the system is coupled to
+     * the heat bath.
+     */
+    void setFriction(double coeff) {
+        friction = coeff;
+    }
+
+
+	  void setProjectionVectorInfo( ProtoMol::EigenvectorInfo* inVectorInfo ){
+		  mProjectionVectorInfo = inVectorInfo;
+	  }
+
+	  ProtoMol::EigenvectorInfo* getProjectionVectorInfo() const{
+		  return mProjectionVectorInfo;
+	  }
+
+    unsigned int getNumProjectionVectors() const{
+      //return mProjectionVectorInfo->myNumEigenvectors;
+      return mProjectionVectorInfo->myNumUsedEigenvectors;
+    }
+
+    double getMinimumLimit() const{
+      return mProjectionVectorInfo->myMinimumLimit *  4.184; //Kcal->KJ
+    }
+
+    bool getProjVecChanged() const{      
+      return mProjectionVectorInfo->myEigVecChanged;
+    }
+
+    double* getProjectionVectors() const{      
+      return mProjectionVectorInfo->myEigenvectors;
+    }
+
+    void setProjVecChanged(bool inChanged){		  
+      mProjectionVectorInfo->myEigVecChanged = inChanged;
+    }
+
+    double getMaxEigenvalue() const{
+        //return max eigenvalue (f^2) in 1/ps from 1/charmm
+      return mProjectionVectorInfo->myMaxEigenvalue * 
+        (1000.0 / 48.88821290839616) * (1000.0 / 48.88821290839616) //(1/charmm to 1/ps)^2
+        * 100.0; //(Angstron to nm [inverse])^2
+    }
+
+    /**
+     * Advance a simulation through time by taking a series of time steps.
+     *
+     * @param steps   the number of time steps to take
+     */
+    void step(int steps);
+    
+    //Minimizer
+    void minimize(int maxsteps);
+    
+protected:
+    /**
+     * This will be called by the OpenMMContext 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.
+     * It will also get called again if the application calls reinitialize() on the OpenMMContext.
+     */
+    void initialize(ContextImpl& context);
+    /**
+     * Get the names of all Kernels used by this Integrator.
+     */
+    std::vector<std::string> getKernelNames();
+private:
+	  ProtoMol::EigenvectorInfo* mProjectionVectorInfo;
+
+    double temperature, friction;
+    ContextImpl* context;
+    Kernel kernel;
+};
+
+} // namespace OpenMM
+
+#endif /*OPENMM_NMLIntegrator_H_*/

plugins/normalModeLangevin/include/NMLKernels.h

+#ifndef NML_KERNELS_H
+#define NML_KERNELS_H
+
+/* -------------------------------------------------------------------------- *
+ *                                   OpenMM                                   *
+ * -------------------------------------------------------------------------- *
+ * This is part of the OpenMM molecular simulation toolkit originating from   *
+ * Simbios, the NIH National Center for Physics-Based Simulation of           *
+ * Biological Structures at Stanford, funded under the NIH Roadmap for        *
+ * Medical Research, grant U54 GM072970. See https://simtk.org.               *
+ *                                                                            *
+ * Portions copyright (c) 2009 Stanford University and the Authors.           *
+ * Authors: Christopher Bruns                                                 *
+ * Contributors:                                                              *
+ *                                                                            *
+ * Permission is hereby granted, free of charge, to any person obtaining a    *
+ * copy of this software and associated documentation files (the "Software"), *
+ * to deal in the Software without restriction, including without limitation  *
+ * the rights to use, copy, modify, merge, publish, distribute, sublicense,   *
+ * and/or sell copies of the Software, and to permit persons to whom the      *
+ * Software is furnished to do so, subject to the following conditions:       *
+ *                                                                            *
+ * The above copyright notice and this permission notice shall be included in *
+ * all copies or substantial portions of the Software.                        *
+ *                                                                            *
+ * THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR *
+ * IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,   *
+ * FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL    *
+ * THE AUTHORS, CONTRIBUTORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM,    *
+ * DAMAGES OR OTHER LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR      *
+ * OTHERWISE, ARISING FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE  *
+ * USE OR OTHER DEALINGS IN THE SOFTWARE.                                     *
+ * -------------------------------------------------------------------------- */
+
+#include "IntegrateNMLStepKernel.h"
+
+#endif /* NML_KERNELS_H */
+

plugins/normalModeLangevin/include/NormalModeLangevin.h

+#ifndef OPENMM_NORMAL_MODE_LANGEVIN_H_
+#define OPENMM_NORMAL_MODE_LANGEVIN_H_
+
+/* -------------------------------------------------------------------------- *
+ *                                   OpenMM                                   *
+ * -------------------------------------------------------------------------- *
+ * This is part of the OpenMM molecular simulation toolkit originating from   *
+ * Simbios, the NIH National Center for Physics-Based Simulation of           *
+ * Biological Structures at Stanford, funded under the NIH Roadmap for        *
+ * Medical Research, grant U54 GM072970. See https://simtk.org.               *
+ *                                                                            *
+ * Portions copyright (c) 2009 Stanford University and the Authors.           *
+ * Authors: Christopher Bruns                                                 *
+ * Contributors:                                                              *
+ *                                                                            *
+ * This program is free software: you can redistribute it and/or modify       *
+ * it under the terms of the GNU Lesser General Public License as published   *
+ * by the Free Software Foundation, either version 3 of the License, or       *
+ * (at your option) any later version.                                        *
+ *                                                                            *
+ * This program is distributed in the hope that it will be useful,            *
+ * but WITHOUT ANY WARRANTY; without even the implied warranty of             *
+ * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the              *
+ * GNU Lesser General Public License for more details.                        *
+ *                                                                            *
+ * You should have received a copy of the GNU Lesser General Public License   *
+ * along with this program.  If not, see <http://www.gnu.org/licenses/>.      *
+ * -------------------------------------------------------------------------- */
+
+
+namespace OpenMM {
+} // namespace OpenMM
+
+#endif /*OPENMM_NORMAL_MODE_LANGEVIN_H_*/
+

plugins/normalModeLangevin/platforms/reference/CMakeLists.txt

+#/* -------------------------------------------------------------------------- *
+# *                                   OpenMM                                   *
+# * -------------------------------------------------------------------------- *
+# * This is part of the OpenMM molecular simulation toolkit originating from   *
+# * Simbios, the NIH National Center for Physics-Based Simulation of           *
+# * Biological Structures at Stanford, funded under the NIH Roadmap for        *
+# * Medical Research, grant U54 GM072970. See https://simtk.org.               *
+# *                                                                            *
+# * Portions copyright (c) 2009 Stanford University and the Authors.           *
+# * Authors: Chris Sweet, Christopher Bruns                                    *
+# * Contributors:                                                              *
+# *                                                                            *
+# * Permission is hereby granted, free of charge, to any person obtaining a    *
+# * copy of this software and associated documentation files (the "Software"), *
+# * to deal in the Software without restriction, including without limitation  *
+# * the rights to use, copy, modify, merge, publish, distribute, sublicense,   *
+# * and/or sell copies of the Software, and to permit persons to whom the      *
+# * Software is furnished to do so, subject to the following conditions:       *
+# *                                                                            *
+# * The above copyright notice and this permission notice shall be included in *
+# * all copies or substantial portions of the Software.                        *
+# *                                                                            *
+# * THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR *
+# * IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,   *
+# * FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL    *
+# * THE AUTHORS, CONTRIBUTORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM,    *
+# * DAMAGES OR OTHER LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR      *
+# * OTHERWISE, ARISING FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE  *
+# * USE OR OTHER DEALINGS IN THE SOFTWARE.                                     *
+# * -------------------------------------------------------------------------- */
+

plugins/normalModeLangevin/platforms/reference/include/ReferenceIntegrateNMLStepKernel.h

+#ifndef REFERENCE_INTEGRATE_NML_STEP_KERNEL_H
+#define REFERENCE_INTEGRATE_NML_STEP_KERNEL_H
+
+
+/* -------------------------------------------------------------------------- *
+ *                                   OpenMM                                   *
+ * -------------------------------------------------------------------------- *
+ * This is part of the OpenMM molecular simulation toolkit originating from   *
+ * Simbios, the NIH National Center for Physics-Based Simulation of           *
+ * Biological Structures at Stanford, funded under the NIH Roadmap for        *
+ * Medical Research, grant U54 GM072970. See https://simtk.org.               *
+ *                                                                            *
+ * Portions copyright (c) 2009 Stanford University and the Authors.           *
+ * Authors: Chris Sweet                                                       *
+ * Contributors: Christopher Bruns                                            *
+ *                                                                            *
+ * Permission is hereby granted, free of charge, to any person obtaining a    *
+ * copy of this software and associated documentation files (the "Software"), *
+ * to deal in the Software without restriction, including without limitation  *
+ * the rights to use, copy, modify, merge, publish, distribute, sublicense,   *
+ * and/or sell copies of the Software, and to permit persons to whom the      *
+ * Software is furnished to do so, subject to the following conditions:       *
+ *                                                                            *
+ * The above copyright notice and this permission notice shall be included in *
+ * all copies or substantial portions of the Software.                        *
+ *                                                                            *
+ * THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR *
+ * IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,   *
+ * FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL    *
+ * THE AUTHORS, CONTRIBUTORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM,    *
+ * DAMAGES OR OTHER LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR      *
+ * OTHERWISE, ARISING FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE  *
+ * USE OR OTHER DEALINGS IN THE SOFTWARE.                                     *
+ * -------------------------------------------------------------------------- */
+
+#include "ReferencePlatform.h"
+#include "ReferenceNMLDynamics.h"
+#include "IntegrateNMLStepKernel.h"
+
+namespace OpenMM {
+
+/**
+ * This kernel is invoked by NMLIntegrator to take one time step.
+ */
+class ReferenceIntegrateNMLStepKernel : public IntegrateNMLStepKernel {
+    public:
+        ReferenceIntegrateNMLStepKernel(std::string name, const Platform& platform, ReferencePlatform::PlatformData& data) : IntegrateNMLStepKernel(name, platform),
+        data(data), dynamics(0), constraints(0), masses(0), constraintDistances(0), constraintIndices(0), projectionVectors(0) {
+        }
+        ~ReferenceIntegrateNMLStepKernel();
+        /**
+         * Initialize the kernel, setting up the particle masses.
+         *
+         * @param system     the System this kernel will be applied to
+         * @param integrator the NMLIntegrator this kernel will be used for
+         */
+        void initialize(const System& system, const NMLIntegrator& integrator);
+        /**
+         * Execute the kernel.
+         *
+         * @param context    the context in which to execute this kernel
+         * @param integrator the NMLIntegrator this kernel is being used for
+         */
+        void execute(ContextImpl& context, const NMLIntegrator& integrator, const double currentPE, const int stepType);
+    private:
+        ReferencePlatform::PlatformData& data;
+        ReferenceNMLDynamics* dynamics;
+        ReferenceConstraintAlgorithm* constraints;
+        RealOpenMM* masses;
+        RealOpenMM* constraintDistances;
+        int** constraintIndices;
+        int numConstraints;
+        double prevTemp, prevFriction, prevStepSize;
+        //double prevTemp, prevFriction, prevErrorTol;
+        RealOpenMM* projectionVectors;
+
+};
+
+} /* namespace OpenMM */
+
+#endif /* REFERENCE_INTEGRATE_NML_STEP_KERNEL_H */
+

plugins/normalModeLangevin/platforms/reference/include/ReferenceNMLDynamics.h

+#ifndef __ReferenceNMLDynamics_H__
+#define __ReferenceNMLDynamics_H__
+
+/* -------------------------------------------------------------------------- *
+ *                                   OpenMM                                   *
+ * -------------------------------------------------------------------------- *
+ * This is part of the OpenMM molecular simulation toolkit originating from   *
+ * Simbios, the NIH National Center for Physics-Based Simulation of           *
+ * Biological Structures at Stanford, funded under the NIH Roadmap for        *
+ * Medical Research, grant U54 GM072970. See https://simtk.org.               *
+ *                                                                            *
+ * Portions copyright (c) 2009 Stanford University and the Authors.           *
+ * Authors: Chris Sweet                                                       *
+ * Contributors: Christopher Bruns, Pande Group                               *
+ *                                                                            *
+ * Permission is hereby granted, free of charge, to any person obtaining a    *
+ * copy of this software and associated documentation files (the "Software"), *
+ * to deal in the Software without restriction, including without limitation  *
+ * the rights to use, copy, modify, merge, publish, distribute, sublicense,   *
+ * and/or sell copies of the Software, and to permit persons to whom the      *
+ * Software is furnished to do so, subject to the following conditions:       *
+ *                                                                            *
+ * The above copyright notice and this permission notice shall be included in *
+ * all copies or substantial portions of the Software.                        *
+ *                                                                            *
+ * THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR *
+ * IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,   *
+ * FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL    *
+ * THE AUTHORS, CONTRIBUTORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM,    *
+ * DAMAGES OR OTHER LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR      *
+ * OTHERWISE, ARISING FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE  *
+ * USE OR OTHER DEALINGS IN THE SOFTWARE.                                     *
+ * -------------------------------------------------------------------------- */
+
+#include "SimTKUtilities/SimTKOpenMMRealType.h"
+#include "SimTKReference/ReferenceDynamics.h"
+
+// ---------------------------------------------------------------------------------------
+
+class ReferenceNMLDynamics : public ReferenceDynamics {
+
+   private:
+
+      //dt, myGamma, fdt, vdt, ndt, sqrtFCoverM from Langevin Leapfrog
+      //enum FixedParameters   { GDT, EPH, EMH, EP, EM, B, C, D, V, X, Yv, Yx, MaxFixedParameters };
+      enum FixedParameters   { DT, GAMMA, FDT, VDT, NDT, SQRTFCOVERM, MaxFixedParameters };
+      enum TwoDArrayIndicies { X2D, V2D, OldV, xPrime2D, vPrime2D, Max2DArrays };
+      enum OneDArrayIndicies { InverseMasses, Max1DArrays };
+
+      RealOpenMM _tau;
+      RealOpenMM _fixedParameters[MaxFixedParameters];
+      RealOpenMM* _projectionVectors;
+      unsigned int _numProjectionVectors;
+      RealOpenMM _minimumLimit, _maxEig;
+      RealOpenMM lastPE, lastSlope;
+
+      /**---------------------------------------------------------------------------------------
+      
+         Set fixed values
+      
+         @return ReferenceDynamics::DefaultReturn
+      
+         --------------------------------------------------------------------------------------- */
+      
+      int _setFixedParameters( void );
+      
+   public:
+
+      /**---------------------------------------------------------------------------------------
+      
+         Constructor
+
+         @param numberOfAtoms  number of atoms
+         @param deltaT         delta t for dynamics
+         @param tau            viscosity
+         @param temperature    temperature
+      
+         --------------------------------------------------------------------------------------- */
+
+       ReferenceNMLDynamics( int numberOfAtoms, RealOpenMM deltaT, RealOpenMM tau, RealOpenMM temperature,
+                   RealOpenMM* projectionVectors, unsigned int numProjectionVectors, 
+                   RealOpenMM minimumLimit, RealOpenMM maxEig
+                    );
+
+      /**---------------------------------------------------------------------------------------
+      
+         Destructor
+      
+         --------------------------------------------------------------------------------------- */
+
+       ~ReferenceNMLDynamics( );
+
+      /**---------------------------------------------------------------------------------------
+      
+         Get tau
+      
+         @return tau
+      
+         --------------------------------------------------------------------------------------- */
+      
+      RealOpenMM getTau( void ) const;
+      
+      /**---------------------------------------------------------------------------------------
+      
+         Get array of fixed parameters indexed by 'FixedParameters' enums
+      
+         @return array
+      
+         --------------------------------------------------------------------------------------- */
+      
+      const RealOpenMM* getFixedParameters( void ) const;
+      
+      /**---------------------------------------------------------------------------------------
+      
+         Print parameters
+      
+         @param message message
+
+         @return ReferenceDynamics::DefaultReturn
+      
+         --------------------------------------------------------------------------------------- */
+      
+      int printParameters( std::stringstream& message ) const;
+      
+      /**---------------------------------------------------------------------------------------
+      
+         Update
+      
+         @param numberOfAtoms       number of atoms
+         @param atomCoordinates     atom coordinates
+         @param velocities          velocities
+         @param forces              forces
+         @param masses              atom masses
+      
+         @return ReferenceDynamics::DefaultReturn
+      
+         --------------------------------------------------------------------------------------- */
+     
+      int update( int numberOfAtoms, RealOpenMM** atomCoordinates,
+                  RealOpenMM** velocities, RealOpenMM** forces, RealOpenMM* masses, const RealOpenMM currentPE, const int stepType );
+     
+    /**---------------------------------------------------------------------------------------
+     
+     Velocity update; Langevin Leapfrog
+     
+     @param numberOfAtoms       number of atoms
+     @param atomCoordinates     atom coordinates
+     @param velocities          velocities
+     @param forces              forces
+     @param masses              atom masses
+     @param inverseMasses       inverse atom masses
+     @param xVector             xVector
+     
+     --------------------------------------------------------------------------------------- */
+    void halfKick( int numberOfAtoms, RealOpenMM** atomCoordinates, RealOpenMM** velocities,
+                    RealOpenMM** forces, RealOpenMM* masses, RealOpenMM* inverseMasses, RealOpenMM** xVector );
+    
+    /**---------------------------------------------------------------------------------------
+     
+     Position update; Langevin Leapfrog
+     
+     @param numberOfAtoms       number of atoms
+     @param atomCoordinates     atom coordinates
+     @param velocities          velocities
+     @param forces              forces
+     @param masses              atom masses
+     @param inverseMasses       inverse atom masses
+     @param xPrime              xPrime
+     
+     --------------------------------------------------------------------------------------- */
+    void drift( int numberOfAtoms, RealOpenMM** atomCoordinates, RealOpenMM** velocities,
+                  RealOpenMM** forces, RealOpenMM* masses, RealOpenMM* inverseMasses,
+                  RealOpenMM** xPrime );
+          
+      /**---------------------------------------------------------------------------------------
+      
+         Write state
+      
+         @param numberOfAtoms       number of atoms
+         @param atomCoordinates     atom coordinates
+         @param velocities          velocities
+         @param forces              forces
+         @param masses              atom masses
+         @param state               0 if initial state; otherwise nonzero
+         @param baseFileName        base file name
+      
+         @return ReferenceDynamics::DefaultReturn
+      
+         --------------------------------------------------------------------------------------- */
+      
+      int writeState( int numberOfAtoms, RealOpenMM** atomCoordinates,
+                      RealOpenMM** velocities, RealOpenMM** forces, RealOpenMM* masses,
+                      int state, const std::string& baseFileName ) const;
+
+      /**---------------------------------------------------------------------------------------
+       Find forces OR positions inside subspace (defined as the span of the 'eigenvectors' Q)
+       Take 'array' as input, 'outArray' as output (may be the same vector).
+       'projectionMode' determines:
+       a) if the projection is for force or positions, bit 1 set=force.
+       b) if the projection is for the sub-space or-complement space, bit 2 set=sub-space.
+       ----------------------------------------------------------------------------------------- */
+      void subspaceProjection(  RealOpenMM** arrayParam, 
+                                RealOpenMM** outArrayParam, 
+                                int numberOfAtoms,
+                                RealOpenMM* masses,
+                                RealOpenMM* inverseMasses,
+                                int projectionMode);
+  
+};
+
+// ---------------------------------------------------------------------------------------
+
+#endif // __ReferenceNMLDynamics_H__

plugins/normalModeLangevin/platforms/reference/include/ReferenceNMLKernelFactory.h

+#ifndef OPENMM_REFERENCE_NMLKERNELFACTORY_H_
+#define OPENMM_REFERENCE_NMLKERNELFACTORY_H_
+
+/* -------------------------------------------------------------------------- *
+ *                                   OpenMM                                   *
+ * -------------------------------------------------------------------------- *
+ * This is part of the OpenMM molecular simulation toolkit originating from   *
+ * Simbios, the NIH National Center for Physics-Based Simulation of           *
+ * Biological Structures at Stanford, funded under the NIH Roadmap for        *
+ * Medical Research, grant U54 GM072970. See https://simtk.org.               *
+ *                                                                            *
+ * Portions copyright (c) 2009 Stanford University and the Authors.           *
+ * Authors: Chris Sweet                                                       *
+ * Contributors: Christopher Bruns                                            *
+ *                                                                            *
+ * Permission is hereby granted, free of charge, to any person obtaining a    *
+ * copy of this software and associated documentation files (the "Software"), *
+ * to deal in the Software without restriction, including without limitation  *
+ * the rights to use, copy, modify, merge, publish, distribute, sublicense,   *
+ * and/or sell copies of the Software, and to permit persons to whom the      *
+ * Software is furnished to do so, subject to the following conditions:       *
+ *                                                                            *
+ * The above copyright notice and this permission notice shall be included in *
+ * all copies or substantial portions of the Software.                        *
+ *                                                                            *
+ * THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR *
+ * IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,   *
+ * FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL    *
+ * THE AUTHORS, CONTRIBUTORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM,    *
+ * DAMAGES OR OTHER LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR      *
+ * OTHERWISE, ARISING FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE  *
+ * USE OR OTHER DEALINGS IN THE SOFTWARE.                                     *
+ * -------------------------------------------------------------------------- */
+
+#include "openmm/KernelFactory.h"
+
+namespace OpenMM {
+
+/**
+ * This KernelFactory creates all kernels for NormalModeLangevin.
+ */
+
+class ReferenceNMLKernelFactory : public KernelFactory {
+public:
+    KernelImpl* createKernelImpl(std::string name, const Platform& platform, ContextImpl& context) const;
+};
+
+} // namespace OpenMM
+
+#endif /*OPENMM_REFERENCE_NMLKERNELFACTORY_H_*/
+

plugins/normalModeLangevin/platforms/reference/src/ReferenceIntegrateNMLStepKernel.cpp

+/* -------------------------------------------------------------------------- *
+ *                                   OpenMM                                   *
+ * -------------------------------------------------------------------------- *
+ * This is part of the OpenMM molecular simulation toolkit originating from   *
+ * Simbios, the NIH National Center for Physics-Based Simulation of           *
+ * Biological Structures at Stanford, funded under the NIH Roadmap for        *
+ * Medical Research, grant U54 GM072970. See https://simtk.org.               *
+ *                                                                            *
+ * Portions copyright (c) 2009 Stanford University and the Authors.           *
+ * Authors: Chris Sweet                                                       *
+ * Contributors: Christopher Bruns                                            *
+ *                                                                            *
+ * Permission is hereby granted, free of charge, to any person obtaining a    *
+ * copy of this software and associated documentation files (the "Software"), *
+ * to deal in the Software without restriction, including without limitation  *
+ * the rights to use, copy, modify, merge, publish, distribute, sublicense,   *
+ * and/or sell copies of the Software, and to permit persons to whom the      *
+ * Software is furnished to do so, subject to the following conditions:       *
+ *                                                                            *
+ * The above copyright notice and this permission notice shall be included in *
+ * all copies or substantial portions of the Software.                        *
+ *                                                                            *
+ * THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR *
+ * IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,   *
+ * FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL    *
+ * THE AUTHORS, CONTRIBUTORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM,    *
+ * DAMAGES OR OTHER LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR      *
+ * OTHERWISE, ARISING FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE  *
+ * USE OR OTHER DEALINGS IN THE SOFTWARE.                                     *
+ * -------------------------------------------------------------------------- */
+
+#include "ReferenceIntegrateNMLStepKernel.h"
+#include "openmm/HarmonicAngleForce.h"
+#include "openmm/internal/ContextImpl.h"
+#include "SimTKReference/ReferenceCCMAAlgorithm.h"
+#include <vector>
+
+using namespace std;
+
+namespace OpenMM {
+
+// These methods lifted from ReferenceKernels.cpp, so NormalModeLangevin could be a plugin
+//   allocateIntArray
+//   disposeIntArray
+//   extractPositions
+//   extractVelocities
+//   extractForces
+//   findAnglesForCCMA
+
+static int** allocateIntArray(int length, int width) {
+    int** array = new int*[length];
+    for (int i = 0; i < length; ++i)
+        array[i] = new int[width];
+    return array;
+}
+
+static void disposeIntArray(int** array, int size) {
+    if (array) {
+        for (int i = 0; i < size; ++i)
+            delete[] array[i];
+        delete[] array;
+    }
+}
+
+static RealOpenMM** extractPositions(ContextImpl& context) {
+    ReferencePlatform::PlatformData* data = reinterpret_cast<ReferencePlatform::PlatformData*>(context.getPlatformData());
+    return (RealOpenMM**) data->positions;
+}
+
+static RealOpenMM** extractVelocities(ContextImpl& context) {
+    ReferencePlatform::PlatformData* data = reinterpret_cast<ReferencePlatform::PlatformData*>(context.getPlatformData());
+    return (RealOpenMM**) data->velocities;
+}
+
+static RealOpenMM** extractForces(ContextImpl& context) {
+    ReferencePlatform::PlatformData* data = reinterpret_cast<ReferencePlatform::PlatformData*>(context.getPlatformData());
+    return (RealOpenMM**) data->forces;
+}
+
+static void findAnglesForCCMA(const System& system, vector<ReferenceCCMAAlgorithm::AngleInfo>& angles) {
+    for (int i = 0; i < system.getNumForces(); i++) {
+        const HarmonicAngleForce* force = dynamic_cast<const HarmonicAngleForce*>(&system.getForce(i));
+        if (force != NULL) {
+            for (int j = 0; j < force->getNumAngles(); j++) {
+                int atom1, atom2, atom3;
+                double angle, k;
+                force->getAngleParameters(j, atom1, atom2, atom3, angle, k);
+                angles.push_back(ReferenceCCMAAlgorithm::AngleInfo(atom1, atom2, atom3, (RealOpenMM)angle));
+            }
+        }
+    }
+}
+
+
+ReferenceIntegrateNMLStepKernel::~ReferenceIntegrateNMLStepKernel() {
+    if (dynamics)
+        delete dynamics;
+    if (constraints)
+        delete constraints;
+    if (masses)
+        delete[] masses;
+    if (constraintIndices)
+        disposeIntArray(constraintIndices, numConstraints);
+    if (constraintDistances)
+        delete[] constraintDistances;
+    if (projectionVectors)
+        delete projectionVectors;
+}
+
+void ReferenceIntegrateNMLStepKernel::initialize(const System& system, const NMLIntegrator& integrator) {
+    int numParticles = system.getNumParticles();
+    masses = new RealOpenMM[numParticles];
+    for (int i = 0; i < numParticles; ++i)
+        masses[i] = static_cast<RealOpenMM>(system.getParticleMass(i));
+    numConstraints = system.getNumConstraints();
+    constraintIndices = allocateIntArray(numConstraints, 2);
+    constraintDistances = new RealOpenMM[numConstraints];
+    for (int i = 0; i < numConstraints; ++i) {
+        int particle1, particle2;
+        double distance;
+        system.getConstraintParameters(i, particle1, particle2, distance);
+        constraintIndices[i][0] = particle1;
+        constraintIndices[i][1] = particle2;
+        constraintDistances[i] = static_cast<RealOpenMM>(distance);
+    }
+    //SimTKOpenMMUtilities::setRandomNumberSeed((unsigned int) integrator.getRandomNumberSeed());
+}
+
+void ReferenceIntegrateNMLStepKernel::execute(ContextImpl& context, const NMLIntegrator& integrator, const double currentPE, const int stepType) {
+
+    double temperature = integrator.getTemperature();
+    double friction = integrator.getFriction();
+    double stepSize = integrator.getStepSize();
+    double* dProjectionVectors = integrator.getProjectionVectors();
+    unsigned int numProjectionVectors = integrator.getNumProjectionVectors();
+    bool projVecChanged = integrator.getProjVecChanged();
+    double minimumLimit = integrator.getMinimumLimit();
+    double maxEig = integrator.getMaxEigenvalue();
+
+
+    RealOpenMM** posData = extractPositions(context);
+    RealOpenMM** velData = extractVelocities(context);
+    RealOpenMM** forceData = extractForces(context);
+
+    //projection vectors
+    if (projectionVectors == 0 || projVecChanged) {
+        std::cout << "Setting vectors " << projVecChanged << std::endl;
+        unsigned int arraySz = numProjectionVectors * context.getSystem().getNumParticles() * 3;
+        if(projectionVectors == 0) {
+            projectionVectors = new RealOpenMM[arraySz];
+        }
+        for(unsigned int i=0; i<arraySz; i++) {
+            projectionVectors[i] = static_cast<RealOpenMM>(dProjectionVectors[i]);
+        }
+
+    }
+
+
+    if (dynamics == 0 ){ //|| temperature != prevTemp || friction != prevFriction || stepSize != prevStepSize) {
+        // Recreate the computation objects with the new parameters.
+
+        if (dynamics) {
+            delete dynamics;
+            delete constraints;
+        }
+        RealOpenMM tau = static_cast<RealOpenMM>( friction == 0.0 ? 0.0 : 1.0/friction );
+
+        dynamics = new ReferenceNMLDynamics(context.getSystem().getNumParticles(),
+                                            static_cast<RealOpenMM>(stepSize),
+                                            static_cast<RealOpenMM>(tau),
+
+                                            static_cast<RealOpenMM>(temperature),
+                                            projectionVectors,
+                                            numProjectionVectors,
+                                            static_cast<RealOpenMM>(minimumLimit),
+                                            static_cast<RealOpenMM>(maxEig));
+
+        vector<ReferenceCCMAAlgorithm::AngleInfo> angles;
+        findAnglesForCCMA(context.getSystem(), angles);
+        constraints = new ReferenceCCMAAlgorithm(context.getSystem().getNumParticles(), numConstraints, constraintIndices, constraintDistances, masses, angles, (RealOpenMM)integrator.getConstraintTolerance());
+        dynamics->setReferenceConstraintAlgorithm(constraints);
+        prevTemp = temperature;
+        prevFriction = friction;
+        prevStepSize = stepSize;
+    }
+    //for minimization quadratic fit we need the potential energy
+    RealOpenMM currentPE2 = static_cast<RealOpenMM>(currentPE);
+    //
+    dynamics->update(context.getSystem().getNumParticles(), posData, velData, forceData, masses, currentPE2, stepType );
+    //update at dynamic step 2
+    if(stepType == 2){
+        data.time += stepSize;
+        data.stepCount++;
+    }
+
+}
+
+} /* namespace OpenMM */
+

plugins/normalModeLangevin/platforms/reference/src/ReferenceNMLDynamics.cpp

+/* -------------------------------------------------------------------------- *
+ *                                   OpenMM                                   *
+ * -------------------------------------------------------------------------- *
+ * This is part of the OpenMM molecular simulation toolkit originating from   *
+ * Simbios, the NIH National Center for Physics-Based Simulation of           *
+ * Biological Structures at Stanford, funded under the NIH Roadmap for        *
+ * Medical Research, grant U54 GM072970. See https://simtk.org.               *
+ *                                                                            *
+ * Portions copyright (c) 2009 Stanford University and the Authors.           *
+ * Authors: Chris Sweet                                                       *
+ * Contributors: Christopher Bruns, Pande Group                               *
+ *                                                                            *
+ * Permission is hereby granted, free of charge, to any person obtaining a    *
+ * copy of this software and associated documentation files (the "Software"), *
+ * to deal in the Software without restriction, including without limitation  *
+ * the rights to use, copy, modify, merge, publish, distribute, sublicense,   *
+ * and/or sell copies of the Software, and to permit persons to whom the      *
+ * Software is furnished to do so, subject to the following conditions:       *
+ *                                                                            *
+ * The above copyright notice and this permission notice shall be included in *
+ * all copies or substantial portions of the Software.                        *
+ *                                                                            *
+ * THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR *
+ * IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,   *
+ * FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL    *
+ * THE AUTHORS, CONTRIBUTORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM,    *
+ * DAMAGES OR OTHER LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR      *
+ * OTHERWISE, ARISING FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE  *
+ * USE OR OTHER DEALINGS IN THE SOFTWARE.                                     *
+ * -------------------------------------------------------------------------- */
+
+#include <string.h>
+#include <sstream>
+
+#include "SimTKUtilities/SimTKOpenMMCommon.h"
+#include "SimTKUtilities/SimTKOpenMMLog.h"
+#include "SimTKUtilities/SimTKOpenMMUtilities.h"
+#include "ReferenceNMLDynamics.h"
+
+/**---------------------------------------------------------------------------------------
+
+   ReferenceNMLDynamics constructor
+
+   @param numberOfAtoms  number of atoms
+   @param deltaT         delta t for dynamics
+   @param tau            viscosity(?)
+   @param temperature    temperature
+
+   --------------------------------------------------------------------------------------- */
+
+ReferenceNMLDynamics::ReferenceNMLDynamics( int numberOfAtoms,
+                                              RealOpenMM deltaT, RealOpenMM tau,
+                                              RealOpenMM temperature,
+                                              RealOpenMM* projectionVectors, 
+                                              unsigned int numProjectionVectors, 
+                                              RealOpenMM minimumLimit, RealOpenMM maxEig  ) : 
+            ReferenceDynamics( numberOfAtoms, deltaT, temperature ), _tau( tau ),
+            _projectionVectors(projectionVectors), _numProjectionVectors(numProjectionVectors), 
+            _minimumLimit(minimumLimit), _maxEig(maxEig)  {
+
+   // ---------------------------------------------------------------------------------------
+
+   static const char* methodName      = "\nReferenceNMLDynamics::ReferenceNMLDynamics";
+
+   static const RealOpenMM zero       =  0.0;
+   static const RealOpenMM one        =  1.0;
+
+   // ---------------------------------------------------------------------------------------
+
+   // insure tau is not zero -- if it is print warning message
+
+   if( _tau == zero ){
+
+      std::stringstream message;
+      message << methodName;
+      message << " input tau value=" << tau << " is invalid -- setting to 1.";
+      SimTKOpenMMLog::printError( message );
+
+      _tau = one;
+     
+   }
+   _setFixedParameters( );
+
+   allocate2DArrays( numberOfAtoms, 3, Max2DArrays );
+   allocate1DArrays( numberOfAtoms, Max1DArrays );
+   
+}
+
+/**---------------------------------------------------------------------------------------
+
+   ReferenceNMLDynamics destructor
+
+   --------------------------------------------------------------------------------------- */
+
+ReferenceNMLDynamics::~ReferenceNMLDynamics( ){
+
+   // ---------------------------------------------------------------------------------------
+
+   // static const char* methodName = "\nReferenceNMLDynamics::~ReferenceNMLDynamics";
+
+   // ---------------------------------------------------------------------------------------
+
+}
+
+/**---------------------------------------------------------------------------------------
+
+   Set fixed parameters
+
+   @return ReferenceDynamics::DefaultReturn
+
+   --------------------------------------------------------------------------------------- */
+
+int ReferenceNMLDynamics::_setFixedParameters( void ){
+
+   // ---------------------------------------------------------------------------------------
+
+   // static const char* methodName = "\nReferenceNMLDynamics::_setFixedParameters";
+
+   static const RealOpenMM zero       =  0.0;
+   static const RealOpenMM one        =  1.0;
+   static const RealOpenMM two        =  2.0;
+   static const RealOpenMM half       =  0.5;
+    static const RealOpenMM eighty       =  80.0;
+
+   // ---------------------------------------------------------------------------------------
+
+    //dt, myGamma, fdt, vdt, ndt, sqrtFCoverM from Langevin Leapfrog
+    
+    _fixedParameters[DT]         = getTimeStep();
+    _fixedParameters[GAMMA]      = (_tau == zero ? eighty : one/_tau);
+    _fixedParameters[FDT]        = ( one - EXP( -half * _fixedParameters[GAMMA] * _fixedParameters[DT] ) ) / _fixedParameters[GAMMA];
+    _fixedParameters[VDT]        = EXP(-half * _fixedParameters[GAMMA] * _fixedParameters[DT]);
+    _fixedParameters[NDT]        = SQRT( ( one - EXP( -_fixedParameters[GAMMA] * _fixedParameters[DT] ) ) / (two * _fixedParameters[GAMMA]) ); 
+    _fixedParameters[SQRTFCOVERM]        = SQRT( two * (RealOpenMM) BOLTZ * getTemperature() * _fixedParameters[GAMMA] );
+    
+   return 0; // ReferenceDynamics::DefaultReturn;
+
+};
+
+/**---------------------------------------------------------------------------------------
+
+   Get tau
+
+   @return tau
+
+   --------------------------------------------------------------------------------------- */
+
+RealOpenMM ReferenceNMLDynamics::getTau( void ) const {
+
+   // ---------------------------------------------------------------------------------------
+
+   // static const char* methodName  = "\nReferenceNMLDynamics::getTau";
+
+   // ---------------------------------------------------------------------------------------
+
+   return _tau;
+}
+
+/**---------------------------------------------------------------------------------------
+
+   Get array of fixed parameters indexed by 'FixedParameters' enums
+
+   @return array
+
+   --------------------------------------------------------------------------------------- */
+   
+const RealOpenMM* ReferenceNMLDynamics::getFixedParameters( void ) const {
+
+   // ---------------------------------------------------------------------------------------
+
+   // static const char* methodName  = "\nReferenceNMLDynamics::getFixedParameters";
+
+   // ---------------------------------------------------------------------------------------
+
+   return _fixedParameters;
+}
+
+/**---------------------------------------------------------------------------------------
+
+   Print parameters
+
+   @param message             message
+
+   @return ReferenceDynamics::DefaultReturn
+
+   --------------------------------------------------------------------------------------- */
+
+int ReferenceNMLDynamics::printParameters( std::stringstream& message ) const {
+
+   // ---------------------------------------------------------------------------------------
+    //dt, myGamma, fdt, vdt, ndt, sqrtFCoverM from Langevin Leapfrog
+
+   static const char* methodName  = "\nReferenceNMLDynamics::printParameters";
+   static const char* parameterNames[MaxFixedParameters] = { "dt", "myGamma", "fdt", "ndt", "sqrtFCoverM" };
+
+   // ---------------------------------------------------------------------------------------
+
+   // print parameters
+
+   ReferenceDynamics::printParameters( message );
+   message << " tau=" << getTau();
+   message << " T=" << getTemperature();
+   int cut = 3;
+   for( int ii = 0; ii < MaxFixedParameters; ii++ ){
+      message << " " << parameterNames[ii] << "=" << _fixedParameters[ii];
+      if( cut++ > 5 ){
+         cut = 0;
+         message << std::endl;
+      }
+   }
+
+   return 0; // ReferenceDynamics::DefaultReturn;
+
+}
+
+/**---------------------------------------------------------------------------------------
+ 
+ Drift update
+ 
+ @param numberOfAtoms       number of atoms
+ @param atomCoordinates     atom coordinates
+ @param velocities          velocities
+ @param forces              forces
+ @param inverseMasses       inverse atom masses
+ @param xPrime              xPrime
+ @param oldVelocities       previous velocities
+ @param xVector             xVector
+ @param vVector             vVector
+ 
+ @return ReferenceDynamics::DefaultReturn
+ 
+ --------------------------------------------------------------------------------------- */
+
+void ReferenceNMLDynamics::halfKick( int numberOfAtoms, RealOpenMM** atomCoordinates,
+                                      RealOpenMM** velocities,
+                                      RealOpenMM** forces, RealOpenMM* masses, 
+                                      RealOpenMM* inverseMasses, RealOpenMM** xVector ){
+    
+    //dt, myGamma, fdt, vdt, ndt, sqrtFCoverM from Langevin Leapfrog
+    static const RealOpenMM dt = getTimeStep(); // in ps
+    static const RealOpenMM myGamma = (_tau == 0.0 ? 80.0 : 1.0/_tau);
+    static const RealOpenMM fdt = ( 1.0 - EXP( -0.5 * myGamma * dt ) ) / myGamma;
+    static const RealOpenMM vdt = EXP(-0.5 * myGamma * dt);
+    static const RealOpenMM ndt = SQRT( ( 1.0 - EXP( -myGamma * dt ) ) / (2.0 * myGamma) ); 
+    static const RealOpenMM sqrtFCoverM = SQRT( 2.0 * (RealOpenMM) BOLTZ * getTemperature() * myGamma );
+    
+    //project forces into sub space. TODO put this somewhere callable from NMLIntegrator, else projecting twice.
+    subspaceProjection(forces, forces, numberOfAtoms, masses, inverseMasses, 0);
+
+    for( int ii = 0; ii < numberOfAtoms; ii++ ){
+        //RealOpenMM sqrtInverseMass = SQRT( inverseMasses[ii] )*fixedParameters[X];
+        for( int jj = 0; jj < 3; jj++ ){
+            xVector[ii][jj] = SimTKOpenMMUtilities::getNormallyDistributedRandomNumber() * SQRT( masses[ii] );
+        }
+    }
+    
+    subspaceProjection(xVector, xVector, numberOfAtoms, masses, inverseMasses, 0);
+    
+    for( int ii = 0; ii < numberOfAtoms; ii++ ){
+        //RealOpenMM sqrtInverseMass = SQRT( inverseMasses[ii] )*fixedParameters[X];
+        for( int jj = 0; jj < 3; jj++ ){
+            xVector[ii][jj] *= inverseMasses[ii] ;
+        }
+    }
+    
+    
+    
+    for( int i = 0; i < numberOfAtoms; i++ ) {
+        // semi-update velocities
+        for( int j = 0; j < 3; j++ ) {
+            velocities[i][j] = velocities[i][j] * vdt
+                                    + forces[i][j] * fdt * inverseMasses[i]
+                                        + xVector[i][j] * sqrtFCoverM * ndt;
+                                        //+ SimTKOpenMMUtilities::getNormallyDistributedRandomNumber() * sqrtFCoverM * ndt;
+        }
+    }
+    
+    //project forces into complement space, put in xPrime
+    subspaceProjection(velocities, velocities, numberOfAtoms, masses, inverseMasses, 1);    //in subspace but position/velocity projection
+
+}
+
+void ReferenceNMLDynamics::drift( int numberOfAtoms, RealOpenMM** atomCoordinates,
+                                   RealOpenMM** velocities,
+                                   RealOpenMM** forces, RealOpenMM* masses, RealOpenMM* inverseMasses,
+                                   RealOpenMM** xPrime ){
+    
+    static const RealOpenMM deltaT = getTimeStep();
+    
+    for( int i = 0; i < numberOfAtoms; i++ ){
+        for( int j = 0; j < 3; j++ ) {
+        
+            xPrime[i][j] = atomCoordinates[i][j] + deltaT * velocities[i][j];
+            
+        }
+    }
+    
+
+}
+    
+/**---------------------------------------------------------------------------------------
+
+   Update -- driver routine for performing stochastic dynamics update of coordinates
+   and velocities
+
+   @param numberOfAtoms       number of atoms
+   @param atomCoordinates     atom coordinates
+   @param velocities          velocities
+   @param forces              forces
+   @param masses              atom masses
+
+   @return ReferenceDynamics::DefaultReturn
+
+   --------------------------------------------------------------------------------------- */
+
+int ReferenceNMLDynamics::update( int numberOfAtoms, RealOpenMM** atomCoordinates,
+                                          RealOpenMM** velocities,
+                                          RealOpenMM** forces, RealOpenMM* masses, 
+                                          const RealOpenMM currentPE, const int stepType ){
+
+   // ---------------------------------------------------------------------------------------
+
+   static const char* methodName      = "\nReferenceNMLDynamics::update";
+
+   static const RealOpenMM zero       =  0.0;
+   static const RealOpenMM one        =  1.0;
+
+   static int debug                   =  0;
+
+   // ---------------------------------------------------------------------------------------
+
+   // get work arrays
+
+   RealOpenMM** xPrime          = get2DArrayAtIndex( xPrime2D );
+   RealOpenMM** oldVelocities   = get2DArrayAtIndex( OldV );
+   RealOpenMM** xVector         = get2DArrayAtIndex( X2D );
+   RealOpenMM** vVector         = get2DArrayAtIndex( V2D );
+
+   RealOpenMM* inverseMasses    = get1DArrayAtIndex( InverseMasses );
+
+   // first-time-through initialization
+
+   if( getTimeStep() == 0 ){
+
+      std::stringstream message;
+      message << methodName;
+      int errors = 0;
+
+      // invert masses
+
+      for( int ii = 0; ii < numberOfAtoms; ii++ ){
+         if( masses[ii] <= zero ){
+            message << "mass at atom index=" << ii << " (" << masses[ii] << ") is <= 0" << std::endl;
+            errors++;
+         } else {
+            inverseMasses[ii] = one/masses[ii];
+         }
+      }
+
+      // set xVector 
+
+      for( int ii = 0; ii < numberOfAtoms; ii++ ){
+         //RealOpenMM sqrtInverseMass = SQRT( inverseMasses[ii] )*fixedParameters[X];
+         for( int jj = 0; jj < 3; jj++ ){
+            xVector[ii][jj] = SimTKOpenMMUtilities::getNormallyDistributedRandomNumber() * SQRT( masses[ii] );
+         }
+      }
+
+      subspaceProjection(xVector, xVector, numberOfAtoms, masses, inverseMasses, 0);
+       
+       for( int ii = 0; ii < numberOfAtoms; ii++ ){
+           //RealOpenMM sqrtInverseMass = SQRT( inverseMasses[ii] )*fixedParameters[X];
+           for( int jj = 0; jj < 3; jj++ ){
+               xVector[ii][jj] *= inverseMasses[ii] ;
+           }
+       }
+       
+      
+      // exit if errors
+
+      if( errors ){
+         SimTKOpenMMLog::printError( message );
+      }
+   }
+
+    std::stringstream messagea;
+    messagea << "Tau " << _tau << "\n";
+    //messagea << "Eigs " << _projectionVectors[0] << "," << _projectionVectors[1] << "," << _projectionVectors[2] << "," << _projectionVectors[3] << ","  <<"\n";
+    //messagea << "Num vect " << _numProjectionVectors << "\n";
+    SimTKOpenMMLog::printMessage( messagea );
+    
+    switch(stepType){
+        case 1: 
+            halfKick( numberOfAtoms, atomCoordinates,
+                           velocities,
+                           forces, masses, inverseMasses, xVector );
+            
+            drift( numberOfAtoms, atomCoordinates,
+                     velocities,
+                     forces, masses, inverseMasses,
+                     xPrime );
+            
+            break;
+            
+        case 2: 
+            halfKick( numberOfAtoms, atomCoordinates,
+                                  velocities,
+                                  forces, masses, inverseMasses, xVector );
+            break;
+            
+        case 3:{
+            //save current PE in case quadratic required
+            lastPE = currentPE;
+            //project forces into complement space, put in xPrime
+            subspaceProjection(forces, xPrime, numberOfAtoms, masses, inverseMasses, 2);
+
+            //just equal to minus dot product of 'proposed position move' and forces (=-\nabla PE)
+            lastSlope = 0.0;
+            for( int ii = 0; ii < numberOfAtoms; ii++ ){
+                for( int jj = 0; jj < 3; jj++ ){
+                    lastSlope -= xPrime[ii][jj] * forces[ii][jj] * inverseMasses[ii]; //posTemp[k].dot((*myForcesP)[k]);
+                }
+            }
+            
+            std::stringstream message1;
+            message1 << "Last slope " << lastSlope  <<"\n";
+            
+            SimTKOpenMMLog::printMessage( message1 );
+            
+            // copy xPrime -> atomCoordinates
+            for( int ii = 0; ii < numberOfAtoms; ii++ ){
+                const RealOpenMM factor = inverseMasses[ii] / _maxEig;
+    
+                xPrime[ii][0] = atomCoordinates[ii][0] + factor * xPrime[ii][0];
+                xPrime[ii][1] = atomCoordinates[ii][1] + factor * xPrime[ii][1];
+                xPrime[ii][2] = atomCoordinates[ii][2] + factor * xPrime[ii][2];
+            }
+            break;
+        }
+        case 4:
+            //project forces into complement space, put in xPrime
+            subspaceProjection(forces, xPrime, numberOfAtoms, masses, inverseMasses, 2);
+            
+            //find slope of PE with /lambda here
+            RealOpenMM lambda = 1.0 / _maxEig;
+            RealOpenMM oldLambda = lambda;
+            std::stringstream message;
+            message << "Lambda "<< lambda << " oldLambda " << oldLambda  <<"\n";
+            
+            SimTKOpenMMLog::printMessage( message );
+            
+            //just equal to minus dot product of 'proposed position move' and forces (=-\nabla PE)
+            //RealOpenMM lambdaSlp = 0.0;
+            //for( int ii = 0; ii < numberOfAtoms; ii++ ){
+            //    for( int jj = 0; jj < 3; jj++ ){
+            //        lambdaSlp -= xPrime[ii][jj] * forces[ii][jj] * inverseMasses[ii]; //posTemp[k].dot((*myForcesP)[k]);
+            //    }
+            //}
+            
+            //solve for minimum for quadratic fit using two PE vales and the slope with /lambda=0
+            RealOpenMM a, b;
+            
+            //a = -(((currentPE - lastPE) / lambda - lambdaSlp) / lambda);
+            //b = lambdaSlp - 2.0 * a * lambda;
+            a = (((currentPE - lastPE) / lambda - lastSlope) / lambda);
+            b = lastSlope;
+            
+            if( a != 0.0){
+                lambda = -b / (2 * a);
+            }else{
+                lambda = oldLambda / 2.0;
+            }
+            
+            //test if lambda negative, if so just use smaller lambda
+            if(lambda <= 0.0) lambda = oldLambda / 2.0;
+            
+            //std::stringstream message;
+            //message << "Lambda "<< lambda << " oldLambda " << oldLambda << " slope " << lambdaSlp << " lastSlope " << lastSlope<<"\n";
+            message << "Lambda "<< lambda << " oldLambda " << oldLambda << " lastSlope " << lastSlope<<"\n";
+            
+            SimTKOpenMMLog::printMessage( message );
+            
+            // copy xPrime -> atomCoordinates
+            for( int ii = 0; ii < numberOfAtoms; ii++ ){
+                const RealOpenMM factor = inverseMasses[ii] * (lambda-oldLambda);
+                
+                xPrime[ii][0] = atomCoordinates[ii][0] + factor * xPrime[ii][0];
+                xPrime[ii][1] = atomCoordinates[ii][1] + factor * xPrime[ii][1];
+                xPrime[ii][2] = atomCoordinates[ii][2] + factor * xPrime[ii][2];
+            }
+            break;
+    }
+    
+  //New NML minimizer
+  
+  std::stringstream message;
+    message << "Array     "<<forces[0][0]<<" " <<forces[0][1]<<" " <<forces[0][2]<<" "<<xPrime[1][0]<<" " <<xPrime[1][1]<<" " <<xPrime[1][2]<<" " <<"\n";
+    RealOpenMM* mydata = &forces[0][0];
+    message << "Lin Array "<<mydata[0]<<" " <<mydata[1]<<" " <<mydata[2]<<" "<<mydata[3]<<" " <<mydata[4]<<" " <<mydata[5]<<" " <<"\n";
+    forces[0][0] = 1.0; forces[0][1] = 1.1; xPrime[1][0] = 2.0;
+    message << "Lin Array "<<mydata[0]<<" " <<mydata[1]<<" " <<mydata[2]<<" "<<mydata[3]<<" " <<mydata[4]<<" " <<mydata[5]<<" " <<"\n";
+    
+    SimTKOpenMMLog::printMessage( message );
+    
+  //End NML Minimizer
+  
+
+   if( debug ){
+      int maxPrint = 5;
+      std::stringstream message;
+      message << methodName << " Post SD2 atoms=" << numberOfAtoms << "\n";
+      RealOpenMM** oldVelocities   = get2DArrayAtIndex( OldV );
+      for( int ii = 0; ii < maxPrint; ii++ ){
+         message << " x[";
+         SimTKOpenMMUtilities::formatRealStringStream( message, atomCoordinates[ii], 3, one );
+         message << "] xp[";
+         SimTKOpenMMUtilities::formatRealStringStream( message, xPrime[ii], 3, one );
+         message << "] v[";
+         SimTKOpenMMUtilities::formatRealStringStream( message, velocities[ii], 3, one );
+         message << "] ov[";
+         SimTKOpenMMUtilities::formatRealStringStream( message, oldVelocities[ii], 3, one );
+         message << "]\n";
+      }
+      SimTKOpenMMLog::printMessage( message );
+   }
+#if 0
+  ReferenceConstraintAlgorithm* referenceConstraintAlgorithm = getReferenceConstraintAlgorithm();
+  if( referenceConstraintAlgorithm ){
+    
+    /*
+     std::stringstream message;
+     message << methodName;
+     message << " calling constrain1\n";
+     SimTKOpenMMLog::printMessage( message );
+     */
+    
+    referenceConstraintAlgorithm->apply( numberOfAtoms, atomCoordinates, xPrime,
+                                        inverseMasses );
+    
+  }
+#endif  
+   // copy xPrime -> atomCoordinates
+
+   for( int ii = 0; ii < numberOfAtoms; ii++ ){
+      atomCoordinates[ii][0] = xPrime[ii][0];
+      atomCoordinates[ii][1] = xPrime[ii][1];
+      atomCoordinates[ii][2] = xPrime[ii][2];
+   }
+
+   incrementTimeStep();
+
+   return 0; // ReferenceDynamics::DefaultReturn;
+
+}
+
+/**---------------------------------------------------------------------------------------
+
+   Write state
+
+   @param numberOfAtoms       number of atoms
+   @param atomCoordinates     atom coordinates
+   @param velocities          velocities
+   @param forces              forces
+   @param masses              atom masses
+   @param state               0 if initial state; otherwise nonzero
+   @param baseFileName        base file name
+
+   @return ReferenceDynamics::DefaultReturn
+
+   --------------------------------------------------------------------------------------- */
+
+int ReferenceNMLDynamics::writeState( int numberOfAtoms, RealOpenMM** atomCoordinates,
+                                             RealOpenMM** velocities,
+                                             RealOpenMM** forces, RealOpenMM* masses,
+                                             int state, const std::string& baseFileName ) const {
+
+   // ---------------------------------------------------------------------------------------
+
+   static const char* methodName      = "\nReferenceNMLDynamics::writeState";
+
+   static const RealOpenMM zero       =  0.0;
+   static const RealOpenMM one        =  1.0;
+   static const int threeI            =  3;
+
+   // ---------------------------------------------------------------------------------------
+
+   std::stringstream stateFileName;
+
+   stateFileName << baseFileName;
+   stateFileName << "_Step" << getTimeStep();
+   // stateFileName << "_State" << state;
+   stateFileName << ".txt";
+
+   // ---------------------------------------------------------------------------------------
+
+   // open file -- return if unsuccessful
+
+   FILE* stateFile = NULL;
+#ifdef WIN32
+   fopen_s( &stateFile, stateFileName.str().c_str(), "w" );
+#else
+   stateFile = fopen( stateFileName.str().c_str(), "w" );
+#endif
+
+   // ---------------------------------------------------------------------------------------
+
+   // diagnostics
+
+   if( stateFile != NULL ){
+      std::stringstream message;
+      message << methodName;
+      message << " Opened file=<" << stateFileName.str() << ">.\n";
+      SimTKOpenMMLog::printMessage( message );
+   } else {
+      std::stringstream message;
+      message << methodName;
+      message << " could not open file=<" << stateFileName.str() << "> -- abort output.\n";
+      SimTKOpenMMLog::printMessage( message );
+      return -1; // ReferenceDynamics::ErrorReturn;
+   }   
+
+   // ---------------------------------------------------------------------------------------
+
+   StringVector scalarNameI;
+   IntVector scalarI;
+
+   StringVector scalarNameR;
+   RealOpenMMVector scalarR;
+
+   StringVector scalarNameR1;
+   RealOpenMMPtrVector scalarR1;
+
+   StringVector scalarNameR2;
+   RealOpenMMPtrPtrVector scalarR2;
+
+   scalarI.push_back( getNumberOfAtoms() );
+   scalarNameI.push_back( "Atoms" );
+
+   scalarI.push_back( getTimeStep() );
+   scalarNameI.push_back( "Timestep" );
+
+   if( state == 0 || state == -1 ){
+
+      scalarR.push_back( getDeltaT() );
+      scalarNameR.push_back( "delta_t" );
+
+      scalarR.push_back( getTemperature() );
+      scalarNameR.push_back( "T" );
+
+      scalarR.push_back( getTau() );
+      scalarNameR.push_back( "Tau" );
+
+      scalarR1.push_back( masses );
+      scalarNameR1.push_back( "mass" );
+
+      scalarR2.push_back( atomCoordinates );
+      scalarNameR2.push_back( "coord" );
+
+      scalarR2.push_back( velocities );
+      scalarNameR2.push_back( "velocities" );
+
+      scalarR2.push_back( forces );
+      scalarNameR2.push_back( "forces" );
+
+      if( state == -1 ){
+
+         RealOpenMM** xPrime          = get2DArrayAtIndex( xPrime2D );
+         RealOpenMM** oldVelocities   = get2DArrayAtIndex( OldV );
+         RealOpenMM** xVector         = get2DArrayAtIndex( X2D );
+         RealOpenMM** vVector         = get2DArrayAtIndex( V2D );
+
+         scalarR2.push_back( xPrime );
+         scalarNameR2.push_back( "xPrime" );
+
+         scalarR2.push_back( oldVelocities);
+         scalarNameR2.push_back( "vold" );
+
+         scalarR2.push_back( xVector );
+         scalarNameR2.push_back( "xVector" );
+
+         scalarR2.push_back( vVector );
+         scalarNameR2.push_back( "vVector" );
+      }
+      
+   } else {
+
+      scalarR2.push_back( atomCoordinates );
+      scalarNameR2.push_back( "coord" );
+
+      scalarR2.push_back( velocities );
+      scalarNameR2.push_back( "velocities" );
+
+   }
+
+   writeStateToFile( stateFile, scalarNameI, scalarI, scalarNameR, scalarR, getNumberOfAtoms(), scalarNameR1, scalarR1, threeI, scalarNameR2, scalarR2 ); 
+
+   (void) fclose( stateFile );
+
+   return 0; // ReferenceDynamics::DefaultReturn;
+
+}
+
+//~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
+// Find forces OR positions inside subspace (defined as the span of the 'eigenvectors' Q)
+// Take 'array' as input, 'outArray' as output (may be the same vector).
+// 'projectionMode' determines:
+// a) if the projection is for force or positions, bit 1 set=force.
+// b) if the projection is for the sub-space or-complement space, bit 2 set=sub-space.
+//~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
+void ReferenceNMLDynamics::subspaceProjection(  RealOpenMM** arrayParam, 
+                                                RealOpenMM** outArrayParam, 
+                                                int numberOfAtoms,
+                                                RealOpenMM* masses,
+                                                RealOpenMM* inverseMasses,
+                                                int projectionMode){
+  
+  //move to linear arrays for Blas etc.
+  RealOpenMM* array = &arrayParam[0][0];
+  RealOpenMM* outArray = &outArrayParam[0][0];
+  
+  //If 'array' and 'outArray are not the same array
+  //copy 'array' into outArray
+  const unsigned int _3N = numberOfAtoms * 3;
+  if (array != outArray) {
+    for (unsigned int i=0; i<_3N; i++)
+      outArray[i] = array[i];
+  }
+  
+  //Temporary array for the 'mode' values, denoted 'c'
+  //the size is equal to the number of eigenvectors
+  RealOpenMM* tmpC = new RealOpenMM[_numProjectionVectors];
+  
+  //~~~~We need to calculate M^{1/2}QQ^TM^{-1/2}force or M^{-1/2}QQ^TM^{1/2}positions~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
+  
+  //First we weight the array by multiplying by 
+  //the square-root of the atomic masses for positions
+  //or the the inverse of the square-root of the atomic masses for forces.
+  //
+  //a'=M^{-1/2}*f for forces, OR a'= M^{1/2}*x for positions
+  //
+  if(!(projectionMode & 0x01)) { //Forces if zero so weight is sqrt(1/m)
+    
+    //array is 3N long, for N atoms
+    for( int i=0; i < numberOfAtoms; i++) {           //N loops
+      RealOpenMM  massWt = SQRT( inverseMasses[i] );
+      for (unsigned int j=0; j<3; j++)                //times 3 loops
+        outArray[i*3 + j] *= massWt;
+    }
+    
+  } else { //positions if first bit set, so weight is sqrt(m)
+    
+    //array is 3N long, for N atoms
+    for( int i=0; i < numberOfAtoms; i++) {           //N loops
+      RealOpenMM  massWt = SQRT( masses[i] );
+      for (unsigned int j=0; j<3; j++)                //times 3 loops
+        outArray[i*3 + j] *= massWt;
+    }
+    
+  }
+  
+  //Project onto mode space by taking the matrix product of 
+  //the transpose of the eigenvectors Q with the array.
+  //
+  //c=Q^T*a', a' from last algorithm step
+  //
+#if defined(HAVE_LAPACK)
+  //~~~~Blas here~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
+  char *transA = "T";							// Transpose, LAPACK checks only first character N/V
+  int m = numberOfAtoms * 3; //eigenvector vector/array length
+  int n = _numProjectionVectors;  //number of eigenvectors
+  int incxy = 1;	//sizes
+  RealOpenMM alpha = 1.0;	RealOpenMM beta = 0.0;
+  
+  //if single precision use Blas 'sgemv'
+  #if RealOpenMMType == 1 
+    //
+    sgemv_ (transA, &m, &n, &alpha, _projectionVectors, &m, outArray, &incxy, &beta, tmpC, &incxy);
+    //
+    
+    //Now find projected force/positions a'' by matrix product with Eigenvectors Q
+    //or the complementary Eigenvectors
+    //a''=Qc
+    char *transB = "N"; /* LAPACK checks only first character N/V */
+    
+    //if projectionMode has bit 2 set then complement space (I-Q^TQ) not (Q^TQ)
+    if(!(projectionMode & 0x02)) {
+      alpha = 1.0;	beta = 0.0;
+    } else {
+      alpha = -1.0;	beta = 1.0;
+    }
+    //
+    sgemv_ (transB, &m, &n, &alpha, _projectionVectors, &m, tmpC, &incxy, &beta, outArray, &incxy);
+    
+    //if double precision use Blas 'dgemv'
+  #else
+    //
+    dgemv_ (transA, &m, &n, &alpha, _projectionVectors), &m, outArray, &incxy, &beta, tmpC, &incxy);
+    //
+    
+    //Now find projected force/positions a'' by matrix product with Eigenvectors Q
+    //a''=Qc
+    char *transB = "N"; /* LAPACK checks only first character N/V */
+    
+    //if projectionMode has bit 2 set then complement space (I-Q^TQ) not (Q^TQ)
+    if(!(projectionMode & 0x02)) {
+      alpha = 1.0;	beta = 0.0;
+    } else {    
+      alpha = -1.0;	beta = 1.0;
+    }
+    //
+    dgemv_ (transB, &m, &n, &alpha, _projectionVectors, &m, tmpC, &incxy, &beta, outArray, &incxy);
+  #endif
+#else
+
+  //~~~~NO Blas here~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
+  //If no Blas is available we need to manually find the product c=A*b
+  //c_i=\sum_{j=1}^n A_{i,j} b_j
+  
+  //c=Q^T*a', a' from last algorithm step
+  //Q is a linear array in column major format
+  //so tmpC_i = \sum_{j=1}^n Q_{j,i} outArray_j
+  //Q_{j,i}=_projectionVectors[j * numberOfAtoms * 3 + i]
+  
+  //const int _3N = numberOfAtoms * 3; //eigenvector vector/array length
+  
+  for(int i=0; i < _numProjectionVectors; i++){   //over all eigenvectors
+    
+    tmpC[i] = 0.0;  //clear
+    
+    for(int j=0; j< _3N; j++){  //over each element in the vector
+      //tmpC[i] += _projectionVectors[j * _3N + i] * outArray[j];
+        tmpC[i] += _projectionVectors[j  + i * _3N] * outArray[j];
+    }
+  }
+  
+  //Now find projected force/positions a'' by matrix product with Eigenvectors Q
+  //a''=Qc
+  //so outArray_i  = \sum_{j=1}^n Q_{i,j} tmpC_i
+  
+  //if projectionMode has bit 2 set then complement space (I-Q^TQ) not (Q^TQ)
+  const bool pMode = !(projectionMode & 0x02);
+  
+  //find product
+  for(int i=0; i< _3N; i++){  //over each element in the vector
+    
+    //if sub-space do Q*c
+    //else do a'-Q(Q^T a') = (I-QQ^T)a'
+      if(pMode){
+          outArray[i] = 0.0; //if not complement
+    
+          for(int j=0; j < _numProjectionVectors; j++){   //over all eigenvectors
+              //outArray[i] += _projectionVectors[i * _3N + j] * tmpC[j];
+              outArray[i] += _projectionVectors[i + j * _3N] * tmpC[j];
+          }
+      }else{
+          for(int j=0; j < _numProjectionVectors; j++){   //over all eigenvectors
+              //outArray[i] -= _projectionVectors[i * _3N + j] * tmpC[j];
+              outArray[i] -= _projectionVectors[i + j * _3N] * tmpC[j];
+          }
+          
+      }
+    
+  }
+  
+  //~~~~End of no-blas~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
+#endif
+    
+  //Finally we weight the array by multiplying by 
+  //the inverse of the square-root of the atomic masses for positions
+  //or the the square-root of the atomic masses for forces.
+  //
+  //a'''=M^{1/2}*a'' or a'''=M^{-1/2}*a''
+  //
+  if(!(projectionMode & 0x01)) { //Forces if zero so weight is sqrt(m)
+    for( int i=0; i < numberOfAtoms; i++) {
+      RealOpenMM  massUnWt = SQRT( masses[i] );
+      for (unsigned int j=0; j<3; j++)
+        outArray[i*3 + j] *= massUnWt;
+    }
+  } else {  //positions if first bit set, so weight is sqrt(1/m)
+    for( int i=0; i < numberOfAtoms; i++) {
+      RealOpenMM  massUnWt = SQRT( inverseMasses[i] );
+      for (unsigned int j=0; j<3; j++)
+        outArray[i*3 + j] *= massUnWt;
+    }
+  }
+  
+  //Clean up
+  
+  //delete temporary array
+  delete [] tmpC;
+  
+}

plugins/normalModeLangevin/platforms/reference/src/ReferenceNMLKernelFactory.cpp

+/* -------------------------------------------------------------------------- *
+ *                                   OpenMM                                   *
+ * -------------------------------------------------------------------------- *
+ * This is part of the OpenMM molecular simulation toolkit originating from   *
+ * Simbios, the NIH National Center for Physics-Based Simulation of           *
+ * Biological Structures at Stanford, funded under the NIH Roadmap for        *
+ * Medical Research, grant U54 GM072970. See https://simtk.org.               *
+ *                                                                            *
+ * Portions copyright (c) 2009 Stanford University and the Authors.           *
+ * Authors: Chris Sweet                                                       *
+ * Contributors: Christopher Bruns                                            *
+ *                                                                            *
+ * Permission is hereby granted, free of charge, to any person obtaining a    *
+ * copy of this software and associated documentation files (the "Software"), *
+ * to deal in the Software without restriction, including without limitation  *
+ * the rights to use, copy, modify, merge, publish, distribute, sublicense,   *
+ * and/or sell copies of the Software, and to permit persons to whom the      *
+ * Software is furnished to do so, subject to the following conditions:       *
+ *                                                                            *
+ * The above copyright notice and this permission notice shall be included in *
+ * all copies or substantial portions of the Software.                        *
+ *                                                                            *
+ * THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR *
+ * IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,   *
+ * FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL    *
+ * THE AUTHORS, CONTRIBUTORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM,    *
+ * DAMAGES OR OTHER LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR      *
+ * OTHERWISE, ARISING FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE  *
+ * USE OR OTHER DEALINGS IN THE SOFTWARE.                                     *
+ * -------------------------------------------------------------------------- */
+
+#include "openmm/OpenMMException.h"
+#include "openmm/internal/ContextImpl.h"
+#include "ReferenceNMLKernelFactory.h"
+#include "ReferenceIntegrateNMLStepKernel.h"
+
+namespace OpenMM {
+
+KernelImpl* ReferenceNMLKernelFactory::createKernelImpl(std::string name, const Platform& platform, ContextImpl& context) const {
+    ReferencePlatform::PlatformData& data = *static_cast<ReferencePlatform::PlatformData*>(context.getPlatformData());
+    if (name == IntegrateNMLStepKernel::Name())
+        return new ReferenceIntegrateNMLStepKernel(name, platform, data);
+    throw OpenMMException((std::string("Tried to create kernel with illegal kernel name '")+name+"'").c_str());
+}
+
+} /* namespace OpenMM */

plugins/normalModeLangevin/src/NMLIntegrator.cpp

+/* -------------------------------------------------------------------------- *
+ *                                   OpenMM                                   *
+ * -------------------------------------------------------------------------- *
+ * This is part of the OpenMM molecular simulation toolkit originating from   *
+ * Simbios, the NIH National Center for Physics-Based Simulation of           *
+ * Biological Structures at Stanford, funded under the NIH Roadmap for        *
+ * Medical Research, grant U54 GM072970. See https://simtk.org.               *
+ *                                                                            *
+ * Portions copyright (c) 2009 Stanford University and the Authors.           *
+ * Authors: Chris Sweet                                                       *
+ * Contributors: Christopher Bruns                                            *
+ *                                                                            *
+ * Permission is hereby granted, free of charge, to any person obtaining a    *
+ * copy of this software and associated documentation files (the "Software"), *
+ * to deal in the Software without restriction, including without limitation  *
+ * the rights to use, copy, modify, merge, publish, distribute, sublicense,   *
+ * and/or sell copies of the Software, and to permit persons to whom the      *
+ * Software is furnished to do so, subject to the following conditions:       *
+ *                                                                            *
+ * The above copyright notice and this permission notice shall be included in *
+ * all copies or substantial portions of the Software.                        *
+ *                                                                            *
+ * THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR *
+ * IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,   *
+ * FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL    *
+ * THE AUTHORS, CONTRIBUTORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM,    *
+ * DAMAGES OR OTHER LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR      *
+ * OTHERWISE, ARISING FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE  *
+ * USE OR OTHER DEALINGS IN THE SOFTWARE.                                     *
+ * -------------------------------------------------------------------------- */
+
+#include "NMLIntegrator.h"
+#include "IntegrateNMLStepKernel.h"
+#include "openmm/Context.h"
+#include "openmm/internal/ContextImpl.h"
+#include "openmm/kernels.h"
+#include <string>
+
+using namespace OpenMM;
+using std::string;
+using std::vector;
+
+NMLIntegrator::NMLIntegrator(double temperature, double frictionCoeff, double stepSize, ProtoMol::EigenvectorInfo* projectionVectorInfo ) {
+    setTemperature(temperature);
+    setFriction(frictionCoeff);
+    setStepSize(stepSize);
+    setProjectionVectorInfo(projectionVectorInfo);
+    setConstraintTolerance(1e-4);
+}
+
+void NMLIntegrator::initialize(ContextImpl& contextRef) {
+    context = &contextRef;
+    kernel = context->getPlatform().createKernel(IntegrateNMLStepKernel::Name(), contextRef);
+    dynamic_cast<IntegrateNMLStepKernel&>(kernel.getImpl()).initialize(contextRef.getSystem(), *this);
+}
+
+vector<string> NMLIntegrator::getKernelNames() {
+    std::vector<std::string> names;
+    names.push_back(IntegrateNMLStepKernel::Name());
+    return names;
+}
+
+void NMLIntegrator::step(int steps) {
+    
+
+    std::cout << "One loop " << std::endl;
+    
+    //updates per loop. Need it here? ####Moved within loop. Single v multiple steps error if not!
+    //context->updateContextState();
+    //context->calcForces();
+    
+    //loop
+    for (int i = 0; i < steps; ++i) {
+
+        //updates per loop. Need it here?
+        context->updateContextState();
+        context->calcForces();
+
+        //do half kick and full drift
+        dynamic_cast<IntegrateNMLStepKernel&>(kernel.getImpl()).execute(*context, *this, 0.0, 1 );    //stepType 1 is halfKick/drift LangevinLeapfrog
+        //in case projection vectors changed, clear flag
+        setProjVecChanged(false);
+
+        //updates per loop, ####done in minimizer
+        //context->updateContextState();
+        //context->calcForces();
+        
+        //minimize compliment space
+        minimize(50);
+
+        //do half kick
+        dynamic_cast<IntegrateNMLStepKernel&>(kernel.getImpl()).execute(*context, *this, 0.0, 2 );    //stepType 1 is halfKick LangevinLeapfrog
+              
+    }
+    
+    //update time
+    context->setTime(context->getTime()+getStepSize() * steps);
+    
+}
+
+void NMLIntegrator::minimize(int maxsteps) {
+    
+    //minimum limit
+    const double minlim = getMinimumLimit();
+    
+    //loop
+    for (int i = 0; i < maxsteps; ++i) {
+        
+        //updates per loop
+        context->updateContextState();
+        context->calcForces();
+        
+        //get initial PE
+        const double initialPE = context->calcPotentialEnergy();
+        
+        //minimize
+        dynamic_cast<IntegrateNMLStepKernel&>(kernel.getImpl()).execute(*context, *this, initialPE, 3 );    //stepType 3 is simple minimizer
+        setProjVecChanged(false);
+        
+        //get PE difference
+        const double currentPE = context->calcPotentialEnergy();
+        const double peDifference =  initialPE - currentPE;
+        
+        //std::cout << "Loop " << i << ", Current " << currentPE << ", Old " << initialPE << ", lim " << minlim << std::endl;
+        //end condition met?
+        if(peDifference < minlim && peDifference >= 0.0){
+            //std::cout << "Breaking at " << i << ", diff " << peDifference << std::endl;
+            break;
+        }
+        
+        //Require quadratic solution
+        if(peDifference < 0.0){
+            //std::cout << "Quadratic at " << i << std::endl;
+            
+            //update for new positions
+            context->updateContextState();
+            context->calcForces();
+            
+            //minimize, uses quadratic soultion as 'stepsize' not forced
+            dynamic_cast<IntegrateNMLStepKernel&>(kernel.getImpl()).execute(*context, *this, currentPE, 4 ); //stepType 4 is quadratic minimizer   
+            setProjVecChanged(false);
+            
+            std::cout << "Quadratic " << i << ", Current " << context->calcPotentialEnergy() << ", Old " << initialPE << ", lim " << minlim << std::endl;
+            
+            const double quadraticPE = context->calcPotentialEnergy();
+            
+            if(quadraticPE > initialPE){
+                break;
+            }
+            
+        }
+        
+        
+    }
+    
+}

plugins/normalModeLangevin/src/NormalModeLangevin.cpp

+/* -------------------------------------------------------------------------- *
+ *                                   OpenMM                                   *
+ * -------------------------------------------------------------------------- *
+ * This is part of the OpenMM molecular simulation toolkit originating from   *
+ * Simbios, the NIH National Center for Physics-Based Simulation of           *
+ * Biological Structures at Stanford, funded under the NIH Roadmap for        *
+ * Medical Research, grant U54 GM072970. See https://simtk.org.               *
+ *                                                                            *
+ * Portions copyright (c) 2009 Stanford University and the Authors.           *
+ * Authors: Chris Sweet                                                       *
+ * Contributors: Christopher Bruns                                            *
+ *                                                                            *
+ * Permission is hereby granted, free of charge, to any person obtaining a    *
+ * copy of this software and associated documentation files (the "Software"), *
+ * to deal in the Software without restriction, including without limitation  *
+ * the rights to use, copy, modify, merge, publish, distribute, sublicense,   *
+ * and/or sell copies of the Software, and to permit persons to whom the      *
+ * Software is furnished to do so, subject to the following conditions:       *
+ *                                                                            *
+ * The above copyright notice and this permission notice shall be included in *
+ * all copies or substantial portions of the Software.                        *
+ *                                                                            *
+ * THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR *
+ * IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,   *
+ * FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL    *
+ * THE AUTHORS, CONTRIBUTORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM,    *
+ * DAMAGES OR OTHER LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR      *
+ * OTHERWISE, ARISING FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE  *
+ * USE OR OTHER DEALINGS IN THE SOFTWARE.                                     *
+ * -------------------------------------------------------------------------- */
+
+#include <iostream>
+#include "OpenMM.h"
+#include "ReferenceNMLKernelFactory.h"
+
+using namespace OpenMM;
+using namespace std;
+
+extern "C" void initOpenMMPlugin() {
+    cout << "Initializing Normal Mode Langevin OpenMM plugin..." << endl;
+    for (int p = 0; p < Platform::getNumPlatforms(); ++p) {
+        Platform& platform = Platform::getPlatform(p);
+        if (platform.getName() == "CudaPlatform")
+            // We don't have a Cuda one yet, but we will soon
+            platform.registerKernelFactory("ReferenceNMLKernelFactory", new ReferenceNMLKernelFactory());
+        else
+            platform.registerKernelFactory("ReferenceNMLKernelFactory", new ReferenceNMLKernelFactory());
+    }
+}

plugins/normalModeLangevin/test/CMakeLists.txt

+
+add_executable(TestNormalModeLangevin TestNormalModeLangevin.cpp)
+target_link_libraries(TestNormalModeLangevin OpenMM)
+add_test(TestNormalModeLangevin ${EXECUTABLE_OUTPUT_PATH}/TestNormalModeLangevin)
+