Skip to content

Commit

Permalink
ENH: Create Topologys from .prmtop files`
Browse files Browse the repository at this point in the history
  • Loading branch information
mattwthompson committed Mar 28, 2023
1 parent 54326e9 commit f16f853
Show file tree
Hide file tree
Showing 6 changed files with 98 additions and 3 deletions.
2 changes: 1 addition & 1 deletion openff/interchange/drivers/openmm.py
Original file line number Diff line number Diff line change
Expand Up @@ -95,7 +95,7 @@ def _get_openmm_energies(
positions: openmm.unit.Quantity,
round_positions: Optional[int],
platform: str,
) -> EnergyReport:
) -> Dict[int, openmm.Force]:
"""Given prepared `openmm` objects, run a single-point energy calculation."""
for index, force in enumerate(system.getForces()):
force.setForceGroup(index)
Expand Down
1 change: 1 addition & 0 deletions openff/interchange/interop/amber/_import/__init__.py
Original file line number Diff line number Diff line change
@@ -0,0 +1 @@
"""Interfaces with Amber."""
91 changes: 91 additions & 0 deletions openff/interchange/interop/amber/_import/_import.py
Original file line number Diff line number Diff line change
@@ -0,0 +1,91 @@
"""Interfaces with Amber."""
from typing import TYPE_CHECKING, Dict, List

if TYPE_CHECKING:
from openff.toolkit import Topology

from openff.interchange import Interchange


def from_prmtop(
file: str,
) -> "Interchange":
"""Import from a prmtop file."""
from openff.interchange import Interchange

interchange = Interchange()

data: Dict[str, List[str]] = dict()

with open(file) as f:
for chunk in f.read().split(r"%FLAG"):
tag, format, *_data = chunk.strip().split()

if tag == "%VERSION":
continue

data[tag] = _data

interchange.topology = _make_topology(data)

return interchange


def _make_topology(data: Dict[str, List[str]]) -> "Topology":
"""Make a topology from the data."""
from openff.toolkit import Topology
from openff.toolkit.topology._mm_molecule import _SimpleMolecule

Topology._add_bond = _add_bond

topology = Topology()

start_index = 0

for molecule_index in range(int(data["POINTERS"][11])):
molecule = _SimpleMolecule()

end_index = start_index + int(data["ATOMS_PER_MOLECULE"][molecule_index])

for atom_index in range(start_index, end_index):
# TODO: Check for isotopes (unsupported) or otherwise botches atomic masses
molecule.add_atom(
atomic_number=int(data["ATOMIC_NUMBER"][atom_index]),
name=data["ATOM_NAME"][atom_index],
)

topology.add_molecule(molecule)

start_index = end_index

bonds: List[str] = data["BONDS_INC_HYDROGEN"] + data["BONDS_WITHOUT_HYDROGEN"]

# third value in each triplet is an index to the bond type
for n1, n2 in zip(bonds[::3], bonds[1::3]):
# See BONDS_INC_HYDROGEN in https://ambermd.org/prmtop.pdf
# For run-time efficiency, the atom indexes are actually indexes into a coordinate array,
# so the actual atom index A is calculated from the coordinate array index N by A = N/3 + 1

a1 = int(int(n1) / 3)
a2 = int(int(n2) / 3)

topology._add_bond(int(a1), int(a2))

return topology


def _add_bond(self, atom1_index: int, atom2_index: int):
atom1 = self.atom(atom1_index)
atom2 = self.atom(atom2_index)

if atom1.molecule is not atom2.molecule:
raise ValueError(
"Cannot add a bond between atoms in different molecules.",
)

molecule = atom1.molecule

molecule.add_bond(
atom1,
atom2,
)
3 changes: 2 additions & 1 deletion openff/interchange/interop/amber/export/__init__.py
Original file line number Diff line number Diff line change
@@ -1 +1,2 @@
"""Exports to Amber files."""
"""Exports to Amber files."""
from openff.interchange.interop.amber.export.export import to_inpcrd, to_prmtop
3 changes: 2 additions & 1 deletion openff/interchange/smirnoff/_gbsa.py
Original file line number Diff line number Diff line change
Expand Up @@ -6,6 +6,7 @@

from openff.interchange.components.potentials import Potential
from openff.interchange.constants import kcal_mol_a2
from openff.interchange.exceptions import InvalidParameterHandlerError
from openff.interchange.smirnoff._base import SMIRNOFFCollection

if TYPE_CHECKING:
Expand Down Expand Up @@ -70,7 +71,7 @@ def store_potentials(self, parameter_handler: GBSAHandler) -> None:
self.potentials[potential_key] = potential

@classmethod
def create( # type: ignore[override]
def create(
cls,
parameter_handler: GBSAHandler,
topology: "Topology",
Expand Down
1 change: 1 addition & 0 deletions setup.cfg
Original file line number Diff line number Diff line change
Expand Up @@ -19,6 +19,7 @@ per-file-ignores =
openff/interchange/components/*:F821
openff/interchange/__init__.py:F401
openff/interchange/tests/data/*:INP001
openff/interchange/interop/amber/export/__init__.py:F401
plugins/*:INP001

[isort]
Expand Down

0 comments on commit f16f853

Please sign in to comment.