coMMParE - Comparing Molecular Mechanics Potential Energies

The goal of this package is to evaluate energies using ParmEd as a middleman and various molecular mechanics engines as the backend for computation. Ideas derive from intermol, parmed, validate.

Some overarching goals:

• Validate FF parameter and molecular model translation to a particular format or data structure
• Assess consistency of MM engines in calculating energies

The broad steps to compare energy via various permutations:

• Input of molecular model (MM input files, Pythonic creation, etc)
  • These inputs get turned into parmed.Structure
• Output of molecular model for energy calculation (MM engine)
  • The parmed.Structure gets converted to its respective MM engine data object/format

While more-sophisticated Python packaging may come, here is a list of possible dependencies, with conda links where available. The packages for "inputs" and "outputs" are not all required, just packages whose conversions/energy evaluations are supported here. The only real, core package necessary is ParmEd, as it is the core data structure used.

Creating and parametrizing molecular models (inputs)

• mbuild (conda)
• foyer (conda)
  • gaff-foyer
• parmed (conda)
• openforcefield (conda)

MM-interconversion tools

• parmed (conda)
• this mbuild PR

MM engines for measuring energy (outputs)

• openmm (conda)
• hoomd (conda)
• gromacs (conda)
  • panedr (conda)

Contributing

There are multiple ways to expand this testing suite - more MM engines or more reference systems.

Contributing new MM engine functionality

• Desired Input: parmed.Structure containing molecular model
• Desired Output: pandas.DataFrame containing energies
• Given a parmed.Structure object, be able to use your MM engine to build your simulation and measure the energy from the intial configuration and molecular model parameters within the parmed.Structure.
  • If input files need to be written, please use the tempfile library to create and dump files to a temp directory (commpare/gromacs utilizes temp directories)
• After building your MM engine's simulation, measure the energy and "canonicalize" the energy terms.
  • Ultimately, we want a pandas.DataFrame whose index is the MM engine with columns: bond, angle, dihedral, nonbond, all. The values should be energies in kJ/mol. The energy breakdown may be further decomposed, but the goal for now is those 4 energy groups. all does not necessarily have to be the sum of bond, angle, dihedral, nonbond if there are force/energy groups unaccounted for and non-canonicalized.
• While a single function gets called in conversion.py, inputs a parmed.Structure (and possibly optional arguments like in commpare/hoomd), and outputs a canonicalized pandas.DataFrame, helper functions are certainly welcomed and encouraged to enhance readability and debugging
• mm_engines.py will need to be updated to detect the new MM engine
• conversion.py will need to be updated to convert and run the parmed.Structure in the MD engine
• Add unit tests in commpare/tests that, given a parmed.Structure, will build the particular MM engine's simulation and return the pandas.DataFrame with canonicalized energy terms.

Contributing new MM tests

• The goal is to test different molecular model/force field terms and observe how well they translate to different MM engines
• Ideally, there should be some directory organization to the MM tests such that each directory is a particular chemical system understood by parmed
• reference_systems.py will need to be updated to detect the locally available reference systems against which to test the local MM engines
• Using pytest, put together some test scripts in commpare/mm_tests that will construct the parmed.Structure from the reference systems, then run commpare.spawn_engine_simulations to "shotgun" a variety of energy calculations

The energies don't agree - what now?

Here's some spots where things could have gone wrong:

• Did ParmEd correctly parse your input into a parmed.Structure? Some functional forms, coefficients, exceptions, scalings, etc. might not be parsed correctly. Check your units and factors of two, ParmEd tries to convert units to Angstrom and kcal/mol
• Did ParmEd or mBuild correctly translate the parmed.Structure into the engine of choice? Again, some functional forms, coefficients, exceptions, scalings, etc. might not be be translated correctly. Check your units and factors of two
• Within coMMParE, are certain simulation/force field treatments being implemented/treated the same across engines?
• Within coMMParE, are the molecular mechanics energies compared "apples to apples" and canonicalized the same? For example, was an engine's 1,4 LJ energy lumped into a dihedral energy, or vice versa?
• Lastly, is the engine calculating these functional forms correctly?

Notes and Warnings

This is NOT meant to be your one-stop-shop package to run molecular simulations. To CORRECTLY run a simulation means using correct simulation parameters, consistent software, and reliable force field parameters. This package assumes certain simulation parameters and certain model implementations that will likely yield PHYSICALLY INACCURATE results. While inaccurate, we aim to CONTROL for some of these choices to validate models and force fields are correctly being translated.