From 7ab3dc3c506e0d40b6bef1889dc12a52980c4e47 Mon Sep 17 00:00:00 2001 From: IAlibay Date: Fri, 5 Dec 2025 18:04:46 +0000 Subject: [PATCH 01/11] Add absolute alchemical settings --- openfe/protocols/openmm_afe/base.py | 2 +- .../openmm_afe/equil_afe_settings.py | 12 ++--- .../openmm_afe/equil_binding_afe_method.py | 8 +-- .../openmm_afe/equil_solvation_afe_method.py | 8 +-- openfe/protocols/openmm_septop/base.py | 2 +- .../openmm_septop/equil_septop_method.py | 8 +-- .../openmm_septop/equil_septop_settings.py | 10 +--- openfe/protocols/openmm_utils/omm_settings.py | 52 +++++++++++++++++++ 8 files changed, 71 insertions(+), 31 deletions(-) diff --git a/openfe/protocols/openmm_afe/base.py b/openfe/protocols/openmm_afe/base.py index 9cafdf47b..574743b08 100644 --- a/openfe/protocols/openmm_afe/base.py +++ b/openfe/protocols/openmm_afe/base.py @@ -337,7 +337,7 @@ def _handle_settings(self) -> dict[str, SettingsBaseModel]: * forcefield_settings : OpenMMSystemGeneratorFFSettings * thermo_settings : ThermoSettings * solvation_settings : BaseSolvationSettings - * alchemical_settings : AlchemicalSettings + * alchemical_settings : AbsoluteAlchemicalSettings * lambda_settings : LambdaSettings * engine_settings : OpenMMEngineSettings * integrator_settings : IntegratorSettings diff --git a/openfe/protocols/openmm_afe/equil_afe_settings.py b/openfe/protocols/openmm_afe/equil_afe_settings.py index e652cfa1b..89e50711b 100644 --- a/openfe/protocols/openmm_afe/equil_afe_settings.py +++ b/openfe/protocols/openmm_afe/equil_afe_settings.py @@ -34,6 +34,7 @@ OpenFFPartialChargeSettings, OpenMMEngineSettings, OpenMMSolvationSettings, + AbsoluteAlchemicalSettings, ) from openfe.protocols.restraint_utils.settings import ( BaseRestraintSettings, @@ -41,13 +42,6 @@ ) -class AlchemicalSettings(SettingsBaseModel): - """Settings for the alchemical protocol - - Empty place holder for right now. - """ - - class LambdaSettings(SettingsBaseModel): """Lambda schedule settings. @@ -210,7 +204,7 @@ def must_be_positive(cls, v): """Settings for solvating the system.""" # Alchemical settings - alchemical_settings: AlchemicalSettings + alchemical_settings: AbsoluteAlchemicalSettings """ Alchemical protocol settings. """ @@ -321,7 +315,7 @@ def must_be_positive(cls, v): """Settings for solvating the system in the complex.""" # Alchemical settings - alchemical_settings: AlchemicalSettings + alchemical_settings: AbsoluteAlchemicalSettings """ Alchemical protocol settings. """ diff --git a/openfe/protocols/openmm_afe/equil_binding_afe_method.py b/openfe/protocols/openmm_afe/equil_binding_afe_method.py index dc12f1a8e..ad704301e 100644 --- a/openfe/protocols/openmm_afe/equil_binding_afe_method.py +++ b/openfe/protocols/openmm_afe/equil_binding_afe_method.py @@ -58,7 +58,7 @@ from openfe.protocols.openmm_afe.equil_afe_settings import ( ABFEPreEquilOutputSettings, AbsoluteBindingSettings, - AlchemicalSettings, + AbsoluteAlchemicalSettings, BoreschRestraintSettings, IntegratorSettings, LambdaSettings, @@ -553,7 +553,7 @@ def _default_settings(cls): temperature=298.15 * offunit.kelvin, pressure=1 * offunit.bar, ), - alchemical_settings=AlchemicalSettings(), + alchemical_settings=AbsoluteAlchemicalSettings(), solvent_lambda_settings=LambdaSettings( lambda_elec=[ 0.0, 0.25, 0.5, 0.75, 1.0, @@ -958,7 +958,7 @@ def _handle_settings(self) -> dict[str, SettingsBaseModel]: * thermo_settings : ThermoSettings * charge_settings : OpenFFPartialChargeSettings * solvation_settings : OpenMMSolvationSettings - * alchemical_settings : AlchemicalSettings + * alchemical_settings : AbsoluteAlchemicalSettings * lambda_settings : LambdaSettings * engine_settings : OpenMMEngineSettings * integrator_settings : IntegratorSettings @@ -1306,7 +1306,7 @@ def _handle_settings(self) -> dict[str, SettingsBaseModel]: * thermo_settings : ThermoSettings * charge_settings : OpenFFPartialChargeSettings * solvation_settings : OpenMMSolvationSettings - * alchemical_settings : AlchemicalSettings + * alchemical_settings : AbsoluteAlchemicalSettings * lambda_settings : LambdaSettings * engine_settings : OpenMMEngineSettings * integrator_settings : IntegratorSettings diff --git a/openfe/protocols/openmm_afe/equil_solvation_afe_method.py b/openfe/protocols/openmm_afe/equil_solvation_afe_method.py index 4ffb94770..ae4de822f 100644 --- a/openfe/protocols/openmm_afe/equil_solvation_afe_method.py +++ b/openfe/protocols/openmm_afe/equil_solvation_afe_method.py @@ -55,7 +55,7 @@ from openfe.due import Doi, due from openfe.protocols.openmm_afe.equil_afe_settings import ( AbsoluteSolvationSettings, - AlchemicalSettings, + AbsoluteAlchemicalSettings, IntegratorSettings, LambdaSettings, MDOutputSettings, @@ -442,7 +442,7 @@ def _default_settings(cls): temperature=298.15 * unit.kelvin, pressure=1 * unit.bar, ), - alchemical_settings=AlchemicalSettings(), + alchemical_settings=AbsoluteAlchemicalSettings(), lambda_settings=LambdaSettings( lambda_elec=[ 0.0, 0.25, 0.5, 0.75, 1.0, 1.0, 1.0, @@ -859,7 +859,7 @@ def _handle_settings(self) -> dict[str, SettingsBaseModel]: * thermo_settings : ThermoSettings * charge_settings : OpenFFPartialChargeSettings * solvation_settings : OpenMMSolvationSettings - * alchemical_settings : AlchemicalSettings + * alchemical_settings : AbsoluteAlchemicalSettings * lambda_settings : LambdaSettings * engine_settings : OpenMMEngineSettings * integrator_settings : IntegratorSettings @@ -934,7 +934,7 @@ def _handle_settings(self) -> dict[str, SettingsBaseModel]: * thermo_settings : ThermoSettings * charge_settings : OpenFFPartialChargeSettings * solvation_settings : OpenMMSolvationSettings - * alchemical_settings : AlchemicalSettings + * alchemical_settings : AbsoluteAlchemicalSettings * lambda_settings : LambdaSettings * engine_settings : OpenMMEngineSettings * integrator_settings : IntegratorSettings diff --git a/openfe/protocols/openmm_septop/base.py b/openfe/protocols/openmm_septop/base.py index 18d0209b8..323a3e6db 100644 --- a/openfe/protocols/openmm_septop/base.py +++ b/openfe/protocols/openmm_septop/base.py @@ -766,7 +766,7 @@ def _handle_settings(self): * forcefield_settings : OpenMMSystemGeneratorFFSettings * thermo_settings : ThermoSettings * solvation_settings : BaseSolvationSettings - * alchemical_settings : AlchemicalSettings + * alchemical_settings : AbsoluteAlchemicalSettings * lambda_settings : LambdaSettings * engine_settings : OpenMMEngineSettings * integrator_settings : IntegratorSettings diff --git a/openfe/protocols/openmm_septop/equil_septop_method.py b/openfe/protocols/openmm_septop/equil_septop_method.py index cefabf359..8391c5ae8 100644 --- a/openfe/protocols/openmm_septop/equil_septop_method.py +++ b/openfe/protocols/openmm_septop/equil_septop_method.py @@ -68,7 +68,7 @@ from openfe.due import Doi, due from openfe.protocols.openmm_septop.equil_septop_settings import ( - AlchemicalSettings, + AbsoluteAlchemicalSettings, IntegratorSettings, LambdaSettings, MDSimulationSettings, @@ -241,7 +241,7 @@ def _handle_settings(self) -> dict[str, SettingsBaseModel]: * thermo_settings : ThermoSettings * charge_settings : OpenFFPartialChargeSettings * solvation_settings : OpenMMSolvationSettings - * alchemical_settings : AlchemicalSettings + * alchemical_settings : AbsoluteAlchemicalSettings * lambda_settings : LambdaSettings * engine_settings : OpenMMEngineSettings * integrator_settings : IntegratorSettings @@ -327,7 +327,7 @@ def _handle_settings(self) -> dict[str, SettingsBaseModel]: * thermo_settings : ThermoSettings * charge_settings : OpenFFPartialChargeSettings * solvation_settings : OpenMMSolvationSettings - * alchemical_settings : AlchemicalSettings + * alchemical_settings : AbsoluteAlchemicalSettings * lambda_settings : LambdaSettings * engine_settings : OpenMMEngineSettings * integrator_settings : IntegratorSettings @@ -894,7 +894,7 @@ def _default_settings(cls): temperature=298.15 * unit.kelvin, pressure=1 * unit.bar, ), - alchemical_settings=AlchemicalSettings(), + alchemical_settings=AbsoluteAlchemicalSettings(), solvent_lambda_settings=LambdaSettings( lambda_elec_A=[ 0.0, diff --git a/openfe/protocols/openmm_septop/equil_septop_settings.py b/openfe/protocols/openmm_septop/equil_septop_settings.py index 6e6522acb..59af68b40 100644 --- a/openfe/protocols/openmm_septop/equil_septop_settings.py +++ b/openfe/protocols/openmm_septop/equil_septop_settings.py @@ -32,17 +32,11 @@ OpenFFPartialChargeSettings, OpenMMEngineSettings, OpenMMSolvationSettings, + AbsoluteAlchemicalSettings, ) from openfe.protocols.restraint_utils.settings import BaseRestraintSettings -class AlchemicalSettings(SettingsBaseModel): - """Settings for the alchemical protocol - - Empty place holder for right now. - """ - - class LambdaSettings(SettingsBaseModel): """Lambda schedule settings. @@ -369,7 +363,7 @@ def must_be_positive(cls, v): """Settings for solvating the complex system.""" # Alchemical settings - alchemical_settings: AlchemicalSettings + alchemical_settings: AbsoluteAlchemicalSettings """ Alchemical protocol settings. """ diff --git a/openfe/protocols/openmm_utils/omm_settings.py b/openfe/protocols/openmm_utils/omm_settings.py index 4e11685c6..a18f6d99f 100644 --- a/openfe/protocols/openmm_utils/omm_settings.py +++ b/openfe/protocols/openmm_utils/omm_settings.py @@ -711,3 +711,55 @@ class MDOutputSettings(OutputSettings): Filename for writing the log of the MD simulation, including timesteps, energies, density, etc. """ + + +class AbsoluteAlchemicalSettings(SettingsBaseModel): + """ + Alchemical settings for Protocols which use the + AbsoluteAlchemicalFactory. + """ + consistent_exceptions: bool = False + """ + If True, the same functional form of the System's + nonbonded method will be used to determine the electrostatic + potential energy contributions of 1,4 exceptions instead + of the classical q1*q2/(4*epsilon*epsilon0*pi*r). Default is False. + """ + disable_alchemical_dispersion_corrrection: bool = False + """ + If True, the long-range dispersion correction will not + be included for the alchemical region, avoiding the need + to recompute the correction. This can improve performance, + at the cost of accuracy. Default is False. + """ + annihilate_sterics: bool = False + """ + If True, sterics (Lennard-Jones) will be annhilated instead + of decoupled. Default is False. + """ + softcore_alpha: float = 0.5 + """ + Alchemical softcore parameter for the Lennard-Jones interactions + (default is 0.5). + + The generalized softcore potential formalism introduced by + Pham and Shirts, J. Chem. Phys. 135, 034114 (2011), equation 13, + is used here. The ``softcore_a``, ``softcore_b``, and + ``softcore_c`` parameters are used alongside ``softcore_alpha`` + to control how the potential is scaled. + """ + softcore_a: float = 1.0 + """ + Scaling constant ``a`` in + Eq. 13 from Pham and Shirts, J. Chem. Phys. 135, 034114 (2011). + """ + softcore_b: float = 1.0 + """ + Scaling constant ``b`` in + Eq. 13 from Pham and Shirts, J. Chem. Phys. 135, 034114 (2011). + """ + softcore_c: float = 6.0 + """ + Scaling constant ``c`` in + Eq. 13 from Pham and Shirts, J. Chem. Phys. 135, 034114 (2011). + """ From 48fe3359a1365630d1f996e76331f0a4125c2467 Mon Sep 17 00:00:00 2001 From: IAlibay Date: Fri, 5 Dec 2025 18:22:40 +0000 Subject: [PATCH 02/11] Expose settings in AFE protocols --- openfe/protocols/openmm_afe/base.py | 34 ++++++++++++++++--- openfe/protocols/openmm_utils/omm_settings.py | 4 +-- 2 files changed, 31 insertions(+), 7 deletions(-) diff --git a/openfe/protocols/openmm_afe/base.py b/openfe/protocols/openmm_afe/base.py index 574743b08..269d3913b 100644 --- a/openfe/protocols/openmm_afe/base.py +++ b/openfe/protocols/openmm_afe/base.py @@ -36,7 +36,7 @@ from openff.units import Quantity, unit from openff.units.openmm import ensure_quantity, from_openmm, to_openmm from openmm import app -from openmm import unit as omm_unit +from openmm import unit as ommunit from openmmforcefields.generators import SystemGenerator from openmmtools import multistate from openmmtools.alchemy import ( @@ -168,10 +168,10 @@ def _pre_equilibrate( self, system: openmm.System, topology: openmm.app.Topology, - positions: omm_unit.Quantity, + positions: ommunit.Quantity, settings: dict[str, SettingsBaseModel], dry: bool, - ) -> tuple[omm_unit.Quantity, omm_unit.Quantity]: + ) -> tuple[ommunit.Quantity, ommunit.Quantity]: """ Run a non-alchemical equilibration to get a stable system. @@ -619,6 +619,7 @@ def _get_alchemical_system( system: openmm.System, comp_resids: dict[Component, npt.NDArray], alchem_comps: dict[str, list[Component]], + alchemical_settings: AbsoluteAlchemicalSettings, ) -> tuple[AbsoluteAlchemicalFactory, openmm.System, list[int]]: """ Get an alchemically modified system and its associated factory @@ -633,6 +634,8 @@ def _get_alchemical_system( A dictionary of residues for each component in the System. alchem_comps : dict[str, list[Component]] A dictionary of alchemical components for each end state. + alchemical_settings : AbsolulteAlchemicalSettings + Settings controlling how the alchemical system is built. Returns ------- @@ -652,9 +655,26 @@ def _get_alchemical_system( alchemical_region = AlchemicalRegion( alchemical_atoms=alchemical_indices, + softcore_alpha=alchemical_settings.softcore_alpha, + annihilate_electrostatics=True, + annihilate_sterics=alchemical_settings.annihilate_sterics, + softcore_a=alchemical_settings.softcore_a, + softcore_b=alchemical_settings.softcore_b, + softcore_c=alchemical_settings.softcore_c, + softcore_beta=0.0, + softcore_d=1.0, + softcore_e=1.0, + softcore_f=2.0, ) - alchemical_factory = AbsoluteAlchemicalFactory() + alchemical_factory = AbsoluteAlchemicalFactory( + consistent_exceptions=alchemical_settings.consistent_exceptions, + switch_width=1.0 * ommunit.angstroms, + alchemical_pme_treatment='exact', + alchemical_rf_treatment='switched', + disable_alchemical_dispersion_correction=alchemical_settings.disable_alchemical_dispersion_correction, + split_alchemical_forces=True + ) alchemical_system = alchemical_factory.create_alchemical_system(system, alchemical_region) return alchemical_factory, alchemical_system, alchemical_indices @@ -1141,7 +1161,11 @@ def run( # 7. Get alchemical system alchem_factory, alchem_system, alchem_indices = self._get_alchemical_system( - omm_topology, restrained_omm_system, comp_resids, alchem_comps + topology=omm_topology, + system=restrained_omm_system, + comp_resids=comp_resids, + alchem_comps=alchem_comps, + alchemical_settings=settings["alchemical_settings"], ) # 8. Get compound and sampler states diff --git a/openfe/protocols/openmm_utils/omm_settings.py b/openfe/protocols/openmm_utils/omm_settings.py index a18f6d99f..643d0f4d9 100644 --- a/openfe/protocols/openmm_utils/omm_settings.py +++ b/openfe/protocols/openmm_utils/omm_settings.py @@ -754,12 +754,12 @@ class AbsoluteAlchemicalSettings(SettingsBaseModel): Eq. 13 from Pham and Shirts, J. Chem. Phys. 135, 034114 (2011). """ softcore_b: float = 1.0 - """ + """ Scaling constant ``b`` in Eq. 13 from Pham and Shirts, J. Chem. Phys. 135, 034114 (2011). """ softcore_c: float = 6.0 - """ + """ Scaling constant ``c`` in Eq. 13 from Pham and Shirts, J. Chem. Phys. 135, 034114 (2011). """ From 89f78147a34311789e697766e193b90a867d09ed Mon Sep 17 00:00:00 2001 From: IAlibay Date: Fri, 5 Dec 2025 18:28:25 +0000 Subject: [PATCH 03/11] expose settings in septop --- openfe/protocols/openmm_afe/base.py | 1 + openfe/protocols/openmm_septop/base.py | 39 +++++++++++++++++-- .../openmm_septop/equil_septop_method.py | 2 + 3 files changed, 38 insertions(+), 4 deletions(-) diff --git a/openfe/protocols/openmm_afe/base.py b/openfe/protocols/openmm_afe/base.py index 269d3913b..af9daca3c 100644 --- a/openfe/protocols/openmm_afe/base.py +++ b/openfe/protocols/openmm_afe/base.py @@ -52,6 +52,7 @@ ) from openfe.protocols.openmm_afe.equil_afe_settings import ( + AbsoluteAlchemicalSettings, BaseSolvationSettings, IntegratorSettings, MultiStateOutputSettings, diff --git a/openfe/protocols/openmm_septop/base.py b/openfe/protocols/openmm_septop/base.py index 323a3e6db..f1a70a4f7 100644 --- a/openfe/protocols/openmm_septop/base.py +++ b/openfe/protocols/openmm_septop/base.py @@ -45,6 +45,7 @@ ) from openfe.protocols.openmm_afe.equil_afe_settings import ( + AbsoluteAlchemicalSettings, BaseSolvationSettings, IntegratorSettings, MultiStateOutputSettings, @@ -285,6 +286,7 @@ def _get_alchemical_system( system: openmm.System, alchem_indices_A: list[int], alchem_indices_B: list[int], + alchemical_settings: AbsoluteAlchemicalSettings, ) -> tuple[AbsoluteAlchemicalFactory, openmm.System]: """ Get an alchemically modified system and its associated factory @@ -309,17 +311,46 @@ def _get_alchemical_system( """ alchemical_factory = AbsoluteAlchemicalFactory( - consistent_exceptions=False, + consistent_exceptions=alchemical_settings.consistent_exceptions, switch_width=1.0 * unit.angstroms, alchemical_pme_treatment="exact", alchemical_rf_treatment="switched", - disable_alchemical_dispersion_correction=False, + disable_alchemical_dispersion_correction=alchemical_settings.disable_alchemical_dispersion_correction, split_alchemical_forces=True, ) + # Alchemical Region for ligand A - alchemical_region_A = AlchemicalRegion(alchemical_atoms=alchem_indices_A, name="A") + alchemical_region_A = AlchemicalRegion( + alchemical_atoms=alchem_indices_A, + name="A", + softcore_alpha=alchemical_settings.softcore_alpha, + annihilate_electrostatics=True, + annihilate_sterics=alchemical_settings.annihilate_sterics, + softcore_a=alchemical_settings.softcore_a, + softcore_b=alchemical_settings.softcore_b, + softcore_c=alchemical_settings.softcore_c, + softcore_beta=0.0, + softcore_d=1.0, + softcore_e=1.0, + softcore_f=2.0, + ) + # Alchemical Region for ligand B - alchemical_region_B = AlchemicalRegion(alchemical_atoms=alchem_indices_B, name="B") + alchemical_region_B = AlchemicalRegion( + alchemical_atoms=alchem_indices_B, + name="B", + softcore_alpha=alchemical_settings.softcore_alpha, + annihilate_electrostatics=True, + annihilate_sterics=alchemical_settings.annihilate_sterics, + softcore_a=alchemical_settings.softcore_a, + softcore_b=alchemical_settings.softcore_b, + softcore_c=alchemical_settings.softcore_c, + softcore_beta=0.0, + softcore_d=1.0, + softcore_e=1.0, + softcore_f=2.0, + ) + alchemical_system = alchemical_factory.create_alchemical_system( system, [alchemical_region_A, alchemical_region_B] ) diff --git a/openfe/protocols/openmm_septop/equil_septop_method.py b/openfe/protocols/openmm_septop/equil_septop_method.py index 8391c5ae8..3d4dee6dd 100644 --- a/openfe/protocols/openmm_septop/equil_septop_method.py +++ b/openfe/protocols/openmm_septop/equil_septop_method.py @@ -1971,6 +1971,7 @@ def run( omm_system_AB, atom_indices_AB_A, atom_indices_AB_B, + settings['alchemical_settings'], ) # 10. Apply Restraints @@ -2254,6 +2255,7 @@ def run( omm_system_AB, atom_indices_AB_A, atom_indices_AB_B, + settings['alchemical_settings'], ) # 8. Apply Restraints From 832b5c10e3081e506f82ac69a90543237f49a4d1 Mon Sep 17 00:00:00 2001 From: IAlibay Date: Fri, 5 Dec 2025 19:35:11 +0000 Subject: [PATCH 04/11] Fix up settings --- openfe/protocols/openmm_afe/base.py | 8 +-- .../openmm_afe/equil_afe_settings.py | 57 ++++++++++++++++++- .../openmm_afe/equil_binding_afe_method.py | 8 +-- .../openmm_afe/equil_solvation_afe_method.py | 8 +-- openfe/protocols/openmm_septop/base.py | 12 ++-- .../openmm_septop/equil_septop_method.py | 8 +-- .../openmm_septop/equil_septop_settings.py | 4 +- openfe/protocols/openmm_utils/omm_settings.py | 54 +----------------- 8 files changed, 80 insertions(+), 79 deletions(-) diff --git a/openfe/protocols/openmm_afe/base.py b/openfe/protocols/openmm_afe/base.py index af9daca3c..b065f6076 100644 --- a/openfe/protocols/openmm_afe/base.py +++ b/openfe/protocols/openmm_afe/base.py @@ -52,7 +52,7 @@ ) from openfe.protocols.openmm_afe.equil_afe_settings import ( - AbsoluteAlchemicalSettings, + AlchemicalSettings, BaseSolvationSettings, IntegratorSettings, MultiStateOutputSettings, @@ -338,7 +338,7 @@ def _handle_settings(self) -> dict[str, SettingsBaseModel]: * forcefield_settings : OpenMMSystemGeneratorFFSettings * thermo_settings : ThermoSettings * solvation_settings : BaseSolvationSettings - * alchemical_settings : AbsoluteAlchemicalSettings + * alchemical_settings : AlchemicalSettings * lambda_settings : LambdaSettings * engine_settings : OpenMMEngineSettings * integrator_settings : IntegratorSettings @@ -620,7 +620,7 @@ def _get_alchemical_system( system: openmm.System, comp_resids: dict[Component, npt.NDArray], alchem_comps: dict[str, list[Component]], - alchemical_settings: AbsoluteAlchemicalSettings, + alchemical_settings: AlchemicalSettings, ) -> tuple[AbsoluteAlchemicalFactory, openmm.System, list[int]]: """ Get an alchemically modified system and its associated factory @@ -635,7 +635,7 @@ def _get_alchemical_system( A dictionary of residues for each component in the System. alchem_comps : dict[str, list[Component]] A dictionary of alchemical components for each end state. - alchemical_settings : AbsolulteAlchemicalSettings + alchemical_settings : AlchemicalSettings Settings controlling how the alchemical system is built. Returns diff --git a/openfe/protocols/openmm_afe/equil_afe_settings.py b/openfe/protocols/openmm_afe/equil_afe_settings.py index 89e50711b..333daa5d3 100644 --- a/openfe/protocols/openmm_afe/equil_afe_settings.py +++ b/openfe/protocols/openmm_afe/equil_afe_settings.py @@ -34,7 +34,6 @@ OpenFFPartialChargeSettings, OpenMMEngineSettings, OpenMMSolvationSettings, - AbsoluteAlchemicalSettings, ) from openfe.protocols.restraint_utils.settings import ( BaseRestraintSettings, @@ -42,6 +41,58 @@ ) +class AlchemicalSettings(SettingsBaseModel): + """ + Alchemical settings for Protocols which use the + AbsoluteAlchemicalFactory. + """ + consistent_exceptions: bool = False + """ + If True, the same functional form of the System's + nonbonded method will be used to determine the electrostatic + potential energy contributions of 1,4 exceptions instead + of the classical q1*q2/(4*epsilon*epsilon0*pi*r). Default is False. + """ + disable_alchemical_dispersion_correction: bool = False + """ + If True, the long-range dispersion correction will not + be included for the alchemical region, avoiding the need + to recompute the correction. This can improve performance, + at the cost of accuracy. Default is False. + """ + annihilate_sterics: bool = False + """ + If True, sterics (Lennard-Jones) will be annhilated instead + of decoupled. Default is False. + """ + softcore_alpha: float = 0.5 + """ + Alchemical softcore parameter for the Lennard-Jones interactions + (default is 0.5). + + The generalized softcore potential formalism introduced by + Pham and Shirts, J. Chem. Phys. 135, 034114 (2011), equation 13, + is used here. The ``softcore_a``, ``softcore_b``, and + ``softcore_c`` parameters are used alongside ``softcore_alpha`` + to control how the potential is scaled. + """ + softcore_a: float = 1.0 + """ + Scaling constant ``a`` in + Eq. 13 from Pham and Shirts, J. Chem. Phys. 135, 034114 (2011). + """ + softcore_b: float = 1.0 + """ + Scaling constant ``b`` in + Eq. 13 from Pham and Shirts, J. Chem. Phys. 135, 034114 (2011). + """ + softcore_c: float = 6.0 + """ + Scaling constant ``c`` in + Eq. 13 from Pham and Shirts, J. Chem. Phys. 135, 034114 (2011). + """ + + class LambdaSettings(SettingsBaseModel): """Lambda schedule settings. @@ -204,7 +255,7 @@ def must_be_positive(cls, v): """Settings for solvating the system.""" # Alchemical settings - alchemical_settings: AbsoluteAlchemicalSettings + alchemical_settings: AlchemicalSettings """ Alchemical protocol settings. """ @@ -315,7 +366,7 @@ def must_be_positive(cls, v): """Settings for solvating the system in the complex.""" # Alchemical settings - alchemical_settings: AbsoluteAlchemicalSettings + alchemical_settings: AlchemicalSettings """ Alchemical protocol settings. """ diff --git a/openfe/protocols/openmm_afe/equil_binding_afe_method.py b/openfe/protocols/openmm_afe/equil_binding_afe_method.py index ad704301e..dc12f1a8e 100644 --- a/openfe/protocols/openmm_afe/equil_binding_afe_method.py +++ b/openfe/protocols/openmm_afe/equil_binding_afe_method.py @@ -58,7 +58,7 @@ from openfe.protocols.openmm_afe.equil_afe_settings import ( ABFEPreEquilOutputSettings, AbsoluteBindingSettings, - AbsoluteAlchemicalSettings, + AlchemicalSettings, BoreschRestraintSettings, IntegratorSettings, LambdaSettings, @@ -553,7 +553,7 @@ def _default_settings(cls): temperature=298.15 * offunit.kelvin, pressure=1 * offunit.bar, ), - alchemical_settings=AbsoluteAlchemicalSettings(), + alchemical_settings=AlchemicalSettings(), solvent_lambda_settings=LambdaSettings( lambda_elec=[ 0.0, 0.25, 0.5, 0.75, 1.0, @@ -958,7 +958,7 @@ def _handle_settings(self) -> dict[str, SettingsBaseModel]: * thermo_settings : ThermoSettings * charge_settings : OpenFFPartialChargeSettings * solvation_settings : OpenMMSolvationSettings - * alchemical_settings : AbsoluteAlchemicalSettings + * alchemical_settings : AlchemicalSettings * lambda_settings : LambdaSettings * engine_settings : OpenMMEngineSettings * integrator_settings : IntegratorSettings @@ -1306,7 +1306,7 @@ def _handle_settings(self) -> dict[str, SettingsBaseModel]: * thermo_settings : ThermoSettings * charge_settings : OpenFFPartialChargeSettings * solvation_settings : OpenMMSolvationSettings - * alchemical_settings : AbsoluteAlchemicalSettings + * alchemical_settings : AlchemicalSettings * lambda_settings : LambdaSettings * engine_settings : OpenMMEngineSettings * integrator_settings : IntegratorSettings diff --git a/openfe/protocols/openmm_afe/equil_solvation_afe_method.py b/openfe/protocols/openmm_afe/equil_solvation_afe_method.py index ae4de822f..4ffb94770 100644 --- a/openfe/protocols/openmm_afe/equil_solvation_afe_method.py +++ b/openfe/protocols/openmm_afe/equil_solvation_afe_method.py @@ -55,7 +55,7 @@ from openfe.due import Doi, due from openfe.protocols.openmm_afe.equil_afe_settings import ( AbsoluteSolvationSettings, - AbsoluteAlchemicalSettings, + AlchemicalSettings, IntegratorSettings, LambdaSettings, MDOutputSettings, @@ -442,7 +442,7 @@ def _default_settings(cls): temperature=298.15 * unit.kelvin, pressure=1 * unit.bar, ), - alchemical_settings=AbsoluteAlchemicalSettings(), + alchemical_settings=AlchemicalSettings(), lambda_settings=LambdaSettings( lambda_elec=[ 0.0, 0.25, 0.5, 0.75, 1.0, 1.0, 1.0, @@ -859,7 +859,7 @@ def _handle_settings(self) -> dict[str, SettingsBaseModel]: * thermo_settings : ThermoSettings * charge_settings : OpenFFPartialChargeSettings * solvation_settings : OpenMMSolvationSettings - * alchemical_settings : AbsoluteAlchemicalSettings + * alchemical_settings : AlchemicalSettings * lambda_settings : LambdaSettings * engine_settings : OpenMMEngineSettings * integrator_settings : IntegratorSettings @@ -934,7 +934,7 @@ def _handle_settings(self) -> dict[str, SettingsBaseModel]: * thermo_settings : ThermoSettings * charge_settings : OpenFFPartialChargeSettings * solvation_settings : OpenMMSolvationSettings - * alchemical_settings : AbsoluteAlchemicalSettings + * alchemical_settings : AlchemicalSettings * lambda_settings : LambdaSettings * engine_settings : OpenMMEngineSettings * integrator_settings : IntegratorSettings diff --git a/openfe/protocols/openmm_septop/base.py b/openfe/protocols/openmm_septop/base.py index f1a70a4f7..2ff4f1af8 100644 --- a/openfe/protocols/openmm_septop/base.py +++ b/openfe/protocols/openmm_septop/base.py @@ -45,7 +45,7 @@ ) from openfe.protocols.openmm_afe.equil_afe_settings import ( - AbsoluteAlchemicalSettings, + AlchemicalSettings, BaseSolvationSettings, IntegratorSettings, MultiStateOutputSettings, @@ -286,7 +286,7 @@ def _get_alchemical_system( system: openmm.System, alchem_indices_A: list[int], alchem_indices_B: list[int], - alchemical_settings: AbsoluteAlchemicalSettings, + alchemical_settings: AlchemicalSettings, ) -> tuple[AbsoluteAlchemicalFactory, openmm.System]: """ Get an alchemically modified system and its associated factory @@ -295,12 +295,14 @@ def _get_alchemical_system( ---------- system : openmm.System System to alchemically modify. - alchem_indices_A: list[int] + alchem_indices_A : list[int] A list of atom indices for the alchemically modified ligand A in the system. - alchem_indices_B: list[int] + alchem_indices_B : list[int] A list of atom indices for the alchemically modified ligand B in the system. + alchemical_settings : AlchemicalSettings + Settings controlling how the alchemical system will be built. Returns ------- @@ -797,7 +799,7 @@ def _handle_settings(self): * forcefield_settings : OpenMMSystemGeneratorFFSettings * thermo_settings : ThermoSettings * solvation_settings : BaseSolvationSettings - * alchemical_settings : AbsoluteAlchemicalSettings + * alchemical_settings : AlchemicalSettings * lambda_settings : LambdaSettings * engine_settings : OpenMMEngineSettings * integrator_settings : IntegratorSettings diff --git a/openfe/protocols/openmm_septop/equil_septop_method.py b/openfe/protocols/openmm_septop/equil_septop_method.py index 3d4dee6dd..bc224d914 100644 --- a/openfe/protocols/openmm_septop/equil_septop_method.py +++ b/openfe/protocols/openmm_septop/equil_septop_method.py @@ -68,7 +68,7 @@ from openfe.due import Doi, due from openfe.protocols.openmm_septop.equil_septop_settings import ( - AbsoluteAlchemicalSettings, + AlchemicalSettings, IntegratorSettings, LambdaSettings, MDSimulationSettings, @@ -241,7 +241,7 @@ def _handle_settings(self) -> dict[str, SettingsBaseModel]: * thermo_settings : ThermoSettings * charge_settings : OpenFFPartialChargeSettings * solvation_settings : OpenMMSolvationSettings - * alchemical_settings : AbsoluteAlchemicalSettings + * alchemical_settings : AlchemicalSettings * lambda_settings : LambdaSettings * engine_settings : OpenMMEngineSettings * integrator_settings : IntegratorSettings @@ -327,7 +327,7 @@ def _handle_settings(self) -> dict[str, SettingsBaseModel]: * thermo_settings : ThermoSettings * charge_settings : OpenFFPartialChargeSettings * solvation_settings : OpenMMSolvationSettings - * alchemical_settings : AbsoluteAlchemicalSettings + * alchemical_settings : AlchemicalSettings * lambda_settings : LambdaSettings * engine_settings : OpenMMEngineSettings * integrator_settings : IntegratorSettings @@ -894,7 +894,7 @@ def _default_settings(cls): temperature=298.15 * unit.kelvin, pressure=1 * unit.bar, ), - alchemical_settings=AbsoluteAlchemicalSettings(), + alchemical_settings=AlchemicalSettings(), solvent_lambda_settings=LambdaSettings( lambda_elec_A=[ 0.0, diff --git a/openfe/protocols/openmm_septop/equil_septop_settings.py b/openfe/protocols/openmm_septop/equil_septop_settings.py index 59af68b40..821717a76 100644 --- a/openfe/protocols/openmm_septop/equil_septop_settings.py +++ b/openfe/protocols/openmm_septop/equil_septop_settings.py @@ -32,8 +32,8 @@ OpenFFPartialChargeSettings, OpenMMEngineSettings, OpenMMSolvationSettings, - AbsoluteAlchemicalSettings, ) +from openfe.protocols.openmm_afe.equil_afe_settings import AlchemicalSettings from openfe.protocols.restraint_utils.settings import BaseRestraintSettings @@ -363,7 +363,7 @@ def must_be_positive(cls, v): """Settings for solvating the complex system.""" # Alchemical settings - alchemical_settings: AbsoluteAlchemicalSettings + alchemical_settings: AlchemicalSettings """ Alchemical protocol settings. """ diff --git a/openfe/protocols/openmm_utils/omm_settings.py b/openfe/protocols/openmm_utils/omm_settings.py index 643d0f4d9..2d7bd91dc 100644 --- a/openfe/protocols/openmm_utils/omm_settings.py +++ b/openfe/protocols/openmm_utils/omm_settings.py @@ -710,56 +710,4 @@ class MDOutputSettings(OutputSettings): """ Filename for writing the log of the MD simulation, including timesteps, energies, density, etc. - """ - - -class AbsoluteAlchemicalSettings(SettingsBaseModel): - """ - Alchemical settings for Protocols which use the - AbsoluteAlchemicalFactory. - """ - consistent_exceptions: bool = False - """ - If True, the same functional form of the System's - nonbonded method will be used to determine the electrostatic - potential energy contributions of 1,4 exceptions instead - of the classical q1*q2/(4*epsilon*epsilon0*pi*r). Default is False. - """ - disable_alchemical_dispersion_corrrection: bool = False - """ - If True, the long-range dispersion correction will not - be included for the alchemical region, avoiding the need - to recompute the correction. This can improve performance, - at the cost of accuracy. Default is False. - """ - annihilate_sterics: bool = False - """ - If True, sterics (Lennard-Jones) will be annhilated instead - of decoupled. Default is False. - """ - softcore_alpha: float = 0.5 - """ - Alchemical softcore parameter for the Lennard-Jones interactions - (default is 0.5). - - The generalized softcore potential formalism introduced by - Pham and Shirts, J. Chem. Phys. 135, 034114 (2011), equation 13, - is used here. The ``softcore_a``, ``softcore_b``, and - ``softcore_c`` parameters are used alongside ``softcore_alpha`` - to control how the potential is scaled. - """ - softcore_a: float = 1.0 - """ - Scaling constant ``a`` in - Eq. 13 from Pham and Shirts, J. Chem. Phys. 135, 034114 (2011). - """ - softcore_b: float = 1.0 - """ - Scaling constant ``b`` in - Eq. 13 from Pham and Shirts, J. Chem. Phys. 135, 034114 (2011). - """ - softcore_c: float = 6.0 - """ - Scaling constant ``c`` in - Eq. 13 from Pham and Shirts, J. Chem. Phys. 135, 034114 (2011). - """ + """ \ No newline at end of file From 910522904eddef358c1739e60faedb1623bd1946 Mon Sep 17 00:00:00 2001 From: "pre-commit-ci[bot]" <66853113+pre-commit-ci[bot]@users.noreply.github.com> Date: Fri, 5 Dec 2025 20:00:13 +0000 Subject: [PATCH 05/11] [pre-commit.ci] auto fixes from pre-commit.com hooks for more information, see https://pre-commit.ci --- openfe/protocols/openmm_afe/base.py | 6 +++--- openfe/protocols/openmm_afe/equil_afe_settings.py | 1 + openfe/protocols/openmm_septop/equil_septop_method.py | 4 ++-- openfe/protocols/openmm_septop/equil_septop_settings.py | 2 +- openfe/protocols/openmm_utils/omm_settings.py | 2 +- 5 files changed, 8 insertions(+), 7 deletions(-) diff --git a/openfe/protocols/openmm_afe/base.py b/openfe/protocols/openmm_afe/base.py index b065f6076..d2ecdcb9c 100644 --- a/openfe/protocols/openmm_afe/base.py +++ b/openfe/protocols/openmm_afe/base.py @@ -671,10 +671,10 @@ def _get_alchemical_system( alchemical_factory = AbsoluteAlchemicalFactory( consistent_exceptions=alchemical_settings.consistent_exceptions, switch_width=1.0 * ommunit.angstroms, - alchemical_pme_treatment='exact', - alchemical_rf_treatment='switched', + alchemical_pme_treatment="exact", + alchemical_rf_treatment="switched", disable_alchemical_dispersion_correction=alchemical_settings.disable_alchemical_dispersion_correction, - split_alchemical_forces=True + split_alchemical_forces=True, ) alchemical_system = alchemical_factory.create_alchemical_system(system, alchemical_region) diff --git a/openfe/protocols/openmm_afe/equil_afe_settings.py b/openfe/protocols/openmm_afe/equil_afe_settings.py index 333daa5d3..390df00f8 100644 --- a/openfe/protocols/openmm_afe/equil_afe_settings.py +++ b/openfe/protocols/openmm_afe/equil_afe_settings.py @@ -46,6 +46,7 @@ class AlchemicalSettings(SettingsBaseModel): Alchemical settings for Protocols which use the AbsoluteAlchemicalFactory. """ + consistent_exceptions: bool = False """ If True, the same functional form of the System's diff --git a/openfe/protocols/openmm_septop/equil_septop_method.py b/openfe/protocols/openmm_septop/equil_septop_method.py index bc224d914..87f7aa7a0 100644 --- a/openfe/protocols/openmm_septop/equil_septop_method.py +++ b/openfe/protocols/openmm_septop/equil_septop_method.py @@ -1971,7 +1971,7 @@ def run( omm_system_AB, atom_indices_AB_A, atom_indices_AB_B, - settings['alchemical_settings'], + settings["alchemical_settings"], ) # 10. Apply Restraints @@ -2255,7 +2255,7 @@ def run( omm_system_AB, atom_indices_AB_A, atom_indices_AB_B, - settings['alchemical_settings'], + settings["alchemical_settings"], ) # 8. Apply Restraints diff --git a/openfe/protocols/openmm_septop/equil_septop_settings.py b/openfe/protocols/openmm_septop/equil_septop_settings.py index 821717a76..e1898c5fa 100644 --- a/openfe/protocols/openmm_septop/equil_septop_settings.py +++ b/openfe/protocols/openmm_septop/equil_septop_settings.py @@ -23,6 +23,7 @@ from openff.units import unit as offunit from pydantic import field_validator +from openfe.protocols.openmm_afe.equil_afe_settings import AlchemicalSettings from openfe.protocols.openmm_utils.omm_settings import ( IntegratorSettings, MDOutputSettings, @@ -33,7 +34,6 @@ OpenMMEngineSettings, OpenMMSolvationSettings, ) -from openfe.protocols.openmm_afe.equil_afe_settings import AlchemicalSettings from openfe.protocols.restraint_utils.settings import BaseRestraintSettings diff --git a/openfe/protocols/openmm_utils/omm_settings.py b/openfe/protocols/openmm_utils/omm_settings.py index 2d7bd91dc..4e11685c6 100644 --- a/openfe/protocols/openmm_utils/omm_settings.py +++ b/openfe/protocols/openmm_utils/omm_settings.py @@ -710,4 +710,4 @@ class MDOutputSettings(OutputSettings): """ Filename for writing the log of the MD simulation, including timesteps, energies, density, etc. - """ \ No newline at end of file + """ From fcf44867a5955b8bf37fc33a546b3a36f4757434 Mon Sep 17 00:00:00 2001 From: IAlibay Date: Sat, 6 Dec 2025 01:06:01 +0000 Subject: [PATCH 06/11] add some tests --- openfe/protocols/openmm_afe/base.py | 2 +- .../openmm_afe/equil_afe_settings.py | 8 - openfe/protocols/openmm_septop/base.py | 2 +- .../openmm_ahfe/test_ahfe_protocol.py | 150 +++++++++++++++++- 4 files changed, 150 insertions(+), 12 deletions(-) diff --git a/openfe/protocols/openmm_afe/base.py b/openfe/protocols/openmm_afe/base.py index d2ecdcb9c..d28a630ca 100644 --- a/openfe/protocols/openmm_afe/base.py +++ b/openfe/protocols/openmm_afe/base.py @@ -669,7 +669,7 @@ def _get_alchemical_system( ) alchemical_factory = AbsoluteAlchemicalFactory( - consistent_exceptions=alchemical_settings.consistent_exceptions, + consistent_exceptions=False, switch_width=1.0 * ommunit.angstroms, alchemical_pme_treatment="exact", alchemical_rf_treatment="switched", diff --git a/openfe/protocols/openmm_afe/equil_afe_settings.py b/openfe/protocols/openmm_afe/equil_afe_settings.py index 390df00f8..778efe145 100644 --- a/openfe/protocols/openmm_afe/equil_afe_settings.py +++ b/openfe/protocols/openmm_afe/equil_afe_settings.py @@ -46,14 +46,6 @@ class AlchemicalSettings(SettingsBaseModel): Alchemical settings for Protocols which use the AbsoluteAlchemicalFactory. """ - - consistent_exceptions: bool = False - """ - If True, the same functional form of the System's - nonbonded method will be used to determine the electrostatic - potential energy contributions of 1,4 exceptions instead - of the classical q1*q2/(4*epsilon*epsilon0*pi*r). Default is False. - """ disable_alchemical_dispersion_correction: bool = False """ If True, the long-range dispersion correction will not diff --git a/openfe/protocols/openmm_septop/base.py b/openfe/protocols/openmm_septop/base.py index 2ff4f1af8..000f0b14e 100644 --- a/openfe/protocols/openmm_septop/base.py +++ b/openfe/protocols/openmm_septop/base.py @@ -313,7 +313,7 @@ def _get_alchemical_system( """ alchemical_factory = AbsoluteAlchemicalFactory( - consistent_exceptions=alchemical_settings.consistent_exceptions, + consistent_exceptions=False, switch_width=1.0 * unit.angstroms, alchemical_pme_treatment="exact", alchemical_rf_treatment="switched", diff --git a/openfe/tests/protocols/openmm_ahfe/test_ahfe_protocol.py b/openfe/tests/protocols/openmm_ahfe/test_ahfe_protocol.py index 33cd36223..44562c914 100644 --- a/openfe/tests/protocols/openmm_ahfe/test_ahfe_protocol.py +++ b/openfe/tests/protocols/openmm_ahfe/test_ahfe_protocol.py @@ -13,7 +13,15 @@ from numpy.testing import assert_allclose from openff.units import unit as offunit from openff.units.openmm import ensure_quantity, from_openmm -from openmm import CustomNonbondedForce, NonbondedForce +from openmm import ( + CustomBondForce, + CustomNonbondedForce, + HarmonicAngleForce, + HarmonicBondForce, + MonteCarloBarostat, + NonbondedForce, + PeriodicTorsionForce, +) from openmmtools.multistate.multistatesampler import MultiStateSampler import openfe @@ -86,6 +94,54 @@ def test_create_independent_repeat_ids(benzene_system): assert len(repeat_ids) == 12 +def _assert_num_forces(system, forcetype, number): + """ + Helper method to check the number of forces of a given + type in a system. + """ + forces = [f for f in system.getForces() if isinstance(f, forcetype)] + assert len(forces) == number + + +def _verify_alchemical_sterics_force_parameters( + force, + long_range=True, + alpha=0.5, + beta=0, + a=1.0, + b=1.0, + c=6.0, + d=1.0, + e=1.0, + f=2.0, +): + assert force.getUseLongRangeCorrection() is long_range + + if force.getNumGlobalParameters() == 8: + shift = 0 + else: + shift = 1 + + # Check the softcore parameters for the sterics forces + assert force.getGlobalParameterName(0 + shift) == 'softcore_alpha' + assert force.getGlobalParameterName(1 + shift) == 'softcore_beta' + assert force.getGlobalParameterName(2 + shift) == 'softcore_a' + assert force.getGlobalParameterName(3 + shift) == 'softcore_b' + assert force.getGlobalParameterName(4 + shift) == 'softcore_c' + assert force.getGlobalParameterName(5 + shift) == 'softcore_d' + assert force.getGlobalParameterName(6 + shift) == 'softcore_e' + assert force.getGlobalParameterName(7 + shift) == 'softcore_f' + + assert force.getGlobalParameterDefaultValue(0 + shift) == alpha + assert force.getGlobalParameterDefaultValue(1 + shift) == beta + assert force.getGlobalParameterDefaultValue(2 + shift) == a + assert force.getGlobalParameterDefaultValue(3 + shift) == b + assert force.getGlobalParameterDefaultValue(4 + shift) == c + assert force.getGlobalParameterDefaultValue(5 + shift) == d + assert force.getGlobalParameterDefaultValue(6 + shift) == e + assert force.getGlobalParameterDefaultValue(7 + shift) == f + + @pytest.mark.parametrize("method", ["repex", "sams", "independent", "InDePeNdENT"]) def test_dry_run_vac_benzene(benzene_system, method, protocol_dry_settings, tmpdir): protocol_dry_settings.vacuum_simulation_settings.sampler_method = method @@ -114,9 +170,99 @@ def test_dry_run_vac_benzene(benzene_system, method, protocol_dry_settings, tmpd assert len(sol_unit) == 1 with tmpdir.as_cwd(): - vac_sampler = vac_unit[0].run(dry=True)["debug"]["sampler"] + debug = vac_unit[0].run(dry=True)["debug"] + vac_sampler = debug["sampler"] assert not vac_sampler.is_periodic + # standard system + system = debug["system"] + assert system.getNumParticles() == 12 + assert len(system.getForces()) == 4 + _assert_num_forces(system, NonbondedForce, 1) + _assert_num_forces(system, HarmonicBondForce, 1) + _assert_num_forces(system, HarmonicAngleForce, 1) + _assert_num_forces(system, PeriodicTorsionForce, 1) + + # alchemical system + alchem_system = debug["alchem_system"] + assert alchem_system.getNumParticles() == 12 + assert len(alchem_system.getForces()) == 12 + _assert_num_forces(alchem_system, NonbondedForce, 1) + _assert_num_forces(alchem_system, CustomNonbondedForce, 4) + _assert_num_forces(alchem_system, CustomBondForce, 4) + _assert_num_forces(alchem_system, HarmonicBondForce, 1) + _assert_num_forces(alchem_system, HarmonicAngleForce, 1) + _assert_num_forces(alchem_system, PeriodicTorsionForce, 1) + + # Check some force contents + stericsf = [ + f for f in alchem_system.getForces() + if isinstance(f, CustomNonbondedForce) and 'U_sterics' in f.getEnergyFunction() + ] + + for force in stericsf: + _verify_alchemical_sterics_force_parameters(force) + + +def test_alchemical_settings_dry_run_vacuum(benzene_system, protocol_dry_settings, tmpdir): + """ + Test non default alchemical settings + """ + protocol_dry_settings.alchemical_settings.softcore_alpha = 0.2 + protocol_dry_settings.alchemical_settings.softcore_a = 2 + protocol_dry_settings.alchemical_settings.softcore_b = 2 + protocol_dry_settings.alchemical_settings.softcore_c = 1 + protocol_dry_settings.alchemical_settings.disable_alchemical_dispersion_correction = True + protocol = openmm_afe.AbsoluteSolvationProtocol(settings=protocol_dry_settings) + + stateA = benzene_system + + stateB = ChemicalSystem({"solvent": SolventComponent()}) + + # Create DAG from protocol, get the vacuum and solvent units + # and eventually dry run the first vacuum unit + dag = protocol.create( + stateA=stateA, + stateB=stateB, + mapping=None, + ) + prot_units = list(dag.protocol_units) + + assert len(prot_units) == 2 + + vac_unit = [u for u in prot_units if isinstance(u, AbsoluteSolvationVacuumUnit)] + sol_unit = [u for u in prot_units if isinstance(u, AbsoluteSolvationSolventUnit)] + + assert len(vac_unit) == 1 + assert len(sol_unit) == 1 + + with tmpdir.as_cwd(): + debug = vac_unit[0].run(dry=True)["debug"] + + alchem_system = debug["alchem_system"] + _assert_num_forces(alchem_system, NonbondedForce, 1) + _assert_num_forces(alchem_system, CustomNonbondedForce, 4) + _assert_num_forces(alchem_system, CustomBondForce, 4) + _assert_num_forces(alchem_system, HarmonicBondForce, 1) + _assert_num_forces(alchem_system, HarmonicAngleForce, 1) + _assert_num_forces(alchem_system, PeriodicTorsionForce, 1) + + # Check some force contents + stericsf = [ + f for f in alchem_system.getForces() + if isinstance(f, CustomNonbondedForce) and 'U_sterics' in f.getEnergyFunction() + ] + + for force in stericsf: + _verify_alchemical_sterics_force_parameters( + force, + long_range=False, + alpha=0.2, + a=2.0, + b=2.0, + c=1.0, + ) + def test_confgen_fail_AFE(benzene_system, protocol_dry_settings, tmpdir): # check system parametrisation works even if confgen fails From 92fce7a2019d43f31c07482a94766bdc5e8f4894 Mon Sep 17 00:00:00 2001 From: IAlibay Date: Sat, 6 Dec 2025 02:22:41 +0000 Subject: [PATCH 07/11] Add septop tests for sterics force settings --- .../openmm_septop/test_septop_protocol.py | 49 ++++++++++++++++++- 1 file changed, 48 insertions(+), 1 deletion(-) diff --git a/openfe/tests/protocols/openmm_septop/test_septop_protocol.py b/openfe/tests/protocols/openmm_septop/test_septop_protocol.py index 25b4f2e85..71e452205 100644 --- a/openfe/tests/protocols/openmm_septop/test_septop_protocol.py +++ b/openfe/tests/protocols/openmm_septop/test_septop_protocol.py @@ -13,11 +13,21 @@ import openmm import openmm.app import openmm.unit +from openmm import ( + CustomBondForce, + CustomNonbondedForce, + HarmonicAngleForce, + HarmonicBondForce, + MonteCarloBarostat, + NonbondedForce, + PeriodicTorsionForce, + CustomCompoundBondForce, + MonteCarloBarostat +) import pytest from numpy.testing import assert_allclose from openff.units import unit as offunit from openff.units.openmm import ensure_quantity, from_openmm -from openmm import CustomNonbondedForce, MonteCarloBarostat, NonbondedForce from openmmtools.alchemy import AbsoluteAlchemicalFactory, AlchemicalRegion from openmmtools.multistate.multistatesampler import MultiStateSampler @@ -38,6 +48,10 @@ from openfe.protocols.openmm_utils import system_validation from openfe.protocols.restraint_utils.geometry.boresch import BoreschRestraintGeometry from openfe.tests.protocols.conftest import compute_energy +from openfe.tests.protocols.openmm_ahfe.test_ahfe_protocol import ( + _assert_num_forces, + _verify_alchemical_sterics_force_parameters, +) E_CHARGE = 1.602176634e-19 * openmm.unit.coulomb EPSILON0 = ( @@ -717,6 +731,22 @@ def test_dry_run_benzene_toluene(benzene_toluene_dag, tmpdir): pdb = md.load_pdb("alchemical_system.pdb") assert pdb.n_atoms == 31 + # Test the solvent system + alchem_system = deserialize(serialized_system) + assert len(alchem_system.getForces()) == 14 + _assert_num_forces(alchem_system, NonbondedForce, 1) + _assert_num_forces(alchem_system, CustomNonbondedForce, 4) + _assert_num_forces(alchem_system, CustomBondForce, 4) + _assert_num_forces(alchem_system, HarmonicBondForce, 2) + _assert_num_forces(alchem_system, HarmonicAngleForce, 1) + _assert_num_forces(alchem_system, PeriodicTorsionForce, 1) + _assert_num_forces(alchem_system, MonteCarloBarostat, 1) + + # Check steric forces + for f in alchem_system.getForces(): + if isinstance(f, CustomNonbondedForce) and 'U_sterics' in f.getEnergyFunction(): + _verify_alchemical_sterics_force_parameters(f) + complex_setup_output = complex_setup_unit[0].run(dry=True) serialized_topology = complex_setup_output["topology"] serialized_system = complex_setup_output["system"] @@ -732,6 +762,23 @@ def test_dry_run_benzene_toluene(benzene_toluene_dag, tmpdir): pdb = md.load_pdb("alchemical_system.pdb") assert pdb.n_atoms == 2687 + # Test the complex system + alchem_system = deserialize(serialized_system) + assert len(alchem_system.getForces()) == 15 + _assert_num_forces(alchem_system, NonbondedForce, 1) + _assert_num_forces(alchem_system, CustomNonbondedForce, 4) + _assert_num_forces(alchem_system, CustomBondForce, 4) + _assert_num_forces(alchem_system, HarmonicBondForce, 1) + _assert_num_forces(alchem_system, HarmonicAngleForce, 1) + _assert_num_forces(alchem_system, PeriodicTorsionForce, 1) + _assert_num_forces(alchem_system, CustomCompoundBondForce, 2) + _assert_num_forces(alchem_system, MonteCarloBarostat, 1) + + # Check steric forces + for f in alchem_system.getForces(): + if isinstance(f, CustomNonbondedForce) and 'U_sterics' in f.getEnergyFunction(): + _verify_alchemical_sterics_force_parameters(f) + @pytest.mark.parametrize("method", ["repex", "sams", "independent"]) def test_dry_run_methods( From 9ac02285a1710e0cca9c9a07bed5741067696771 Mon Sep 17 00:00:00 2001 From: "pre-commit-ci[bot]" <66853113+pre-commit-ci[bot]@users.noreply.github.com> Date: Sat, 6 Dec 2025 02:25:27 +0000 Subject: [PATCH 08/11] [pre-commit.ci] auto fixes from pre-commit.com hooks for more information, see https://pre-commit.ci --- .../openmm_afe/equil_afe_settings.py | 1 + .../openmm_ahfe/test_ahfe_protocol.py | 26 ++++++++++--------- .../openmm_septop/test_septop_protocol.py | 15 +++++------ 3 files changed, 22 insertions(+), 20 deletions(-) diff --git a/openfe/protocols/openmm_afe/equil_afe_settings.py b/openfe/protocols/openmm_afe/equil_afe_settings.py index 778efe145..afb74528c 100644 --- a/openfe/protocols/openmm_afe/equil_afe_settings.py +++ b/openfe/protocols/openmm_afe/equil_afe_settings.py @@ -46,6 +46,7 @@ class AlchemicalSettings(SettingsBaseModel): Alchemical settings for Protocols which use the AbsoluteAlchemicalFactory. """ + disable_alchemical_dispersion_correction: bool = False """ If True, the long-range dispersion correction will not diff --git a/openfe/tests/protocols/openmm_ahfe/test_ahfe_protocol.py b/openfe/tests/protocols/openmm_ahfe/test_ahfe_protocol.py index 44562c914..4f24435ab 100644 --- a/openfe/tests/protocols/openmm_ahfe/test_ahfe_protocol.py +++ b/openfe/tests/protocols/openmm_ahfe/test_ahfe_protocol.py @@ -123,14 +123,14 @@ def _verify_alchemical_sterics_force_parameters( shift = 1 # Check the softcore parameters for the sterics forces - assert force.getGlobalParameterName(0 + shift) == 'softcore_alpha' - assert force.getGlobalParameterName(1 + shift) == 'softcore_beta' - assert force.getGlobalParameterName(2 + shift) == 'softcore_a' - assert force.getGlobalParameterName(3 + shift) == 'softcore_b' - assert force.getGlobalParameterName(4 + shift) == 'softcore_c' - assert force.getGlobalParameterName(5 + shift) == 'softcore_d' - assert force.getGlobalParameterName(6 + shift) == 'softcore_e' - assert force.getGlobalParameterName(7 + shift) == 'softcore_f' + assert force.getGlobalParameterName(0 + shift) == "softcore_alpha" + assert force.getGlobalParameterName(1 + shift) == "softcore_beta" + assert force.getGlobalParameterName(2 + shift) == "softcore_a" + assert force.getGlobalParameterName(3 + shift) == "softcore_b" + assert force.getGlobalParameterName(4 + shift) == "softcore_c" + assert force.getGlobalParameterName(5 + shift) == "softcore_d" + assert force.getGlobalParameterName(6 + shift) == "softcore_e" + assert force.getGlobalParameterName(7 + shift) == "softcore_f" assert force.getGlobalParameterDefaultValue(0 + shift) == alpha assert force.getGlobalParameterDefaultValue(1 + shift) == beta @@ -196,8 +196,9 @@ def test_dry_run_vac_benzene(benzene_system, method, protocol_dry_settings, tmpd # Check some force contents stericsf = [ - f for f in alchem_system.getForces() - if isinstance(f, CustomNonbondedForce) and 'U_sterics' in f.getEnergyFunction() + f + for f in alchem_system.getForces() + if isinstance(f, CustomNonbondedForce) and "U_sterics" in f.getEnergyFunction() ] for force in stericsf: @@ -249,8 +250,9 @@ def test_alchemical_settings_dry_run_vacuum(benzene_system, protocol_dry_setting # Check some force contents stericsf = [ - f for f in alchem_system.getForces() - if isinstance(f, CustomNonbondedForce) and 'U_sterics' in f.getEnergyFunction() + f + for f in alchem_system.getForces() + if isinstance(f, CustomNonbondedForce) and "U_sterics" in f.getEnergyFunction() ] for force in stericsf: diff --git a/openfe/tests/protocols/openmm_septop/test_septop_protocol.py b/openfe/tests/protocols/openmm_septop/test_septop_protocol.py index 71e452205..83090368f 100644 --- a/openfe/tests/protocols/openmm_septop/test_septop_protocol.py +++ b/openfe/tests/protocols/openmm_septop/test_septop_protocol.py @@ -13,21 +13,20 @@ import openmm import openmm.app import openmm.unit +import pytest +from numpy.testing import assert_allclose +from openff.units import unit as offunit +from openff.units.openmm import ensure_quantity, from_openmm from openmm import ( CustomBondForce, + CustomCompoundBondForce, CustomNonbondedForce, HarmonicAngleForce, HarmonicBondForce, MonteCarloBarostat, NonbondedForce, PeriodicTorsionForce, - CustomCompoundBondForce, - MonteCarloBarostat ) -import pytest -from numpy.testing import assert_allclose -from openff.units import unit as offunit -from openff.units.openmm import ensure_quantity, from_openmm from openmmtools.alchemy import AbsoluteAlchemicalFactory, AlchemicalRegion from openmmtools.multistate.multistatesampler import MultiStateSampler @@ -744,7 +743,7 @@ def test_dry_run_benzene_toluene(benzene_toluene_dag, tmpdir): # Check steric forces for f in alchem_system.getForces(): - if isinstance(f, CustomNonbondedForce) and 'U_sterics' in f.getEnergyFunction(): + if isinstance(f, CustomNonbondedForce) and "U_sterics" in f.getEnergyFunction(): _verify_alchemical_sterics_force_parameters(f) complex_setup_output = complex_setup_unit[0].run(dry=True) @@ -776,7 +775,7 @@ def test_dry_run_benzene_toluene(benzene_toluene_dag, tmpdir): # Check steric forces for f in alchem_system.getForces(): - if isinstance(f, CustomNonbondedForce) and 'U_sterics' in f.getEnergyFunction(): + if isinstance(f, CustomNonbondedForce) and "U_sterics" in f.getEnergyFunction(): _verify_alchemical_sterics_force_parameters(f) From 60e5323d43086712549b0c069289f82b0afb5972 Mon Sep 17 00:00:00 2001 From: IAlibay Date: Sat, 6 Dec 2025 02:28:39 +0000 Subject: [PATCH 09/11] Add news item --- news/absolute_settings.rst | 27 +++++++++++++++++++++++++++ 1 file changed, 27 insertions(+) create mode 100644 news/absolute_settings.rst diff --git a/news/absolute_settings.rst b/news/absolute_settings.rst new file mode 100644 index 000000000..dc4d715c1 --- /dev/null +++ b/news/absolute_settings.rst @@ -0,0 +1,27 @@ +**Added:** + +* New options have been added to the ``AlchemicalSettings`` + of the ``SepTopProtocol``, ``AbsoluteSolvationProtocol`` + and ``AbsoluteBindingProtocol``. Notably, these options allow users to + control the softcore parameters as well as the use of + long range dispersion corrections. + +**Changed:** + +* + +**Deprecated:** + +* + +**Removed:** + +* + +**Fixed:** + +* + +**Security:** + +* From ee1733d141fbf510c9ac7bdfd5ff509b9818cf7c Mon Sep 17 00:00:00 2001 From: IAlibay Date: Mon, 8 Dec 2025 13:38:56 +0000 Subject: [PATCH 10/11] add an extra parameterize option --- .../openmm_ahfe/test_ahfe_protocol.py | 29 ++++++++++++------- 1 file changed, 18 insertions(+), 11 deletions(-) diff --git a/openfe/tests/protocols/openmm_ahfe/test_ahfe_protocol.py b/openfe/tests/protocols/openmm_ahfe/test_ahfe_protocol.py index 4f24435ab..64cf623bd 100644 --- a/openfe/tests/protocols/openmm_ahfe/test_ahfe_protocol.py +++ b/openfe/tests/protocols/openmm_ahfe/test_ahfe_protocol.py @@ -205,15 +205,22 @@ def test_dry_run_vac_benzene(benzene_system, method, protocol_dry_settings, tmpd _verify_alchemical_sterics_force_parameters(force) -def test_alchemical_settings_dry_run_vacuum(benzene_system, protocol_dry_settings, tmpdir): +@pytest.mark.parametrize('alpha, a, b, c, correction', [ + [0.2, 2, 2, 1, True], + [0.35, 2.2, 1.5, 0, False], +]) +def test_alchemical_settings_dry_run_vacuum( + alpha, a, b, c, correction, + benzene_system, protocol_dry_settings, tmpdir +): """ Test non default alchemical settings """ - protocol_dry_settings.alchemical_settings.softcore_alpha = 0.2 - protocol_dry_settings.alchemical_settings.softcore_a = 2 - protocol_dry_settings.alchemical_settings.softcore_b = 2 - protocol_dry_settings.alchemical_settings.softcore_c = 1 - protocol_dry_settings.alchemical_settings.disable_alchemical_dispersion_correction = True + protocol_dry_settings.alchemical_settings.softcore_alpha = alpha + protocol_dry_settings.alchemical_settings.softcore_a = a + protocol_dry_settings.alchemical_settings.softcore_b = b + protocol_dry_settings.alchemical_settings.softcore_c = c + protocol_dry_settings.alchemical_settings.disable_alchemical_dispersion_correction = correction protocol = openmm_afe.AbsoluteSolvationProtocol(settings=protocol_dry_settings) stateA = benzene_system @@ -258,11 +265,11 @@ def test_alchemical_settings_dry_run_vacuum(benzene_system, protocol_dry_setting for force in stericsf: _verify_alchemical_sterics_force_parameters( force, - long_range=False, - alpha=0.2, - a=2.0, - b=2.0, - c=1.0, + long_range=not correction, + alpha=alpha, + a=a, + b=b, + c=c, ) From 098c4801d977b2f370d028da921b9c665bc6a413 Mon Sep 17 00:00:00 2001 From: "pre-commit-ci[bot]" <66853113+pre-commit-ci[bot]@users.noreply.github.com> Date: Mon, 8 Dec 2025 13:39:29 +0000 Subject: [PATCH 11/11] [pre-commit.ci] auto fixes from pre-commit.com hooks for more information, see https://pre-commit.ci --- .../protocols/openmm_ahfe/test_ahfe_protocol.py | 14 ++++++++------ 1 file changed, 8 insertions(+), 6 deletions(-) diff --git a/openfe/tests/protocols/openmm_ahfe/test_ahfe_protocol.py b/openfe/tests/protocols/openmm_ahfe/test_ahfe_protocol.py index 64cf623bd..206ff98c7 100644 --- a/openfe/tests/protocols/openmm_ahfe/test_ahfe_protocol.py +++ b/openfe/tests/protocols/openmm_ahfe/test_ahfe_protocol.py @@ -205,13 +205,15 @@ def test_dry_run_vac_benzene(benzene_system, method, protocol_dry_settings, tmpd _verify_alchemical_sterics_force_parameters(force) -@pytest.mark.parametrize('alpha, a, b, c, correction', [ - [0.2, 2, 2, 1, True], - [0.35, 2.2, 1.5, 0, False], -]) +@pytest.mark.parametrize( + "alpha, a, b, c, correction", + [ + [0.2, 2, 2, 1, True], + [0.35, 2.2, 1.5, 0, False], + ], +) def test_alchemical_settings_dry_run_vacuum( - alpha, a, b, c, correction, - benzene_system, protocol_dry_settings, tmpdir + alpha, a, b, c, correction, benzene_system, protocol_dry_settings, tmpdir ): """ Test non default alchemical settings