From 3cbe2f3ed7e18a6dfb153ae4cab3a3fa6b82b0de Mon Sep 17 00:00:00 2001 From: "Matthew W. Thompson" Date: Thu, 7 Oct 2021 14:54:48 -0500 Subject: [PATCH 01/24] Add Interchange + GROMACS to existing GROMACS example --- .../convert_to_amber_gromacs.ipynb | 144 ++++++++++++++++-- 1 file changed, 133 insertions(+), 11 deletions(-) diff --git a/examples/using_smirnoff_in_amber_or_gromacs/convert_to_amber_gromacs.ipynb b/examples/using_smirnoff_in_amber_or_gromacs/convert_to_amber_gromacs.ipynb index 1bff64170..bf4296572 100644 --- a/examples/using_smirnoff_in_amber_or_gromacs/convert_to_amber_gromacs.ipynb +++ b/examples/using_smirnoff_in_amber_or_gromacs/convert_to_amber_gromacs.ipynb @@ -84,7 +84,7 @@ "cell_type": "markdown", "metadata": {}, "source": [ - "### Convert OpenMM System to AMBER and GROMACS files\n", + "### Convert OpenMM System to AMBER and GROMACS files using ParmEd\n", "\n", "\n", "First, we convert the OpenMM `System` into a ParmEd `Structure`. We'll use the atom positions in the PDB to create the coordinate files.\n", @@ -227,16 +227,16 @@ "System loaded from AMBER files:\n", "-------------------------------\n", "{'HarmonicAngleForce': Quantity(value=155.70130920410156, unit=kilojoule/mole),\n", - " 'HarmonicBondForce': Quantity(value=0.4403935372829437, unit=kilojoule/mole),\n", - " 'NonbondedForce': Quantity(value=2.8259036956507657, unit=kilojoule/mole),\n", - " 'PeriodicTorsionForce': Quantity(value=19.563232421875, unit=kilojoule/mole)}\n", + " 'HarmonicBondForce': Quantity(value=0.440397173166275, unit=kilojoule/mole),\n", + " 'NonbondedForce': Quantity(value=-5.549589580022797, unit=kilojoule/mole),\n", + " 'PeriodicTorsionForce': Quantity(value=19.563228607177734, unit=kilojoule/mole)}\n", "\n", "Original OpenMM System:\n", "-----------------------\n", - "{'HarmonicAngleForce': Quantity(value=155.70130920410156, unit=kilojoule/mole),\n", - " 'HarmonicBondForce': Quantity(value=0.44038787484169006, unit=kilojoule/mole),\n", - " 'NonbondedForce': Quantity(value=2.8251933539722245, unit=kilojoule/mole),\n", - " 'PeriodicTorsionForce': Quantity(value=19.563230514526367, unit=kilojoule/mole)}\n" + "{'HarmonicAngleForce': Quantity(value=155.7012939453125, unit=kilojoule/mole),\n", + " 'HarmonicBondForce': Quantity(value=0.44039154052734375, unit=kilojoule/mole),\n", + " 'NonbondedForce': Quantity(value=-5.55016902838662, unit=kilojoule/mole),\n", + " 'PeriodicTorsionForce': Quantity(value=19.563228607177734, unit=kilojoule/mole)}\n" ] } ], @@ -253,17 +253,139 @@ "pprint(omm_energies)" ] }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "### Export to GROMACS files using Interchange\n", + "\n", + "First, we construct an `Interchange` object." + ] + }, + { + "cell_type": "code", + "execution_count": 9, + "metadata": {}, + "outputs": [], + "source": [ + "from openff.interchange.components.interchange import Interchange" + ] + }, { "cell_type": "code", - "execution_count": null, + "execution_count": 10, "metadata": {}, "outputs": [], - "source": [] + "source": [ + "interchange = Interchange.from_smirnoff(force_field=forcefield, topology=off_topology)\n", + "interchange.positions = pdbfile.positions" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "This object can be exported to a variety of formats, including GROMACS, OpenMM, and LAMMPS" + ] + }, + { + "cell_type": "code", + "execution_count": 11, + "metadata": {}, + "outputs": [], + "source": [ + "interchange.to_gro(\"system.gro\")\n", + "interchange.to_top(\"system.top\")\n", + "openmm_system = interchange.to_openmm(combine_nonbonded_forces=True)" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "The Interchange project provides a `drivers` module with functions that compute single-point energies of `Interchange` objects with a variety of simulation engines. We can use it to validate that the exports to each engine produce files and/or objects that produce sufficiently similar single-point energies. Compare these values to those in cell 8 above." + ] + }, + { + "cell_type": "code", + "execution_count": 12, + "metadata": {}, + "outputs": [], + "source": [ + "from openff.interchange.drivers import get_gromacs_energies, get_openmm_energies" + ] + }, + { + "cell_type": "code", + "execution_count": 13, + "metadata": {}, + "outputs": [ + { + "data": { + "text/plain": [ + "{'Bond': 0.440394788980484 ,\n", + " 'Angle': 155.7017059326172 ,\n", + " 'Torsion': 19.56324005126953 ,\n", + " 'vdW': 9.66697609052062 ,\n", + " 'Electrostatics': -15.235475540161133 }" + ] + }, + "execution_count": 13, + "metadata": {}, + "output_type": "execute_result" + } + ], + "source": [ + "get_gromacs_energies(interchange).energies" + ] + }, + { + "cell_type": "code", + "execution_count": 14, + "metadata": { + "scrolled": true + }, + "outputs": [ + { + "data": { + "text/plain": [ + "{'Bond': 0.44039154052734375 ,\n", + " 'Angle': 155.7012939453125 ,\n", + " 'Torsion': 19.563228607177734 ,\n", + " 'Nonbonded': -5.5570980420777545 }" + ] + }, + "execution_count": 14, + "metadata": {}, + "output_type": "execute_result" + } + ], + "source": [ + "get_openmm_energies(interchange).energies" + ] + }, + { + "cell_type": "code", + "execution_count": 15, + "metadata": {}, + "outputs": [], + "source": [ + "# TODO: Add get_summary_data when it makes it into a release (v0.1.4, maybe?)" + ] + }, + { + "cell_type": "code", + "execution_count": 16, + "metadata": {}, + "outputs": [], + "source": [ + "# TODO: Add Amber exports when prmtop writer is more trustworthy" + ] } ], "metadata": { "kernelspec": { - "display_name": "Python 3", + "display_name": "Python 3 (ipykernel)", "language": "python", "name": "python3" }, From 325a35b6ff72ec43a179a03e100558ca3792b5f0 Mon Sep 17 00:00:00 2001 From: "Matthew W. Thompson" Date: Wed, 27 Oct 2021 09:13:13 -0500 Subject: [PATCH 02/24] Overhaul Amber + GROMACS example --- .../convert_to_amber_gromacs.ipynb | 407 ----------------- .../smirnoff_amber_gromacs.ipynb | 409 ++++++++++++++++++ 2 files changed, 409 insertions(+), 407 deletions(-) delete mode 100644 examples/using_smirnoff_in_amber_or_gromacs/convert_to_amber_gromacs.ipynb create mode 100644 examples/using_smirnoff_in_amber_or_gromacs/smirnoff_amber_gromacs.ipynb diff --git a/examples/using_smirnoff_in_amber_or_gromacs/convert_to_amber_gromacs.ipynb b/examples/using_smirnoff_in_amber_or_gromacs/convert_to_amber_gromacs.ipynb deleted file mode 100644 index bf4296572..000000000 --- a/examples/using_smirnoff_in_amber_or_gromacs/convert_to_amber_gromacs.ipynb +++ /dev/null @@ -1,407 +0,0 @@ -{ - "cells": [ - { - "cell_type": "markdown", - "metadata": {}, - "source": [ - "## Convert Open Forcefield System to AMBER and GROMACS input files\n", - "\n", - "The Open Forcefield Toolkit can create parametrized `System` objects that can be natively simulated with OpenMM. This example shows how you can convert such an OpenMM `System` into AMBER prmtop/inpcrd and GROMACS top/gro input files through the ParmEd library." - ] - }, - { - "cell_type": "markdown", - "metadata": {}, - "source": [ - "### Create an OpenMM System\n", - "\n", - "We start by loading a PDB file containing one copy of ethanol and cyclohexane. Our goal is to create an OFF `Topology` object describing this system that we can parametrize with the SMIRNOFF-format \"Sage\" force field.\n", - "\n", - "The two `Molecule` objects created from the SMILES strings can contain information such as partial charges and stereochemistry that is not included in an OpenMM topology. In this example, partial charges are not explicitly given, and `ForceField` will assign AM1/BCC charges as specified by the \"Sage\" force field. Note that the Open Force Field Toolkit produces deterministic partial charges that do not depend on the input conformation of parameterized molecules. See the [FAQ](https://open-forcefield-toolkit.readthedocs.io/en/latest/faq.html#the-partial-charges-generated-by-the-toolkit-don-t-seem-to-depend-on-the-molecule-s-conformation-is-this-a-bug) for more information." - ] - }, - { - "cell_type": "code", - "execution_count": 1, - "metadata": {}, - "outputs": [ - { - "name": "stderr", - "output_type": "stream", - "text": [ - "Warning: Unable to load toolkit 'OpenEye Toolkit'. The Open Force Field Toolkit does not require the OpenEye Toolkits, and can use RDKit/AmberTools instead. However, if you have a valid license for the OpenEye Toolkits, consider installing them for faster performance and additional file format support: https://docs.eyesopen.com/toolkits/python/quickstart-python/linuxosx.html OpenEye offers free Toolkit licenses for academics: https://www.eyesopen.com/academic-licensing\n" - ] - } - ], - "source": [ - "try:\n", - " from openmm.app import PDBFile\n", - "except ImportError:\n", - " from simtk.openmm.app import PDBFile\n", - "\n", - "from openff.toolkit.topology import Molecule, Topology\n", - "from openff.toolkit.typing.engines.smirnoff import ForceField" - ] - }, - { - "cell_type": "code", - "execution_count": 2, - "metadata": {}, - "outputs": [], - "source": [ - "ethanol = Molecule.from_smiles(\"CCO\")\n", - "cyclohexane = Molecule.from_smiles(\"C1CCCCC1\")\n", - "\n", - "# Obtain the OpenMM Topology object from the PDB file.\n", - "pdbfile = PDBFile(\"1_cyclohexane_1_ethanol.pdb\")\n", - "omm_topology = pdbfile.topology\n", - "\n", - "# Create the Open Forcefield Topology.\n", - "off_topology = Topology.from_openmm(\n", - " omm_topology, unique_molecules=[ethanol, cyclohexane]\n", - ")" - ] - }, - { - "cell_type": "markdown", - "metadata": {}, - "source": [ - "Now we parametrize the OFF `Topology` to create an OpenMM `System`. Since ParmEd will run with the `constraints=HBonds` keyword later, we use the _unconstrained_ version of the Sage force field here." - ] - }, - { - "cell_type": "code", - "execution_count": 3, - "metadata": {}, - "outputs": [], - "source": [ - "# Load the \"Sage\" force field.\n", - "forcefield = ForceField(\"openff_unconstrained-2.0.0.offxml\")\n", - "omm_system = forcefield.create_openmm_system(off_topology)" - ] - }, - { - "cell_type": "markdown", - "metadata": {}, - "source": [ - "### Convert OpenMM System to AMBER and GROMACS files using ParmEd\n", - "\n", - "\n", - "First, we convert the OpenMM `System` into a ParmEd `Structure`. We'll use the atom positions in the PDB to create the coordinate files.\n", - "\n", - "
\n", - " Warning: ParmEd's Structure model is inspired by AMBER, and some information in an OpenMM System are not directly translatable into a Structure. In particular, as of today (4/2/2019), long-range interaction treatment method (e.g., PME, CutoffPeriodic) and parameters (e.g., cutoff and cutoff switching distance, PME error tolerance) are known to be lost during the conversion.\n", - "
" - ] - }, - { - "cell_type": "code", - "execution_count": 4, - "metadata": {}, - "outputs": [], - "source": [ - "import parmed\n", - "\n", - "# Convert OpenMM System to a ParmEd structure.\n", - "parmed_structure = parmed.openmm.load_topology(\n", - " omm_topology, omm_system, pdbfile.positions\n", - ")" - ] - }, - { - "cell_type": "markdown", - "metadata": {}, - "source": [ - "We can then use ParmEd to convert an OpenMM `System` to prmtop/inpcrd or top/gro files that can be read by AMBER and GROMACS respectively. ParmEd is capable of converting parametrized files to other formats as well. For further information, see ParmEd's documentation: https://parmed.github.io/ParmEd/html/readwrite.html ." - ] - }, - { - "cell_type": "code", - "execution_count": 5, - "metadata": {}, - "outputs": [], - "source": [ - "# Export AMBER files.\n", - "parmed_structure.save(\"system.prmtop\", overwrite=True)\n", - "parmed_structure.save(\"system.inpcrd\", overwrite=True)\n", - "\n", - "# Export GROMACS files.\n", - "parmed_structure.save(\"system.top\", overwrite=True)\n", - "parmed_structure.save(\"system.gro\", overwrite=True)" - ] - }, - { - "cell_type": "markdown", - "metadata": {}, - "source": [ - "### Validate the conversion\n", - "\n", - "ParmEd is generally a reliable and robust library, but we can easily check that everything went as expected during the conversion by loading the exported files into an OpenMM `System` and comparing it with the original. Note that you'll have to specify the correct nonbonded method and cutoff settings for the energy comparison to make sense since they are not included in the AMBER prmtop (or GROMACS top/gro) files." - ] - }, - { - "cell_type": "code", - "execution_count": 6, - "metadata": {}, - "outputs": [ - { - "name": "stdout", - "output_type": "stream", - "text": [ - "0.9 nm\n", - "False\n", - "True\n" - ] - } - ], - "source": [ - "try:\n", - " import openmm\n", - "except ImportError:\n", - " from simtk import openmm\n", - "\n", - "for force in omm_system.getForces():\n", - " if isinstance(force, openmm.NonbondedForce):\n", - " break\n", - "print(force.getCutoffDistance())\n", - "print(force.getUseSwitchingFunction())\n", - "print(force.getNonbondedMethod() == openmm.NonbondedForce.PME)" - ] - }, - { - "cell_type": "code", - "execution_count": 7, - "metadata": {}, - "outputs": [], - "source": [ - "try:\n", - " from openmm import unit\n", - " from openmm.app import PME, HBonds\n", - "except ImportError:\n", - " from simtk import unit\n", - " from simtk.openmm.app import PME, HBonds\n", - "\n", - "from openff.toolkit.tests.utils import (\n", - " compare_system_energies,\n", - " compare_system_parameters,\n", - ")\n", - "\n", - "# Load the prmtop/inpcrd files into a ParmEd Structure.as an OpenMM System object.\n", - "amber_structure = parmed.load_file(\"system.prmtop\", \"system.inpcrd\")\n", - "\n", - "# Convert the Structure to an OpenMM System. Note that by\n", - "# default ParmEd will add a CMMotionRemover force to the\n", - "# System, and won't constrain the hydrogen bonds.\n", - "amber_system = amber_structure.createSystem(\n", - " nonbondedMethod=PME,\n", - " nonbondedCutoff=9.0 * unit.angstrom,\n", - " switchDistance=0.0 * unit.angstrom,\n", - " constraints=HBonds,\n", - " removeCMMotion=False,\n", - ")\n", - "\n", - "# Compare the parameters of the original and converted Systems.\n", - "# This raises FailedParameterComparisonError if the comparison fails.\n", - "compare_system_parameters(omm_system, amber_system)\n", - "\n", - "# Compare the energies by force.\n", - "# This raises FailedEnergyComparisonError if the comparison fails.\n", - "amber_energies, omm_energies = compare_system_energies(\n", - " amber_system,\n", - " omm_system,\n", - " amber_structure.positions,\n", - " amber_structure.box_vectors,\n", - " rtol=1e-3,\n", - ")" - ] - }, - { - "cell_type": "code", - "execution_count": 8, - "metadata": {}, - "outputs": [ - { - "name": "stdout", - "output_type": "stream", - "text": [ - "System loaded from AMBER files:\n", - "-------------------------------\n", - "{'HarmonicAngleForce': Quantity(value=155.70130920410156, unit=kilojoule/mole),\n", - " 'HarmonicBondForce': Quantity(value=0.440397173166275, unit=kilojoule/mole),\n", - " 'NonbondedForce': Quantity(value=-5.549589580022797, unit=kilojoule/mole),\n", - " 'PeriodicTorsionForce': Quantity(value=19.563228607177734, unit=kilojoule/mole)}\n", - "\n", - "Original OpenMM System:\n", - "-----------------------\n", - "{'HarmonicAngleForce': Quantity(value=155.7012939453125, unit=kilojoule/mole),\n", - " 'HarmonicBondForce': Quantity(value=0.44039154052734375, unit=kilojoule/mole),\n", - " 'NonbondedForce': Quantity(value=-5.55016902838662, unit=kilojoule/mole),\n", - " 'PeriodicTorsionForce': Quantity(value=19.563228607177734, unit=kilojoule/mole)}\n" - ] - } - ], - "source": [ - "# Pretty-print the energies by component.\n", - "from pprint import pprint\n", - "\n", - "print(\"System loaded from AMBER files:\")\n", - "print(\"-------------------------------\")\n", - "pprint(amber_energies)\n", - "\n", - "print(\"\\nOriginal OpenMM System:\")\n", - "print(\"-----------------------\")\n", - "pprint(omm_energies)" - ] - }, - { - "cell_type": "markdown", - "metadata": {}, - "source": [ - "### Export to GROMACS files using Interchange\n", - "\n", - "First, we construct an `Interchange` object." - ] - }, - { - "cell_type": "code", - "execution_count": 9, - "metadata": {}, - "outputs": [], - "source": [ - "from openff.interchange.components.interchange import Interchange" - ] - }, - { - "cell_type": "code", - "execution_count": 10, - "metadata": {}, - "outputs": [], - "source": [ - "interchange = Interchange.from_smirnoff(force_field=forcefield, topology=off_topology)\n", - "interchange.positions = pdbfile.positions" - ] - }, - { - "cell_type": "markdown", - "metadata": {}, - "source": [ - "This object can be exported to a variety of formats, including GROMACS, OpenMM, and LAMMPS" - ] - }, - { - "cell_type": "code", - "execution_count": 11, - "metadata": {}, - "outputs": [], - "source": [ - "interchange.to_gro(\"system.gro\")\n", - "interchange.to_top(\"system.top\")\n", - "openmm_system = interchange.to_openmm(combine_nonbonded_forces=True)" - ] - }, - { - "cell_type": "markdown", - "metadata": {}, - "source": [ - "The Interchange project provides a `drivers` module with functions that compute single-point energies of `Interchange` objects with a variety of simulation engines. We can use it to validate that the exports to each engine produce files and/or objects that produce sufficiently similar single-point energies. Compare these values to those in cell 8 above." - ] - }, - { - "cell_type": "code", - "execution_count": 12, - "metadata": {}, - "outputs": [], - "source": [ - "from openff.interchange.drivers import get_gromacs_energies, get_openmm_energies" - ] - }, - { - "cell_type": "code", - "execution_count": 13, - "metadata": {}, - "outputs": [ - { - "data": { - "text/plain": [ - "{'Bond': 0.440394788980484 ,\n", - " 'Angle': 155.7017059326172 ,\n", - " 'Torsion': 19.56324005126953 ,\n", - " 'vdW': 9.66697609052062 ,\n", - " 'Electrostatics': -15.235475540161133 }" - ] - }, - "execution_count": 13, - "metadata": {}, - "output_type": "execute_result" - } - ], - "source": [ - "get_gromacs_energies(interchange).energies" - ] - }, - { - "cell_type": "code", - "execution_count": 14, - "metadata": { - "scrolled": true - }, - "outputs": [ - { - "data": { - "text/plain": [ - "{'Bond': 0.44039154052734375 ,\n", - " 'Angle': 155.7012939453125 ,\n", - " 'Torsion': 19.563228607177734 ,\n", - " 'Nonbonded': -5.5570980420777545 }" - ] - }, - "execution_count": 14, - "metadata": {}, - "output_type": "execute_result" - } - ], - "source": [ - "get_openmm_energies(interchange).energies" - ] - }, - { - "cell_type": "code", - "execution_count": 15, - "metadata": {}, - "outputs": [], - "source": [ - "# TODO: Add get_summary_data when it makes it into a release (v0.1.4, maybe?)" - ] - }, - { - "cell_type": "code", - "execution_count": 16, - "metadata": {}, - "outputs": [], - "source": [ - "# TODO: Add Amber exports when prmtop writer is more trustworthy" - ] - } - ], - "metadata": { - "kernelspec": { - "display_name": "Python 3 (ipykernel)", - "language": "python", - "name": "python3" - }, - "language_info": { - "codemirror_mode": { - "name": "ipython", - "version": 3 - }, - "file_extension": ".py", - "mimetype": "text/x-python", - "name": "python", - "nbconvert_exporter": "python", - "pygments_lexer": "ipython3", - "version": "3.8.6" - } - }, - "nbformat": 4, - "nbformat_minor": 4 -} diff --git a/examples/using_smirnoff_in_amber_or_gromacs/smirnoff_amber_gromacs.ipynb b/examples/using_smirnoff_in_amber_or_gromacs/smirnoff_amber_gromacs.ipynb new file mode 100644 index 000000000..5b68f685e --- /dev/null +++ b/examples/using_smirnoff_in_amber_or_gromacs/smirnoff_amber_gromacs.ipynb @@ -0,0 +1,409 @@ +{ + "cells": [ + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "## Using OpenFF force fields in Amber and GROMACS\n", + "\n", + "The Open Forcefield Toolkit can create parametrized `openmm.System` objects that can be natively simulated with OpenMM. This example shows the Interchange project can enable parallel workflows using Amber and GROMACS." + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "### Preparing an OpenFF Topology\n", + "\n", + "We start by loading a PDB file containing one copy of ethanol and cyclohexane. Our goal is to create an OpenFF `Topology` object describing this system that we can parametrize with the SMIRNOFF-format \"Sage\" force field.\n", + "\n", + "The two `Molecule` objects created from the SMILES strings can contain information such as partial charges and stereochemistry that is not included in an OpenMM topology. In this example, partial charges are not explicitly given, and `ForceField` will assign AM1/BCC charges as specified by the \"Sage\" force field. Note that the OpenFF Toolkit produces partial charges that do not depend on the input conformation of parameterized molecules. See the [FAQ](https://open-forcefield-toolkit.readthedocs.io/en/latest/faq.html#the-partial-charges-generated-by-the-toolkit-don-t-seem-to-depend-on-the-molecule-s-conformation-is-this-a-bug) for more information." + ] + }, + { + "cell_type": "code", + "execution_count": 1, + "metadata": {}, + "outputs": [ + { + "name": "stderr", + "output_type": "stream", + "text": [ + "Warning: importing 'simtk.openmm' is deprecated. Import 'openmm' instead.\n" + ] + } + ], + "source": [ + "try:\n", + " from openmm import app\n", + "except ImportError:\n", + " from simtk.openmm import app\n", + "\n", + "from openff.toolkit.topology import Molecule, Topology\n", + "from openff.toolkit.typing.engines.smirnoff import ForceField" + ] + }, + { + "cell_type": "code", + "execution_count": 2, + "metadata": {}, + "outputs": [ + { + "data": { + "text/plain": [ + "" + ] + }, + "execution_count": 2, + "metadata": {}, + "output_type": "execute_result" + } + ], + "source": [ + "ethanol = Molecule.from_smiles(\"CCO\")\n", + "cyclohexane = Molecule.from_smiles(\"C1CCCCC1\")\n", + "\n", + "# Load the PDB file using OpenMM and save the OpenMM Topology\n", + "pdbfile = app.PDBFile(\"1_cyclohexane_1_ethanol.pdb\")\n", + "omm_topology = pdbfile.topology\n", + "\n", + "# Create the OpenFF Topology.\n", + "off_topology = Topology.from_openmm(\n", + " omm_topology, unique_molecules=[ethanol, cyclohexane]\n", + ")\n", + "off_topology" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "### Preparing an OpenFF ForceField\n", + "\n", + "Once the `ForceField` class is imported, the only decision to make is which force field to use. An exhaustive list of force fields released by the Open Force Field Initiative can be found [here](from openff.toolkit.typing.engines.smirnoff import ForceField\n", + ").\n", + "\n", + "Here we will use force field from the \"Sage\" line." + ] + }, + { + "cell_type": "code", + "execution_count": 3, + "metadata": {}, + "outputs": [ + { + "data": { + "text/plain": [ + "" + ] + }, + "execution_count": 3, + "metadata": {}, + "output_type": "execute_result" + } + ], + "source": [ + "forcefield = ForceField(\"openff-2.0.0.offxml\")\n", + "forcefield" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "### Preparing an OpenMM System\n", + "\n", + "Once a force field and topology have been loaded, an `openmm.System` can be generated natively with the OpenFF Toolkit." + ] + }, + { + "cell_type": "code", + "execution_count": 4, + "metadata": {}, + "outputs": [ + { + "data": { + "text/plain": [ + " >" + ] + }, + "execution_count": 4, + "metadata": {}, + "output_type": "execute_result" + } + ], + "source": [ + "omm_system = forcefield.create_openmm_system(off_topology)\n", + "omm_system" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "### Preparing an Interchange object\n", + "\n", + "To exports to engines other than OpenMM, we will make use of the [Interchange](https://openff-interchange.readthedocs.io/) project. There is a high-level `Interchange.from_smirnoff` function that consumes OpenFF Toolkit and ForceField objects and produces an `Interchange` object which can then be exported to formats understood by other molecular simulation engines. This extra step is needed to provide a clean interface between _applied_ parameters and engines. Note also that this step does not require an OpenMM System to be generated; `ForceField.create_openmm_system` does not need to be called to use Amber and GROMACS." + ] + }, + { + "cell_type": "code", + "execution_count": 5, + "metadata": {}, + "outputs": [ + { + "data": { + "text/plain": [ + "Interchange with 27 atoms, periodic topology" + ] + }, + "execution_count": 5, + "metadata": {}, + "output_type": "execute_result" + } + ], + "source": [ + "from openff.interchange.components.interchange import Interchange\n", + "\n", + "interchange = Interchange.from_smirnoff(\n", + " force_field=forcefield,\n", + " topology=off_topology,\n", + ")\n", + "interchange.positions = pdbfile.positions\n", + "interchange" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "### Exporting to Amber and GROMACS files\n", + "\n", + "Once an `Interchange` object has been constructed, its API can be used to export to files understood by GROMACS, Amber, and more." + ] + }, + { + "cell_type": "code", + "execution_count": 6, + "metadata": {}, + "outputs": [], + "source": [ + "# Export AMBER files.\n", + "interchange.to_prmtop(\"system.prmtop\")\n", + "interchange.to_inpcrd(\"system.inpcrd\")\n", + "\n", + "# Export GROMACS files.\n", + "interchange.to_top(\"system.top\")\n", + "interchange.to_gro(\"system.gro\")" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "### Validate the conversion\n", + "\n", + "The Interchange project includes functions that take in an `Interchange` object and call out to simulation engines to run single-point energy calculations (with no minimization or dynamics) for the purpose of validating the export layer with each engine. Under the hood, each of these functions calls API points like those used above while converting to files understood by each engine. These rely on having each engine installed and accessible in the current `$PATH`." + ] + }, + { + "cell_type": "code", + "execution_count": 7, + "metadata": {}, + "outputs": [], + "source": [ + "from openff.interchange.drivers import (\n", + " get_all_energies,\n", + " get_amber_energies,\n", + " get_gromacs_energies,\n", + " get_openmm_energies,\n", + ")\n", + "from openff.interchange.drivers.all import get_summary_data" + ] + }, + { + "cell_type": "code", + "execution_count": 8, + "metadata": {}, + "outputs": [ + { + "data": { + "text/plain": [ + "{'Bond': 0.29528921842575073 ,\n", + " 'Angle': 155.7017059326172 ,\n", + " 'Torsion': 19.56324005126953 ,\n", + " 'vdW': 9.66697609052062 ,\n", + " 'Electrostatics': -15.235475540161133 }" + ] + }, + "execution_count": 8, + "metadata": {}, + "output_type": "execute_result" + } + ], + "source": [ + "get_gromacs_energies(interchange).energies" + ] + }, + { + "cell_type": "code", + "execution_count": 9, + "metadata": { + "scrolled": true + }, + "outputs": [ + { + "data": { + "text/plain": [ + "{'Bond': 0.2952878475189209 ,\n", + " 'Angle': 155.7012939453125 ,\n", + " 'Torsion': 19.563228607177734 ,\n", + " 'vdW': 9.6586229795419 ,\n", + " 'Electrostatics': -15.215721021619657 }" + ] + }, + "execution_count": 9, + "metadata": {}, + "output_type": "execute_result" + } + ], + "source": [ + "get_openmm_energies(interchange).energies" + ] + }, + { + "cell_type": "code", + "execution_count": 10, + "metadata": { + "scrolled": false + }, + "outputs": [ + { + "data": { + "text/plain": [ + "{'Bond': 0.0741 ,\n", + " 'Angle': 37.2023 ,\n", + " 'Torsion': 4.6778 ,\n", + " 'vdW': 9.7231976 ,\n", + " 'Electrostatics': -15.212187200000006 }" + ] + }, + "execution_count": 10, + "metadata": {}, + "output_type": "execute_result" + } + ], + "source": [ + "get_amber_energies(interchange).energies" + ] + }, + { + "cell_type": "code", + "execution_count": 11, + "metadata": {}, + "outputs": [ + { + "data": { + "text/html": [ + "
\n", + "\n", + "\n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + " \n", + "
BondAngleTorsionvdWElectrostatics
OpenMM0.295288155.70129419.5632299.658623-15.215721
GROMACS0.295289155.70170619.5632409.666976-15.235476
LAMMPS0.440418155.70134319.5632339.630999-15.219631
Amber0.310034155.65442319.5719159.723198-15.212187
\n", + "
" + ], + "text/plain": [ + " Bond Angle Torsion vdW Electrostatics\n", + "OpenMM 0.295288 155.701294 19.563229 9.658623 -15.215721\n", + "GROMACS 0.295289 155.701706 19.563240 9.666976 -15.235476\n", + "LAMMPS 0.440418 155.701343 19.563233 9.630999 -15.219631\n", + "Amber 0.310034 155.654423 19.571915 9.723198 -15.212187" + ] + }, + "execution_count": 11, + "metadata": {}, + "output_type": "execute_result" + } + ], + "source": [ + "get_summary_data(interchange)" + ] + } + ], + "metadata": { + "kernelspec": { + "display_name": "Python 3 (ipykernel)", + "language": "python", + "name": "python3" + }, + "language_info": { + "codemirror_mode": { + "name": "ipython", + "version": 3 + }, + "file_extension": ".py", + "mimetype": "text/x-python", + "name": "python", + "nbconvert_exporter": "python", + "pygments_lexer": "ipython3", + "version": "3.8.12" + } + }, + "nbformat": 4, + "nbformat_minor": 4 +} From a964ae68da0ae58995156e2849037794767e19ae Mon Sep 17 00:00:00 2001 From: "Matthew W. Thompson" Date: Wed, 27 Oct 2021 09:16:41 -0500 Subject: [PATCH 03/24] Update conda env used by binder --- environment.yml | 4 ++++ 1 file changed, 4 insertions(+) diff --git a/environment.yml b/environment.yml index dd63dcf0a..9fe2d0b82 100644 --- a/environment.yml +++ b/environment.yml @@ -3,6 +3,7 @@ name: binder-env channels: - conda-forge - defaults + - bioconda dependencies: # Let conda pull in the toolkit and its deps - openff-toolkit @@ -13,3 +14,6 @@ dependencies: - qcportal - qcengine - openmmforcefields + - openff-interchange >=0.1.2 + - gromacs + - lammps From f594cd4e47ec5f12f9c63fc27e46cf7d22195fd1 Mon Sep 17 00:00:00 2001 From: "Matthew W. Thompson" Date: Wed, 27 Oct 2021 12:40:45 -0500 Subject: [PATCH 04/24] Refactor protein-ligand examples to not use ParmEd --- .../BRD4_inhibitor_benchmark.ipynb | 208 ++++-------------- .../README.md | 8 +- .../toluene_in_T4_lysozyme.ipynb | 182 ++++----------- 3 files changed, 86 insertions(+), 312 deletions(-) diff --git a/examples/using_smirnoff_with_amber_protein_forcefield/BRD4_inhibitor_benchmark.ipynb b/examples/using_smirnoff_with_amber_protein_forcefield/BRD4_inhibitor_benchmark.ipynb index a26ffabd0..f8955e98f 100644 --- a/examples/using_smirnoff_with_amber_protein_forcefield/BRD4_inhibitor_benchmark.ipynb +++ b/examples/using_smirnoff_with_amber_protein_forcefield/BRD4_inhibitor_benchmark.ipynb @@ -2,18 +2,18 @@ "cells": [ { "cell_type": "markdown", - "metadata": {}, + "metadata": { + "tags": [] + }, "source": [ "## Using SMIRNOFF with Amber on BRD4:inhibitor complexes: Exporting parameterized complexes to Amber, Gromacs, and CHARMM\n", "\n", - "This example applies SMIRNOFF-format parameters to BRD4 inhibitors from the [living review on binding free energy benchmark systems](https://www.annualreviews.org/doi/abs/10.1146/annurev-biophys-070816-033654) by Mobley and Gilson. The BRD4 system comes from the [accompanying GitHub repository](https://github.com/MobleyLab/benchmarksets/tree/master/input_files/BRD4).\n", - "\n", - "This example uses [ParmEd](http://parmed.github.io) to combine a protein parameterized with Amber with a ligand parameterized with SMIRNOFF. This example is meant to illustrate how to apply parameters to a single ligand, but it's also easy to process many ligands." + "This example applies SMIRNOFF-format parameters to BRD4 inhibitors from the [living review on binding free energy benchmark systems](https://www.annualreviews.org/doi/abs/10.1146/annurev-biophys-070816-033654) by Mobley and Gilson. The BRD4 system comes from the [accompanying GitHub repository](https://github.com/MobleyLab/benchmarksets/tree/master/input_files/BRD4)." ] }, { "cell_type": "code", - "execution_count": 1, + "execution_count": null, "metadata": {}, "outputs": [], "source": [ @@ -34,227 +34,103 @@ " open(filename, \"w\").write(r.text)" ] }, - { - "cell_type": "markdown", - "metadata": {}, - "source": [ - "### Parametrize a molecule with the SMIRNOFF-format \"Sage\" force field\n", - "\n", - "First, we parametrize the ligand with the Sage force field." - ] - }, - { - "cell_type": "markdown", - "metadata": {}, - "source": [ - "We use both a PDB file and an SDF file for the ligand." - ] - }, - { - "cell_type": "code", - "execution_count": 2, - "metadata": {}, - "outputs": [ - { - "name": "stderr", - "output_type": "stream", - "text": [ - "Warning: Unable to load toolkit 'OpenEye Toolkit'. The Open Force Field Toolkit does not require the OpenEye Toolkits, and can use RDKit/AmberTools instead. However, if you have a valid license for the OpenEye Toolkits, consider installing them for faster performance and additional file format support: https://docs.eyesopen.com/toolkits/python/quickstart-python/linuxosx.html OpenEye offers free Toolkit licenses for academics: https://www.eyesopen.com/academic-licensing\n" - ] - } - ], - "source": [ - "# Create an OpenFF Molecule object from the ligand SDF fiel\n", - "from openff.toolkit.topology import Molecule\n", - "\n", - "ligand_off_molecule = Molecule(\"ligand.sdf\")\n", - "\n", - "# Load the SMIRNOFF-format Sage force field\n", - "from openff.toolkit.typing.engines.smirnoff import ForceField\n", - "\n", - "force_field = ForceField(\"openff_unconstrained-2.0.0.offxml\")\n", - "\n", - "# Parametrize the ligand molecule by creating a Topology object from it\n", - "ligand_system = force_field.create_openmm_system(ligand_off_molecule.to_topology())" - ] - }, - { - "cell_type": "markdown", - "metadata": {}, - "source": [ - "and we convert the OpenMM `System` to a ParmEd `Structure` that we'll be able to mix with the protein\n", - "\n", - "
\n", - " Warning: ParmEd's Structure model is inspired by AMBER. Some information in an OpenMM System are not directly translatable into a Structure. In particular, long-range interaction treatment method (e.g., PME, CutoffPeriodic) and parameters (e.g., cutoff and cutoff switching distance, PME error tolerance) are known to be lost during the conversion.\n", - "
" - ] - }, { "cell_type": "code", - "execution_count": 3, + "execution_count": null, "metadata": {}, "outputs": [], "source": [ - "try:\n", - " from openmm import XmlSerializer, app, unit\n", - " from openmm.app import HBonds, NoCutoff, PDBFile\n", - "except ImportError:\n", - " from simtk import unit\n", - " from simtk.openmm import XmlSerializer, app\n", - " from simtk.openmm.app import HBonds, NoCutoff, PDBFile" + "from openff.interchange.components.interchange import Interchange\n", + "\n", + "from openff.toolkit.topology import Molecule\n", + "from openff.toolkit.typing.engines.smirnoff import ForceField" ] }, { "cell_type": "code", - "execution_count": 4, + "execution_count": null, "metadata": {}, "outputs": [], "source": [ - "# Read in the coordinates of the ligand from the PDB file\n", - "ligand_pdbfile = PDBFile(\"ligand.pdb\")\n", - "\n", - "# Convert OpenMM System object containing ligand parameters into a ParmEd Structure.\n", - "import parmed\n", - "\n", - "ligand_structure = parmed.openmm.load_topology(\n", - " ligand_pdbfile.topology, ligand_system, xyz=ligand_pdbfile.positions\n", - ")" - ] - }, - { - "cell_type": "markdown", - "metadata": {}, - "source": [ - "### Create a ParmEd `Structure` of an AMBER-parametrized receptor\n", - "\n", - "We have to create a ParmEd `Structure` containing positions and parameters of the receptor (BRD4) that we will then combine with the positions and parameters in the ligand `Structure` we created above. \n", + "ligand_molecule = Molecule(\"ligand.sdf\")\n", + "sage = ForceField(\"openff-2.0.0.offxml\")\n", "\n", - "First, we use OpenMM to assign Amber14 parameters using OpenMM, but you can also create a `Structure` [from Amber `prmtop` and `inpcrd` files](http://parmed.github.io/ParmEd/html/structure.html).\n", + "ligand = Interchange.from_smirnoff(\n", + " force_field=sage, topology=ligand_molecule.to_topology()\n", + ")\n", "\n", - "
\n", - " Note: If you already have AMBER (prmtop/inpcrd), GROMACS (top/gro), or any other file specifying the protein parameters supported by ParmEd, you can simply load the files directly into a Structure using ParmEd's functionalities. See https://parmed.github.io/ParmEd/html/readwrite.html .\n", - "
" - ] - }, - { - "cell_type": "code", - "execution_count": 5, - "metadata": {}, - "outputs": [], - "source": [ - "receptor_pdbfile = PDBFile(\"receptor.pdb\")\n", + "receptor_molecule = Molecule.from_pdb(\"receptor.pdb\")\n", "\n", - "# Load the AMBER protein force field through OpenMM.\n", - "omm_forcefield = app.ForceField(\"amber14-all.xml\")\n", + "ff14sb = ForceField(\"ff14sb_off_impropers_0.0.1.offxml\")\n", + "ff14sb.deregister_parameter_handler(\"ImproperTorsions\")\n", "\n", - "# Parameterize the protein.\n", - "receptor_system = omm_forcefield.createSystem(receptor_pdbfile.topology)\n", + "receptor = Interchange.from_smirnoff(\n", + " force_field=ff14sb, topology=receptor_molecule.to_topology()\n", + ")\n", "\n", - "# Convert the protein System into a ParmEd Structure.\n", - "receptor_structure = parmed.openmm.load_topology(\n", - " receptor_pdbfile.topology, receptor_system, xyz=receptor_pdbfile.positions\n", - ")" - ] - }, - { - "cell_type": "markdown", - "metadata": {}, - "source": [ - "### Combine receptor and ligand structures\n", + "complex_system = receptor + ligand\n", "\n", - "We can then merge the receptor and ligand `Structure` objects into a single `Structure` containing both positions and parameters." - ] - }, - { - "cell_type": "code", - "execution_count": 6, - "metadata": {}, - "outputs": [], - "source": [ - "complex_structure = receptor_structure + ligand_structure" + "# Make sure ...\n", + "# complex.box = receptor.box\n", + "# complex.positions = np.vstack([receptor.positions, ligand.positions])" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ - "### Export to OpenMM\n", - "\n", - "Once we have the `Structure` for the complex, containing both positions and parameters, we can chose to create an OpenMM `System` object that we can simulate using [OpenMM](http://openmm.org):" - ] - }, - { - "cell_type": "code", - "execution_count": 7, - "metadata": {}, - "outputs": [], - "source": [ - "# Convert the Structure to an OpenMM System in vacuum.\n", - "complex_system = complex_structure.createSystem(\n", - " nonbondedMethod=NoCutoff,\n", - " nonbondedCutoff=9.0 * unit.angstrom,\n", - " constraints=HBonds,\n", - " removeCMMotion=False,\n", - ")" + "### Export to OpenMM" ] }, { "cell_type": "code", - "execution_count": 8, + "execution_count": null, "metadata": {}, "outputs": [], "source": [ - "# Export the System to an OpenMM System XML and PDB file\n", - "complex_structure.save(\"complex.pdb\", overwrite=True)\n", - "\n", - "with open(\"complex.xml\", \"w\") as f:\n", - " f.write(XmlSerializer.serialize(complex_system))" + "complex_system.to_openmm()" ] }, { "cell_type": "markdown", - "metadata": {}, + "metadata": { + "tags": [] + }, "source": [ - "### Export to Amber\n", - "\n", - "We can also export the `Structure` to Amber `prmtop` and `inpcrd` files using [ParmEd](http://parmed.github.io/ParmEd/html/amber.html):" + "### Export to Amber" ] }, { "cell_type": "code", - "execution_count": 9, + "execution_count": null, "metadata": {}, "outputs": [], "source": [ - "# Export the Structure to AMBER files\n", - "complex_structure.save(\"complex.prmtop\", overwrite=True)\n", - "complex_structure.save(\"complex.inpcrd\", overwrite=True)" + "complex_system.to_inpcrd(\"complex.inpcrd\")\n", + "complex_system.to_prmtop(\"complex.prmtop\")" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ - "### Export to gromacs\n", - "\n", - "We can export the `Structure` to gromacs `gro` and `top` files using [ParmEd](http://parmed.github.io/ParmEd/html/gromacs.html):" + "### Export to GROMACS" ] }, { "cell_type": "code", - "execution_count": 10, + "execution_count": null, "metadata": {}, "outputs": [], "source": [ - "# Export the Structure to Gromacs files\n", - "complex_structure.save(\"complex.gro\", overwrite=True)\n", - "complex_structure.save(\"complex.top\", overwrite=True)" + "complex_structure.to_gro(\"complex.gro\")\n", + "complex_structure.to_top(\"complex.top\")" ] } ], "metadata": { "kernelspec": { - "display_name": "Python 3", + "display_name": "Python 3 (ipykernel)", "language": "python", "name": "python3" }, @@ -268,7 +144,7 @@ "name": "python", "nbconvert_exporter": "python", "pygments_lexer": "ipython3", - "version": "3.7.10" + "version": "3.9.7" } }, "nbformat": 4, diff --git a/examples/using_smirnoff_with_amber_protein_forcefield/README.md b/examples/using_smirnoff_with_amber_protein_forcefield/README.md index 1fda09eb2..a22d271b8 100644 --- a/examples/using_smirnoff_with_amber_protein_forcefield/README.md +++ b/examples/using_smirnoff_with_amber_protein_forcefield/README.md @@ -1,17 +1,15 @@ -## Combining a SMIRNOFF parameterized small molecule with an Amber parameterized protein using ParmEd +## Combining a small molecule with Sage and a protein parameterized with an Amber protein force field [![Binder](https://mybinder.org/badge_logo.svg)](https://mybinder.org/v2/gh/openforcefield/openff-toolkit/stable?filepath=examples%2Fusing_smirnoff_with_amber_protein_forcefield) -These examples illustrate how the [ParmEd](http://parmed.github.io/ParmEd/html/index.html) utility can be used to merge a small molecule parameterized by SMIRNOFF with a traditionally parameterized protein (or other biopolymer) to create a fully parameterized protein-ligand system. - ### BRD4:inhibitor complex [![Binder](https://mybinder.org/badge_logo.svg)](https://mybinder.org/v2/gh/openforcefield/openff-toolkit/stable?filepath=examples%2Fusing_smirnoff_with_amber_protein_forcefield%2FBRD4_inhibitor_benchmark.ipynb) -[`BRD4_inhibitor_benchmark.ipynb`](BRD4_inhibitor_benchmark.ipynb) contains an example illustrating applying SMIRNOFF and Amber14 parameters to a BRD4:inhibitor complex taken from the [free energy benchmark systems living review](https://www.annualreviews.org/doi/abs/10.1146/annurev-biophys-070816-033654) [GitHub repo](https://github.com/MobleyLab/benchmarksets/tree/master/input_files/BRD4). +[`BRD4_inhibitor_benchmark.ipynb`](BRD4_inhibitor_benchmark.ipynb) contains an example illustrating applying Sage and ff14SB parameters to a BRD4:inhibitor complex taken from the [free energy benchmark systems living review](https://www.annualreviews.org/doi/abs/10.1146/annurev-biophys-070816-033654) [GitHub repo](https://github.com/MobleyLab/benchmarksets/tree/master/input_files/BRD4). ### Toluene in complex with T4 lysozyme L99A in TIP3P-FB water [![Binder](https://mybinder.org/badge_logo.svg)](https://mybinder.org/v2/gh/openforcefield/openff-toolkit/stable?filepath=examples%2Fusing_smirnoff_with_amber_protein_forcefield%2Ftoluene_in_T4_lysozyme.ipynb) -[`toluene_in_T4_lysozyme.ipynb`](toluene_in_T4_lysozyme.ipynb) contains an example illustrating applying SMIRNOFF and Amber99SB-ILDN parameters to toluene complexed with T4 lysozyme L99A in a TIP3P-FB water box. +[`toluene_in_T4_lysozyme.ipynb`](toluene_in_T4_lysozyme.ipynb) contains an example illustrating applying Sage and ff14SB parameters to toluene complexed with T4 lysozyme L99A in a TIP3P water box. diff --git a/examples/using_smirnoff_with_amber_protein_forcefield/toluene_in_T4_lysozyme.ipynb b/examples/using_smirnoff_with_amber_protein_forcefield/toluene_in_T4_lysozyme.ipynb index 8427d10b4..2881831dc 100644 --- a/examples/using_smirnoff_with_amber_protein_forcefield/toluene_in_T4_lysozyme.ipynb +++ b/examples/using_smirnoff_with_amber_protein_forcefield/toluene_in_T4_lysozyme.ipynb @@ -4,11 +4,7 @@ "cell_type": "markdown", "metadata": {}, "source": [ - "## Create a system mixing SMIRNOFF and non-SMIRNOFF-formatted force fields\n", - "\n", - "This example shows how to create a receptor-ligand OpenMM `System` where the ligand (toluene) is parametrized with a SMIRNOFF force field and the protein (T4 Lysozyme) solvated in water is assigned AMBER and TIP3P-FB parameters through the ParmEd library.\n", - "\n", - "We'll need two PDB files. One for the ligand in vacuum, and one for the solvated protein without ligand. The coordinates of the protein-ligand complex will be determined by these PDB files, so their positions needs to be consistent with the ligand being positioned in the binding pocket if this is the desired initial configuration of your simulation." + "## Create a system mixing SMIRNOFF and non-SMIRNOFF-formatted force fields" ] }, { @@ -22,186 +18,90 @@ }, { "cell_type": "code", - "execution_count": 1, + "execution_count": null, "metadata": {}, - "outputs": [ - { - "name": "stderr", - "output_type": "stream", - "text": [ - "Warning: Unable to load toolkit 'OpenEye Toolkit'. The Open Force Field Toolkit does not require the OpenEye Toolkits, and can use RDKit/AmberTools instead. However, if you have a valid license for the OpenEye Toolkits, consider installing them for faster performance and additional file format support: https://docs.eyesopen.com/toolkits/python/quickstart-python/linuxosx.html OpenEye offers free Toolkit licenses for academics: https://www.eyesopen.com/academic-licensing\n" - ] - } - ], + "outputs": [], "source": [ "try:\n", - " from openmm import app, unit\n", - " from openmm.app import HBonds, NoCutoff, PDBFile\n", + " from openmm import app\n", + " from openmm import unit as openmm_unit\n", "except ImportError:\n", - " from simtk import unit\n", + " from simtk import unit as openmm_unit\n", " from simtk.openmm import app\n", - " from simtk.openmm.app import PDBFile, HBonds, NoCutoff\n", - "\n", - "from openff.toolkit.topology import Molecule, Topology\n", - "from openff.toolkit.typing.engines.smirnoff import ForceField\n", - "from openff.toolkit.utils import get_data_file_path" - ] - }, - { - "cell_type": "code", - "execution_count": 2, - "metadata": {}, - "outputs": [], - "source": [ - "# Create an OpenFF Topology of toluene from a pdb file.\n", - "toluene_pdb_file_path = get_data_file_path(\"molecules/toluene.pdb\")\n", - "toluene_pdbfile = PDBFile(toluene_pdb_file_path)\n", - "toluene = Molecule.from_smiles(\"Cc1ccccc1\")\n", - "off_topology = Topology.from_openmm(\n", - " openmm_topology=toluene_pdbfile.topology, unique_molecules=[toluene]\n", - ")\n", - "\n", - "# Load the Sage force field from disk.\n", - "force_field = ForceField(\"openff_unconstrained-2.0.0.offxml\")\n", "\n", - "# Parametrize the toluene molecule.\n", - "toluene_system = force_field.create_openmm_system(off_topology)" - ] - }, - { - "cell_type": "markdown", - "metadata": {}, - "source": [ - "and we convert the OpenMM `System` to a ParmEd `Structure` that we'll be able to mix with the protein.\n", - "\n", - "
\n", - " Warning: ParmEd's Structure model is inspired by AMBER. Some information in an OpenMM System are not directly translatable into a Structure. In particular, long-range interaction treatment method (e.g., PME, CutoffPeriodic) and parameters (e.g., cutoff and cutoff switching distance, PME error tolerance) are known to be lost during the conversion.\n", - "
" - ] - }, - { - "cell_type": "code", - "execution_count": 3, - "metadata": {}, - "outputs": [], - "source": [ - "import parmed\n", - "\n", - "# Convert OpenMM System into a ParmEd Structure.\n", - "toluene_structure = parmed.openmm.load_topology(\n", - " toluene_pdbfile.topology, toluene_system, xyz=toluene_pdbfile.positions\n", - ")" - ] - }, - { - "cell_type": "markdown", - "metadata": {}, - "source": [ - "### Create a ParmEd `Structure` of an AMBER-parametrized receptor in TIP3P-FB water\n", + "import numpy as np\n", + "from openff.interchange.components.interchange import Interchange\n", + "from openff.units import unit\n", "\n", - "We have to create a ParmEd `Structure` of the receptor (T4 Lysozyme) to combine to the toluene `Structure`. Here we assign `amber99sbildn` to a T4 Lysozyme receptor solvated in TIP3P-FB water using OpenMM.\n", + "from openff.toolkit.topology import Molecule\n", + "from openff.toolkit.typing.engines.smirnoff import ForceField\n", + "from openff.toolkit.utils import get_data_file_path\n", "\n", - "First, we load the parameters and the PDB file including positions for the protein, water, and ion atoms.\n", "\n", - "
\n", - " Note: If you already have AMBER (prmtop/inpcrd), GROMACS (top/gro), or any other file supported by ParmEd specifying the parameters for the solvated protein, you can simply load the files directly into a Structure using ParmEd's functionalities. See https://parmed.github.io/ParmEd/html/readwrite.html .\n", - "
" + "def converter(x):\n", + " return x.value_in_unit(openmm_unit.nanometer)" ] }, { "cell_type": "code", - "execution_count": 4, - "metadata": {}, + "execution_count": null, + "metadata": { + "tags": [] + }, "outputs": [], "source": [ - "# Load the AMBER protein force field parameters through OpenMM.\n", - "omm_forcefield = app.ForceField(\"amber99sbildn.xml\", \"tip3pfb.xml\")\n", - "\n", - "# Load the solvated receptor from a PDB file.\n", - "t4_pdb_file_path = get_data_file_path(\"systems/test_systems/T4_lysozyme_water_ions.pdb\")\n", - "t4_pdb_file = PDBFile(t4_pdb_file_path)\n", - "\n", - "# Obtain the updated OpenMM Topology and positions.\n", - "omm_topology = t4_pdb_file.getTopology()\n", - "positions = t4_pdb_file.getPositions()" - ] - }, - { - "cell_type": "markdown", - "metadata": {}, - "source": [ - "We then create a parametrized OpenMM `System` and convert it to a `Structure`.\n", + "toluene_positions = app.PDBFile(get_data_file_path(\"molecules/toluene.pdb\")).positions\n", + "toluene_topology = Molecule.from_smiles(\"Cc1ccccc1\").to_topology()\n", "\n", - "Note the `rigidWater=False` argument in `ForceField.createSystem()`. This is necessary to work around a problem araising with ParmEd in reading the parameters of constrained bonds from an OpenMM `System` (see https://github.com/openforcefield/openff-toolkit/issues/259 for more details). We'll re-add the hydrogen bond constraints when we'll create the `System` for the complex.\n", + "sage = ForceField(\"openff-2.0.0.offxml\")\n", "\n", - "
\n", - " Note: If you don't solvate the system or if you load it directly from AMBER, GROMACS, or other files directly to ParmEd, you won't need extra precautions.\n", - "
" + "toluene = Interchange.from_smirnoff(force_field=sage, topology=toluene_topology)\n", + "toluene.positions = np.frompyfunc(converter, 1, 1)(toluene_positions).astype(float)" ] }, { "cell_type": "code", - "execution_count": 5, + "execution_count": null, "metadata": {}, "outputs": [], "source": [ - "# Parameterize the protein.\n", - "t4_system = omm_forcefield.createSystem(omm_topology, rigidWater=False)\n", + "ff14sb = ForceField(\"ff14sb_off_impropers_0.0.1.offxml\")\n", + "ff14sb.deregister_parameter_handler(\"ImproperTorsions\")\n", "\n", - "# Convert the protein System into a ParmEd Structure.\n", - "t4_structure = parmed.openmm.load_topology(omm_topology, t4_system, xyz=positions)" - ] - }, - { - "cell_type": "markdown", - "metadata": {}, - "source": [ - "### Combine receptor and ligand structures\n", - "\n", - "We can then merge the receptor and ligand `Structure` objects to form the complex. Note that the coordinates of protein and ligand in this example are determined by the PDB file, and they are already consistent with the ligand being positioned in the binding pocket.\n", + "# TODO: This PDB file also has water an ions, which are not yet intelligently parsed here.\n", + "# Unclear if we should have a Topology.from_pdb, just drop the waters, hack through OpenMM ...\n", + "t4_lysozyme_molecule = Molecule.from_pdb(\n", + " get_data_file_path(\"systems/test_systems/T4_lysozyme_water_ions.pdb\")\n", + ")\n", "\n", - "
\n", - " Note: If you want to include water molecules in the binding site, you will have to be careful to place them so that they won't create steric clashes once the ligand is inserted.\n", - "
" + "t4_lysozyme = Interchange.from_smirnoff(\n", + " force_field=ff14sb, topology=t4_lysozyme_molecule.to_topoplogy()\n", + ")\n", + "t4_lysozyme.positions = t4_lysozyme_molecule.conformers[0]" ] }, { "cell_type": "code", - "execution_count": 6, + "execution_count": null, "metadata": {}, "outputs": [], "source": [ - "complex_structure = toluene_structure + t4_structure" - ] - }, - { - "cell_type": "markdown", - "metadata": {}, - "source": [ - "### Convert back the structure into an OpenMM System\n", - "\n", - "Once we have the `Structure` of the complex, we can chose to create a `System` object that we can simulate with OpenMM." + "complex_system = t4_lysozyme + toluene" ] }, { "cell_type": "code", - "execution_count": 7, + "execution_count": null, "metadata": {}, "outputs": [], "source": [ - "# Convert the Structure to an OpenMM System in vacuum.\n", - "complex_system = complex_structure.createSystem(\n", - " nonbondedMethod=NoCutoff,\n", - " nonbondedCutoff=9.0 * unit.angstrom,\n", - " constraints=HBonds,\n", - " removeCMMotion=False,\n", - ")" + "complex_system.to_openmm()" ] } ], "metadata": { "kernelspec": { - "display_name": "Python 3", + "display_name": "Python 3 (ipykernel)", "language": "python", "name": "python3" }, @@ -215,7 +115,7 @@ "name": "python", "nbconvert_exporter": "python", "pygments_lexer": "ipython3", - "version": "3.7.10" + "version": "3.9.7" } }, "nbformat": 4, From dd28d089b678e375dcdf394bf4995315f4bbee2e Mon Sep 17 00:00:00 2001 From: "Matthew W. Thompson" Date: Wed, 27 Oct 2021 15:42:50 -0500 Subject: [PATCH 05/24] Add back ParmEd-based example --- .../convert_to_amber_gromacs.ipynb | 285 ++++++++++++++++++ 1 file changed, 285 insertions(+) create mode 100644 examples/using_smirnoff_in_amber_or_gromacs/convert_to_amber_gromacs.ipynb diff --git a/examples/using_smirnoff_in_amber_or_gromacs/convert_to_amber_gromacs.ipynb b/examples/using_smirnoff_in_amber_or_gromacs/convert_to_amber_gromacs.ipynb new file mode 100644 index 000000000..1bff64170 --- /dev/null +++ b/examples/using_smirnoff_in_amber_or_gromacs/convert_to_amber_gromacs.ipynb @@ -0,0 +1,285 @@ +{ + "cells": [ + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "## Convert Open Forcefield System to AMBER and GROMACS input files\n", + "\n", + "The Open Forcefield Toolkit can create parametrized `System` objects that can be natively simulated with OpenMM. This example shows how you can convert such an OpenMM `System` into AMBER prmtop/inpcrd and GROMACS top/gro input files through the ParmEd library." + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "### Create an OpenMM System\n", + "\n", + "We start by loading a PDB file containing one copy of ethanol and cyclohexane. Our goal is to create an OFF `Topology` object describing this system that we can parametrize with the SMIRNOFF-format \"Sage\" force field.\n", + "\n", + "The two `Molecule` objects created from the SMILES strings can contain information such as partial charges and stereochemistry that is not included in an OpenMM topology. In this example, partial charges are not explicitly given, and `ForceField` will assign AM1/BCC charges as specified by the \"Sage\" force field. Note that the Open Force Field Toolkit produces deterministic partial charges that do not depend on the input conformation of parameterized molecules. See the [FAQ](https://open-forcefield-toolkit.readthedocs.io/en/latest/faq.html#the-partial-charges-generated-by-the-toolkit-don-t-seem-to-depend-on-the-molecule-s-conformation-is-this-a-bug) for more information." + ] + }, + { + "cell_type": "code", + "execution_count": 1, + "metadata": {}, + "outputs": [ + { + "name": "stderr", + "output_type": "stream", + "text": [ + "Warning: Unable to load toolkit 'OpenEye Toolkit'. The Open Force Field Toolkit does not require the OpenEye Toolkits, and can use RDKit/AmberTools instead. However, if you have a valid license for the OpenEye Toolkits, consider installing them for faster performance and additional file format support: https://docs.eyesopen.com/toolkits/python/quickstart-python/linuxosx.html OpenEye offers free Toolkit licenses for academics: https://www.eyesopen.com/academic-licensing\n" + ] + } + ], + "source": [ + "try:\n", + " from openmm.app import PDBFile\n", + "except ImportError:\n", + " from simtk.openmm.app import PDBFile\n", + "\n", + "from openff.toolkit.topology import Molecule, Topology\n", + "from openff.toolkit.typing.engines.smirnoff import ForceField" + ] + }, + { + "cell_type": "code", + "execution_count": 2, + "metadata": {}, + "outputs": [], + "source": [ + "ethanol = Molecule.from_smiles(\"CCO\")\n", + "cyclohexane = Molecule.from_smiles(\"C1CCCCC1\")\n", + "\n", + "# Obtain the OpenMM Topology object from the PDB file.\n", + "pdbfile = PDBFile(\"1_cyclohexane_1_ethanol.pdb\")\n", + "omm_topology = pdbfile.topology\n", + "\n", + "# Create the Open Forcefield Topology.\n", + "off_topology = Topology.from_openmm(\n", + " omm_topology, unique_molecules=[ethanol, cyclohexane]\n", + ")" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "Now we parametrize the OFF `Topology` to create an OpenMM `System`. Since ParmEd will run with the `constraints=HBonds` keyword later, we use the _unconstrained_ version of the Sage force field here." + ] + }, + { + "cell_type": "code", + "execution_count": 3, + "metadata": {}, + "outputs": [], + "source": [ + "# Load the \"Sage\" force field.\n", + "forcefield = ForceField(\"openff_unconstrained-2.0.0.offxml\")\n", + "omm_system = forcefield.create_openmm_system(off_topology)" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "### Convert OpenMM System to AMBER and GROMACS files\n", + "\n", + "\n", + "First, we convert the OpenMM `System` into a ParmEd `Structure`. We'll use the atom positions in the PDB to create the coordinate files.\n", + "\n", + "
\n", + " Warning: ParmEd's Structure model is inspired by AMBER, and some information in an OpenMM System are not directly translatable into a Structure. In particular, as of today (4/2/2019), long-range interaction treatment method (e.g., PME, CutoffPeriodic) and parameters (e.g., cutoff and cutoff switching distance, PME error tolerance) are known to be lost during the conversion.\n", + "
" + ] + }, + { + "cell_type": "code", + "execution_count": 4, + "metadata": {}, + "outputs": [], + "source": [ + "import parmed\n", + "\n", + "# Convert OpenMM System to a ParmEd structure.\n", + "parmed_structure = parmed.openmm.load_topology(\n", + " omm_topology, omm_system, pdbfile.positions\n", + ")" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "We can then use ParmEd to convert an OpenMM `System` to prmtop/inpcrd or top/gro files that can be read by AMBER and GROMACS respectively. ParmEd is capable of converting parametrized files to other formats as well. For further information, see ParmEd's documentation: https://parmed.github.io/ParmEd/html/readwrite.html ." + ] + }, + { + "cell_type": "code", + "execution_count": 5, + "metadata": {}, + "outputs": [], + "source": [ + "# Export AMBER files.\n", + "parmed_structure.save(\"system.prmtop\", overwrite=True)\n", + "parmed_structure.save(\"system.inpcrd\", overwrite=True)\n", + "\n", + "# Export GROMACS files.\n", + "parmed_structure.save(\"system.top\", overwrite=True)\n", + "parmed_structure.save(\"system.gro\", overwrite=True)" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "### Validate the conversion\n", + "\n", + "ParmEd is generally a reliable and robust library, but we can easily check that everything went as expected during the conversion by loading the exported files into an OpenMM `System` and comparing it with the original. Note that you'll have to specify the correct nonbonded method and cutoff settings for the energy comparison to make sense since they are not included in the AMBER prmtop (or GROMACS top/gro) files." + ] + }, + { + "cell_type": "code", + "execution_count": 6, + "metadata": {}, + "outputs": [ + { + "name": "stdout", + "output_type": "stream", + "text": [ + "0.9 nm\n", + "False\n", + "True\n" + ] + } + ], + "source": [ + "try:\n", + " import openmm\n", + "except ImportError:\n", + " from simtk import openmm\n", + "\n", + "for force in omm_system.getForces():\n", + " if isinstance(force, openmm.NonbondedForce):\n", + " break\n", + "print(force.getCutoffDistance())\n", + "print(force.getUseSwitchingFunction())\n", + "print(force.getNonbondedMethod() == openmm.NonbondedForce.PME)" + ] + }, + { + "cell_type": "code", + "execution_count": 7, + "metadata": {}, + "outputs": [], + "source": [ + "try:\n", + " from openmm import unit\n", + " from openmm.app import PME, HBonds\n", + "except ImportError:\n", + " from simtk import unit\n", + " from simtk.openmm.app import PME, HBonds\n", + "\n", + "from openff.toolkit.tests.utils import (\n", + " compare_system_energies,\n", + " compare_system_parameters,\n", + ")\n", + "\n", + "# Load the prmtop/inpcrd files into a ParmEd Structure.as an OpenMM System object.\n", + "amber_structure = parmed.load_file(\"system.prmtop\", \"system.inpcrd\")\n", + "\n", + "# Convert the Structure to an OpenMM System. Note that by\n", + "# default ParmEd will add a CMMotionRemover force to the\n", + "# System, and won't constrain the hydrogen bonds.\n", + "amber_system = amber_structure.createSystem(\n", + " nonbondedMethod=PME,\n", + " nonbondedCutoff=9.0 * unit.angstrom,\n", + " switchDistance=0.0 * unit.angstrom,\n", + " constraints=HBonds,\n", + " removeCMMotion=False,\n", + ")\n", + "\n", + "# Compare the parameters of the original and converted Systems.\n", + "# This raises FailedParameterComparisonError if the comparison fails.\n", + "compare_system_parameters(omm_system, amber_system)\n", + "\n", + "# Compare the energies by force.\n", + "# This raises FailedEnergyComparisonError if the comparison fails.\n", + "amber_energies, omm_energies = compare_system_energies(\n", + " amber_system,\n", + " omm_system,\n", + " amber_structure.positions,\n", + " amber_structure.box_vectors,\n", + " rtol=1e-3,\n", + ")" + ] + }, + { + "cell_type": "code", + "execution_count": 8, + "metadata": {}, + "outputs": [ + { + "name": "stdout", + "output_type": "stream", + "text": [ + "System loaded from AMBER files:\n", + "-------------------------------\n", + "{'HarmonicAngleForce': Quantity(value=155.70130920410156, unit=kilojoule/mole),\n", + " 'HarmonicBondForce': Quantity(value=0.4403935372829437, unit=kilojoule/mole),\n", + " 'NonbondedForce': Quantity(value=2.8259036956507657, unit=kilojoule/mole),\n", + " 'PeriodicTorsionForce': Quantity(value=19.563232421875, unit=kilojoule/mole)}\n", + "\n", + "Original OpenMM System:\n", + "-----------------------\n", + "{'HarmonicAngleForce': Quantity(value=155.70130920410156, unit=kilojoule/mole),\n", + " 'HarmonicBondForce': Quantity(value=0.44038787484169006, unit=kilojoule/mole),\n", + " 'NonbondedForce': Quantity(value=2.8251933539722245, unit=kilojoule/mole),\n", + " 'PeriodicTorsionForce': Quantity(value=19.563230514526367, unit=kilojoule/mole)}\n" + ] + } + ], + "source": [ + "# Pretty-print the energies by component.\n", + "from pprint import pprint\n", + "\n", + "print(\"System loaded from AMBER files:\")\n", + "print(\"-------------------------------\")\n", + "pprint(amber_energies)\n", + "\n", + "print(\"\\nOriginal OpenMM System:\")\n", + "print(\"-----------------------\")\n", + "pprint(omm_energies)" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [] + } + ], + "metadata": { + "kernelspec": { + "display_name": "Python 3", + "language": "python", + "name": "python3" + }, + "language_info": { + "codemirror_mode": { + "name": "ipython", + "version": 3 + }, + "file_extension": ".py", + "mimetype": "text/x-python", + "name": "python", + "nbconvert_exporter": "python", + "pygments_lexer": "ipython3", + "version": "3.8.6" + } + }, + "nbformat": 4, + "nbformat_minor": 4 +} From 72d8352b65d0f4e4e00045dea28ceaf61f66f3e9 Mon Sep 17 00:00:00 2001 From: "Matthew W. Thompson" Date: Wed, 27 Oct 2021 15:55:27 -0500 Subject: [PATCH 06/24] Fix other files --- .../smirnoff_amber_gromacs.ipynb | 29 ++- .../BRD4_inhibitor_benchmark.ipynb | 208 ++++++++++++++---- .../toluene_in_T4_lysozyme.ipynb | 182 +++++++++++---- 3 files changed, 320 insertions(+), 99 deletions(-) diff --git a/examples/using_smirnoff_in_amber_or_gromacs/smirnoff_amber_gromacs.ipynb b/examples/using_smirnoff_in_amber_or_gromacs/smirnoff_amber_gromacs.ipynb index 5b68f685e..4bac9f501 100644 --- a/examples/using_smirnoff_in_amber_or_gromacs/smirnoff_amber_gromacs.ipynb +++ b/examples/using_smirnoff_in_amber_or_gromacs/smirnoff_amber_gromacs.ipynb @@ -24,15 +24,7 @@ "cell_type": "code", "execution_count": 1, "metadata": {}, - "outputs": [ - { - "name": "stderr", - "output_type": "stream", - "text": [ - "Warning: importing 'simtk.openmm' is deprecated. Import 'openmm' instead.\n" - ] - } - ], + "outputs": [], "source": [ "try:\n", " from openmm import app\n", @@ -51,7 +43,7 @@ { "data": { "text/plain": [ - "" + "" ] }, "execution_count": 2, @@ -94,7 +86,7 @@ { "data": { "text/plain": [ - "" + "" ] }, "execution_count": 3, @@ -124,7 +116,7 @@ { "data": { "text/plain": [ - " >" + " >" ] }, "execution_count": 4, @@ -151,6 +143,13 @@ "execution_count": 5, "metadata": {}, "outputs": [ + { + "name": "stderr", + "output_type": "stream", + "text": [ + "Warning: importing 'simtk.openmm' is deprecated. Import 'openmm' instead.\n" + ] + }, { "data": { "text/plain": [ @@ -233,7 +232,7 @@ " 'Angle': 155.7017059326172 ,\n", " 'Torsion': 19.56324005126953 ,\n", " 'vdW': 9.66697609052062 ,\n", - " 'Electrostatics': -15.235475540161133 }" + " 'Electrostatics': -15.235476016998291 }" ] }, "execution_count": 8, @@ -274,9 +273,7 @@ { "cell_type": "code", "execution_count": 10, - "metadata": { - "scrolled": false - }, + "metadata": {}, "outputs": [ { "data": { diff --git a/examples/using_smirnoff_with_amber_protein_forcefield/BRD4_inhibitor_benchmark.ipynb b/examples/using_smirnoff_with_amber_protein_forcefield/BRD4_inhibitor_benchmark.ipynb index f8955e98f..a26ffabd0 100644 --- a/examples/using_smirnoff_with_amber_protein_forcefield/BRD4_inhibitor_benchmark.ipynb +++ b/examples/using_smirnoff_with_amber_protein_forcefield/BRD4_inhibitor_benchmark.ipynb @@ -2,18 +2,18 @@ "cells": [ { "cell_type": "markdown", - "metadata": { - "tags": [] - }, + "metadata": {}, "source": [ "## Using SMIRNOFF with Amber on BRD4:inhibitor complexes: Exporting parameterized complexes to Amber, Gromacs, and CHARMM\n", "\n", - "This example applies SMIRNOFF-format parameters to BRD4 inhibitors from the [living review on binding free energy benchmark systems](https://www.annualreviews.org/doi/abs/10.1146/annurev-biophys-070816-033654) by Mobley and Gilson. The BRD4 system comes from the [accompanying GitHub repository](https://github.com/MobleyLab/benchmarksets/tree/master/input_files/BRD4)." + "This example applies SMIRNOFF-format parameters to BRD4 inhibitors from the [living review on binding free energy benchmark systems](https://www.annualreviews.org/doi/abs/10.1146/annurev-biophys-070816-033654) by Mobley and Gilson. The BRD4 system comes from the [accompanying GitHub repository](https://github.com/MobleyLab/benchmarksets/tree/master/input_files/BRD4).\n", + "\n", + "This example uses [ParmEd](http://parmed.github.io) to combine a protein parameterized with Amber with a ligand parameterized with SMIRNOFF. This example is meant to illustrate how to apply parameters to a single ligand, but it's also easy to process many ligands." ] }, { "cell_type": "code", - "execution_count": null, + "execution_count": 1, "metadata": {}, "outputs": [], "source": [ @@ -35,102 +35,226 @@ ] }, { - "cell_type": "code", - "execution_count": null, + "cell_type": "markdown", "metadata": {}, - "outputs": [], "source": [ - "from openff.interchange.components.interchange import Interchange\n", + "### Parametrize a molecule with the SMIRNOFF-format \"Sage\" force field\n", "\n", + "First, we parametrize the ligand with the Sage force field." + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "We use both a PDB file and an SDF file for the ligand." + ] + }, + { + "cell_type": "code", + "execution_count": 2, + "metadata": {}, + "outputs": [ + { + "name": "stderr", + "output_type": "stream", + "text": [ + "Warning: Unable to load toolkit 'OpenEye Toolkit'. The Open Force Field Toolkit does not require the OpenEye Toolkits, and can use RDKit/AmberTools instead. However, if you have a valid license for the OpenEye Toolkits, consider installing them for faster performance and additional file format support: https://docs.eyesopen.com/toolkits/python/quickstart-python/linuxosx.html OpenEye offers free Toolkit licenses for academics: https://www.eyesopen.com/academic-licensing\n" + ] + } + ], + "source": [ + "# Create an OpenFF Molecule object from the ligand SDF fiel\n", "from openff.toolkit.topology import Molecule\n", - "from openff.toolkit.typing.engines.smirnoff import ForceField" + "\n", + "ligand_off_molecule = Molecule(\"ligand.sdf\")\n", + "\n", + "# Load the SMIRNOFF-format Sage force field\n", + "from openff.toolkit.typing.engines.smirnoff import ForceField\n", + "\n", + "force_field = ForceField(\"openff_unconstrained-2.0.0.offxml\")\n", + "\n", + "# Parametrize the ligand molecule by creating a Topology object from it\n", + "ligand_system = force_field.create_openmm_system(ligand_off_molecule.to_topology())" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "and we convert the OpenMM `System` to a ParmEd `Structure` that we'll be able to mix with the protein\n", + "\n", + "
\n", + " Warning: ParmEd's Structure model is inspired by AMBER. Some information in an OpenMM System are not directly translatable into a Structure. In particular, long-range interaction treatment method (e.g., PME, CutoffPeriodic) and parameters (e.g., cutoff and cutoff switching distance, PME error tolerance) are known to be lost during the conversion.\n", + "
" + ] + }, + { + "cell_type": "code", + "execution_count": 3, + "metadata": {}, + "outputs": [], + "source": [ + "try:\n", + " from openmm import XmlSerializer, app, unit\n", + " from openmm.app import HBonds, NoCutoff, PDBFile\n", + "except ImportError:\n", + " from simtk import unit\n", + " from simtk.openmm import XmlSerializer, app\n", + " from simtk.openmm.app import HBonds, NoCutoff, PDBFile" ] }, { "cell_type": "code", - "execution_count": null, + "execution_count": 4, "metadata": {}, "outputs": [], "source": [ - "ligand_molecule = Molecule(\"ligand.sdf\")\n", - "sage = ForceField(\"openff-2.0.0.offxml\")\n", + "# Read in the coordinates of the ligand from the PDB file\n", + "ligand_pdbfile = PDBFile(\"ligand.pdb\")\n", "\n", - "ligand = Interchange.from_smirnoff(\n", - " force_field=sage, topology=ligand_molecule.to_topology()\n", - ")\n", + "# Convert OpenMM System object containing ligand parameters into a ParmEd Structure.\n", + "import parmed\n", + "\n", + "ligand_structure = parmed.openmm.load_topology(\n", + " ligand_pdbfile.topology, ligand_system, xyz=ligand_pdbfile.positions\n", + ")" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "### Create a ParmEd `Structure` of an AMBER-parametrized receptor\n", "\n", - "receptor_molecule = Molecule.from_pdb(\"receptor.pdb\")\n", + "We have to create a ParmEd `Structure` containing positions and parameters of the receptor (BRD4) that we will then combine with the positions and parameters in the ligand `Structure` we created above. \n", "\n", - "ff14sb = ForceField(\"ff14sb_off_impropers_0.0.1.offxml\")\n", - "ff14sb.deregister_parameter_handler(\"ImproperTorsions\")\n", + "First, we use OpenMM to assign Amber14 parameters using OpenMM, but you can also create a `Structure` [from Amber `prmtop` and `inpcrd` files](http://parmed.github.io/ParmEd/html/structure.html).\n", "\n", - "receptor = Interchange.from_smirnoff(\n", - " force_field=ff14sb, topology=receptor_molecule.to_topology()\n", - ")\n", + "
\n", + " Note: If you already have AMBER (prmtop/inpcrd), GROMACS (top/gro), or any other file specifying the protein parameters supported by ParmEd, you can simply load the files directly into a Structure using ParmEd's functionalities. See https://parmed.github.io/ParmEd/html/readwrite.html .\n", + "
" + ] + }, + { + "cell_type": "code", + "execution_count": 5, + "metadata": {}, + "outputs": [], + "source": [ + "receptor_pdbfile = PDBFile(\"receptor.pdb\")\n", + "\n", + "# Load the AMBER protein force field through OpenMM.\n", + "omm_forcefield = app.ForceField(\"amber14-all.xml\")\n", "\n", - "complex_system = receptor + ligand\n", + "# Parameterize the protein.\n", + "receptor_system = omm_forcefield.createSystem(receptor_pdbfile.topology)\n", "\n", - "# Make sure ...\n", - "# complex.box = receptor.box\n", - "# complex.positions = np.vstack([receptor.positions, ligand.positions])" + "# Convert the protein System into a ParmEd Structure.\n", + "receptor_structure = parmed.openmm.load_topology(\n", + " receptor_pdbfile.topology, receptor_system, xyz=receptor_pdbfile.positions\n", + ")" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ - "### Export to OpenMM" + "### Combine receptor and ligand structures\n", + "\n", + "We can then merge the receptor and ligand `Structure` objects into a single `Structure` containing both positions and parameters." ] }, { "cell_type": "code", - "execution_count": null, + "execution_count": 6, "metadata": {}, "outputs": [], "source": [ - "complex_system.to_openmm()" + "complex_structure = receptor_structure + ligand_structure" ] }, { "cell_type": "markdown", - "metadata": { - "tags": [] - }, + "metadata": {}, + "source": [ + "### Export to OpenMM\n", + "\n", + "Once we have the `Structure` for the complex, containing both positions and parameters, we can chose to create an OpenMM `System` object that we can simulate using [OpenMM](http://openmm.org):" + ] + }, + { + "cell_type": "code", + "execution_count": 7, + "metadata": {}, + "outputs": [], + "source": [ + "# Convert the Structure to an OpenMM System in vacuum.\n", + "complex_system = complex_structure.createSystem(\n", + " nonbondedMethod=NoCutoff,\n", + " nonbondedCutoff=9.0 * unit.angstrom,\n", + " constraints=HBonds,\n", + " removeCMMotion=False,\n", + ")" + ] + }, + { + "cell_type": "code", + "execution_count": 8, + "metadata": {}, + "outputs": [], + "source": [ + "# Export the System to an OpenMM System XML and PDB file\n", + "complex_structure.save(\"complex.pdb\", overwrite=True)\n", + "\n", + "with open(\"complex.xml\", \"w\") as f:\n", + " f.write(XmlSerializer.serialize(complex_system))" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, "source": [ - "### Export to Amber" + "### Export to Amber\n", + "\n", + "We can also export the `Structure` to Amber `prmtop` and `inpcrd` files using [ParmEd](http://parmed.github.io/ParmEd/html/amber.html):" ] }, { "cell_type": "code", - "execution_count": null, + "execution_count": 9, "metadata": {}, "outputs": [], "source": [ - "complex_system.to_inpcrd(\"complex.inpcrd\")\n", - "complex_system.to_prmtop(\"complex.prmtop\")" + "# Export the Structure to AMBER files\n", + "complex_structure.save(\"complex.prmtop\", overwrite=True)\n", + "complex_structure.save(\"complex.inpcrd\", overwrite=True)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ - "### Export to GROMACS" + "### Export to gromacs\n", + "\n", + "We can export the `Structure` to gromacs `gro` and `top` files using [ParmEd](http://parmed.github.io/ParmEd/html/gromacs.html):" ] }, { "cell_type": "code", - "execution_count": null, + "execution_count": 10, "metadata": {}, "outputs": [], "source": [ - "complex_structure.to_gro(\"complex.gro\")\n", - "complex_structure.to_top(\"complex.top\")" + "# Export the Structure to Gromacs files\n", + "complex_structure.save(\"complex.gro\", overwrite=True)\n", + "complex_structure.save(\"complex.top\", overwrite=True)" ] } ], "metadata": { "kernelspec": { - "display_name": "Python 3 (ipykernel)", + "display_name": "Python 3", "language": "python", "name": "python3" }, @@ -144,7 +268,7 @@ "name": "python", "nbconvert_exporter": "python", "pygments_lexer": "ipython3", - "version": "3.9.7" + "version": "3.7.10" } }, "nbformat": 4, diff --git a/examples/using_smirnoff_with_amber_protein_forcefield/toluene_in_T4_lysozyme.ipynb b/examples/using_smirnoff_with_amber_protein_forcefield/toluene_in_T4_lysozyme.ipynb index 2881831dc..8427d10b4 100644 --- a/examples/using_smirnoff_with_amber_protein_forcefield/toluene_in_T4_lysozyme.ipynb +++ b/examples/using_smirnoff_with_amber_protein_forcefield/toluene_in_T4_lysozyme.ipynb @@ -4,7 +4,11 @@ "cell_type": "markdown", "metadata": {}, "source": [ - "## Create a system mixing SMIRNOFF and non-SMIRNOFF-formatted force fields" + "## Create a system mixing SMIRNOFF and non-SMIRNOFF-formatted force fields\n", + "\n", + "This example shows how to create a receptor-ligand OpenMM `System` where the ligand (toluene) is parametrized with a SMIRNOFF force field and the protein (T4 Lysozyme) solvated in water is assigned AMBER and TIP3P-FB parameters through the ParmEd library.\n", + "\n", + "We'll need two PDB files. One for the ligand in vacuum, and one for the solvated protein without ligand. The coordinates of the protein-ligand complex will be determined by these PDB files, so their positions needs to be consistent with the ligand being positioned in the binding pocket if this is the desired initial configuration of your simulation." ] }, { @@ -18,90 +22,186 @@ }, { "cell_type": "code", - "execution_count": null, + "execution_count": 1, "metadata": {}, - "outputs": [], + "outputs": [ + { + "name": "stderr", + "output_type": "stream", + "text": [ + "Warning: Unable to load toolkit 'OpenEye Toolkit'. The Open Force Field Toolkit does not require the OpenEye Toolkits, and can use RDKit/AmberTools instead. However, if you have a valid license for the OpenEye Toolkits, consider installing them for faster performance and additional file format support: https://docs.eyesopen.com/toolkits/python/quickstart-python/linuxosx.html OpenEye offers free Toolkit licenses for academics: https://www.eyesopen.com/academic-licensing\n" + ] + } + ], "source": [ "try:\n", - " from openmm import app\n", - " from openmm import unit as openmm_unit\n", + " from openmm import app, unit\n", + " from openmm.app import HBonds, NoCutoff, PDBFile\n", "except ImportError:\n", - " from simtk import unit as openmm_unit\n", + " from simtk import unit\n", " from simtk.openmm import app\n", + " from simtk.openmm.app import PDBFile, HBonds, NoCutoff\n", "\n", - "import numpy as np\n", - "from openff.interchange.components.interchange import Interchange\n", - "from openff.units import unit\n", - "\n", - "from openff.toolkit.topology import Molecule\n", + "from openff.toolkit.topology import Molecule, Topology\n", "from openff.toolkit.typing.engines.smirnoff import ForceField\n", - "from openff.toolkit.utils import get_data_file_path\n", + "from openff.toolkit.utils import get_data_file_path" + ] + }, + { + "cell_type": "code", + "execution_count": 2, + "metadata": {}, + "outputs": [], + "source": [ + "# Create an OpenFF Topology of toluene from a pdb file.\n", + "toluene_pdb_file_path = get_data_file_path(\"molecules/toluene.pdb\")\n", + "toluene_pdbfile = PDBFile(toluene_pdb_file_path)\n", + "toluene = Molecule.from_smiles(\"Cc1ccccc1\")\n", + "off_topology = Topology.from_openmm(\n", + " openmm_topology=toluene_pdbfile.topology, unique_molecules=[toluene]\n", + ")\n", "\n", + "# Load the Sage force field from disk.\n", + "force_field = ForceField(\"openff_unconstrained-2.0.0.offxml\")\n", "\n", - "def converter(x):\n", - " return x.value_in_unit(openmm_unit.nanometer)" + "# Parametrize the toluene molecule.\n", + "toluene_system = force_field.create_openmm_system(off_topology)" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "and we convert the OpenMM `System` to a ParmEd `Structure` that we'll be able to mix with the protein.\n", + "\n", + "
\n", + " Warning: ParmEd's Structure model is inspired by AMBER. Some information in an OpenMM System are not directly translatable into a Structure. In particular, long-range interaction treatment method (e.g., PME, CutoffPeriodic) and parameters (e.g., cutoff and cutoff switching distance, PME error tolerance) are known to be lost during the conversion.\n", + "
" ] }, { "cell_type": "code", - "execution_count": null, - "metadata": { - "tags": [] - }, + "execution_count": 3, + "metadata": {}, "outputs": [], "source": [ - "toluene_positions = app.PDBFile(get_data_file_path(\"molecules/toluene.pdb\")).positions\n", - "toluene_topology = Molecule.from_smiles(\"Cc1ccccc1\").to_topology()\n", + "import parmed\n", "\n", - "sage = ForceField(\"openff-2.0.0.offxml\")\n", + "# Convert OpenMM System into a ParmEd Structure.\n", + "toluene_structure = parmed.openmm.load_topology(\n", + " toluene_pdbfile.topology, toluene_system, xyz=toluene_pdbfile.positions\n", + ")" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "### Create a ParmEd `Structure` of an AMBER-parametrized receptor in TIP3P-FB water\n", + "\n", + "We have to create a ParmEd `Structure` of the receptor (T4 Lysozyme) to combine to the toluene `Structure`. Here we assign `amber99sbildn` to a T4 Lysozyme receptor solvated in TIP3P-FB water using OpenMM.\n", "\n", - "toluene = Interchange.from_smirnoff(force_field=sage, topology=toluene_topology)\n", - "toluene.positions = np.frompyfunc(converter, 1, 1)(toluene_positions).astype(float)" + "First, we load the parameters and the PDB file including positions for the protein, water, and ion atoms.\n", + "\n", + "
\n", + " Note: If you already have AMBER (prmtop/inpcrd), GROMACS (top/gro), or any other file supported by ParmEd specifying the parameters for the solvated protein, you can simply load the files directly into a Structure using ParmEd's functionalities. See https://parmed.github.io/ParmEd/html/readwrite.html .\n", + "
" ] }, { "cell_type": "code", - "execution_count": null, + "execution_count": 4, "metadata": {}, "outputs": [], "source": [ - "ff14sb = ForceField(\"ff14sb_off_impropers_0.0.1.offxml\")\n", - "ff14sb.deregister_parameter_handler(\"ImproperTorsions\")\n", + "# Load the AMBER protein force field parameters through OpenMM.\n", + "omm_forcefield = app.ForceField(\"amber99sbildn.xml\", \"tip3pfb.xml\")\n", "\n", - "# TODO: This PDB file also has water an ions, which are not yet intelligently parsed here.\n", - "# Unclear if we should have a Topology.from_pdb, just drop the waters, hack through OpenMM ...\n", - "t4_lysozyme_molecule = Molecule.from_pdb(\n", - " get_data_file_path(\"systems/test_systems/T4_lysozyme_water_ions.pdb\")\n", - ")\n", + "# Load the solvated receptor from a PDB file.\n", + "t4_pdb_file_path = get_data_file_path(\"systems/test_systems/T4_lysozyme_water_ions.pdb\")\n", + "t4_pdb_file = PDBFile(t4_pdb_file_path)\n", "\n", - "t4_lysozyme = Interchange.from_smirnoff(\n", - " force_field=ff14sb, topology=t4_lysozyme_molecule.to_topoplogy()\n", - ")\n", - "t4_lysozyme.positions = t4_lysozyme_molecule.conformers[0]" + "# Obtain the updated OpenMM Topology and positions.\n", + "omm_topology = t4_pdb_file.getTopology()\n", + "positions = t4_pdb_file.getPositions()" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "We then create a parametrized OpenMM `System` and convert it to a `Structure`.\n", + "\n", + "Note the `rigidWater=False` argument in `ForceField.createSystem()`. This is necessary to work around a problem araising with ParmEd in reading the parameters of constrained bonds from an OpenMM `System` (see https://github.com/openforcefield/openff-toolkit/issues/259 for more details). We'll re-add the hydrogen bond constraints when we'll create the `System` for the complex.\n", + "\n", + "
\n", + " Note: If you don't solvate the system or if you load it directly from AMBER, GROMACS, or other files directly to ParmEd, you won't need extra precautions.\n", + "
" ] }, { "cell_type": "code", - "execution_count": null, + "execution_count": 5, "metadata": {}, "outputs": [], "source": [ - "complex_system = t4_lysozyme + toluene" + "# Parameterize the protein.\n", + "t4_system = omm_forcefield.createSystem(omm_topology, rigidWater=False)\n", + "\n", + "# Convert the protein System into a ParmEd Structure.\n", + "t4_structure = parmed.openmm.load_topology(omm_topology, t4_system, xyz=positions)" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "### Combine receptor and ligand structures\n", + "\n", + "We can then merge the receptor and ligand `Structure` objects to form the complex. Note that the coordinates of protein and ligand in this example are determined by the PDB file, and they are already consistent with the ligand being positioned in the binding pocket.\n", + "\n", + "
\n", + " Note: If you want to include water molecules in the binding site, you will have to be careful to place them so that they won't create steric clashes once the ligand is inserted.\n", + "
" + ] + }, + { + "cell_type": "code", + "execution_count": 6, + "metadata": {}, + "outputs": [], + "source": [ + "complex_structure = toluene_structure + t4_structure" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "### Convert back the structure into an OpenMM System\n", + "\n", + "Once we have the `Structure` of the complex, we can chose to create a `System` object that we can simulate with OpenMM." ] }, { "cell_type": "code", - "execution_count": null, + "execution_count": 7, "metadata": {}, "outputs": [], "source": [ - "complex_system.to_openmm()" + "# Convert the Structure to an OpenMM System in vacuum.\n", + "complex_system = complex_structure.createSystem(\n", + " nonbondedMethod=NoCutoff,\n", + " nonbondedCutoff=9.0 * unit.angstrom,\n", + " constraints=HBonds,\n", + " removeCMMotion=False,\n", + ")" ] } ], "metadata": { "kernelspec": { - "display_name": "Python 3 (ipykernel)", + "display_name": "Python 3", "language": "python", "name": "python3" }, @@ -115,7 +215,7 @@ "name": "python", "nbconvert_exporter": "python", "pygments_lexer": "ipython3", - "version": "3.9.7" + "version": "3.7.10" } }, "nbformat": 4, From fdf19b9d36adca612f94a271c117fb4f75faf684 Mon Sep 17 00:00:00 2001 From: "Matthew W. Thompson" Date: Thu, 28 Oct 2021 10:35:58 -0500 Subject: [PATCH 07/24] Rename notebooks --- .../{convert_to_amber_gromacs.ipynb => convert_with_parmed.ipynb} | 0 ...smirnoff_amber_gromacs.ipynb => export_with_interchange.ipynb} | 0 2 files changed, 0 insertions(+), 0 deletions(-) rename examples/using_smirnoff_in_amber_or_gromacs/{convert_to_amber_gromacs.ipynb => convert_with_parmed.ipynb} (100%) rename examples/using_smirnoff_in_amber_or_gromacs/{smirnoff_amber_gromacs.ipynb => export_with_interchange.ipynb} (100%) diff --git a/examples/using_smirnoff_in_amber_or_gromacs/convert_to_amber_gromacs.ipynb b/examples/using_smirnoff_in_amber_or_gromacs/convert_with_parmed.ipynb similarity index 100% rename from examples/using_smirnoff_in_amber_or_gromacs/convert_to_amber_gromacs.ipynb rename to examples/using_smirnoff_in_amber_or_gromacs/convert_with_parmed.ipynb diff --git a/examples/using_smirnoff_in_amber_or_gromacs/smirnoff_amber_gromacs.ipynb b/examples/using_smirnoff_in_amber_or_gromacs/export_with_interchange.ipynb similarity index 100% rename from examples/using_smirnoff_in_amber_or_gromacs/smirnoff_amber_gromacs.ipynb rename to examples/using_smirnoff_in_amber_or_gromacs/export_with_interchange.ipynb From 7b353f92de8248afcd97b4160fe7f4000245e00c Mon Sep 17 00:00:00 2001 From: "Matthew W. Thompson" Date: Thu, 28 Oct 2021 10:47:22 -0500 Subject: [PATCH 08/24] Add header notes to each notebook --- .../convert_with_parmed.ipynb | 20 ++++++++++++++++--- .../export_with_interchange.ipynb | 9 +++++++++ 2 files changed, 26 insertions(+), 3 deletions(-) diff --git a/examples/using_smirnoff_in_amber_or_gromacs/convert_with_parmed.ipynb b/examples/using_smirnoff_in_amber_or_gromacs/convert_with_parmed.ipynb index 1bff64170..6e1e19167 100644 --- a/examples/using_smirnoff_in_amber_or_gromacs/convert_with_parmed.ipynb +++ b/examples/using_smirnoff_in_amber_or_gromacs/convert_with_parmed.ipynb @@ -3,6 +3,17 @@ { "cell_type": "markdown", "metadata": {}, + "source": [ + "
\n", + " Note: \n", + " Please take a look at the other notebook in the folder (export_with_interchange.ipynb), showcasing how to use OpenFF Interchange. Interchange is well-validated for small molecule interoperability with several formats and is in the process of replacing ParmEd in OpenFF workflows.
" + ] + }, + { + "cell_type": "markdown", + "metadata": { + "tags": [] + }, "source": [ "## Convert Open Forcefield System to AMBER and GROMACS input files\n", "\n", @@ -82,7 +93,10 @@ }, { "cell_type": "markdown", - "metadata": {}, + "metadata": { + "jp-MarkdownHeadingCollapsed": true, + "tags": [] + }, "source": [ "### Convert OpenMM System to AMBER and GROMACS files\n", "\n", @@ -263,7 +277,7 @@ ], "metadata": { "kernelspec": { - "display_name": "Python 3", + "display_name": "Python 3 (ipykernel)", "language": "python", "name": "python3" }, @@ -277,7 +291,7 @@ "name": "python", "nbconvert_exporter": "python", "pygments_lexer": "ipython3", - "version": "3.8.6" + "version": "3.8.12" } }, "nbformat": 4, diff --git a/examples/using_smirnoff_in_amber_or_gromacs/export_with_interchange.ipynb b/examples/using_smirnoff_in_amber_or_gromacs/export_with_interchange.ipynb index 4bac9f501..fa3877bfc 100644 --- a/examples/using_smirnoff_in_amber_or_gromacs/export_with_interchange.ipynb +++ b/examples/using_smirnoff_in_amber_or_gromacs/export_with_interchange.ipynb @@ -1,5 +1,14 @@ { "cells": [ + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "
\n", + " Note: \n", + " Interchange is in the process of replacing ParmEd in many workflows, but it still in an alpha testing phase. Our internal tests indicate it is reliable for many small-molecule systems, but it is not yet reliable for complex, multi-component systems and there are likely still rough edges throughout. Feedback is welcome on the Interchange issue tracker.
" + ] + }, { "cell_type": "markdown", "metadata": {}, From 6df34eb61646689506c4c104348ae1afc512daad Mon Sep 17 00:00:00 2001 From: Jeff Wagner Date: Mon, 1 Nov 2021 14:27:57 -0700 Subject: [PATCH 09/24] add new example deps to example yaml --- examples/environment.yaml | 3 +++ 1 file changed, 3 insertions(+) diff --git a/examples/environment.yaml b/examples/environment.yaml index f59f25bbd..12b2decae 100644 --- a/examples/environment.yaml +++ b/examples/environment.yaml @@ -17,3 +17,6 @@ dependencies: - numpy - openmm - openff-forcefields >=2.0.0 + - openff-interchange >=0.1.2 + - gromacs + - lammps From 70cd1424f423ba17b9cd545802b65d9ee95cfc67 Mon Sep 17 00:00:00 2001 From: "Matthew W. Thompson" Date: Tue, 2 Nov 2021 12:03:43 -0500 Subject: [PATCH 10/24] Skip Interchange example when AmberTools is not installed --- .github/workflows/examples.yml | 1 + 1 file changed, 1 insertion(+) diff --git a/.github/workflows/examples.yml b/.github/workflows/examples.yml index 25242579f..07bae52da 100644 --- a/.github/workflows/examples.yml +++ b/.github/workflows/examples.yml @@ -117,6 +117,7 @@ jobs: if [[ "$RDKIT" == false ]]; then PYTEST_ARGS+=" --ignore=examples/check_dataset_parameter_coverage" PYTEST_ARGS+=" --ignore=examples/QCArchive_interface" + PYTEST_ARGS+=" --ignore=examples/using_smirnoff_in_amber_or_gromacs" fi pytest $PYTEST_ARGS openff/toolkit/tests/test_examples.py From cc4b1e04ad5d1c3f23a2f99ef7e6f676c5cad56f Mon Sep 17 00:00:00 2001 From: "Matthew W. Thompson" Date: Tue, 2 Nov 2021 15:00:19 -0500 Subject: [PATCH 11/24] Do not assume LAMMPS and GROMACS are installed --- .github/workflows/examples.yml | 2 + .../export_with_interchange.ipynb | 279 +++++------------- 2 files changed, 75 insertions(+), 206 deletions(-) diff --git a/.github/workflows/examples.yml b/.github/workflows/examples.yml index 07bae52da..448dea493 100644 --- a/.github/workflows/examples.yml +++ b/.github/workflows/examples.yml @@ -101,6 +101,8 @@ jobs: - name: Install package shell: bash -l {0} run: | + pip install git+https://github.com/openforcefield/openff-interchange.git@skip-drivers + # Remove the packaged openff-toolkit, installed as a dependency of openmmforcefields conda remove --force openff-toolkit-base python setup.py develop --no-deps diff --git a/examples/using_smirnoff_in_amber_or_gromacs/export_with_interchange.ipynb b/examples/using_smirnoff_in_amber_or_gromacs/export_with_interchange.ipynb index fa3877bfc..a40520b95 100644 --- a/examples/using_smirnoff_in_amber_or_gromacs/export_with_interchange.ipynb +++ b/examples/using_smirnoff_in_amber_or_gromacs/export_with_interchange.ipynb @@ -31,7 +31,7 @@ }, { "cell_type": "code", - "execution_count": 1, + "execution_count": null, "metadata": {}, "outputs": [], "source": [ @@ -46,20 +46,9 @@ }, { "cell_type": "code", - "execution_count": 2, + "execution_count": null, "metadata": {}, - "outputs": [ - { - "data": { - "text/plain": [ - "" - ] - }, - "execution_count": 2, - "metadata": {}, - "output_type": "execute_result" - } - ], + "outputs": [], "source": [ "ethanol = Molecule.from_smiles(\"CCO\")\n", "cyclohexane = Molecule.from_smiles(\"C1CCCCC1\")\n", @@ -89,20 +78,9 @@ }, { "cell_type": "code", - "execution_count": 3, + "execution_count": null, "metadata": {}, - "outputs": [ - { - "data": { - "text/plain": [ - "" - ] - }, - "execution_count": 3, - "metadata": {}, - "output_type": "execute_result" - } - ], + "outputs": [], "source": [ "forcefield = ForceField(\"openff-2.0.0.offxml\")\n", "forcefield" @@ -119,20 +97,9 @@ }, { "cell_type": "code", - "execution_count": 4, + "execution_count": null, "metadata": {}, - "outputs": [ - { - "data": { - "text/plain": [ - " >" - ] - }, - "execution_count": 4, - "metadata": {}, - "output_type": "execute_result" - } - ], + "outputs": [], "source": [ "omm_system = forcefield.create_openmm_system(off_topology)\n", "omm_system" @@ -149,27 +116,9 @@ }, { "cell_type": "code", - "execution_count": 5, + "execution_count": null, "metadata": {}, - "outputs": [ - { - "name": "stderr", - "output_type": "stream", - "text": [ - "Warning: importing 'simtk.openmm' is deprecated. Import 'openmm' instead.\n" - ] - }, - { - "data": { - "text/plain": [ - "Interchange with 27 atoms, periodic topology" - ] - }, - "execution_count": 5, - "metadata": {}, - "output_type": "execute_result" - } - ], + "outputs": [], "source": [ "from openff.interchange.components.interchange import Interchange\n", "\n", @@ -192,7 +141,7 @@ }, { "cell_type": "code", - "execution_count": 6, + "execution_count": null, "metadata": {}, "outputs": [], "source": [ @@ -209,184 +158,102 @@ "cell_type": "markdown", "metadata": {}, "source": [ - "### Validate the conversion\n", + "### Validating the conversion to Amber files\n", "\n", "The Interchange project includes functions that take in an `Interchange` object and call out to simulation engines to run single-point energy calculations (with no minimization or dynamics) for the purpose of validating the export layer with each engine. Under the hood, each of these functions calls API points like those used above while converting to files understood by each engine. These rely on having each engine installed and accessible in the current `$PATH`." ] }, { "cell_type": "code", - "execution_count": 7, + "execution_count": null, "metadata": {}, "outputs": [], "source": [ - "from openff.interchange.drivers import (\n", - " get_all_energies,\n", - " get_amber_energies,\n", - " get_gromacs_energies,\n", - " get_openmm_energies,\n", - ")\n", - "from openff.interchange.drivers.all import get_summary_data" + "from openff.interchange.drivers import get_amber_energies, get_openmm_energies" ] }, { "cell_type": "code", - "execution_count": 8, - "metadata": {}, - "outputs": [ - { - "data": { - "text/plain": [ - "{'Bond': 0.29528921842575073 ,\n", - " 'Angle': 155.7017059326172 ,\n", - " 'Torsion': 19.56324005126953 ,\n", - " 'vdW': 9.66697609052062 ,\n", - " 'Electrostatics': -15.235476016998291 }" - ] - }, - "execution_count": 8, - "metadata": {}, - "output_type": "execute_result" - } - ], + "execution_count": null, + "metadata": { + "scrolled": true + }, + "outputs": [], "source": [ - "get_gromacs_energies(interchange).energies" + "openmm_energies = get_openmm_energies(interchange)\n", + "openmm_energies.energies" ] }, { "cell_type": "code", - "execution_count": 9, + "execution_count": null, "metadata": { "scrolled": true }, - "outputs": [ - { - "data": { - "text/plain": [ - "{'Bond': 0.2952878475189209 ,\n", - " 'Angle': 155.7012939453125 ,\n", - " 'Torsion': 19.563228607177734 ,\n", - " 'vdW': 9.6586229795419 ,\n", - " 'Electrostatics': -15.215721021619657 }" - ] - }, - "execution_count": 9, - "metadata": {}, - "output_type": "execute_result" - } - ], + "outputs": [], "source": [ - "get_openmm_energies(interchange).energies" + "amber_energies = get_amber_energies(interchange)\n", + "amber_energies.energies" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "### Appendix: Validating the conversion to GROMACS and LAMMPS files\n", + "\n", + "If GROMACS and/or LAMMPS are installed on your machine, the same comparisons can also take place with those engines. They are available via `conda` by running a command like:\n", + "\n", + "```conda install gromacs lammps -c conda-forge -c bioconda```" ] }, { "cell_type": "code", - "execution_count": 10, + "execution_count": null, "metadata": {}, - "outputs": [ - { - "data": { - "text/plain": [ - "{'Bond': 0.0741 ,\n", - " 'Angle': 37.2023 ,\n", - " 'Torsion': 4.6778 ,\n", - " 'vdW': 9.7231976 ,\n", - " 'Electrostatics': -15.212187200000006 }" - ] - }, - "execution_count": 10, - "metadata": {}, - "output_type": "execute_result" - } - ], + "outputs": [], "source": [ - "get_amber_energies(interchange).energies" + "from distutils.spawn import find_executable\n", + "from pprint import pprint\n", + "\n", + "from openff.interchange.drivers import get_gromacs_energies, get_lammps_energies" ] }, { "cell_type": "code", - "execution_count": 11, + "execution_count": null, "metadata": {}, - "outputs": [ - { - "data": { - "text/html": [ - "
\n", - "\n", - "\n", - " \n", - " \n", - " \n", - " \n", - " \n", - " \n", - " \n", - " \n", - " \n", - " \n", - " \n", - " \n", - " \n", - " \n", - " \n", - " \n", - " \n", - " \n", - " \n", - " \n", - " \n", - " \n", - " \n", - " \n", - " \n", - " \n", - " \n", - " \n", - " \n", - " \n", - " \n", - " \n", - " \n", - " \n", - " \n", - " \n", - " \n", - " \n", - " \n", - " \n", - " \n", - " \n", - " \n", - " \n", - "
BondAngleTorsionvdWElectrostatics
OpenMM0.295288155.70129419.5632299.658623-15.215721
GROMACS0.295289155.70170619.5632409.666976-15.235476
LAMMPS0.440418155.70134319.5632339.630999-15.219631
Amber0.310034155.65442319.5719159.723198-15.212187
\n", - "
" - ], - "text/plain": [ - " Bond Angle Torsion vdW Electrostatics\n", - "OpenMM 0.295288 155.701294 19.563229 9.658623 -15.215721\n", - "GROMACS 0.295289 155.701706 19.563240 9.666976 -15.235476\n", - "LAMMPS 0.440418 155.701343 19.563233 9.630999 -15.219631\n", - "Amber 0.310034 155.654423 19.571915 9.723198 -15.212187" - ] - }, - "execution_count": 11, - "metadata": {}, - "output_type": "execute_result" - } - ], + "outputs": [], "source": [ + "if find_executable(\"lmp_serial\"):\n", + " pprint(get_lammps_energies(interchange).energies)" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "if find_executable(\"gmx\"):\n", + " pprint(get_gromacs_energies(interchange).energies)" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "Finally, there is a helper function `get_summary_data` that will attempt to run drivers of each engine. A summary reported is prepared as a Pandas `DataFrame`." + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "from openff.interchange.drivers.all import get_summary_data\n", + "\n", "get_summary_data(interchange)" ] } @@ -407,7 +274,7 @@ "name": "python", "nbconvert_exporter": "python", "pygments_lexer": "ipython3", - "version": "3.8.12" + "version": "3.9.7" } }, "nbformat": 4, From 9b9d823ed6f2ee1971aafb724975fee5b1429549 Mon Sep 17 00:00:00 2001 From: "Matthew W. Thompson" Date: Mon, 8 Nov 2021 14:22:34 -0600 Subject: [PATCH 12/24] Fix how GROMACS/Amber notebook is detected --- .github/workflows/examples.yml | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/.github/workflows/examples.yml b/.github/workflows/examples.yml index 448dea493..da3314475 100644 --- a/.github/workflows/examples.yml +++ b/.github/workflows/examples.yml @@ -119,7 +119,6 @@ jobs: if [[ "$RDKIT" == false ]]; then PYTEST_ARGS+=" --ignore=examples/check_dataset_parameter_coverage" PYTEST_ARGS+=" --ignore=examples/QCArchive_interface" - PYTEST_ARGS+=" --ignore=examples/using_smirnoff_in_amber_or_gromacs" fi pytest $PYTEST_ARGS openff/toolkit/tests/test_examples.py @@ -130,6 +129,7 @@ jobs: NB_ARGS+=" --ignore=examples/QCArchive_interface" NB_ARGS+=" --ignore=examples/check_dataset_parameter_coverage" NB_ARGS+=" --ignore=examples/conformer_energies" + NB_ARGS+=" --ignore=examples/using_smirnoff_in_amber_or_gromacs" fi pytest $NB_ARGS examples From 23d674f9ee49b348c348440b2e91ec0e00c6fb46 Mon Sep 17 00:00:00 2001 From: "Matthew W. Thompson" Date: Tue, 9 Nov 2021 10:00:29 -0600 Subject: [PATCH 13/24] Debug files --- .../export_with_interchange.ipynb | 36 +++++++++++++++++++ 1 file changed, 36 insertions(+) diff --git a/examples/using_smirnoff_in_amber_or_gromacs/export_with_interchange.ipynb b/examples/using_smirnoff_in_amber_or_gromacs/export_with_interchange.ipynb index a40520b95..ca0cf3d71 100644 --- a/examples/using_smirnoff_in_amber_or_gromacs/export_with_interchange.ipynb +++ b/examples/using_smirnoff_in_amber_or_gromacs/export_with_interchange.ipynb @@ -163,6 +163,42 @@ "The Interchange project includes functions that take in an `Interchange` object and call out to simulation engines to run single-point energy calculations (with no minimization or dynamics) for the purpose of validating the export layer with each engine. Under the hood, each of these functions calls API points like those used above while converting to files understood by each engine. These rely on having each engine installed and accessible in the current `$PATH`." ] }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "!cat system.prmtop" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "!cat system.inpcrd" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "!cat system.top" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "!cat system.gro" + ] + }, { "cell_type": "code", "execution_count": null, From 5b26465e326d45cb8f93f0adf558049382edfe5a Mon Sep 17 00:00:00 2001 From: "Matthew W. Thompson" Date: Tue, 9 Nov 2021 10:02:59 -0600 Subject: [PATCH 14/24] Move debug so it might get printed --- .../export_with_interchange.ipynb | 40 ++----------------- 1 file changed, 4 insertions(+), 36 deletions(-) diff --git a/examples/using_smirnoff_in_amber_or_gromacs/export_with_interchange.ipynb b/examples/using_smirnoff_in_amber_or_gromacs/export_with_interchange.ipynb index ca0cf3d71..7136f877e 100644 --- a/examples/using_smirnoff_in_amber_or_gromacs/export_with_interchange.ipynb +++ b/examples/using_smirnoff_in_amber_or_gromacs/export_with_interchange.ipynb @@ -163,42 +163,6 @@ "The Interchange project includes functions that take in an `Interchange` object and call out to simulation engines to run single-point energy calculations (with no minimization or dynamics) for the purpose of validating the export layer with each engine. Under the hood, each of these functions calls API points like those used above while converting to files understood by each engine. These rely on having each engine installed and accessible in the current `$PATH`." ] }, - { - "cell_type": "code", - "execution_count": null, - "metadata": {}, - "outputs": [], - "source": [ - "!cat system.prmtop" - ] - }, - { - "cell_type": "code", - "execution_count": null, - "metadata": {}, - "outputs": [], - "source": [ - "!cat system.inpcrd" - ] - }, - { - "cell_type": "code", - "execution_count": null, - "metadata": {}, - "outputs": [], - "source": [ - "!cat system.top" - ] - }, - { - "cell_type": "code", - "execution_count": null, - "metadata": {}, - "outputs": [], - "source": [ - "!cat system.gro" - ] - }, { "cell_type": "code", "execution_count": null, @@ -228,6 +192,10 @@ }, "outputs": [], "source": [ + "!cat system.inpcrd\n", + "!cat system.prmtop\n", + "!cat system.top\n", + "!cat system.gro\n", "amber_energies = get_amber_energies(interchange)\n", "amber_energies.energies" ] From eca429d2321fd0cef56269f5b051a98728686098 Mon Sep 17 00:00:00 2001 From: "Matthew W. Thompson" Date: Tue, 9 Nov 2021 13:20:59 -0600 Subject: [PATCH 15/24] Split out Amber/GROMACS example into earlier step --- .github/workflows/examples.yml | 7 +++++++ 1 file changed, 7 insertions(+) diff --git a/.github/workflows/examples.yml b/.github/workflows/examples.yml index b8962db69..4f81ccff7 100644 --- a/.github/workflows/examples.yml +++ b/.github/workflows/examples.yml @@ -122,6 +122,13 @@ jobs: fi pytest $PYTEST_ARGS openff/toolkit/tests/test_examples.py + - name: Run only the Amber/GROMACS notebooks + shell: bash -l {0} + run: | + if [[ "$RDKIT" == true ]]; then + pytest $NB_ARGS examples/using_smirnoff_in_amber_or_gromacs/ + fi + - name: Run example notebooks shell: bash -l {0} run: | From 511a445dc1a09f299f319412e06ff279a286b2ec Mon Sep 17 00:00:00 2001 From: "Matthew W. Thompson" Date: Tue, 9 Nov 2021 14:27:22 -0600 Subject: [PATCH 16/24] Run examples with Python 3.9 --- .github/workflows/examples.yml | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/.github/workflows/examples.yml b/.github/workflows/examples.yml index 4f81ccff7..79e281949 100644 --- a/.github/workflows/examples.yml +++ b/.github/workflows/examples.yml @@ -27,11 +27,11 @@ jobs: - ubuntu-latest - macos-latest cfg: - - python-version: 3.7 + - python-version: 3.9 rdkit: "true" openeye: "false" - - python-version: 3.7 + - python-version: 3.9 rdkit: "false" openeye: "true" From 618bc7b24c3ebde2d875938a0da5e65c72f440ef Mon Sep 17 00:00:00 2001 From: Jeff Wagner Date: Tue, 9 Nov 2021 17:55:25 -0800 Subject: [PATCH 17/24] try not using mamba in case its causing trouble with ambertools --- .github/workflows/examples.yml | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/.github/workflows/examples.yml b/.github/workflows/examples.yml index 79e281949..7ac4890b8 100644 --- a/.github/workflows/examples.yml +++ b/.github/workflows/examples.yml @@ -61,7 +61,7 @@ jobs: auto-activate-base: false miniforge-version: latest miniforge-variant: Mambaforge - use-mamba: true + use-mamba: false - uses: conda-incubator/setup-miniconda@v2.1.1 name: Install with OpenEye toolkits if: ${{ matrix.cfg.rdkit == 'FALSE' && matrix.cfg.openeye == 'TRUE' }} From db9f1869dcc93532f0f9bd73e22cbc548f97c76a Mon Sep 17 00:00:00 2001 From: Jeff Wagner Date: Wed, 10 Nov 2021 03:02:47 +0000 Subject: [PATCH 18/24] re enable mamba, force ambertools downgrade --- .github/workflows/examples.yml | 2 +- examples/environment.yaml | 1 + 2 files changed, 2 insertions(+), 1 deletion(-) diff --git a/.github/workflows/examples.yml b/.github/workflows/examples.yml index 7ac4890b8..79e281949 100644 --- a/.github/workflows/examples.yml +++ b/.github/workflows/examples.yml @@ -61,7 +61,7 @@ jobs: auto-activate-base: false miniforge-version: latest miniforge-variant: Mambaforge - use-mamba: false + use-mamba: true - uses: conda-incubator/setup-miniconda@v2.1.1 name: Install with OpenEye toolkits if: ${{ matrix.cfg.rdkit == 'FALSE' && matrix.cfg.openeye == 'TRUE' }} diff --git a/examples/environment.yaml b/examples/environment.yaml index 12b2decae..c7c1dfaa4 100644 --- a/examples/environment.yaml +++ b/examples/environment.yaml @@ -20,3 +20,4 @@ dependencies: - openff-interchange >=0.1.2 - gromacs - lammps + - ambertools =21.3 From 637eae99dd3b8f66ee76d30da0911d5c42857a4b Mon Sep 17 00:00:00 2001 From: "Matthew W. Thompson" Date: Wed, 10 Nov 2021 10:10:19 -0600 Subject: [PATCH 19/24] Try running it as a test script --- .../convert_with_parmed.ipynb | 7 +- .../export_with_interchange.py | 172 ++++++++++++++++++ 2 files changed, 177 insertions(+), 2 deletions(-) create mode 100644 examples/using_smirnoff_in_amber_or_gromacs/export_with_interchange.py diff --git a/examples/using_smirnoff_in_amber_or_gromacs/convert_with_parmed.ipynb b/examples/using_smirnoff_in_amber_or_gromacs/convert_with_parmed.ipynb index 6e1e19167..b30c6377e 100644 --- a/examples/using_smirnoff_in_amber_or_gromacs/convert_with_parmed.ipynb +++ b/examples/using_smirnoff_in_amber_or_gromacs/convert_with_parmed.ipynb @@ -51,7 +51,8 @@ " from simtk.openmm.app import PDBFile\n", "\n", "from openff.toolkit.topology import Molecule, Topology\n", - "from openff.toolkit.typing.engines.smirnoff import ForceField" + "from openff.toolkit.typing.engines.smirnoff import ForceField\n", + "from openff.toolkit.utils import get_data_file_path" ] }, { @@ -64,7 +65,9 @@ "cyclohexane = Molecule.from_smiles(\"C1CCCCC1\")\n", "\n", "# Obtain the OpenMM Topology object from the PDB file.\n", - "pdbfile = PDBFile(\"1_cyclohexane_1_ethanol.pdb\")\n", + "pdbfile = PDBFile(\n", + " get_data_file_path(\"systems/test_systems/1_cyclohexane_1_ethanol.pdb\")\n", + ")\n", "omm_topology = pdbfile.topology\n", "\n", "# Create the Open Forcefield Topology.\n", diff --git a/examples/using_smirnoff_in_amber_or_gromacs/export_with_interchange.py b/examples/using_smirnoff_in_amber_or_gromacs/export_with_interchange.py new file mode 100644 index 000000000..e1ea62259 --- /dev/null +++ b/examples/using_smirnoff_in_amber_or_gromacs/export_with_interchange.py @@ -0,0 +1,172 @@ +#!/usr/bin/env python +# coding: utf-8 + +#
+# Note: +# Interchange is in the process of replacing ParmEd in many workflows, but it still in an alpha testing phase. Our internal tests indicate it is reliable for many small-molecule systems, but it is not yet reliable for complex, multi-component systems and there are likely still rough edges throughout. Feedback is welcome on the Interchange issue tracker.
+ +# ## Using OpenFF force fields in Amber and GROMACS +# +# The Open Forcefield Toolkit can create parametrized `openmm.System` objects that can be natively simulated with OpenMM. This example shows the Interchange project can enable parallel workflows using Amber and GROMACS. + +# ### Preparing an OpenFF Topology +# +# We start by loading a PDB file containing one copy of ethanol and cyclohexane. Our goal is to create an OpenFF `Topology` object describing this system that we can parametrize with the SMIRNOFF-format "Sage" force field. +# +# The two `Molecule` objects created from the SMILES strings can contain information such as partial charges and stereochemistry that is not included in an OpenMM topology. In this example, partial charges are not explicitly given, and `ForceField` will assign AM1/BCC charges as specified by the "Sage" force field. Note that the OpenFF Toolkit produces partial charges that do not depend on the input conformation of parameterized molecules. See the [FAQ](https://open-forcefield-toolkit.readthedocs.io/en/latest/faq.html#the-partial-charges-generated-by-the-toolkit-don-t-seem-to-depend-on-the-molecule-s-conformation-is-this-a-bug) for more information. + +# In[ ]: + + +try: + from openmm import app +except ImportError: + from simtk.openmm import app + +from openff.toolkit.utils import get_data_file_path +from openff.toolkit.topology import Molecule, Topology +from openff.toolkit.typing.engines.smirnoff import ForceField + + +# In[ ]: + + +ethanol = Molecule.from_smiles("CCO") +cyclohexane = Molecule.from_smiles("C1CCCCC1") + +# Load the PDB file using OpenMM and save the OpenMM Topology +pdbfile = app.PDBFile( + get_data_file_path("systems/test_systems/1_cyclohexane_1_ethanol.pdb") +) +omm_topology = pdbfile.topology + +# Create the OpenFF Topology. +off_topology = Topology.from_openmm( + omm_topology, unique_molecules=[ethanol, cyclohexane] +) +off_topology + + +# ### Preparing an OpenFF ForceField +# +# Once the `ForceField` class is imported, the only decision to make is which force field to use. An exhaustive list of force fields released by the Open Force Field Initiative can be found [here](from openff.toolkit.typing.engines.smirnoff import ForceField +# ). +# +# Here we will use force field from the "Sage" line. + +# In[ ]: + + +forcefield = ForceField("openff-2.0.0.offxml") +forcefield + + +# ### Preparing an OpenMM System +# +# Once a force field and topology have been loaded, an `openmm.System` can be generated natively with the OpenFF Toolkit. + +# In[ ]: + + +omm_system = forcefield.create_openmm_system(off_topology) +omm_system + + +# ### Preparing an Interchange object +# +# To exports to engines other than OpenMM, we will make use of the [Interchange](https://openff-interchange.readthedocs.io/) project. There is a high-level `Interchange.from_smirnoff` function that consumes OpenFF Toolkit and ForceField objects and produces an `Interchange` object which can then be exported to formats understood by other molecular simulation engines. This extra step is needed to provide a clean interface between _applied_ parameters and engines. Note also that this step does not require an OpenMM System to be generated; `ForceField.create_openmm_system` does not need to be called to use Amber and GROMACS. + +# In[ ]: + + +from openff.interchange.components.interchange import Interchange + +interchange = Interchange.from_smirnoff( + force_field=forcefield, + topology=off_topology, +) +interchange.positions = pdbfile.positions +interchange + + +# ### Exporting to Amber and GROMACS files +# +# Once an `Interchange` object has been constructed, its API can be used to export to files understood by GROMACS, Amber, and more. + +# In[ ]: + + +# Export AMBER files. +interchange.to_prmtop("system.prmtop") +interchange.to_inpcrd("system.inpcrd") + +# Export GROMACS files. +interchange.to_top("system.top") +interchange.to_gro("system.gro") + + +# ### Validating the conversion to Amber files +# +# The Interchange project includes functions that take in an `Interchange` object and call out to simulation engines to run single-point energy calculations (with no minimization or dynamics) for the purpose of validating the export layer with each engine. Under the hood, each of these functions calls API points like those used above while converting to files understood by each engine. These rely on having each engine installed and accessible in the current `$PATH`. + +# In[ ]: + + +from openff.interchange.drivers import get_amber_energies, get_openmm_energies + + +# In[ ]: + + +openmm_energies = get_openmm_energies(interchange) +openmm_energies.energies + + +# In[ ]: + + +# get_ipython().system('cat system.inpcrd') +# get_ipython().system('cat system.prmtop') +# get_ipython().system('cat system.top') +# get_ipython().system('cat system.gro') +amber_energies = get_amber_energies(interchange) +amber_energies.energies + + +# ### Appendix: Validating the conversion to GROMACS and LAMMPS files +# +# If GROMACS and/or LAMMPS are installed on your machine, the same comparisons can also take place with those engines. They are available via `conda` by running a command like: +# +# ```conda install gromacs lammps -c conda-forge -c bioconda``` + +# In[ ]: + + +from distutils.spawn import find_executable +from pprint import pprint + +from openff.interchange.drivers import get_gromacs_energies, get_lammps_energies + + +# In[ ]: + + +if find_executable("lmp_serial"): + pprint(get_lammps_energies(interchange).energies) + + +# In[ ]: + + +if find_executable("gmx"): + pprint(get_gromacs_energies(interchange).energies) + + +# Finally, there is a helper function `get_summary_data` that will attempt to run drivers of each engine. A summary reported is prepared as a Pandas `DataFrame`. + +# In[ ]: + + +from openff.interchange.drivers.all import get_summary_data + +get_summary_data(interchange) From f7e2d3b4a32c48fc9678796906dd5cea59f2bdc5 Mon Sep 17 00:00:00 2001 From: "Matthew W. Thompson" Date: Wed, 10 Nov 2021 10:15:32 -0600 Subject: [PATCH 20/24] Unpin AmberTools --- examples/environment.yaml | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/examples/environment.yaml b/examples/environment.yaml index c7c1dfaa4..923bc43de 100644 --- a/examples/environment.yaml +++ b/examples/environment.yaml @@ -20,4 +20,4 @@ dependencies: - openff-interchange >=0.1.2 - gromacs - lammps - - ambertools =21.3 + - ambertools From 64c4af22c80d67d721a99127416d2e3ef3db1c12 Mon Sep 17 00:00:00 2001 From: "Matthew W. Thompson" Date: Wed, 10 Nov 2021 15:05:31 -0600 Subject: [PATCH 21/24] Revert some debugging in workflow --- .github/workflows/examples.yml | 11 ++--------- 1 file changed, 2 insertions(+), 9 deletions(-) diff --git a/.github/workflows/examples.yml b/.github/workflows/examples.yml index 79e281949..b8962db69 100644 --- a/.github/workflows/examples.yml +++ b/.github/workflows/examples.yml @@ -27,11 +27,11 @@ jobs: - ubuntu-latest - macos-latest cfg: - - python-version: 3.9 + - python-version: 3.7 rdkit: "true" openeye: "false" - - python-version: 3.9 + - python-version: 3.7 rdkit: "false" openeye: "true" @@ -122,13 +122,6 @@ jobs: fi pytest $PYTEST_ARGS openff/toolkit/tests/test_examples.py - - name: Run only the Amber/GROMACS notebooks - shell: bash -l {0} - run: | - if [[ "$RDKIT" == true ]]; then - pytest $NB_ARGS examples/using_smirnoff_in_amber_or_gromacs/ - fi - - name: Run example notebooks shell: bash -l {0} run: | From 7219b0b1f14a4cc843256aa4571f875402f8250d Mon Sep 17 00:00:00 2001 From: "Matthew W. Thompson" Date: Wed, 10 Nov 2021 15:14:34 -0600 Subject: [PATCH 22/24] Remove script used for debugging --- .../export_with_interchange.py | 172 ------------------ 1 file changed, 172 deletions(-) delete mode 100644 examples/using_smirnoff_in_amber_or_gromacs/export_with_interchange.py diff --git a/examples/using_smirnoff_in_amber_or_gromacs/export_with_interchange.py b/examples/using_smirnoff_in_amber_or_gromacs/export_with_interchange.py deleted file mode 100644 index e1ea62259..000000000 --- a/examples/using_smirnoff_in_amber_or_gromacs/export_with_interchange.py +++ /dev/null @@ -1,172 +0,0 @@ -#!/usr/bin/env python -# coding: utf-8 - -#
-# Note: -# Interchange is in the process of replacing ParmEd in many workflows, but it still in an alpha testing phase. Our internal tests indicate it is reliable for many small-molecule systems, but it is not yet reliable for complex, multi-component systems and there are likely still rough edges throughout. Feedback is welcome on the Interchange issue tracker.
- -# ## Using OpenFF force fields in Amber and GROMACS -# -# The Open Forcefield Toolkit can create parametrized `openmm.System` objects that can be natively simulated with OpenMM. This example shows the Interchange project can enable parallel workflows using Amber and GROMACS. - -# ### Preparing an OpenFF Topology -# -# We start by loading a PDB file containing one copy of ethanol and cyclohexane. Our goal is to create an OpenFF `Topology` object describing this system that we can parametrize with the SMIRNOFF-format "Sage" force field. -# -# The two `Molecule` objects created from the SMILES strings can contain information such as partial charges and stereochemistry that is not included in an OpenMM topology. In this example, partial charges are not explicitly given, and `ForceField` will assign AM1/BCC charges as specified by the "Sage" force field. Note that the OpenFF Toolkit produces partial charges that do not depend on the input conformation of parameterized molecules. See the [FAQ](https://open-forcefield-toolkit.readthedocs.io/en/latest/faq.html#the-partial-charges-generated-by-the-toolkit-don-t-seem-to-depend-on-the-molecule-s-conformation-is-this-a-bug) for more information. - -# In[ ]: - - -try: - from openmm import app -except ImportError: - from simtk.openmm import app - -from openff.toolkit.utils import get_data_file_path -from openff.toolkit.topology import Molecule, Topology -from openff.toolkit.typing.engines.smirnoff import ForceField - - -# In[ ]: - - -ethanol = Molecule.from_smiles("CCO") -cyclohexane = Molecule.from_smiles("C1CCCCC1") - -# Load the PDB file using OpenMM and save the OpenMM Topology -pdbfile = app.PDBFile( - get_data_file_path("systems/test_systems/1_cyclohexane_1_ethanol.pdb") -) -omm_topology = pdbfile.topology - -# Create the OpenFF Topology. -off_topology = Topology.from_openmm( - omm_topology, unique_molecules=[ethanol, cyclohexane] -) -off_topology - - -# ### Preparing an OpenFF ForceField -# -# Once the `ForceField` class is imported, the only decision to make is which force field to use. An exhaustive list of force fields released by the Open Force Field Initiative can be found [here](from openff.toolkit.typing.engines.smirnoff import ForceField -# ). -# -# Here we will use force field from the "Sage" line. - -# In[ ]: - - -forcefield = ForceField("openff-2.0.0.offxml") -forcefield - - -# ### Preparing an OpenMM System -# -# Once a force field and topology have been loaded, an `openmm.System` can be generated natively with the OpenFF Toolkit. - -# In[ ]: - - -omm_system = forcefield.create_openmm_system(off_topology) -omm_system - - -# ### Preparing an Interchange object -# -# To exports to engines other than OpenMM, we will make use of the [Interchange](https://openff-interchange.readthedocs.io/) project. There is a high-level `Interchange.from_smirnoff` function that consumes OpenFF Toolkit and ForceField objects and produces an `Interchange` object which can then be exported to formats understood by other molecular simulation engines. This extra step is needed to provide a clean interface between _applied_ parameters and engines. Note also that this step does not require an OpenMM System to be generated; `ForceField.create_openmm_system` does not need to be called to use Amber and GROMACS. - -# In[ ]: - - -from openff.interchange.components.interchange import Interchange - -interchange = Interchange.from_smirnoff( - force_field=forcefield, - topology=off_topology, -) -interchange.positions = pdbfile.positions -interchange - - -# ### Exporting to Amber and GROMACS files -# -# Once an `Interchange` object has been constructed, its API can be used to export to files understood by GROMACS, Amber, and more. - -# In[ ]: - - -# Export AMBER files. -interchange.to_prmtop("system.prmtop") -interchange.to_inpcrd("system.inpcrd") - -# Export GROMACS files. -interchange.to_top("system.top") -interchange.to_gro("system.gro") - - -# ### Validating the conversion to Amber files -# -# The Interchange project includes functions that take in an `Interchange` object and call out to simulation engines to run single-point energy calculations (with no minimization or dynamics) for the purpose of validating the export layer with each engine. Under the hood, each of these functions calls API points like those used above while converting to files understood by each engine. These rely on having each engine installed and accessible in the current `$PATH`. - -# In[ ]: - - -from openff.interchange.drivers import get_amber_energies, get_openmm_energies - - -# In[ ]: - - -openmm_energies = get_openmm_energies(interchange) -openmm_energies.energies - - -# In[ ]: - - -# get_ipython().system('cat system.inpcrd') -# get_ipython().system('cat system.prmtop') -# get_ipython().system('cat system.top') -# get_ipython().system('cat system.gro') -amber_energies = get_amber_energies(interchange) -amber_energies.energies - - -# ### Appendix: Validating the conversion to GROMACS and LAMMPS files -# -# If GROMACS and/or LAMMPS are installed on your machine, the same comparisons can also take place with those engines. They are available via `conda` by running a command like: -# -# ```conda install gromacs lammps -c conda-forge -c bioconda``` - -# In[ ]: - - -from distutils.spawn import find_executable -from pprint import pprint - -from openff.interchange.drivers import get_gromacs_energies, get_lammps_energies - - -# In[ ]: - - -if find_executable("lmp_serial"): - pprint(get_lammps_energies(interchange).energies) - - -# In[ ]: - - -if find_executable("gmx"): - pprint(get_gromacs_energies(interchange).energies) - - -# Finally, there is a helper function `get_summary_data` that will attempt to run drivers of each engine. A summary reported is prepared as a Pandas `DataFrame`. - -# In[ ]: - - -from openff.interchange.drivers.all import get_summary_data - -get_summary_data(interchange) From 5363f563d8056843af52176e896f0b4ce86b0708 Mon Sep 17 00:00:00 2001 From: "Matthew W. Thompson" Date: Fri, 12 Nov 2021 18:07:13 -0600 Subject: [PATCH 23/24] Pin to interchange v0.1.3+ --- examples/environment.yaml | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/examples/environment.yaml b/examples/environment.yaml index 923bc43de..9c1354a9c 100644 --- a/examples/environment.yaml +++ b/examples/environment.yaml @@ -17,7 +17,7 @@ dependencies: - numpy - openmm - openff-forcefields >=2.0.0 - - openff-interchange >=0.1.2 + - openff-interchange >=0.1.3 - gromacs - lammps - ambertools From f41cee62d955c1d8afd1e3393264bc92fa6f325d Mon Sep 17 00:00:00 2001 From: "Matthew W. Thompson" Date: Tue, 16 Nov 2021 10:07:46 -0600 Subject: [PATCH 24/24] Update CHANGELOG --- docs/releasehistory.md | 5 +++++ 1 file changed, 5 insertions(+) diff --git a/docs/releasehistory.md b/docs/releasehistory.md index 04102e1f3..1393475e8 100644 --- a/docs/releasehistory.md +++ b/docs/releasehistory.md @@ -8,6 +8,11 @@ Releases follow the `major.minor.micro` scheme recommended by [PEP440](https://w ## Current Development +### Examples added + +- [PR #1113](https://github.com/openforcefield/openff-toolkit/pull/1113): Updates the Amber/GROMACS + example to use Interchange. + ## 0.10.1 Minor feature and bugfix release