diff --git a/news/abfe_host_smc.rst b/news/abfe_host_smc.rst new file mode 100644 index 000000000..b7c31a871 --- /dev/null +++ b/news/abfe_host_smc.rst @@ -0,0 +1,25 @@ +**Added:** + +* Add support for SmallMoleculeComponents to act + as host in the ``AbsoluteBindingProtocol`` + (`PR 2020 `_). + +**Changed:** + +* + +**Deprecated:** + +* + +**Removed:** + +* + +**Fixed:** + +* + +**Security:** + +* diff --git a/news/abfe_user_dfn_restraints.rst b/news/abfe_user_dfn_restraints.rst new file mode 100644 index 000000000..65c1d3b75 --- /dev/null +++ b/news/abfe_user_dfn_restraints.rst @@ -0,0 +1,24 @@ +**Added:** + +* Added support for user-defined Boresch restraints in the + ABFE Protocol (`PR #2019 `_). + +**Changed:** + +* + +**Deprecated:** + +* + +**Removed:** + +* + +**Fixed:** + +* + +**Security:** + +* diff --git a/src/openfe/protocols/openmm_afe/abfe_units.py b/src/openfe/protocols/openmm_afe/abfe_units.py index 343cbd2e8..3969db0cd 100644 --- a/src/openfe/protocols/openmm_afe/abfe_units.py +++ b/src/openfe/protocols/openmm_afe/abfe_units.py @@ -26,7 +26,7 @@ from rdkit import Chem from openfe.protocols.openmm_afe.equil_afe_settings import ( - BoreschRestraintSettings, + ABFEBoreschRestraintSettings, SettingsBaseModel, ) from openfe.protocols.openmm_utils import system_validation @@ -218,7 +218,7 @@ def _get_boresch_restraint( guest_atom_ids: list[int], host_atom_ids: list[int], temperature: Quantity, - settings: BoreschRestraintSettings, + settings: ABFEBoreschRestraintSettings, ) -> tuple[BoreschRestraintGeometry, BoreschRestraint]: """ Get a Boresch-like restraint Geometry and OpenMM restraint force @@ -236,7 +236,7 @@ def _get_boresch_restraint( A list of atom indices defining the host molecules in the universe. temperature : openff.units.Quantity The temperature of the simulation where the restraint will be added. - settings : BoreschRestraintSettings + settings : ABFEBoreschRestraintSettings Settings on how the Boresch-like restraint should be defined. Returns @@ -254,6 +254,12 @@ def _get_boresch_restraint( guest_rdmol=guest_rdmol, guest_idxs=guest_atom_ids, host_idxs=host_atom_ids, + guest_restraint_atoms_idxs=list(settings.guest_restraint_ids) + if settings.guest_restraint_ids is not None + else None, + host_restraint_atoms_idxs=list(settings.host_restraint_ids) + if settings.host_restraint_ids is not None + else None, host_selection=settings.host_selection, anchor_finding_strategy=settings.anchor_finding_strategy, dssp_filter=settings.dssp_filter, @@ -355,7 +361,7 @@ def _add_restraints( self.shared_basepath / settings["equil_output_settings"].production_trajectory_filename, ) - if isinstance(settings["restraint_settings"], BoreschRestraintSettings): + if isinstance(settings["restraint_settings"], ABFEBoreschRestraintSettings): rest_geom, restraint = self._get_boresch_restraint( univ, guest_rdmol, diff --git a/src/openfe/protocols/openmm_afe/equil_afe_settings.py b/src/openfe/protocols/openmm_afe/equil_afe_settings.py index 0941b6784..5853b4f76 100644 --- a/src/openfe/protocols/openmm_afe/equil_afe_settings.py +++ b/src/openfe/protocols/openmm_afe/equil_afe_settings.py @@ -213,6 +213,29 @@ def must_be_all(cls, v): return v +class ABFEBoreschRestraintSettings(BoreschRestraintSettings): + host_restraint_ids: tuple[int, int, int] | None = None + """ + The indices of the host component atoms to restrain. + The entries define the H0, H1, and H2 atoms in order. + If defined, these will override any automatic selection. + """ + guest_restraint_ids: tuple[int, int, int] | None = None + """ + The indices of the guest component atoms to restraint. + The entries define the G0, G1, and G2 atoms in order. + If defined, these will override any automatic selection. + """ + + @field_validator("guest_restraint_ids", "host_restraint_ids") + def positive_idxs_three_tuple(cls, v): + if v is not None: + if any([i < 0 for i in v]): + errmsg = "``guest_atoms`` and ``host_atoms`` cannot have negative indices." + raise ValueError(errmsg) + return v + + # This subclasses from SettingsBaseModel as it has vacuum_forcefield and # solvent_forcefield fields, not just a single forcefield_settings field class AbsoluteSolvationSettings(SettingsBaseModel): diff --git a/src/openfe/protocols/openmm_afe/equil_binding_afe_method.py b/src/openfe/protocols/openmm_afe/equil_binding_afe_method.py index 06dcbaed0..d7b14cea2 100644 --- a/src/openfe/protocols/openmm_afe/equil_binding_afe_method.py +++ b/src/openfe/protocols/openmm_afe/equil_binding_afe_method.py @@ -43,10 +43,10 @@ from openfe.due import Doi, due from openfe.protocols.openmm_afe.equil_afe_settings import ( + ABFEBoreschRestraintSettings, ABFEPreEquilOutputSettings, AbsoluteBindingSettings, AlchemicalSettings, - BoreschRestraintSettings, IntegratorSettings, LambdaSettings, MDSimulationSettings, @@ -174,7 +174,7 @@ def _default_settings(cls): engine_settings=OpenMMEngineSettings(), solvent_integrator_settings=IntegratorSettings(), complex_integrator_settings=IntegratorSettings(), - restraint_settings=BoreschRestraintSettings(), + restraint_settings=ABFEBoreschRestraintSettings(), solvent_equil_simulation_settings=MDSimulationSettings( equilibration_length_nvt=0.1 * offunit.nanosecond, equilibration_length=0.2 * offunit.nanosecond, @@ -278,8 +278,12 @@ def _validate_endstates( If the alchemical species is charged. """ if not (stateA.contains(ProteinComponent) and stateB.contains(ProteinComponent)): - errmsg = "No ProteinComponent found" - raise ValueError(errmsg) + # Check if there is a suitable SmallMoleculeComponent that could + # be the host molecule instead. + smcs = stateA.get_components_of_type(SmallMoleculeComponent) + if all(smc in stateA.component_diff(stateB)[0] for smc in smcs): + errmsg = "No suitable host molecule found" + raise ValueError(errmsg) if not (stateA.contains(SolventComponent) and stateB.contains(SolventComponent)): errmsg = "No SolventComponent found" diff --git a/src/openfe/protocols/restraint_utils/settings.py b/src/openfe/protocols/restraint_utils/settings.py index 231b6c141..8ab8dcd86 100644 --- a/src/openfe/protocols/restraint_utils/settings.py +++ b/src/openfe/protocols/restraint_utils/settings.py @@ -174,17 +174,6 @@ class BoreschRestraintSettings(BaseRestraintSettings): Boresch-like restraint search parameter. The maximum distance between any host atom and the guest G0 atom. Must be in units compatible with nanometer. """ - # TODO: re-enable this (Issue #1555) - # host_atoms: Optional[list[int]] = None - # """ - # The indices of the host component atoms to restrain. - # If defined, these will override any automatic selection. - # """ - # guest_atoms: Optional[list[int]] = None - # """ - # The indices of the guest component atoms to restraint. - # If defined, these will override any automatic selection. - # """ anchor_finding_strategy: Literal["multi-residue", "bonded"] = "bonded" """ The Boresch atom picking strategy to use. @@ -193,11 +182,3 @@ class BoreschRestraintSettings(BaseRestraintSettings): * `bonded`: pick host atoms that are bonded to each other. * `multi-residue`: pick host atoms which can span multiple residues. """ - - -# @field_validator("guest_atoms", "host_atoms") -# def positive_idxs_list(cls, v): -# if v is not None and any([i < 0 for i in v]): -# errmsg = "negative indices passed" -# raise ValueError(errmsg) -# return v diff --git a/src/openfe/tests/conftest.py b/src/openfe/tests/conftest.py index fcd41412f..3e9046b12 100644 --- a/src/openfe/tests/conftest.py +++ b/src/openfe/tests/conftest.py @@ -233,6 +233,17 @@ def benzene_modifications(): return files +@pytest.fixture(scope="session") +def benzene_modifications_am1bcc(): + files = {} + with resources.as_file(resources.files("openfe.tests.data")) as d: + fn = str(d / "benzene_modifications_am1bcc.sdf") + supp = Chem.SDMolSupplier(str(fn), removeHs=False) + for rdmol in supp: + files[rdmol.GetProp("_Name")] = SmallMoleculeComponent(rdmol) + return files + + @pytest.fixture(scope="session") def charged_benzene_modifications(): files = {} diff --git a/src/openfe/tests/data/benzene_modifications_am1bcc.sdf b/src/openfe/tests/data/benzene_modifications_am1bcc.sdf new file mode 100644 index 000000000..206c82d61 --- /dev/null +++ b/src/openfe/tests/data/benzene_modifications_am1bcc.sdf @@ -0,0 +1,295 @@ +benzene + RDKit 3D + + 12 12 0 0 0 0 0 0 0 0999 V2000 + 25.9780 5.3270 4.7790 C 0 0 0 0 0 0 0 0 0 0 0 0 + 26.3950 5.0740 3.4990 C 0 0 0 0 0 0 0 0 0 0 0 0 + 27.3400 5.8600 2.9020 C 0 0 0 0 0 0 0 0 0 0 0 0 + 27.8370 6.9210 3.5690 C 0 0 0 0 0 0 0 0 0 0 0 0 + 27.4200 7.1960 4.8560 C 0 0 0 0 0 0 0 0 0 0 0 0 + 26.4980 6.3790 5.4690 C 0 0 0 0 0 0 0 0 0 0 0 0 + 25.2298 4.6859 5.2451 H 0 0 0 0 0 0 0 0 0 0 0 0 + 25.9676 4.2351 2.9497 H 0 0 0 0 0 0 0 0 0 0 0 0 + 27.6890 5.6311 1.8951 H 0 0 0 0 0 0 0 0 0 0 0 0 + 28.5730 7.5660 3.0889 H 0 0 0 0 0 0 0 0 0 0 0 0 + 27.8209 8.0598 5.3863 H 0 0 0 0 0 0 0 0 0 0 0 0 + 26.1874 6.5720 6.4958 H 0 0 0 0 0 0 0 0 0 0 0 0 + 1 2 2 0 + 1 6 1 0 + 1 7 1 0 + 2 3 1 0 + 2 8 1 0 + 3 4 2 0 + 3 9 1 0 + 4 5 1 0 + 4 10 1 0 + 5 6 2 0 + 5 11 1 0 + 6 12 1 0 +M END + +> +benzene + +> +-0.13 -0.13 -0.13 -0.13 -0.13 -0.13 0.13 0.13 0.13 0.13 0.13 0.13 + +$$$$ +toluene + RDKit 3D + + 15 15 0 0 0 0 0 0 0 0999 V2000 + 28.9072 8.7434 5.1220 H 0 0 0 0 0 0 0 0 0 0 0 0 + 28.1966 8.1433 6.6393 H 0 0 0 0 0 0 0 0 0 0 0 0 + 27.9864 8.4164 5.6052 C 0 0 0 0 0 0 0 0 0 0 0 0 + 27.2579 9.2269 5.5838 H 0 0 0 0 0 0 0 0 0 0 0 0 + 25.9780 5.3270 4.7790 C 0 0 0 0 0 0 0 0 0 0 0 0 + 26.3950 5.0740 3.4990 C 0 0 0 0 0 0 0 0 0 0 0 0 + 27.3400 5.8600 2.9020 C 0 0 0 0 0 0 0 0 0 0 0 0 + 27.8370 6.9210 3.5690 C 0 0 0 0 0 0 0 0 0 0 0 0 + 27.4200 7.1960 4.8560 C 0 0 0 0 0 0 0 0 0 0 0 0 + 26.4980 6.3790 5.4690 C 0 0 0 0 0 0 0 0 0 0 0 0 + 25.2298 4.6859 5.2451 H 0 0 0 0 0 0 0 0 0 0 0 0 + 25.9676 4.2351 2.9497 H 0 0 0 0 0 0 0 0 0 0 0 0 + 27.6890 5.6311 1.8951 H 0 0 0 0 0 0 0 0 0 0 0 0 + 28.5730 7.5660 3.0889 H 0 0 0 0 0 0 0 0 0 0 0 0 + 26.1874 6.5720 6.4958 H 0 0 0 0 0 0 0 0 0 0 0 0 + 1 3 1 0 + 2 3 1 0 + 3 4 1 0 + 5 6 2 0 + 5 10 1 0 + 5 11 1 0 + 6 7 1 0 + 6 12 1 0 + 7 8 2 0 + 7 13 1 0 + 8 9 1 0 + 8 14 1 0 + 3 9 1 0 + 9 10 2 0 + 10 15 1 0 +M END + +> +toluene + +> +0.044033066666666669 0.044033066666666669 -0.053799933333333334 0.044033066666666669 -0.12699993333333334 -0.13499993333333335 -0.12699993333333334 -0.13099993333333335 +-0.077299933333333321 -0.13099993333333335 0.13000006666666666 0.13000006666666666 0.13000006666666666 0.13000006666666666 0.13000006666666666 + +$$$$ +phenol + RDKit 3D + + 13 13 0 0 0 0 0 0 0 0999 V2000 + 25.9780 5.3270 4.7790 C 0 0 0 0 0 0 0 0 0 0 0 0 + 26.3950 5.0740 3.4990 C 0 0 0 0 0 0 0 0 0 0 0 0 + 27.3400 5.8600 2.9020 C 0 0 0 0 0 0 0 0 0 0 0 0 + 27.8370 6.9210 3.5690 C 0 0 0 0 0 0 0 0 0 0 0 0 + 27.4200 7.1960 4.8560 C 0 0 0 0 0 0 0 0 0 0 0 0 + 26.4980 6.3790 5.4690 C 0 0 0 0 0 0 0 0 0 0 0 0 + 25.2298 4.6859 5.2451 H 0 0 0 0 0 0 0 0 0 0 0 0 + 25.9676 4.2351 2.9497 H 0 0 0 0 0 0 0 0 0 0 0 0 + 27.6890 5.6311 1.8951 H 0 0 0 0 0 0 0 0 0 0 0 0 + 28.5730 7.5660 3.0889 H 0 0 0 0 0 0 0 0 0 0 0 0 + 28.1311 8.0887 6.4624 H 0 0 0 0 0 0 0 0 0 0 0 0 + 26.1874 6.5720 6.4958 H 0 0 0 0 0 0 0 0 0 0 0 0 + 27.9460 8.3293 5.5517 O 0 0 0 0 0 0 0 0 0 0 0 0 + 1 2 2 0 + 1 6 1 0 + 1 7 1 0 + 2 3 1 0 + 2 8 1 0 + 3 4 2 0 + 3 9 1 0 + 4 5 1 0 + 4 10 1 0 + 5 6 2 0 + 5 13 1 0 + 6 12 1 0 + 11 13 1 0 +M END + +> +phenol + +> +-0.094423076923076915 -0.16592307692307692 -0.094423076923076915 -0.18492307692307691 0.12317692307692309 -0.18492307692307691 0.13307692307692309 0.13307692307692309 0.13307692307692309 +0.14157692307692307 0.41807692307692307 0.14157692307692307 -0.4990230769230769 + +$$$$ +benzonitrile + RDKit 3D + + 13 13 0 0 0 0 0 0 0 0999 V2000 + 28.5559 9.5700 6.2831 N 0 0 0 0 0 0 0 0 0 0 0 0 + 27.9981 8.4043 5.5824 C 0 0 0 0 0 0 0 0 0 0 0 0 + 25.9780 5.3270 4.7790 C 0 0 0 0 0 0 0 0 0 0 0 0 + 26.3950 5.0740 3.4990 C 0 0 0 0 0 0 0 0 0 0 0 0 + 27.3400 5.8600 2.9020 C 0 0 0 0 0 0 0 0 0 0 0 0 + 27.8370 6.9210 3.5690 C 0 0 0 0 0 0 0 0 0 0 0 0 + 27.4200 7.1960 4.8560 C 0 0 0 0 0 0 0 0 0 0 0 0 + 26.4980 6.3790 5.4690 C 0 0 0 0 0 0 0 0 0 0 0 0 + 25.2298 4.6859 5.2451 H 0 0 0 0 0 0 0 0 0 0 0 0 + 25.9676 4.2351 2.9497 H 0 0 0 0 0 0 0 0 0 0 0 0 + 27.6890 5.6311 1.8951 H 0 0 0 0 0 0 0 0 0 0 0 0 + 28.5730 7.5660 3.0889 H 0 0 0 0 0 0 0 0 0 0 0 0 + 26.1874 6.5720 6.4958 H 0 0 0 0 0 0 0 0 0 0 0 0 + 1 2 3 0 + 3 4 2 0 + 3 8 1 0 + 3 9 1 0 + 4 5 1 0 + 4 10 1 0 + 5 6 2 0 + 5 11 1 0 + 6 7 1 0 + 6 12 1 0 + 2 7 1 0 + 7 8 2 0 + 8 13 1 0 +M END + +> +benzonitrile + +> +-0.36380000000000001 0.23380000000000001 -0.13500000000000001 -0.107 -0.13500000000000001 -0.090999999999999998 -0.019000000000000003 -0.090999999999999998 0.14000000000000001 +0.13800000000000001 0.14000000000000001 0.14499999999999999 0.14499999999999999 + +$$$$ +anisole + RDKit 3D + + 16 16 0 0 0 0 0 0 0 0999 V2000 + 29.2873 8.8784 4.9226 C 0 0 0 0 0 0 0 0 0 0 0 0 + 29.5502 9.7990 5.4437 H 0 0 0 0 0 0 0 0 0 0 0 0 + 28.0548 8.3459 5.4720 O 0 0 0 0 0 0 0 0 0 0 0 0 + 30.0866 8.1484 5.0502 H 0 0 0 0 0 0 0 0 0 0 0 0 + 29.1525 9.0868 3.8612 H 0 0 0 0 0 0 0 0 0 0 0 0 + 25.9780 5.3270 4.7790 C 0 0 0 0 0 0 0 0 0 0 0 0 + 26.3950 5.0740 3.4990 C 0 0 0 0 0 0 0 0 0 0 0 0 + 27.3400 5.8600 2.9020 C 0 0 0 0 0 0 0 0 0 0 0 0 + 27.8370 6.9210 3.5690 C 0 0 0 0 0 0 0 0 0 0 0 0 + 27.4200 7.1960 4.8560 C 0 0 0 0 0 0 0 0 0 0 0 0 + 26.4980 6.3790 5.4690 C 0 0 0 0 0 0 0 0 0 0 0 0 + 25.2298 4.6859 5.2451 H 0 0 0 0 0 0 0 0 0 0 0 0 + 25.9676 4.2351 2.9497 H 0 0 0 0 0 0 0 0 0 0 0 0 + 27.6890 5.6311 1.8951 H 0 0 0 0 0 0 0 0 0 0 0 0 + 28.5730 7.5660 3.0889 H 0 0 0 0 0 0 0 0 0 0 0 0 + 26.1874 6.5720 6.4958 H 0 0 0 0 0 0 0 0 0 0 0 0 + 1 2 1 0 + 1 4 1 0 + 1 5 1 0 + 1 3 1 0 + 6 7 2 0 + 6 11 1 0 + 6 12 1 0 + 7 8 1 0 + 7 13 1 0 + 8 9 2 0 + 8 14 1 0 + 9 10 1 0 + 9 15 1 0 + 3 10 1 0 + 10 11 2 0 + 11 16 1 0 +M END + +> +anisole + +> +0.113825 0.043825000000000003 -0.32877500000000004 0.043825000000000003 0.043825000000000003 -0.097875000000000004 -0.16487500000000002 -0.097875000000000004 -0.17987500000000001 0.123225 +-0.17987500000000001 0.13212499999999999 0.13212499999999999 0.13212499999999999 0.14212499999999997 0.14212499999999997 + +$$$$ +benzaldehyde + RDKit 3D + + 14 14 0 0 0 0 0 0 0 0999 V2000 + 29.2079 8.8492 4.9632 O 0 0 0 0 0 0 0 0 0 0 0 0 + 27.5482 8.8691 6.4597 H 0 0 0 0 0 0 0 0 0 0 0 0 + 27.9981 8.4043 5.5824 C 0 0 0 0 0 0 0 0 0 0 0 0 + 25.9780 5.3270 4.7790 C 0 0 0 0 0 0 0 0 0 0 0 0 + 26.3950 5.0740 3.4990 C 0 0 0 0 0 0 0 0 0 0 0 0 + 27.3400 5.8600 2.9020 C 0 0 0 0 0 0 0 0 0 0 0 0 + 27.8370 6.9210 3.5690 C 0 0 0 0 0 0 0 0 0 0 0 0 + 27.4200 7.1960 4.8560 C 0 0 0 0 0 0 0 0 0 0 0 0 + 26.4980 6.3790 5.4690 C 0 0 0 0 0 0 0 0 0 0 0 0 + 25.2298 4.6859 5.2451 H 0 0 0 0 0 0 0 0 0 0 0 0 + 25.9676 4.2351 2.9497 H 0 0 0 0 0 0 0 0 0 0 0 0 + 27.6890 5.6311 1.8951 H 0 0 0 0 0 0 0 0 0 0 0 0 + 28.5730 7.5660 3.0889 H 0 0 0 0 0 0 0 0 0 0 0 0 + 26.1874 6.5720 6.4958 H 0 0 0 0 0 0 0 0 0 0 0 0 + 1 3 2 0 + 2 3 1 0 + 4 5 2 0 + 4 9 1 0 + 4 10 1 0 + 5 6 1 0 + 5 11 1 0 + 6 7 2 0 + 6 12 1 0 + 7 8 1 0 + 7 13 1 0 + 3 8 1 0 + 8 9 2 0 + 9 14 1 0 +M END + +> +benzaldehyde + +> +-0.52817142857142862 -0.0028714285714285795 0.5754285714285714 -0.14507142857142857 -0.098071428571428587 -0.14507142857142857 -0.078071428571428583 -0.19767142857142858 +-0.078071428571428583 0.13742857142857143 0.13492857142857143 0.13742857142857143 0.14392857142857141 0.14392857142857141 + +$$$$ +styrene + RDKit 3D + + 16 16 0 0 0 0 0 0 0 0999 V2000 + 29.2873 8.8784 4.9226 C 0 0 0 0 0 0 0 0 0 0 0 0 + 29.6609 8.3486 4.0463 H 0 0 0 0 0 0 0 0 0 0 0 0 + 29.8344 9.7353 5.3157 H 0 0 0 0 0 0 0 0 0 0 0 0 + 27.5365 8.8812 6.4825 H 0 0 0 0 0 0 0 0 0 0 0 0 + 27.9864 8.4164 5.6052 C 0 0 0 0 0 0 0 0 0 0 0 0 + 25.9780 5.3270 4.7790 C 0 0 0 0 0 0 0 0 0 0 0 0 + 26.3950 5.0740 3.4990 C 0 0 0 0 0 0 0 0 0 0 0 0 + 27.3400 5.8600 2.9020 C 0 0 0 0 0 0 0 0 0 0 0 0 + 27.8370 6.9210 3.5690 C 0 0 0 0 0 0 0 0 0 0 0 0 + 27.4200 7.1960 4.8560 C 0 0 0 0 0 0 0 0 0 0 0 0 + 26.4980 6.3790 5.4690 C 0 0 0 0 0 0 0 0 0 0 0 0 + 25.2298 4.6859 5.2451 H 0 0 0 0 0 0 0 0 0 0 0 0 + 25.9676 4.2351 2.9497 H 0 0 0 0 0 0 0 0 0 0 0 0 + 27.6890 5.6311 1.8951 H 0 0 0 0 0 0 0 0 0 0 0 0 + 28.5730 7.5660 3.0889 H 0 0 0 0 0 0 0 0 0 0 0 0 + 26.1874 6.5720 6.4958 H 0 0 0 0 0 0 0 0 0 0 0 0 + 1 2 1 0 + 1 3 1 0 + 1 5 2 0 + 4 5 1 0 + 6 7 2 0 + 6 11 1 0 + 6 12 1 0 + 7 8 1 0 + 7 13 1 0 + 8 9 2 0 + 8 14 1 0 + 9 10 1 0 + 9 15 1 0 + 5 10 1 0 + 10 11 2 0 + 11 16 1 0 +M END + +> +styrene + +> +-0.20899999999999999 0.11349999999999999 0.11349999999999999 0.123 -0.1132 -0.13100000000000001 -0.127 -0.13100000000000001 -0.11849999999999999 -0.057800000000000004 -0.11849999999999999 +0.13100000000000001 0.13100000000000001 0.13100000000000001 0.13150000000000001 0.13150000000000001 + +$$$$ diff --git a/src/openfe/tests/protocols/openmm_abfe/test_abfe_protocol.py b/src/openfe/tests/protocols/openmm_abfe/test_abfe_protocol.py index 6ba27ca2b..29139c285 100644 --- a/src/openfe/tests/protocols/openmm_abfe/test_abfe_protocol.py +++ b/src/openfe/tests/protocols/openmm_abfe/test_abfe_protocol.py @@ -41,6 +41,7 @@ from openfe.protocols.openmm_afe import ( AbsoluteBindingProtocol, ) +from openfe.protocols.restraint_utils.geometry import BoreschRestraintGeometry from .utils import UNIT_TYPES, _get_units @@ -614,6 +615,51 @@ def settings(self): return s +def test_user_restraint(benzene_modifications_am1bcc, T4_protein_component, tmp_path): + s = openmm_afe.AbsoluteBindingProtocol.default_settings() + s.protocol_repeats = 1 + s.engine_settings.compute_platform = "cpu" + s.restraint_settings.guest_restraint_ids = [0, 1, 2] + # Ca and C from VAL 87, and N from TYR 88 + s.restraint_settings.host_restraint_ids = [1383, 1384, 1398] + + protocol = openmm_afe.AbsoluteBindingProtocol(settings=s) + + stateA = gufe.ChemicalSystem( + { + "protein": T4_protein_component, + "benzene": benzene_modifications_am1bcc["benzene"], + "solvent": gufe.SolventComponent(), + } + ) + + stateB = gufe.ChemicalSystem( + { + "protein": T4_protein_component, + "solvent": gufe.SolventComponent(), + } + ) + + dag = protocol.create(stateA=stateA, stateB=stateB, mapping=None) + + complex_setup_units = _get_units(dag.protocol_units, UNIT_TYPES["complex"]["setup"]) + + results = complex_setup_units[0].run( + dry=True, scratch_basepath=tmp_path, shared_basepath=tmp_path + ) + + geom = BoreschRestraintGeometry.model_validate(results["restraint_geometry"]) + # This should be C1, C2, and C3 on the benzene + assert geom.guest_atoms == [2613, 2614, 2615] + assert geom.host_atoms == [1383, 1384, 1398] + assert pytest.approx(geom.r_aA0.to("nanometer").m, rel=1e-4) == 0.510798 + assert pytest.approx(geom.theta_A0.to("radians").m, rel=1e-4) == 1.20278 + assert pytest.approx(geom.theta_B0.to("radians").m, rel=1e-4) == 1.25705 + assert pytest.approx(geom.phi_A0.to("radians").m, rel=1e-4) == 0.86035 + assert pytest.approx(geom.phi_B0.to("radians").m, rel=1e-4) == 1.59444 + assert pytest.approx(geom.phi_C0.to("radians").m, rel=1e-4) == 2.92365 + + def test_user_charges(benzene_modifications, T4_protein_component, tmp_path): s = openmm_afe.AbsoluteBindingProtocol.default_settings() s.protocol_repeats = 1 diff --git a/src/openfe/tests/protocols/openmm_abfe/test_abfe_settings.py b/src/openfe/tests/protocols/openmm_abfe/test_abfe_settings.py index fad97e6ad..013c775ff 100644 --- a/src/openfe/tests/protocols/openmm_abfe/test_abfe_settings.py +++ b/src/openfe/tests/protocols/openmm_abfe/test_abfe_settings.py @@ -7,6 +7,7 @@ AbsoluteBindingProtocol, AbsoluteBindingSettings, ) +from openfe.protocols.openmm_afe.equil_afe_settings import ABFEBoreschRestraintSettings @pytest.fixture() @@ -70,6 +71,15 @@ def test_monotonic_lambda_windows(val, default_settings): lambda_settings.lambda_restraints = val["restraints"] +@pytest.mark.parametrize("parameter", ["host_restraint_ids", "guest_restraint_ids"]) +def test_boresch_restraints_restraint_negative_ids(parameter): + setting = ABFEBoreschRestraintSettings() + + errmsg = "``guest_atoms`` and ``host_atoms`` cannot have negative indices." + with pytest.raises(ValueError, match=errmsg): + setattr(setting, parameter, [1, 2, -3]) + + def test_equil_not_all_complex(default_settings): with pytest.raises(ValueError, match="output_indices must be all"): default_settings.complex_equil_output_settings.output_indices = "not water" diff --git a/src/openfe/tests/protocols/openmm_abfe/test_abfe_validation.py b/src/openfe/tests/protocols/openmm_abfe/test_abfe_validation.py index 378702434..e8f7108c3 100644 --- a/src/openfe/tests/protocols/openmm_abfe/test_abfe_validation.py +++ b/src/openfe/tests/protocols/openmm_abfe/test_abfe_validation.py @@ -121,16 +121,39 @@ def test_validate_no_protcomp( stateB = ChemicalSystem( { - "benzene": benzene_modifications["benzene"], "solvent": SolventComponent(), } ) - errmsg = "No ProteinComponent found" + errmsg = "No suitable host molecule found" with pytest.raises(ValueError, match=errmsg): AbsoluteBindingProtocol._validate_endstates(stateA, stateB) +def test_validate_smc_host( + benzene_modifications, +): + """ + Should pass if there's another smc present that could be a host. + """ + stateA = ChemicalSystem( + { + "benzene": benzene_modifications["benzene"], + "toluene": benzene_modifications["toluene"], + "solvent": SolventComponent(), + } + ) + + stateB = ChemicalSystem( + { + "toluene": benzene_modifications["toluene"], + "solvent": SolventComponent(), + } + ) + + AbsoluteBindingProtocol._validate_endstates(stateA, stateB) + + def test_validate_endstates_nosolvcomp_stateA(benzene_modifications, T4_protein_component): stateA = ChemicalSystem( {