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

Respect mapped SMILES #73

Merged
merged 2 commits into from
Oct 27, 2024
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
7 changes: 4 additions & 3 deletions absolv/runner.py
Original file line number Diff line number Diff line change
@@ -1,4 +1,5 @@
"""Run calculations defined by a config."""

import pathlib
import tempfile
import typing
Expand Down Expand Up @@ -219,7 +220,7 @@ def _run_phase_hremd(
n_steps_per_cycle=protocol.production_protocol.n_steps_per_cycle,
n_cycles=protocol.production_protocol.n_cycles,
trajectory_interval=protocol.production_protocol.trajectory_interval,
trajectory_enforce_pbc=protocol.production_protocol.trajectory_enforce_pbc
trajectory_enforce_pbc=protocol.production_protocol.trajectory_enforce_pbc,
)

with tempfile.TemporaryDirectory() as tmp_dir:
Expand Down Expand Up @@ -254,7 +255,7 @@ def run_eq(
config.temperature,
prepared_system_a,
platform,
None if output_dir is None else output_dir / "solvant-a",
None if output_dir is None else output_dir / "solvent-a",
)

dg_a, dg_a_std = results_a["ddG_kcal_mol"], results_a["ddG_error_kcal_mol"]
Expand All @@ -265,7 +266,7 @@ def run_eq(
config.temperature,
prepared_system_b,
platform,
None if output_dir is None else output_dir / "solvant-b",
None if output_dir is None else output_dir / "solvent-b",
)

dg_b, dg_b_std = results_b["ddG_kcal_mol"], results_b["ddG_error_kcal_mol"]
Expand Down
18 changes: 17 additions & 1 deletion absolv/setup.py
Original file line number Diff line number Diff line change
@@ -1,4 +1,5 @@
"""Setup systems ready for calculations."""

import errno
import functools
import logging
Expand All @@ -9,6 +10,7 @@

import numpy
import openff.toolkit
import openff.toolkit.utils.exceptions
import openff.utilities
import openmm

Expand Down Expand Up @@ -96,6 +98,17 @@ def _generate_input_file(
)


def _molecule_from_smiles(smiles: str) -> openff.toolkit.Molecule:
"""Create a molecule from a SMILES string."""

try:
molecule = openff.toolkit.Molecule.from_mapped_smiles(smiles)
except (openff.toolkit.utils.exceptions.SmilesParsingError, ValueError):
molecule = openff.toolkit.Molecule.from_smiles(smiles)

return molecule


def setup_system(
components: list[tuple[str, int]],
box_target_density: openmm.unit.Quantity = 0.95 * _G_PER_ML,
Expand All @@ -110,6 +123,9 @@ def setup_system(
``smiles_i`` is the SMILES representation of component `i` and
``count_i`` is the number of corresponding instances of that component
to create.

If any SMILES is fully indexed an attempt will be made to create components
with the same atom ordering as the indices.
box_target_density: Target mass density when approximating the box size for the
final system with units compatible with g / mL.
box_scale_factor: The amount to scale the approximate box size by.
Expand Down Expand Up @@ -139,7 +155,7 @@ def setup_system(
if smiles in molecules:
continue

molecule = openff.toolkit.Molecule.from_smiles(smiles)
molecule = _molecule_from_smiles(smiles)
molecule.generate_conformers(n_conformers=1)
molecule.name = f"component-{len(molecules)}.xyz"
molecules[smiles] = molecule
Expand Down
18 changes: 18 additions & 0 deletions absolv/tests/test_setup.py
Original file line number Diff line number Diff line change
Expand Up @@ -7,6 +7,7 @@
from absolv.setup import (
_approximate_box_size_by_density,
_generate_input_file,
_molecule_from_smiles,
setup_system,
)

Expand Down Expand Up @@ -53,6 +54,23 @@ def test_generate_input_file(monkeypatch):
assert actual_input_file == expected_input_file


@pytest.mark.parametrize(
"smiles, expected_order",
[
("O", [8, 1, 1]),
("[O:1]([H:2])[H:3]", [8, 1, 1]),
("[H:2][O:1][H:3]", [8, 1, 1]),
("[O:2]([H:1])[H:3]", [1, 8, 1]),
("[O:2]([H])[H:3]", [8, 1, 1]),
]
)
def test_molecule_from_smiles(smiles, expected_order):
molecule = _molecule_from_smiles(smiles)

assert [atom.atomic_number for atom in molecule.atoms] == expected_order



def test_setup_system_packmol_missing(monkeypatch):
monkeypatch.setattr(shutil, "which", lambda *_: None)

Expand Down
Loading