Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Allow parsing of offxml strings #288

Merged
merged 5 commits into from
Jul 27, 2023
Merged
Show file tree
Hide file tree
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
19 changes: 12 additions & 7 deletions openmmforcefields/generators/template_generators.py
Original file line number Diff line number Diff line change
Expand Up @@ -1289,14 +1289,19 @@ def __init__(self, molecules=None, cache=None, forcefield=None):

# Create ForceField object
import openff.toolkit.typing.engines.smirnoff
try:
filename = forcefield
if not filename.endswith('.offxml'):
filename += '.offxml'

# check for an installed force field
available_force_fields = openff.toolkit.typing.engines.smirnoff.get_available_force_fields()
if (filename := forcefield + ".offxml") in available_force_fields or (filename := forcefield) in available_force_fields:
self._smirnoff_forcefield = openff.toolkit.typing.engines.smirnoff.ForceField(filename)
except Exception as e:
_logger.error(e)
raise ValueError(f"Can't find specified SMIRNOFF force field ({forcefield}) in install paths") from e

# just try parsing the input and let openff handle the error
else:
try:
self._smirnoff_forcefield = openff.toolkit.typing.engines.smirnoff.ForceField(forcefield)
except Exception as e:
_logger.error(e)
raise ValueError(f"Can't find specified SMIRNOFF force field ({forcefield}) in install paths or parse the input as a string.") from e

# Delete constraints, if present
if 'Constraints' in self._smirnoff_forcefield._parameter_handlers:
Expand Down
37 changes: 37 additions & 0 deletions openmmforcefields/tests/test_template_generators.py
Original file line number Diff line number Diff line change
@@ -1,13 +1,16 @@
import copy
import logging
import os
import pytest
import tempfile
import unittest
import pytest

import numpy as np
import openmm
from openff.toolkit.topology import Molecule
from openff.toolkit.typing.engines.smirnoff import ForceField as OFFForceField
from openff.units import unit as OFFUnit
from openmm.app import PME, ForceField, Modeller, NoCutoff, PDBFile

from openmmforcefields.generators import (
Expand Down Expand Up @@ -855,6 +858,40 @@ def test_version(self):
assert generator.smirnoff_filename.endswith(forcefield + '.offxml')
assert os.path.exists(generator.smirnoff_filename)

def test_bespoke_force_field(self):
"""
Make sure a molecule can be parameterised using a bespoke force field passed as a string to
the template generator.
"""

custom_sage = OFFForceField("openff-2.0.0.offxml")
# Create a simple molecule with one bond type
ethane = Molecule.from_smiles("C")
# Label ethane to get the bond type (not hard coded incase this changes in future)
bond_parameter = custom_sage.label_molecules(ethane.to_topology())[0]["Bonds"][(0, 1)]
# Edit the bond parameter
bonds = custom_sage.get_parameter_handler("Bonds")
new_parameter = bonds[bond_parameter.smirks]
new_parameter.length = 2 * OFFUnit.angstrom

# Use the custom sage passed as string to build a template and an openmm system
generator = SMIRNOFFTemplateGenerator(molecules=ethane, forcefield=custom_sage.to_string())

# Create a ForceField
openmm_forcefield = openmm.app.ForceField()
# Register the template generator
openmm_forcefield.registerTemplateGenerator(generator.generator)
# Use OpenMM app to generate the system
openmm_system = openmm_forcefield.createSystem(ethane.to_topology().to_openmm(), removeCMMotion=False,
nonbondedMethod=NoCutoff)

# Check the bond length has been updated
forces = {force.__class__.__name__: force for force in openmm_system.getForces()}
bond_force = forces["HarmonicBondForce"]
for i in range(bond_force.getNumBonds()):
_, _, length, _ = bond_force.getBondParameters(i)
assert pytest.approx(length.value_in_unit(openmm.unit.angstrom)) == 2


@pytest.mark.espaloma
class TestEspalomaTemplateGenerator(TestGAFFTemplateGenerator):
Expand Down