Merge pull request #837 from peastman/rpmd

Created a barostat for use with RPMD

plugins/rpmd/openmmapi/include/openmm/RPMDMonteCarloBarostat.h

+#ifndef OPENMM_RPMDMONTECARLOBAROSTAT_H_
+#define OPENMM_RPMDMONTECARLOBAROSTAT_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) 2010-2015 Stanford University and the Authors.      *
+ * Authors: Peter Eastman                                                     *
+ * 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 "openmm/Force.h"
+#include <string>
+#include "internal/windowsExportRpmd.h"
+
+namespace OpenMM {
+
+/**
+ * This class is very similar to MonteCarloBarostat, but it is specifically designed for use
+ * with RPMDIntegrator.  For each trial move, it scales all copies of the system by the same
+ * amount, then accepts or rejects the move based on the change to the total energy of the
+ * ring polymer (as returned by the integrator's getTotalEnergy() method).
+ */
+
+class OPENMM_EXPORT_RPMD RPMDMonteCarloBarostat : public Force {
+public:
+    /**
+     * This is the name of the parameter which stores the current pressure acting on
+     * the system (in bar).
+     */
+    static const std::string& Pressure() {
+        static const std::string key = "RPMDMonteCarloPressure";
+        return key;
+    }
+    /**
+     * Create a MonteCarloBarostat.
+     *
+     * @param defaultPressure   the default pressure acting on the system (in bar)
+     * @param frequency         the frequency at which Monte Carlo pressure changes should be attempted (in time steps)
+     */
+    RPMDMonteCarloBarostat(double defaultPressure, int frequency = 25);
+    /**
+     * Get the default pressure acting on the system (in bar).
+     *
+     * @return the default pressure acting on the system, measured in bar.
+     */
+    double getDefaultPressure() const {
+        return defaultPressure;
+    }
+    /**
+     * Set the default pressure acting on the system.  This will affect any new Contexts you create,
+     * but not ones that already exist.
+     *
+     * @param pressure   the default pressure acting on the system, measured in bar.
+     */
+    void setDefaultPressure(double pressure) {
+        defaultPressure = pressure;
+    }
+    /**
+     * Get the frequency (in time steps) at which Monte Carlo pressure changes should be attempted.  If this is set to
+     * 0, the barostat is disabled.
+     */
+    int getFrequency() const {
+        return frequency;
+    }
+    /**
+     * Set the frequency (in time steps) at which Monte Carlo pressure changes should be attempted.  If this is set to
+     * 0, the barostat is disabled.
+     */
+    void setFrequency(int freq) {
+        frequency = freq;
+    }
+    /**
+     * Get the random number seed.  See setRandomNumberSeed() for details.
+     */
+    int getRandomNumberSeed() const {
+        return randomNumberSeed;
+    }
+    /**
+     * Set the random number seed.  It is guaranteed that if two simulations are run
+     * with different random number seeds, the sequence of Monte Carlo steps will be different.  On
+     * the other hand, no guarantees are made about the behavior of simulations that use the same seed.
+     * In particular, Platforms are permitted to use non-deterministic algorithms which produce different
+     * results on successive runs, even if those runs were initialized identically.
+     *
+     * If seed is set to 0 (which is the default value assigned), a unique seed is chosen when a Context
+     * is created from this Force. This is done to ensure that each Context receives unique random seeds
+     * without you needing to set them explicitly.
+     */
+    void setRandomNumberSeed(int seed) {
+        randomNumberSeed = seed;
+    }
+    /**
+     * Returns whether or not this force makes use of periodic boundary
+     * conditions.
+     *
+     * @returns true if force uses PBC and false otherwise
+     */
+    bool usesPeriodicBoundaryConditions() const {
+        return true;
+    }
+protected:
+    ForceImpl* createImpl() const;
+private:
+    double defaultPressure;
+    int frequency, randomNumberSeed;
+};
+
+} // namespace OpenMM
+
+#endif /*OPENMM_RPMDMONTECARLOBAROSTAT_H_*/

plugins/rpmd/openmmapi/include/openmm/RPMDUpdater.h

+#ifndef OPENMM_RPMDUPDATER_H_
+#define OPENMM_RPMDUPDATER_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) 2015 Stanford University and the Authors.           *
+ * Authors: Peter Eastman                                                     *
+ * 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 "openmm/internal/ForceImpl.h"
+#include "internal/windowsExportRpmd.h"
+
+namespace OpenMM {
+
+/**
+ * This ForceImpl subclass is for forces specifically designed to work with RPMDIntegrator.
+ * It adds an updateRPMDState() method that is invoked at the start of each time step
+ * (in contrast to updateContextState(), which gets invoked many times per time step, once
+ * for each copy).
+ */
+
+class OPENMM_EXPORT_RPMD RPMDUpdater : public ForceImpl {
+public:
+    virtual void updateRPMDState(ContextImpl& context) = 0;
+};
+
+} // namespace OpenMM
+
+#endif /*OPENMM_RPMDUPDATER_H_*/

plugins/rpmd/openmmapi/include/openmm/internal/RPMDMonteCarloBarostatImpl.h

+#ifndef OPENMM_RPMDMONTECARLOBAROSTATIMPL_H_
+#define OPENMM_RPMDMONTECARLOBAROSTATIMPL_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) 2010-2015 Stanford University and the Authors.      *
+ * Authors: Peter Eastman                                                     *
+ * 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 "openmm/RPMDMonteCarloBarostat.h"
+#include "openmm/RPMDUpdater.h"
+#include "openmm/Kernel.h"
+#include "openmm/Vec3.h"
+#include "sfmt/SFMT.h"
+#include <string>
+#include <vector>
+
+namespace OpenMM {
+
+/**
+ * This is the internal implementation of RPMDMonteCarloBarostat.
+ */
+
+class RPMDMonteCarloBarostatImpl : public RPMDUpdater {
+public:
+    RPMDMonteCarloBarostatImpl(const RPMDMonteCarloBarostat& owner);
+    void initialize(ContextImpl& context);
+    const RPMDMonteCarloBarostat& getOwner() const {
+        return owner;
+    }
+    void updateRPMDState(ContextImpl& context);
+    void updateContextState(ContextImpl& context) {
+        // This is unused, since the updating is done in updateRPMDState().
+    }
+    double calcForcesAndEnergy(ContextImpl& context, bool includeForces, bool includeEnergy, int groups) {
+        // This force doesn't apply forces to particles.
+        return 0.0;
+    }
+    std::map<std::string, double> getDefaultParameters();
+    std::vector<std::string> getKernelNames();
+private:
+    const RPMDMonteCarloBarostat& owner;
+    int step, numAttempted, numAccepted;
+    double volumeScale;
+    OpenMM_SFMT::SFMT random;
+    std::vector<std::vector<Vec3> > savedPositions;
+    Kernel kernel;
+};
+
+} // namespace OpenMM
+
+#endif /*OPENMM_RPMDMONTECARLOBAROSTATIMPL_H_*/

plugins/rpmd/openmmapi/src/RPMDIntegrator.cpp

@@ -34,6 +34,7 @@
 #include "openmm/OpenMMException.h"
 #include "openmm/internal/ContextImpl.h"
 #include "openmm/RpmdKernels.h"
+#include "openmm/RPMDUpdater.h"
 #include "SimTKOpenMMRealType.h"
 #include <cmath>
 #include <string>
@@ -190,6 +191,12 @@ void RPMDIntegrator::step(int steps) {
         context->getOwner().setPositions(p);
         isFirstStep = false;
     }
+    vector<ForceImpl*>& forceImpls = context->getForceImpls();
+    for (int i = 0; i < (int) forceImpls.size(); i++) {
+        RPMDUpdater* updater = dynamic_cast<RPMDUpdater*>(forceImpls[i]);
+        if (updater != NULL)
+            updater->updateRPMDState(*context);
+    }
     for (int i = 0; i < steps; ++i) {
         kernel.getAs<IntegrateRPMDStepKernel>().execute(*context, *this, forcesAreValid);
         forcesAreValid = true;

plugins/rpmd/openmmapi/src/RPMDMonteCarloBarostat.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) 2010-2015 Stanford University and the Authors.      *
+ * Authors: Peter Eastman                                                     *
+ * 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 "openmm/RPMDMonteCarloBarostat.h"
+#include "openmm/internal/RPMDMonteCarloBarostatImpl.h"
+
+using namespace OpenMM;
+
+RPMDMonteCarloBarostat::RPMDMonteCarloBarostat(double defaultPressure, int frequency) :
+        defaultPressure(defaultPressure), frequency(frequency) {
+    setRandomNumberSeed(0);
+}
+
+ForceImpl* RPMDMonteCarloBarostat::createImpl() const {
+    return new RPMDMonteCarloBarostatImpl(*this);
+}

plugins/rpmd/openmmapi/src/RPMDMonteCarloBarostatImpl.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) 2010-2015 Stanford University and the Authors.      *
+ * Authors: Peter Eastman                                                     *
+ * 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 "openmm/internal/RPMDMonteCarloBarostatImpl.h"
+#include "openmm/internal/ContextImpl.h"
+#include "openmm/internal/OSRngSeed.h"
+#include "openmm/Context.h"
+#include "openmm/kernels.h"
+#include "openmm/OpenMMException.h"
+#include "openmm/RPMDIntegrator.h"
+#include <cmath>
+#include <vector>
+#include <algorithm>
+
+using namespace OpenMM;
+using namespace OpenMM_SFMT;
+using std::vector;
+
+const float BOLTZMANN = 1.380658e-23f; // (J/K)
+const float AVOGADRO = 6.0221367e23f;
+const float RGAS = BOLTZMANN*AVOGADRO; // (J/(mol K))
+const float BOLTZ = RGAS/1000;         // (kJ/(mol K))
+
+RPMDMonteCarloBarostatImpl::RPMDMonteCarloBarostatImpl(const RPMDMonteCarloBarostat& owner) : owner(owner), step(0) {
+}
+
+void RPMDMonteCarloBarostatImpl::initialize(ContextImpl& context) {
+    RPMDIntegrator* integrator = dynamic_cast<RPMDIntegrator*>(&context.getIntegrator());
+    if (integrator == NULL)
+        throw OpenMMException("RPMDMonteCarloBarostat must be used with an RPMDIntegrator");;
+    if (!integrator->getApplyThermostat())
+        throw OpenMMException("RPMDMonteCarloBarostat requires the integrator's thermostat to be enabled");;
+    kernel = context.getPlatform().createKernel(ApplyMonteCarloBarostatKernel::Name(), context);
+    kernel.getAs<ApplyMonteCarloBarostatKernel>().initialize(context.getSystem(), owner);
+    savedPositions.resize(integrator->getNumCopies());
+    Vec3 box[3];
+    context.getPeriodicBoxVectors(box[0], box[1], box[2]);
+    double volume = box[0][0]*box[1][1]*box[2][2];
+    volumeScale = 0.01*volume;
+    numAttempted = 0;
+    numAccepted = 0;
+    int randSeed = owner.getRandomNumberSeed();
+    // A random seed of 0 means use a unique one
+    if (randSeed == 0) randSeed = osrngseed();
+    init_gen_rand(randSeed, random);
+}
+
+void RPMDMonteCarloBarostatImpl::updateRPMDState(ContextImpl& context) {
+    if (++step < owner.getFrequency() || owner.getFrequency() == 0)
+        return;
+    step = 0;
+
+    // Compute the current potential energy.
+
+    RPMDIntegrator& integrator = dynamic_cast<RPMDIntegrator&>(context.getIntegrator());
+
+    // Record the initial positions and energy
+
+    double initialEnergy = 0;
+    int numCopies = integrator.getNumCopies();
+    for (int i = 0; i < numCopies; i++) {
+        State state = integrator.getState(i, State::Positions | State::Energy);
+        savedPositions[i] = state.getPositions();
+        initialEnergy += state.getPotentialEnergy();
+    }
+
+    // Compute the centroid.
+
+    int numParticles = context.getSystem().getNumParticles();
+    vector<Vec3> centroid(numParticles, Vec3());
+    for (int i = 0; i < numParticles; i++) {
+        for (int j = 0; j < numCopies; j++)
+            centroid[i] += savedPositions[j][i];
+        centroid[i] *= 1.0/numCopies;
+    }
+
+    // Modify the periodic box size and scale the coordinates of the centroid.
+
+    Vec3 box[3];
+    context.getPeriodicBoxVectors(box[0], box[1], box[2]);
+    double volume = box[0][0]*box[1][1]*box[2][2];
+    double deltaVolume = volumeScale*2*(genrand_real2(random)-0.5);
+    double newVolume = volume+deltaVolume;
+    double lengthScale = std::pow(newVolume/volume, 1.0/3.0);
+    context.setPositions(centroid);
+    kernel.getAs<ApplyMonteCarloBarostatKernel>().scaleCoordinates(context, lengthScale, lengthScale, lengthScale);
+    context.getOwner().setPeriodicBoxVectors(box[0]*lengthScale, box[1]*lengthScale, box[2]*lengthScale);
+    State scaledState = context.getOwner().getState(State::Positions);
+
+    // Now apply the same offset to all the copies.
+
+    vector<Vec3> delta(numParticles);
+    for (int i = 0; i < numParticles; i++)
+        delta[i] = scaledState.getPositions()[i]-centroid[i];
+    double finalEnergy = 0;
+    vector<Vec3> positions(numParticles);
+    for (int copy = 0; copy < numCopies; copy++) {
+        for (int i = 0; i < numParticles; i++)
+            positions[i] = savedPositions[copy][i]+delta[i];
+        integrator.setPositions(copy, positions);
+        finalEnergy += integrator.getState(copy, State::Energy).getPotentialEnergy();
+    }
+
+    // Compute the energy of the modified system.
+
+    double pressure = context.getParameter(RPMDMonteCarloBarostat::Pressure())*(AVOGADRO*1e-25);
+    double kT = BOLTZ*integrator.getTemperature();
+    double w = (finalEnergy-initialEnergy)/numCopies + pressure*deltaVolume - context.getMolecules().size()*kT*std::log(newVolume/volume);
+    if (w > 0 && genrand_real2(random) > std::exp(-w/kT)) {
+        // Reject the step.
+
+        for (int copy = 0; copy < numCopies; copy++)
+            integrator.setPositions(copy, savedPositions[copy]);
+        context.getOwner().setPeriodicBoxVectors(box[0], box[1], box[2]);
+        volume = newVolume;
+    }
+    else
+        numAccepted++;
+    numAttempted++;
+    if (numAttempted >= 10) {
+        if (numAccepted < 0.25*numAttempted) {
+            volumeScale /= 1.1;
+            numAttempted = 0;
+            numAccepted = 0;
+        }
+        else if (numAccepted > 0.75*numAttempted) {
+            volumeScale = std::min(volumeScale*1.1, volume*0.3);
+            numAttempted = 0;
+            numAccepted = 0;
+        }
+    }
+}
+
+std::map<std::string, double> RPMDMonteCarloBarostatImpl::getDefaultParameters() {
+    std::map<std::string, double> parameters;
+    parameters[RPMDMonteCarloBarostat::Pressure()] = getOwner().getDefaultPressure();
+    return parameters;
+}
+
+std::vector<std::string> RPMDMonteCarloBarostatImpl::getKernelNames() {
+    std::vector<std::string> names;
+    names.push_back(ApplyMonteCarloBarostatKernel::Name());
+    return names;
+}

plugins/rpmd/platforms/cuda/src/CudaRpmdKernels.cpp

@@ -257,13 +257,8 @@ void CudaIntegrateRPMDStepKernel::computeForces(ContextImpl& context) {
         context.updateContextState();
         Vec3 finalBox[3];
         context.getPeriodicBoxVectors(finalBox[0], finalBox[1], finalBox[2]);
-        if (initialBox[0] != finalBox[0] || initialBox[1] != finalBox[1] || initialBox[2] != finalBox[2]) {
-            // A barostat was applied during updateContextState().  Adjust the particle positions in all the
-            // other copies to match this one.
-            
-            void* args[] = {&positions->getDevicePointer(), &cu.getPosq().getDevicePointer(), &cu.getAtomIndexArray().getDevicePointer(), &i};
-            cu.executeKernel(translateKernel, args, cu.getNumAtoms());
-        }
+        if (initialBox[0] != finalBox[0] || initialBox[1] != finalBox[1] || initialBox[2] != finalBox[2])
+            throw OpenMMException("Standard barostats cannot be used with RPMDIntegrator.  Use RPMDMonteCarloBarostat instead.");
         context.calcForcesAndEnergy(true, false, groupsNotContracted);
         void* copyFromContextArgs[] = {&cu.getForce().getDevicePointer(), &forces->getDevicePointer(), &cu.getVelm().getDevicePointer(),
                 &velocities->getDevicePointer(), &cu.getPosq().getDevicePointer(), &positions->getDevicePointer(), &cu.getAtomIndexArray().getDevicePointer(), &i};

plugins/rpmd/platforms/cuda/tests/TestCudaRpmd.cpp

@@ -38,11 +38,11 @@
 #include "openmm/Context.h"
 #include "openmm/CustomNonbondedForce.h"
 #include "openmm/HarmonicBondForce.h"
-#include "openmm/MonteCarloBarostat.h"
 #include "openmm/NonbondedForce.h"
 #include "openmm/Platform.h"
 #include "openmm/System.h"
 #include "openmm/RPMDIntegrator.h"
+#include "openmm/RPMDMonteCarloBarostat.h"
 #include "openmm/VirtualSite.h"
 #include "SimTKOpenMMUtilities.h"
 #include "sfmt/SFMT.h"
@@ -371,7 +371,7 @@ void testContractions() {
     system.addForce(bonds);
     NonbondedForce* nonbonded = new NonbondedForce();
     nonbonded->setCutoffDistance(cutoff);
-    nonbonded->setNonbondedMethod(NonbondedForce::PME);
+    nonbonded->setNonbondedMethod(NonbondedForce::CutoffPeriodic);
     nonbonded->setForceGroup(1);
     nonbonded->setReciprocalSpaceForceGroup(2);
     system.addForce(nonbonded);
@@ -491,11 +491,11 @@ void testWithBarostat() {
     system.addForce(bonds);
     NonbondedForce* nonbonded = new NonbondedForce();
     nonbonded->setCutoffDistance(cutoff);
-    nonbonded->setNonbondedMethod(NonbondedForce::PME);
+    nonbonded->setNonbondedMethod(NonbondedForce::CutoffPeriodic);
     nonbonded->setForceGroup(1);
     nonbonded->setReciprocalSpaceForceGroup(2);
     system.addForce(nonbonded);
-    system.addForce(new MonteCarloBarostat(0.5, temperature));
+    system.addForce(new RPMDMonteCarloBarostat(0.5, 10));
 
     // Create a cloud of molecules.
 

plugins/rpmd/platforms/opencl/src/OpenCLRpmdKernels.cpp

@@ -273,13 +273,8 @@ void OpenCLIntegrateRPMDStepKernel::computeForces(ContextImpl& context) {
         context.updateContextState();
         Vec3 finalBox[3];
         context.getPeriodicBoxVectors(finalBox[0], finalBox[1], finalBox[2]);
-        if (initialBox[0] != finalBox[0] || initialBox[1] != finalBox[1] || initialBox[2] != finalBox[2]) {
-            // A barostat was applied during updateContextState().  Adjust the particle positions in all the
-            // other copies to match this one.
-            
-            translateKernel.setArg<cl_int>(3, i);
-            cl.executeKernel(translateKernel, cl.getNumAtoms());
-        }
+        if (initialBox[0] != finalBox[0] || initialBox[1] != finalBox[1] || initialBox[2] != finalBox[2])
+            throw OpenMMException("Standard barostats cannot be used with RPMDIntegrator.  Use RPMDMonteCarloBarostat instead.");
         context.calcForcesAndEnergy(true, false, groupsNotContracted);
         copyFromContextKernel.setArg<cl_int>(7, i);
         cl.executeKernel(copyFromContextKernel, cl.getNumAtoms());

plugins/rpmd/platforms/opencl/tests/TestOpenCLRpmd.cpp

@@ -38,11 +38,11 @@
 #include "openmm/Context.h"
 #include "openmm/CustomNonbondedForce.h"
 #include "openmm/HarmonicBondForce.h"
-#include "openmm/MonteCarloBarostat.h"
 #include "openmm/NonbondedForce.h"
 #include "openmm/Platform.h"
 #include "openmm/System.h"
 #include "openmm/RPMDIntegrator.h"
+#include "openmm/RPMDMonteCarloBarostat.h"
 #include "openmm/VirtualSite.h"
 #include "SimTKOpenMMUtilities.h"
 #include "sfmt/SFMT.h"
@@ -372,7 +372,7 @@ void testContractions() {
     system.addForce(bonds);
     NonbondedForce* nonbonded = new NonbondedForce();
     nonbonded->setCutoffDistance(cutoff);
-    nonbonded->setNonbondedMethod(NonbondedForce::PME);
+    nonbonded->setNonbondedMethod(NonbondedForce::CutoffPeriodic);
     nonbonded->setForceGroup(1);
     nonbonded->setReciprocalSpaceForceGroup(2);
     system.addForce(nonbonded);
@@ -492,11 +492,11 @@ void testWithBarostat() {
     system.addForce(bonds);
     NonbondedForce* nonbonded = new NonbondedForce();
     nonbonded->setCutoffDistance(cutoff);
-    nonbonded->setNonbondedMethod(NonbondedForce::PME);
+    nonbonded->setNonbondedMethod(NonbondedForce::CutoffPeriodic);
     nonbonded->setForceGroup(1);
     nonbonded->setReciprocalSpaceForceGroup(2);
     system.addForce(nonbonded);
-    system.addForce(new MonteCarloBarostat(0.5, temperature));
+    system.addForce(new RPMDMonteCarloBarostat(0.5, 10));
 
     // Create a cloud of molecules.
 

plugins/rpmd/platforms/reference/src/ReferenceRpmdKernels.cpp

@@ -280,17 +280,8 @@ void ReferenceIntegrateRPMDStepKernel::computeForces(ContextImpl& context, const
         context.updateContextState();
         Vec3 finalBox[3];
         context.getPeriodicBoxVectors(finalBox[0], finalBox[1], finalBox[2]);
-        if (initialBox[0] != finalBox[0] || initialBox[1] != finalBox[1] || initialBox[2] != finalBox[2]) {
-            // A barostat was applied during updateContextState().  Adjust the particle positions in all the
-            // other copies to match this one.
-            
-            for (int j = 0; j < numParticles; j++) {
-                Vec3 delta = pos[j]-positions[i][j];
-                for (int k = 0; k < totalCopies; k++)
-                    if (k != i)
-                        positions[k][j] += delta;
-            }
-        }
+        if (initialBox[0] != finalBox[0] || initialBox[1] != finalBox[1] || initialBox[2] != finalBox[2])
+            throw OpenMMException("Standard barostats cannot be used with RPMDIntegrator.  Use RPMDMonteCarloBarostat instead.");
         positions[i] = pos;
         velocities[i] = vel;
         context.calcForcesAndEnergy(true, false, groupsNotContracted);

plugins/rpmd/platforms/reference/tests/TestReferenceRpmd.cpp

@@ -37,11 +37,11 @@
 #include "openmm/CMMotionRemover.h"
 #include "openmm/Context.h"
 #include "openmm/HarmonicBondForce.h"
-#include "openmm/MonteCarloBarostat.h"
 #include "openmm/NonbondedForce.h"
 #include "openmm/Platform.h"
 #include "openmm/System.h"
 #include "openmm/RPMDIntegrator.h"
+#include "openmm/RPMDMonteCarloBarostat.h"
 #include "openmm/VirtualSite.h"
 #include "SimTKOpenMMUtilities.h"
 #include "sfmt/SFMT.h"
@@ -255,7 +255,7 @@ void testContractions() {
     system.addForce(bonds);
     NonbondedForce* nonbonded = new NonbondedForce();
     nonbonded->setCutoffDistance(cutoff);
-    nonbonded->setNonbondedMethod(NonbondedForce::PME);
+    nonbonded->setNonbondedMethod(NonbondedForce::CutoffPeriodic);
     nonbonded->setForceGroup(1);
     nonbonded->setReciprocalSpaceForceGroup(2);
     system.addForce(nonbonded);
@@ -375,11 +375,11 @@ void testWithBarostat() {
     system.addForce(bonds);
     NonbondedForce* nonbonded = new NonbondedForce();
     nonbonded->setCutoffDistance(cutoff);
-    nonbonded->setNonbondedMethod(NonbondedForce::PME);
+    nonbonded->setNonbondedMethod(NonbondedForce::CutoffPeriodic);
     nonbonded->setForceGroup(1);
     nonbonded->setReciprocalSpaceForceGroup(2);
     system.addForce(nonbonded);
-    system.addForce(new MonteCarloBarostat(0.5, temperature));
+    system.addForce(new RPMDMonteCarloBarostat(0.5, 10));
 
     // Create a cloud of molecules.
 

wrappers/python/src/swig_doxygen/OpenMM.i

@@ -48,6 +48,7 @@ namespace std {
 #include "OpenMM.h"
 #include "OpenMMAmoeba.h"
 #include "openmm/RPMDIntegrator.h"
+#include "openmm/RPMDMonteCarloBarostat.h"
 #include "OpenMMDrude.h"
 #include "openmm/serialization/SerializationNode.h"
 #include "openmm/serialization/SerializationProxy.h"

wrappers/python/src/swig_doxygen/doxygen/Doxyfile.in

@@ -621,6 +621,7 @@ EXCLUDE_PATTERNS       = */tests/* \
                          *amoebaKernels.h \
                          *DrudeKernels.h \
                          *RpmdKernels.h \
+                         *RPMDUpdater.h \
                          *OpenMMFortranModule.f90 \
                          *OpenMMCWrapper.h