search

openmm / openmm-plumed

OpenMM plugin to interface with PLUMED
55 stars 23 forks source link

Implement setMasses and getMasses #30

Closed raimis closed 4 years ago

raimis commented 4 years ago

PLUMED allows to set particle via its API (https://www.plumed.org/doc-v2.6/developer-doc/html/_how_to_plumed_your_m_d.html):

plumed_cmd(plumedmain, "setMasses", &masses);

This PR implements a corresponding functionality for OpenMM-PLUMED API:

void PlumedForce::setMasses(const std::vector<double>& masses);
const std::vector<double>& PlumedForce::getMasses() const;

By default, OpenMM-PLUMED uses the masses from OpenMM::System. If the hydrogen mass repartitioning is used, these masses are nonphysical and shouldn't be used to compute CVs.

peastman commented 4 years ago

This makes me very nervous. It opens the door to a lot of problems!

Whether the masses are nonphysical is beside the point. All MD simulations are nonphysical in a lot of ways. But they're still well defined models, and they produce useful results in many cases despite their nonphysical features.

Mass repartitioning is an approximation. If it's good enough for your purposes then you can use it, and if it isn't then you don't, but it needs to be used consistently. If you use repartitioned masses in some parts of the calculation but not in others, then you no longer have a well defined model. You're likely to get all sorts of subtle but vitally important artifacts. For example, equations of motion that don't conserve energy, which invalidates the assumptions behind most sampling methods.

raimis commented 4 years ago
  1. PLUMED don't use these masses for the equations of motion. The integration is done by the MD engine and PLUMED just returns additional forces. So for a typical use (i.e. restraints or metadynamics), the conservation of energy won't be affected.
  2. PLUMED uses these masses for the definition of collective variables. For example to compute the centers of mass, but this is an arbitrary choice. The centers could be defined using any other weights (e.g. geometric centers) and this won't make the model more correct or incorrect, in general.
  3. And, yes, PLUMED is flexible enough, so probably there is a way to abuse these masses, but even without them there are many other ways to make simulations wrong.

For example, a metadynamics simulation where a CV needs a COM and the integration is done without the mass repartitioning. Assuming, it is valid to use the mass repartitioning in case of an MD simulation, it should be a case for the metadynamics simulation too. And the results shouldn't change significantly, but it won't be the case, if the the integrator setup affects the CV definition.

So the problem I want to solve with this PR is to decouple the integrator setup in OpenMM and the collective variable definitions in PLUMED.

raimis commented 4 years ago

@peastman this is ready!

peastman commented 4 years ago

It still makes me nervous, but I guess as long as PLUMED never tries to do anything with a momentum or kinetic energy it will be ok.