diff --git a/devtools/data/gen-serialized-results.py b/devtools/data/gen-serialized-results.py index d11d57307..ca7ff2200 100644 --- a/devtools/data/gen-serialized-results.py +++ b/devtools/data/gen-serialized-results.py @@ -1,3 +1,14 @@ +""" +Dev script to generate some result jsons that are used for testing + +Generates +- AHFEProtocol_json_results.gz + - used in afe_solvation_json fixture +- RHFEProtocol_json_results.gz + - used in rfe_transformation_json fixture +- MDProtocol_json_results.gz + - used in md_json fixture +""" import gzip import json import logging @@ -60,7 +71,7 @@ def generate_md_json(smc): settings.simulation_settings.equilibration_length_nvt = 0.01 * unit.nanosecond settings.simulation_settings.equilibration_length = 0.01 * unit.nanosecond settings.simulation_settings.production_length = 0.01 * unit.nanosecond - settings.system_settings.nonbonded_method = "nocutoff" + settings.forcefield_settings.nonbonded_method = "nocutoff" protocol = PlainMDProtocol(settings=settings) system = openfe.ChemicalSystem({"ligand": smc}) dag = protocol.create(stateA=system, stateB=system, mapping=None) @@ -80,9 +91,11 @@ def generate_ahfe_json(smc): settings.lambda_settings.lambda_vdw = [0.0, 0.0, 0.0, 0.0, 0.0, 0.12, 0.24, 0.36, 0.48, 0.6, 0.7, 0.77, 0.85, 1.0] - settings.alchemsampler_settings.n_repeats = 3 - settings.alchemsampler_settings.n_replicas = 14 - settings.alchemsampler_settings.online_analysis_target_error = 0.2 * unit.boltzmann_constant * unit.kelvin + settings.protocol_repeats = 3 + settings.solvent_simulation_settings.n_replicas = 14 + settings.vacuum_simulation_settings.n_replicas = 14 + settings.solvent_simulation_settings.early_termination_target_error = 0.12 * unit.kilocalorie_per_mole + settings.vacuum_simulation_settings.early_termination_target_error = 0.12 * unit.kilocalorie_per_mole settings.vacuum_engine_settings.compute_platform = 'CPU' settings.solvent_engine_settings.compute_platform = 'CUDA' @@ -103,7 +116,7 @@ def generate_rfe_json(smcA, smcB): settings = RelativeHybridTopologyProtocol.default_settings() settings.simulation_settings.equilibration_length = 10 * unit.picosecond settings.simulation_settings.production_length = 250 * unit.picosecond - settings.system_settings.nonbonded_method = "nocutoff" + settings.forcefield_settings.nonbonded_method = "nocutoff" protocol = RelativeHybridTopologyProtocol(settings=settings) a_smcB = align_mol_shape(smcB, ref_mol=smcA) diff --git a/docs/reference/api/openmm_rfe.rst b/docs/reference/api/openmm_rfe.rst index 1f395847d..21f340a8c 100644 --- a/docs/reference/api/openmm_rfe.rst +++ b/docs/reference/api/openmm_rfe.rst @@ -58,17 +58,6 @@ Below are the settings which can be tweaked in the protocol. The default setting :inherited-members: SettingsBaseModel :member-order: bysource -.. autopydantic_model:: AlchemicalSamplerSettings - :model-show-json: False - :model-show-field-summary: False - :model-show-config-member: False - :model-show-config-summary: False - :model-show-validator-members: False - :model-show-validator-summary: False - :field-list-validators: False - :inherited-members: SettingsBaseModel - :member-order: bysource - .. autopydantic_model:: AlchemicalSettings :model-show-json: False :model-show-field-summary: False @@ -113,7 +102,7 @@ Below are the settings which can be tweaked in the protocol. The default setting :inherited-members: SettingsBaseModel :member-order: bysource -.. autopydantic_model:: SimulationSettings +.. autopydantic_model:: MultiStateSimulationSettings :model-show-json: False :model-show-field-summary: False :model-show-config-member: False @@ -134,14 +123,3 @@ Below are the settings which can be tweaked in the protocol. The default setting :field-list-validators: False :inherited-members: SettingsBaseModel :member-order: bysource - -.. autopydantic_model:: SystemSettings - :model-show-json: False - :model-show-field-summary: False - :model-show-config-member: False - :model-show-config-summary: False - :model-show-validator-members: False - :model-show-validator-summary: False - :field-list-validators: False - :inherited-members: SettingsBaseModel - :member-order: bysource diff --git a/docs/reference/api/openmm_solvation_afe.rst b/docs/reference/api/openmm_solvation_afe.rst index 81d0877a0..b4d4706ff 100644 --- a/docs/reference/api/openmm_solvation_afe.rst +++ b/docs/reference/api/openmm_solvation_afe.rst @@ -61,18 +61,6 @@ Below are the settings which can be tweaked in the protocol. The default setting :member-order: bysource :noindex: -.. autopydantic_model:: AlchemicalSamplerSettings - :model-show-json: False - :model-show-field-summary: False - :model-show-config-member: False - :model-show-config-summary: False - :model-show-validator-members: False - :model-show-validator-summary: False - :field-list-validators: False - :inherited-members: SettingsBaseModel - :member-order: bysource - :noindex: - .. autopydantic_model:: AlchemicalSettings :model-show-json: False :model-show-field-summary: False @@ -121,7 +109,7 @@ Below are the settings which can be tweaked in the protocol. The default setting :member-order: bysource :noindex: -.. autopydantic_model:: SimulationSettings +.. autopydantic_model:: MultiStateSimulationSettings :model-show-json: False :model-show-field-summary: False :model-show-config-member: False @@ -144,15 +132,3 @@ Below are the settings which can be tweaked in the protocol. The default setting :inherited-members: SettingsBaseModel :member-order: bysource :noindex: - -.. autopydantic_model:: SystemSettings - :model-show-json: False - :model-show-field-summary: False - :model-show-config-member: False - :model-show-config-summary: False - :model-show-validator-members: False - :model-show-validator-summary: False - :field-list-validators: False - :inherited-members: SettingsBaseModel - :member-order: bysource - :noindex: diff --git a/openfe/protocols/openmm_afe/base.py b/openfe/protocols/openmm_afe/base.py index 1a963ba9c..de46195e9 100644 --- a/openfe/protocols/openmm_afe/base.py +++ b/openfe/protocols/openmm_afe/base.py @@ -52,8 +52,9 @@ ) from openfe.protocols.openmm_afe.equil_afe_settings import ( SolvationSettings, - AlchemicalSamplerSettings, OpenMMEngineSettings, - IntegratorSettings, SimulationSettings, LambdaSettings, + MultiStateSimulationSettings, OpenMMEngineSettings, + IntegratorSettings, LambdaSettings, OutputSettings, + ThermoSettings, ) from openfe.protocols.openmm_rfe._rfe_utils import compute from ..openmm_utils import ( @@ -230,14 +231,13 @@ def _handle_settings(self): Get a dictionary with the following entries: * forcefield_settings : OpenMMSystemGeneratorFFSettings * thermo_settings : ThermoSettings - * system_settings : SystemSettings * solvation_settings : SolvationSettings * alchemical_settings : AlchemicalSettings * lambda_settings : LambdaSettings - * sampler_settings : AlchemicalSamplerSettings * engine_settings : OpenMMEngineSettings * integrator_settings : IntegratorSettings - * simulation_settings : SimulationSettings + * simulation_settings : MultiStateSimulationSettings + * output_settings: OutputSettings Settings may change depending on what type of simulation you are running. Cherry pick them and return them to be available later on. @@ -270,15 +270,14 @@ def _get_system_generator( system_generator : openmmforcefields.generator.SystemGenerator System Generator to parameterise this unit. """ - ffcache = settings['simulation_settings'].forcefield_cache + ffcache = settings['output_settings'].forcefield_cache if ffcache is not None: ffcache = self.shared_basepath / ffcache system_generator = system_creation.get_system_generator( forcefield_settings=settings['forcefield_settings'], - thermo_settings=settings['thermo_settings'], integrator_settings=settings['integrator_settings'], - system_settings=settings['system_settings'], + thermo_settings=settings['thermo_settings'], cache=ffcache, has_solvent=solvent_comp is not None, ) @@ -537,7 +536,7 @@ def _get_reporter( self, topology: app.Topology, positions: openmm.unit.Quantity, - simulation_settings: SimulationSettings, + output_settings: OutputSettings, ) -> multistate.MultiStateReporter: """ Get a MultistateReporter for the simulation you are running. @@ -546,8 +545,8 @@ def _get_reporter( ---------- topology : app.Topology A Topology of the system being created. - simulation_settings : SimulationSettings - Settings for the simulation. + output_settings: OutputSettings + Output settings for the simulations Returns ------- @@ -557,16 +556,16 @@ def _get_reporter( mdt_top = mdt.Topology.from_openmm(topology) selection_indices = mdt_top.select( - simulation_settings.output_indices + output_settings.output_indices ) - nc = self.shared_basepath / simulation_settings.output_filename - chk = simulation_settings.checkpoint_storage + nc = self.shared_basepath / output_settings.output_filename + chk = output_settings.checkpoint_storage_filename reporter = multistate.MultiStateReporter( storage=nc, analysis_particle_indices=selection_indices, - checkpoint_interval=simulation_settings.checkpoint_interval.m, + checkpoint_interval=output_settings.checkpoint_interval.m, checkpoint_storage=chk, ) @@ -577,7 +576,7 @@ def _get_reporter( mdt_top.subset(selection_indices), ) traj.save_pdb( - self.shared_basepath / simulation_settings.output_structure + self.shared_basepath / output_settings.output_structure ) return reporter @@ -614,9 +613,10 @@ def _get_ctx_caches( return energy_context_cache, sampler_context_cache + @staticmethod def _get_integrator( - self, - integrator_settings: IntegratorSettings + integrator_settings: IntegratorSettings, + simulation_settings: MultiStateSimulationSettings ) -> openmmtools.mcmc.LangevinDynamicsMove: """ Return a LangevinDynamicsMove integrator @@ -624,16 +624,21 @@ def _get_integrator( Parameters ---------- integrator_settings : IntegratorSettings + simulation_settings : MultiStateSimulationSettings Returns ------- integrator : openmmtools.mcmc.LangevinDynamicsMove A configured integrator object. """ + steps_per_iteration = settings_validation.convert_steps_per_iteration( + simulation_settings, integrator_settings + ) + integrator = openmmtools.mcmc.LangevinDynamicsMove( timestep=to_openmm(integrator_settings.timestep), - collision_rate=to_openmm(integrator_settings.collision_rate), - n_steps=integrator_settings.n_steps.m, + collision_rate=to_openmm(integrator_settings.langevin_collision_rate), + n_steps=steps_per_iteration, reassign_velocities=integrator_settings.reassign_velocities, n_restart_attempts=integrator_settings.n_restart_attempts, constraint_tolerance=integrator_settings.constraint_tolerance, @@ -641,11 +646,12 @@ def _get_integrator( return integrator + @staticmethod def _get_sampler( - self, integrator: openmmtools.mcmc.LangevinDynamicsMove, reporter: openmmtools.multistate.MultiStateReporter, - sampler_settings: AlchemicalSamplerSettings, + simulation_settings: MultiStateSimulationSettings, + thermo_settings: ThermoSettings, cmp_states: list[ThermodynamicState], sampler_states: list[SamplerState], energy_context_cache: openmmtools.cache.ContextCache, @@ -660,8 +666,10 @@ def _get_sampler( The simulation integrator. reporter : openmmtools.multistate.MultiStateReporter The reporter to hook up to the sampler. - sampler_settings : AlchemicalSamplerSettings + simulation_settings : MultiStateSimulationSettings Settings for the alchemical sampler. + thermo_settings : ThermoSettings + Thermodynamic settings cmp_states : list[ThermodynamicState] A list of thermodynamic states to sample. sampler_states : list[SamplerState] @@ -676,30 +684,37 @@ def _get_sampler( sampler : multistate.MultistateSampler A sampler configured for the chosen sampling method. """ + rta_its, rta_min_its = settings_validation.convert_real_time_analysis_iterations( + simulation_settings=simulation_settings, + ) + et_target_err = settings_validation.convert_target_error_from_kcal_per_mole_to_kT( + thermo_settings.temperature, + simulation_settings.early_termination_target_error, + ) # Select the right sampler # Note: doesn't need else, settings already validates choices - if sampler_settings.sampler_method.lower() == "repex": + if simulation_settings.sampler_method.lower() == "repex": sampler = multistate.ReplicaExchangeSampler( mcmc_moves=integrator, - online_analysis_interval=sampler_settings.online_analysis_interval, - online_analysis_target_error=sampler_settings.online_analysis_target_error.m, - online_analysis_minimum_iterations=sampler_settings.online_analysis_minimum_iterations + online_analysis_interval=rta_its, + online_analysis_target_error=et_target_err, + online_analysis_minimum_iterations=rta_min_its ) - elif sampler_settings.sampler_method.lower() == "sams": + elif simulation_settings.sampler_method.lower() == "sams": sampler = multistate.SAMSSampler( mcmc_moves=integrator, - online_analysis_interval=sampler_settings.online_analysis_interval, - online_analysis_minimum_iterations=sampler_settings.online_analysis_minimum_iterations, - flatness_criteria=sampler_settings.flatness_criteria, - gamma0=sampler_settings.gamma0, + online_analysis_interval=rta_its, + online_analysis_minimum_iterations=rta_min_its, + flatness_criteria=simulation_settings.sams_flatness_criteria, + gamma0=simulation_settings.sams_gamma0, ) - elif sampler_settings.sampler_method.lower() == 'independent': + elif simulation_settings.sampler_method.lower() == 'independent': sampler = multistate.MultiStateSampler( mcmc_moves=integrator, - online_analysis_interval=sampler_settings.online_analysis_interval, - online_analysis_target_error=sampler_settings.online_analysis_target_error.m, - online_analysis_minimum_iterations=sampler_settings.online_analysis_minimum_iterations + online_analysis_interval=rta_its, + online_analysis_target_error=et_target_err, + online_analysis_minimum_iterations=rta_min_its, ) sampler.create( @@ -741,7 +756,10 @@ def _run_simulation( if not a dry run. """ # Get the relevant simulation steps - mc_steps = settings['integrator_settings'].n_steps.m + mc_steps = settings_validation.convert_steps_per_iteration( + simulation_settings=settings['simulation_settings'], + integrator_settings=settings['integrator_settings'], + ) equil_steps = settings_validation.get_simsteps( sim_length=settings['simulation_settings'].equilibration_length, @@ -780,7 +798,7 @@ def _run_simulation( analyzer = multistate_analysis.MultistateEquilFEAnalysis( reporter, - sampling_method=settings['sampler_settings'].sampler_method.lower(), + sampling_method=settings['simulation_settings'].sampler_method.lower(), result_units=unit.kilocalorie_per_mole ) analyzer.plot(filepath=self.shared_basepath, filename_prefix="") @@ -793,8 +811,8 @@ def _run_simulation( reporter.close() # clean up the reporter file - fns = [self.shared_basepath / settings['simulation_settings'].output_filename, - self.shared_basepath / settings['simulation_settings'].checkpoint_storage] + fns = [self.shared_basepath / settings['output_settings'].output_filename, + self.shared_basepath / settings['output_settings'].checkpoint_storage_filename] for fn in fns: os.remove(fn) @@ -878,7 +896,7 @@ def run(self, dry=False, verbose=True, # 11. Create the multistate reporter & create PDB reporter = self._get_reporter( omm_topology, positions, - settings['simulation_settings'], + settings['output_settings'], ) # Wrap in try/finally to avoid memory leak issues @@ -889,11 +907,15 @@ def run(self, dry=False, verbose=True, ) # 13. Get integrator - integrator = self._get_integrator(settings['integrator_settings']) + integrator = self._get_integrator( + settings['integrator_settings'], + settings['simulation_settings'], + ) # 14. Get sampler sampler = self._get_sampler( - integrator, reporter, settings['sampler_settings'], + integrator, reporter, settings['simulation_settings'], + settings['thermo_settings'], cmp_states, sampler_states, energy_ctx_cache, sampler_ctx_cache ) @@ -925,8 +947,8 @@ def run(self, dry=False, verbose=True, del integrator, sampler if not dry: - nc = self.shared_basepath / settings['simulation_settings'].output_filename - chk = settings['simulation_settings'].checkpoint_storage + nc = self.shared_basepath / settings['output_settings'].output_filename + chk = settings['output_settings'].checkpoint_storage_filename return { 'nc': nc, 'last_checkpoint': chk, diff --git a/openfe/protocols/openmm_afe/equil_afe_settings.py b/openfe/protocols/openmm_afe/equil_afe_settings.py index 482dc2056..32c519dd2 100644 --- a/openfe/protocols/openmm_afe/equil_afe_settings.py +++ b/openfe/protocols/openmm_afe/equil_afe_settings.py @@ -16,22 +16,19 @@ """ from gufe.settings import ( - Settings, SettingsBaseModel, OpenMMSystemGeneratorFFSettings, ThermoSettings, ) from openfe.protocols.openmm_utils.omm_settings import ( - SystemSettings, + MultiStateSimulationSettings, SolvationSettings, - AlchemicalSamplerSettings, OpenMMEngineSettings, IntegratorSettings, - SimulationSettings + OutputSettings, ) import numpy as np - try: from pydantic.v1 import validator except ImportError: @@ -52,12 +49,13 @@ class LambdaSettings(SettingsBaseModel): transformation. Notes - -------- + ----- * In all cases a lambda value of 0 defines a fully interacting state A and - a non-interacting state B, whilst a value of 1 defines a fully interacting - state B and a non-interacting state A. + a non-interacting state B, whilst a value of 1 defines a fully interacting + state B and a non-interacting state A. * ``lambda_elec``, `lambda_vdw``, and ``lambda_restraints`` must all be of - the same length, defining all the windows of the transformation. + the same length, defining all the windows of the transformation. + """ lambda_elec: list[float] = [ 0.0, 0.25, 0.5, 0.75, 1.0, 1.0, 1.0, 1.0, 1.0, 1.0, @@ -91,7 +89,8 @@ class LambdaSettings(SettingsBaseModel): def must_be_between_0_and_1(cls, v): for window in v: if not 0 <= window <= 1: - errmsg = "Lambda windows must be between 0 and 1." + errmsg = ("Lambda windows must be between 0 and 1, got a" + f" window with value {window}.") raise ValueError(errmsg) return v @@ -101,13 +100,15 @@ def must_be_monotonic(cls, v): difference = np.diff(v) if not all(i >= 0. for i in difference): - errmsg = "The lambda schedule is not monotonic." + errmsg = f"The lambda schedule is not monotonic, got schedule {v}." raise ValueError(errmsg) return v -class AbsoluteSolvationSettings(Settings): +# This subclasses from SettingsBaseModel as it has vacuum_forcefield and +# solvent_forcefield fields, not just a single forcefield_settings field +class AbsoluteSolvationSettings(SettingsBaseModel): """ Configuration object for ``AbsoluteSolvationProtocol``. @@ -115,26 +116,27 @@ class AbsoluteSolvationSettings(Settings): -------- openfe.protocols.openmm_afe.AbsoluteSolvationProtocol """ - class Config: - arbitrary_types_allowed = True + protocol_repeats: int + """ + The number of completely independent repeats of the entire sampling + process. The mean of the repeats defines the final estimate of FE + difference, while the variance between repeats is used as the uncertainty. + """ + + @validator('protocol_repeats') + def must_be_positive(cls, v): + if v <= 0: + errmsg = f"protocol_repeats must be a positive value, got {v}." + raise ValueError(errmsg) + return v # Inherited things - forcefield_settings: OpenMMSystemGeneratorFFSettings + solvent_forcefield_settings: OpenMMSystemGeneratorFFSettings + vacuum_forcefield_settings: OpenMMSystemGeneratorFFSettings """Parameters to set up the force field with OpenMM Force Fields""" thermo_settings: ThermoSettings """Settings for thermodynamic parameters""" - # Things for creating the systems - vacuum_system_settings: SystemSettings - """ - Simulation system settings including the - long-range non-bonded methods for the vacuum transformation. - """ - solvent_system_settings: SystemSettings - """ - Simulation system settings including the - long-range non-bonded methods for the solvent transformation. - """ solvation_settings: SolvationSettings """Settings for solvating the system.""" @@ -148,11 +150,6 @@ class Config: Settings for controlling the lambda schedule for the different components (vdw, elec, restraints). """ - alchemsampler_settings: AlchemicalSamplerSettings - """ - Settings for controlling how we sample alchemical space, including the - number of repeats. - """ # MD Engine things vacuum_engine_settings: OpenMMEngineSettings @@ -174,13 +171,21 @@ class Config: """ # Simulation run settings - vacuum_simulation_settings: SimulationSettings + vacuum_simulation_settings: MultiStateSimulationSettings + """ + Simulation control settings, including simulation lengths + for the vacuum transformation. + """ + solvent_simulation_settings: MultiStateSimulationSettings + """ + Simulation control settings, including simulation lengths + for the solvent transformation. + """ + vacuum_output_settings: OutputSettings """ - Simulation control settings, including simulation lengths and - record-keeping for the vacuum transformation. + Simulation output settings for the vacuum transformation. """ - solvent_simulation_settings: SimulationSettings + solvent_output_settings: OutputSettings """ - Simulation control settings, including simulation lengths and - record-keeping for the solvent transformation. + Simulation output settings for the solvent transformation. """ diff --git a/openfe/protocols/openmm_afe/equil_solvation_afe_method.py b/openfe/protocols/openmm_afe/equil_solvation_afe_method.py index cecba32ac..99b1f9e3d 100644 --- a/openfe/protocols/openmm_afe/equil_solvation_afe_method.py +++ b/openfe/protocols/openmm_afe/equil_solvation_afe_method.py @@ -34,13 +34,12 @@ from collections import defaultdict import gufe from gufe.components import Component -from openff.toolkit.topology import Molecule as OFFMolecule import itertools import numpy as np import numpy.typing as npt from openff.units import unit from openmmtools import multistate -from typing import Dict, Optional, Union +from typing import Optional, Union from typing import Any, Iterable import uuid @@ -50,10 +49,10 @@ ProteinComponent, SolventComponent ) from openfe.protocols.openmm_afe.equil_afe_settings import ( - AbsoluteSolvationSettings, SystemSettings, + AbsoluteSolvationSettings, SolvationSettings, AlchemicalSettings, LambdaSettings, - AlchemicalSamplerSettings, OpenMMEngineSettings, - IntegratorSettings, SimulationSettings, + MultiStateSimulationSettings, OpenMMEngineSettings, + IntegratorSettings, OutputSettings, SettingsBaseModel, ) from ..openmm_utils import system_validation, settings_validation @@ -397,13 +396,15 @@ def _default_settings(cls): a set of default settings """ return AbsoluteSolvationSettings( - forcefield_settings=settings.OpenMMSystemGeneratorFFSettings(), + protocol_repeats=3, + solvent_forcefield_settings=settings.OpenMMSystemGeneratorFFSettings(), + vacuum_forcefield_settings=settings.OpenMMSystemGeneratorFFSettings( + nonbonded_method='nocutoff', + ), thermo_settings=settings.ThermoSettings( temperature=298.15 * unit.kelvin, pressure=1 * unit.bar, ), - solvent_system_settings=SystemSettings(), - vacuum_system_settings=SystemSettings(nonbonded_method='nocutoff'), alchemical_settings=AlchemicalSettings(), lambda_settings=LambdaSettings( lambda_elec=[ @@ -414,24 +415,27 @@ def _default_settings(cls): 0.0, 0.0, 0.0, 0.0, 0.0, 0.12, 0.24, 0.36, 0.48, 0.6, 0.7, 0.77, 0.85, 1.0], ), - alchemsampler_settings=AlchemicalSamplerSettings( - n_replicas=14, - ), solvation_settings=SolvationSettings(), vacuum_engine_settings=OpenMMEngineSettings(), solvent_engine_settings=OpenMMEngineSettings(), integrator_settings=IntegratorSettings(), - solvent_simulation_settings=SimulationSettings( + solvent_simulation_settings=MultiStateSimulationSettings( + n_replicas=14, equilibration_length=1.0 * unit.nanosecond, production_length=10.0 * unit.nanosecond, + ), + solvent_output_settings=OutputSettings( output_filename='solvent.nc', - checkpoint_storage='solvent_checkpoint.nc', + checkpoint_storage_filename='solvent_checkpoint.nc', ), - vacuum_simulation_settings=SimulationSettings( + vacuum_simulation_settings=MultiStateSimulationSettings( + n_replicas=14, equilibration_length=0.5 * unit.nanosecond, production_length=2.0 * unit.nanosecond, + ), + vacuum_output_settings=OutputSettings( output_filename='vacuum.nc', - checkpoint_storage='vacuum_checkpoint.nc' + checkpoint_storage_filename='vacuum_checkpoint.nc' ), ) @@ -530,15 +534,17 @@ def _validate_alchemical_components( @staticmethod def _validate_lambda_schedule( lambda_settings: LambdaSettings, - alchemsampler_settings: AlchemicalSamplerSettings, + simulation_settings: MultiStateSimulationSettings, ) -> None: """ Checks that the lambda schedule is set up correctly. Parameters ---------- - settings : AbsoluteSolvationSettings - Settings object. + lambda_settings : LambdaSettings + the lambda schedule Settings + simulation_settings : MultiStateSimulationSettings + the settings for either the vacuum or solvent phase Raises ------ @@ -553,7 +559,7 @@ def _validate_lambda_schedule( lambda_elec = lambda_settings.lambda_elec lambda_vdw = lambda_settings.lambda_vdw lambda_restraints = lambda_settings.lambda_restraints - n_replicas = alchemsampler_settings.n_replicas + n_replicas = simulation_settings.n_replicas # Ensure that all lambda components have equal amount of windows lambda_components = [lambda_vdw, lambda_elec] @@ -596,7 +602,7 @@ def _create( self, stateA: ChemicalSystem, stateB: ChemicalSystem, - mapping: Optional[Dict[str, gufe.ComponentMapping]] = None, + mapping: Optional[dict[str, gufe.ComponentMapping]] = None, extends: Optional[gufe.ProtocolDAGResult] = None, ) -> list[gufe.ProtocolUnit]: # TODO: extensions @@ -612,11 +618,13 @@ def _create( # Validate the lambda schedule self._validate_lambda_schedule(self.settings.lambda_settings, - self.settings.alchemsampler_settings) + self.settings.solvent_simulation_settings) + self._validate_lambda_schedule(self.settings.lambda_settings, + self.settings.vacuum_simulation_settings) # Check nonbond & solvent compatibility - solv_nonbonded_method = self.settings.solvent_system_settings.nonbonded_method - vac_nonbonded_method = self.settings.vacuum_system_settings.nonbonded_method + solv_nonbonded_method = self.settings.solvent_forcefield_settings.nonbonded_method + vac_nonbonded_method = self.settings.vacuum_forcefield_settings.nonbonded_method # Use the more complete system validation solvent checks system_validation.validate_solvent(stateA, solv_nonbonded_method) # Gas phase is always gas phase @@ -640,7 +648,7 @@ def _create( name=(f"Absolute Solvation, {alchname} solvent leg: " f"repeat {i} generation 0"), ) - for i in range(self.settings.alchemsampler_settings.n_repeats) + for i in range(self.settings.protocol_repeats) ] vacuum_units = [ @@ -654,14 +662,14 @@ def _create( name=(f"Absolute Solvation, {alchname} vacuum leg: " f"repeat {i} generation 0"), ) - for i in range(self.settings.alchemsampler_settings.n_repeats) + for i in range(self.settings.protocol_repeats) ] return solvent_units + vacuum_units def _gather( self, protocol_dag_results: Iterable[gufe.ProtocolDAGResult] - ) -> Dict[str, Dict[str, Any]]: + ) -> dict[str, dict[str, Any]]: # result units will have a repeat_id and generation # first group according to repeat_id unsorted_solvent_repeats = defaultdict(list) @@ -731,28 +739,26 @@ def _handle_settings(self) -> dict[str, SettingsBaseModel]: A dictionary with the following entries: * forcefield_settings : OpenMMSystemGeneratorFFSettings * thermo_settings : ThermoSettings - * system_settings : SystemSettings * solvation_settings : SolvationSettings * alchemical_settings : AlchemicalSettings * lambda_settings : LambdaSettings - * sampler_settings : AlchemicalSamplerSettings * engine_settings : OpenMMEngineSettings * integrator_settings : IntegratorSettings * simulation_settings : SimulationSettings + * output_settings: OutputSettings """ prot_settings = self._inputs['settings'] settings = {} - settings['forcefield_settings'] = prot_settings.forcefield_settings + settings['forcefield_settings'] = prot_settings.vacuum_forcefield_settings settings['thermo_settings'] = prot_settings.thermo_settings - settings['system_settings'] = prot_settings.vacuum_system_settings settings['solvation_settings'] = prot_settings.solvation_settings settings['alchemical_settings'] = prot_settings.alchemical_settings settings['lambda_settings'] = prot_settings.lambda_settings - settings['sampler_settings'] = prot_settings.alchemsampler_settings settings['engine_settings'] = prot_settings.vacuum_engine_settings settings['integrator_settings'] = prot_settings.integrator_settings settings['simulation_settings'] = prot_settings.vacuum_simulation_settings + settings['output_settings'] = prot_settings.vacuum_output_settings settings_validation.validate_timestep( settings['forcefield_settings'].hydrogen_mass, @@ -763,7 +769,7 @@ def _handle_settings(self) -> dict[str, SettingsBaseModel]: def _execute( self, ctx: gufe.Context, **kwargs, - ) -> Dict[str, Any]: + ) -> dict[str, Any]: log_system_probe(logging.INFO, paths=[ctx.scratch]) with without_oechem_backend(): @@ -817,28 +823,26 @@ def _handle_settings(self) -> dict[str, SettingsBaseModel]: A dictionary with the following entries: * forcefield_settings : OpenMMSystemGeneratorFFSettings * thermo_settings : ThermoSettings - * system_settings : SystemSettings * solvation_settings : SolvationSettings * alchemical_settings : AlchemicalSettings * lambda_settings : LambdaSettings - * sampler_settings : AlchemicalSamplerSettings * engine_settings : OpenMMEngineSettings * integrator_settings : IntegratorSettings - * simulation_settings : SimulationSettings + * simulation_settings : MultiStateSimulationSettings + * output_settings: OutputSettings """ prot_settings = self._inputs['settings'] settings = {} - settings['forcefield_settings'] = prot_settings.forcefield_settings + settings['forcefield_settings'] = prot_settings.solvent_forcefield_settings settings['thermo_settings'] = prot_settings.thermo_settings - settings['system_settings'] = prot_settings.solvent_system_settings settings['solvation_settings'] = prot_settings.solvation_settings settings['alchemical_settings'] = prot_settings.alchemical_settings settings['lambda_settings'] = prot_settings.lambda_settings - settings['sampler_settings'] = prot_settings.alchemsampler_settings settings['engine_settings'] = prot_settings.solvent_engine_settings settings['integrator_settings'] = prot_settings.integrator_settings settings['simulation_settings'] = prot_settings.solvent_simulation_settings + settings['output_settings'] = prot_settings.solvent_output_settings settings_validation.validate_timestep( settings['forcefield_settings'].hydrogen_mass, @@ -849,7 +853,7 @@ def _handle_settings(self) -> dict[str, SettingsBaseModel]: def _execute( self, ctx: gufe.Context, **kwargs, - ) -> Dict[str, Any]: + ) -> dict[str, Any]: log_system_probe(logging.INFO, paths=[ctx.scratch]) with without_oechem_backend(): diff --git a/openfe/protocols/openmm_md/plain_md_methods.py b/openfe/protocols/openmm_md/plain_md_methods.py index da097f972..988c22008 100644 --- a/openfe/protocols/openmm_md/plain_md_methods.py +++ b/openfe/protocols/openmm_md/plain_md_methods.py @@ -33,10 +33,9 @@ ProteinComponent, SolventComponent ) from openfe.protocols.openmm_md.plain_md_settings import ( - PlainMDProtocolSettings, SystemSettings, + PlainMDProtocolSettings, SolvationSettings, OpenMMEngineSettings, - IntegratorSettings, SimulationSettingsMD, - RepeatSettings, + IntegratorSettings, MDSimulationSettings, MDOutputSettings, ) from openff.toolkit.topology import Molecule as OFFMolecule @@ -122,16 +121,16 @@ def _default_settings(cls): temperature=298.15 * unit.kelvin, pressure=1 * unit.bar, ), - system_settings=SystemSettings(), solvation_settings=SolvationSettings(), engine_settings=OpenMMEngineSettings(), integrator_settings=IntegratorSettings(), - simulation_settings=SimulationSettingsMD( + simulation_settings=MDSimulationSettings( equilibration_length_nvt=0.1 * unit.nanosecond, equilibration_length=1.0 * unit.nanosecond, production_length=5.0 * unit.nanosecond, ), - repeat_settings=RepeatSettings(), + output_settings=MDOutputSettings(), + protocol_repeats=1, ) def _create( @@ -146,7 +145,7 @@ def _create( raise NotImplementedError("Can't extend simulations yet") # Validate solvent component - nonbond = self.settings.system_settings.nonbonded_method + nonbond = self.settings.forcefield_settings.nonbonded_method system_validation.validate_solvent(stateA, nonbond) # Validate protein component @@ -168,7 +167,7 @@ def _create( system_name += f" {comp_type}:{comp_name}" # our DAG has no dependencies, so just list units - n_repeats = self.settings.repeat_settings.n_repeats + n_repeats = self.settings.protocol_repeats units = [PlainMDProtocolUnit( stateA=stateA, settings=self.settings, @@ -246,7 +245,8 @@ def __init__(self, *, @staticmethod def _run_MD(simulation: openmm.app.Simulation, positions: omm_unit.Quantity, - simulation_settings: SimulationSettingsMD, + simulation_settings: MDSimulationSettings, + output_settings: MDOutputSettings, temperature: settings.ThermoSettings.temperature, barostat_frequency: IntegratorSettings.barostat_frequency, equil_steps_nvt: int, @@ -267,6 +267,8 @@ def _run_MD(simulation: openmm.app.Simulation, Initial positions for the system. simulation_settings : SimulationSettingsMD Settings for MD simulation + output_settings: OutputSettingsMD + Settings for output of MD simulation temperature: settings.ThermoSettings.temperature temperature setting barostat_frequency: IntegratorSettings.barostat_frequency @@ -300,7 +302,7 @@ def _run_MD(simulation: openmm.app.Simulation, # Get the sub selection of the system to save coords for selection_indices = mdtraj.Topology.from_openmm( - simulation.topology).select(simulation_settings.output_indices) + simulation.topology).select(output_settings.output_indices) positions = to_openmm(from_openmm( simulation.context.getState(getPositions=True, @@ -314,7 +316,7 @@ def _run_MD(simulation: openmm.app.Simulation, ) traj.save_pdb( - shared_basepath / simulation_settings.minimized_structure + shared_basepath / output_settings.minimized_structure ) # equilibrate # NVT equilibration @@ -348,7 +350,7 @@ def _run_MD(simulation: openmm.app.Simulation, mdtraj_top.subset(selection_indices), ) traj.save_pdb( - shared_basepath / simulation_settings.equil_NVT_structure + shared_basepath / output_settings.equil_NVT_structure ) # NPT equilibration @@ -380,7 +382,7 @@ def _run_MD(simulation: openmm.app.Simulation, mdtraj_top.subset(selection_indices), ) traj.save_pdb( - shared_basepath / simulation_settings.equil_NPT_structure + shared_basepath / output_settings.equil_NPT_structure ) # production @@ -389,15 +391,15 @@ def _run_MD(simulation: openmm.app.Simulation, # Setup the reporters simulation.reporters.append(XTCReporter( - file=str(shared_basepath / simulation_settings.production_trajectory_filename), - reportInterval=simulation_settings.trajectory_write_interval.m, + file=str(shared_basepath / output_settings.production_trajectory_filename), + reportInterval=output_settings.trajectory_write_interval.m, atomSubset=selection_indices)) simulation.reporters.append(openmm.app.CheckpointReporter( - file=str(shared_basepath / simulation_settings.checkpoint_storage_filename), - reportInterval=simulation_settings.checkpoint_interval.m)) + file=str(shared_basepath / output_settings.checkpoint_storage_filename), + reportInterval=output_settings.checkpoint_interval.m)) simulation.reporters.append(openmm.app.StateDataReporter( - str(shared_basepath / simulation_settings.log_output), - simulation_settings.checkpoint_interval.m, + str(shared_basepath / output_settings.log_output), + output_settings.checkpoint_interval.m, step=True, time=True, potentialEnergy=True, @@ -463,13 +465,11 @@ def run(self, *, dry=False, verbose=True, protocol_settings.forcefield_settings thermo_settings: settings.ThermoSettings = \ protocol_settings.thermo_settings - system_settings: SystemSettings = protocol_settings.system_settings solvation_settings: SolvationSettings = \ protocol_settings.solvation_settings - sim_settings: SimulationSettingsMD = \ - protocol_settings.simulation_settings + sim_settings: MDSimulationSettings = protocol_settings.simulation_settings + output_settings: MDOutputSettings = protocol_settings.output_settings timestep = protocol_settings.integrator_settings.timestep - mc_steps = protocol_settings.integrator_settings.n_steps.m integrator_settings = protocol_settings.integrator_settings # is the timestep good for the mass? @@ -479,31 +479,30 @@ def run(self, *, dry=False, verbose=True, equil_steps_nvt = settings_validation.get_simsteps( sim_length=sim_settings.equilibration_length_nvt, - timestep=timestep, mc_steps=mc_steps + timestep=timestep, mc_steps=1, ) equil_steps_npt = settings_validation.get_simsteps( sim_length=sim_settings.equilibration_length, - timestep=timestep, mc_steps=mc_steps + timestep=timestep, mc_steps=1, ) prod_steps = settings_validation.get_simsteps( sim_length=sim_settings.production_length, - timestep=timestep, mc_steps=mc_steps + timestep=timestep, mc_steps=1, ) solvent_comp, protein_comp, small_mols = system_validation.get_components(stateA) # 1. Create stateA system # a. get a system generator - if sim_settings.forcefield_cache is not None: - ffcache = shared_basepath / sim_settings.forcefield_cache + if output_settings.forcefield_cache is not None: + ffcache = shared_basepath / output_settings.forcefield_cache else: ffcache = None system_generator = system_creation.get_system_generator( forcefield_settings=forcefield_settings, - thermo_settings=thermo_settings, integrator_settings=integrator_settings, - system_settings=system_settings, + thermo_settings=thermo_settings, cache=ffcache, has_solvent=solvent_comp is not None, ) @@ -548,7 +547,7 @@ def run(self, *, dry=False, verbose=True, ) # f. Save pdb of entire system - with open(shared_basepath / sim_settings.preminimized_structure, "w") as f: + with open(shared_basepath / output_settings.preminimized_structure, "w") as f: openmm.app.PDBFile.writeFile( stateA_topology, stateA_positions, file=f, keepIds=True ) @@ -561,7 +560,7 @@ def run(self, *, dry=False, verbose=True, # 11. Set the integrator integrator = openmm.LangevinMiddleIntegrator( to_openmm(thermo_settings.temperature), - to_openmm(integrator_settings.collision_rate), + to_openmm(integrator_settings.langevin_collision_rate), to_openmm(timestep), ) @@ -578,6 +577,7 @@ def run(self, *, dry=False, verbose=True, self._run_MD(simulation, stateA_positions, sim_settings, + output_settings, thermo_settings.temperature, integrator_settings.barostat_frequency, equil_steps_nvt, @@ -593,12 +593,12 @@ def run(self, *, dry=False, verbose=True, if not dry: # pragma: no-cover return { - 'system_pdb': shared_basepath / sim_settings.preminimized_structure, - 'minimized_pdb': shared_basepath / sim_settings.minimized_structure, - 'nvt_equil_pdb': shared_basepath / sim_settings.equil_NVT_structure, - 'npt_equil_pdb': shared_basepath / sim_settings.equil_NPT_structure, - 'nc': shared_basepath / sim_settings.production_trajectory_filename, - 'last_checkpoint': shared_basepath / sim_settings.checkpoint_storage_filename, + 'system_pdb': shared_basepath / output_settings.preminimized_structure, + 'minimized_pdb': shared_basepath / output_settings.minimized_structure, + 'nvt_equil_pdb': shared_basepath / output_settings.equil_NVT_structure, + 'npt_equil_pdb': shared_basepath / output_settings.equil_NPT_structure, + 'nc': shared_basepath / output_settings.production_trajectory_filename, + 'last_checkpoint': shared_basepath / output_settings.checkpoint_storage_filename, } else: return {'debug': {'system': stateA_system}} diff --git a/openfe/protocols/openmm_md/plain_md_settings.py b/openfe/protocols/openmm_md/plain_md_settings.py index 09b4b3b4e..8e9cd3612 100644 --- a/openfe/protocols/openmm_md/plain_md_settings.py +++ b/openfe/protocols/openmm_md/plain_md_settings.py @@ -8,28 +8,34 @@ """ from openfe.protocols.openmm_utils.omm_settings import ( - Settings, SystemSettings, - SolvationSettings, OpenMMEngineSettings, SimulationSettingsMD, - IntegratorSettings + Settings, + SolvationSettings, OpenMMEngineSettings, MDSimulationSettings, + IntegratorSettings, MDOutputSettings, ) from gufe.settings import SettingsBaseModel +try: + from pydantic.v1 import validator +except ImportError: + from pydantic import validator # type: ignore[assignment] -class RepeatSettings(SettingsBaseModel): - """Settings for how many independent MD runs to perform.""" +class PlainMDProtocolSettings(Settings): + class Config: + arbitrary_types_allowed = True - n_repeats: int = 1 + protocol_repeats: int """ - Number of independent repeats to run. Default 1 + Number of independent MD runs to perform. """ - -class PlainMDProtocolSettings(Settings): - class Config: - arbitrary_types_allowed = True + @validator('protocol_repeats') + def must_be_positive(cls, v): + if v <= 0: + errmsg = f"protocol_repeats must be a positive value, got {v}." + raise ValueError(errmsg) + return v # Things for creating the systems - system_settings: SystemSettings solvation_settings: SolvationSettings # MD Engine things @@ -39,7 +45,7 @@ class Config: integrator_settings: IntegratorSettings # Simulation run settings - simulation_settings: SimulationSettingsMD + simulation_settings: MDSimulationSettings - # Setting number of repeats of md simulation - repeat_settings: RepeatSettings + # Simulations output settings + output_settings: MDOutputSettings diff --git a/openfe/protocols/openmm_rfe/_rfe_utils/lambdaprotocol.py b/openfe/protocols/openmm_rfe/_rfe_utils/lambdaprotocol.py index 71ba31285..a8e764dd2 100644 --- a/openfe/protocols/openmm_rfe/_rfe_utils/lambdaprotocol.py +++ b/openfe/protocols/openmm_rfe/_rfe_utils/lambdaprotocol.py @@ -9,7 +9,7 @@ class LambdaProtocol(object): - """Protocols for perturbing each of the compent energy terms in alchemical + """Protocols for perturbing each of the component energy terms in alchemical free energy simulations. TODO diff --git a/openfe/protocols/openmm_rfe/equil_rfe_methods.py b/openfe/protocols/openmm_rfe/equil_rfe_methods.py index 15c822227..b275ac72b 100644 --- a/openfe/protocols/openmm_rfe/equil_rfe_methods.py +++ b/openfe/protocols/openmm_rfe/equil_rfe_methods.py @@ -52,10 +52,10 @@ ) from .equil_rfe_settings import ( - RelativeHybridTopologyProtocolSettings, SystemSettings, + RelativeHybridTopologyProtocolSettings, SolvationSettings, AlchemicalSettings, LambdaSettings, - AlchemicalSamplerSettings, OpenMMEngineSettings, - IntegratorSettings, SimulationSettings + MultiStateSimulationSettings, OpenMMEngineSettings, + IntegratorSettings, OutputSettings, ) from ..openmm_utils import ( system_validation, settings_validation, system_creation, @@ -448,22 +448,22 @@ def _default_settings(cls): a set of default settings """ return RelativeHybridTopologyProtocolSettings( + protocol_repeats=3, forcefield_settings=settings.OpenMMSystemGeneratorFFSettings(), thermo_settings=settings.ThermoSettings( temperature=298.15 * unit.kelvin, pressure=1 * unit.bar, ), - system_settings=SystemSettings(), solvation_settings=SolvationSettings(), - alchemical_settings=AlchemicalSettings(), + alchemical_settings=AlchemicalSettings(softcore_LJ='gapsys'), lambda_settings=LambdaSettings(), - alchemical_sampler_settings=AlchemicalSamplerSettings(), - engine_settings=OpenMMEngineSettings(), - integrator_settings=IntegratorSettings(), - simulation_settings=SimulationSettings( + simulation_settings=MultiStateSimulationSettings( equilibration_length=1.0 * unit.nanosecond, production_length=5.0 * unit.nanosecond, - ) + ), + engine_settings=OpenMMEngineSettings(), + integrator_settings=IntegratorSettings(), + output_settings=OutputSettings(), ) def _create( @@ -487,7 +487,7 @@ def _create( ligandmapping = list(mapping.values())[0] # type: ignore # Validate solvent component - nonbond = self.settings.system_settings.nonbonded_method + nonbond = self.settings.forcefield_settings.nonbonded_method system_validation.validate_solvent(stateA, nonbond) # Validate protein component @@ -497,7 +497,7 @@ def _create( Anames = ','.join(c.name for c in alchem_comps['stateA']) Bnames = ','.join(c.name for c in alchem_comps['stateB']) # our DAG has no dependencies, so just list units - n_repeats = self.settings.alchemical_sampler_settings.n_repeats + n_repeats = self.settings.protocol_repeats units = [RelativeHybridTopologyProtocolUnit( stateA=stateA, stateB=stateB, ligandmapping=ligandmapping, settings=self.settings, @@ -628,10 +628,9 @@ def run(self, *, dry=False, verbose=True, thermo_settings: settings.ThermoSettings = protocol_settings.thermo_settings alchem_settings: AlchemicalSettings = protocol_settings.alchemical_settings lambda_settings: LambdaSettings = protocol_settings.lambda_settings - system_settings: SystemSettings = protocol_settings.system_settings solvation_settings: SolvationSettings = protocol_settings.solvation_settings - sampler_settings: AlchemicalSamplerSettings = protocol_settings.alchemical_sampler_settings - sim_settings: SimulationSettings = protocol_settings.simulation_settings + sampler_settings: MultiStateSimulationSettings = protocol_settings.simulation_settings + output_settings: OutputSettings = protocol_settings.output_settings integrator_settings: IntegratorSettings = protocol_settings.integrator_settings # is the timestep good for the mass? @@ -639,15 +638,22 @@ def run(self, *, dry=False, verbose=True, forcefield_settings.hydrogen_mass, integrator_settings.timestep ) + # TODO: Also validate various conversions? + # Convert various time based inputs to steps/iterations + steps_per_iteration = settings_validation.convert_steps_per_iteration( + simulation_settings=sampler_settings, + integrator_settings=integrator_settings, + ) + equil_steps = settings_validation.get_simsteps( - sim_length=sim_settings.equilibration_length, + sim_length=sampler_settings.equilibration_length, timestep=integrator_settings.timestep, - mc_steps=integrator_settings.n_steps.m, + mc_steps=steps_per_iteration, ) prod_steps = settings_validation.get_simsteps( - sim_length=sim_settings.production_length, + sim_length=sampler_settings.production_length, timestep=integrator_settings.timestep, - mc_steps=integrator_settings.n_steps.m, + mc_steps=steps_per_iteration, ) solvent_comp, protein_comp, small_mols = system_validation.get_components(stateA) @@ -656,23 +662,22 @@ def run(self, *, dry=False, verbose=True, # and check if the charge correction used is appropriate charge_difference = _get_alchemical_charge_difference( mapping, - system_settings.nonbonded_method, + forcefield_settings.nonbonded_method, alchem_settings.explicit_charge_correction, solvent_comp, ) # 1. Create stateA system # a. get a system generator - if sim_settings.forcefield_cache is not None: - ffcache = shared_basepath / sim_settings.forcefield_cache + if output_settings.forcefield_cache is not None: + ffcache = shared_basepath / output_settings.forcefield_cache else: ffcache = None system_generator = system_creation.get_system_generator( forcefield_settings=forcefield_settings, - thermo_settings=thermo_settings, integrator_settings=integrator_settings, - system_settings=system_settings, + thermo_settings=thermo_settings, cache=ffcache, has_solvent=solvent_comp is not None, ) @@ -781,6 +786,12 @@ def run(self, *, dry=False, verbose=True, ) # 3. Create the hybrid topology + # a. Get softcore potential settings + if alchem_settings.softcore_LJ.lower() == 'gapsys': + softcore_LJ_v2 = True + elif alchem_settings.softcore_LJ.lower() == 'beutler': + softcore_LJ_v2 = False + # b. Get hybrid topology factory hybrid_factory = _rfe_utils.relative.HybridTopologyFactory( stateA_system, stateA_positions, stateA_topology, stateB_system, stateB_positions, stateB_topology, @@ -788,9 +799,9 @@ def run(self, *, dry=False, verbose=True, old_to_new_core_atom_map=ligand_mappings['old_to_new_core_atom_map'], use_dispersion_correction=alchem_settings.use_dispersion_correction, softcore_alpha=alchem_settings.softcore_alpha, - softcore_LJ_v2=alchem_settings.softcore_LJ_v2, + softcore_LJ_v2=softcore_LJ_v2, softcore_LJ_v2_alpha=alchem_settings.softcore_alpha, - interpolate_old_and_new_14s=alchem_settings.interpolate_old_and_new_14s, + interpolate_old_and_new_14s=alchem_settings.turn_off_core_unique_exceptions, ) # 4. Create lambda schedule @@ -812,16 +823,20 @@ def run(self, *, dry=False, verbose=True, # 9. Create the multistate reporter # Get the sub selection of the system to print coords for selection_indices = hybrid_factory.hybrid_topology.select( - sim_settings.output_indices + output_settings.output_indices ) # a. Create the multistate reporter - nc = shared_basepath / sim_settings.output_filename - chk = sim_settings.checkpoint_storage + # convert checkpoint_interval from time to steps + checkpoint_fs = output_settings.checkpoint_interval.to(unit.femtosecond).m + ts_fs = integrator_settings.timestep.to(unit.femtosecond).m + checkpoint_int = int(round(checkpoint_fs / ts_fs)) + nc = shared_basepath / output_settings.output_filename + chk = output_settings.checkpoint_storage_filename reporter = multistate.MultiStateReporter( storage=nc, analysis_particle_indices=selection_indices, - checkpoint_interval=sim_settings.checkpoint_interval.m, + checkpoint_interval=checkpoint_int, checkpoint_storage=chk, ) @@ -837,7 +852,7 @@ def run(self, *, dry=False, verbose=True, hybrid_factory.hybrid_positions[selection_indices, :], hybrid_factory.hybrid_topology.subset(selection_indices), ).save_pdb( - shared_basepath / sim_settings.output_structure, + shared_basepath / output_settings.output_structure, bfactors=bfactors, ) @@ -859,8 +874,8 @@ def run(self, *, dry=False, verbose=True, # b. create langevin integrator integrator = openmmtools.mcmc.LangevinDynamicsMove( timestep=to_openmm(integrator_settings.timestep), - collision_rate=to_openmm(integrator_settings.collision_rate), - n_steps=integrator_settings.n_steps.m, + collision_rate=to_openmm(integrator_settings.langevin_collision_rate), + n_steps=steps_per_iteration, reassign_velocities=integrator_settings.reassign_velocities, n_restart_attempts=integrator_settings.n_restart_attempts, constraint_tolerance=integrator_settings.constraint_tolerance, @@ -868,30 +883,39 @@ def run(self, *, dry=False, verbose=True, # 12. Create sampler self.logger.info("Creating and setting up the sampler") + rta_its, rta_min_its = settings_validation.convert_real_time_analysis_iterations( + simulation_settings=sampler_settings, + ) + # convert early_termination_target_error from kcal/mol to kT + early_termination_target_error = settings_validation.convert_target_error_from_kcal_per_mole_to_kT( + thermo_settings.temperature, + sampler_settings.early_termination_target_error, + ) + if sampler_settings.sampler_method.lower() == "repex": sampler = _rfe_utils.multistate.HybridRepexSampler( mcmc_moves=integrator, hybrid_factory=hybrid_factory, - online_analysis_interval=sampler_settings.online_analysis_interval, - online_analysis_target_error=sampler_settings.online_analysis_target_error.m, - online_analysis_minimum_iterations=sampler_settings.online_analysis_minimum_iterations + online_analysis_interval=rta_its, + online_analysis_target_error=early_termination_target_error, + online_analysis_minimum_iterations=rta_min_its, ) elif sampler_settings.sampler_method.lower() == "sams": sampler = _rfe_utils.multistate.HybridSAMSSampler( mcmc_moves=integrator, hybrid_factory=hybrid_factory, - online_analysis_interval=sampler_settings.online_analysis_interval, - online_analysis_minimum_iterations=sampler_settings.online_analysis_minimum_iterations, - flatness_criteria=sampler_settings.flatness_criteria, - gamma0=sampler_settings.gamma0, + online_analysis_interval=rta_its, + online_analysis_minimum_iterations=rta_min_its, + flatness_criteria=sampler_settings.sams_flatness_criteria, + gamma0=sampler_settings.sams_gamma0, ) elif sampler_settings.sampler_method.lower() == 'independent': sampler = _rfe_utils.multistate.HybridMultiStateSampler( mcmc_moves=integrator, hybrid_factory=hybrid_factory, - online_analysis_interval=sampler_settings.online_analysis_interval, - online_analysis_target_error=sampler_settings.online_analysis_target_error.m, - online_analysis_minimum_iterations=sampler_settings.online_analysis_minimum_iterations + online_analysis_interval=rta_its, + online_analysis_target_error=early_termination_target_error, + online_analysis_minimum_iterations=rta_min_its, ) else: @@ -902,7 +926,7 @@ def run(self, *, dry=False, verbose=True, reporter=reporter, lambda_protocol=lambdas, temperature=to_openmm(thermo_settings.temperature), - endstates=alchem_settings.unsampled_endstates, + endstates=alchem_settings.endstate_dispersion_correction, minimization_platform=platform.getName(), ) @@ -924,14 +948,14 @@ def run(self, *, dry=False, verbose=True, if verbose: self.logger.info("Running minimization") - sampler.minimize(max_iterations=sim_settings.minimization_steps) + sampler.minimize(max_iterations=sampler_settings.minimization_steps) # equilibrate if verbose: self.logger.info("Running equilibration phase") sampler.equilibrate( - int(equil_steps / integrator_settings.n_steps.m) # type: ignore + int(equil_steps / steps_per_iteration) # type: ignore ) # production @@ -939,7 +963,7 @@ def run(self, *, dry=False, verbose=True, self.logger.info("Running production phase") sampler.extend( - int(prod_steps / integrator_settings.n_steps.m) # type: ignore + int(prod_steps / steps_per_iteration) # type: ignore ) self.logger.info("Production phase complete") @@ -957,8 +981,8 @@ def run(self, *, dry=False, verbose=True, else: # clean up the reporter file - fns = [shared_basepath / sim_settings.output_filename, - shared_basepath / sim_settings.checkpoint_storage] + fns = [shared_basepath / output_settings.output_filename, + shared_basepath / output_settings.checkpoint_storage_filename] for fn in fns: os.remove(fn) finally: diff --git a/openfe/protocols/openmm_rfe/equil_rfe_settings.py b/openfe/protocols/openmm_rfe/equil_rfe_settings.py index e71fbf18d..c4d0b11e8 100644 --- a/openfe/protocols/openmm_rfe/equil_rfe_settings.py +++ b/openfe/protocols/openmm_rfe/equil_rfe_settings.py @@ -8,10 +8,9 @@ """ from __future__ import annotations -from typing import Optional +from typing import Optional, Literal from openff.units import unit from openff.models.types import FloatQuantity -import os from gufe.settings import ( Settings, @@ -20,8 +19,9 @@ ThermoSettings, ) from openfe.protocols.openmm_utils.omm_settings import ( - SystemSettings, SolvationSettings, AlchemicalSamplerSettings, - OpenMMEngineSettings, IntegratorSettings, SimulationSettings + SolvationSettings, MultiStateSimulationSettings, + OpenMMEngineSettings, IntegratorSettings, + OutputSettings, ) try: @@ -59,7 +59,7 @@ class Config: This describes the creation of the hybrid system. """ - unsampled_endstates = False + endstate_dispersion_correction = False """ Whether to have extra unsampled endstate windows for long range correction. Default False. @@ -71,18 +71,19 @@ class Config: Whether to use dispersion correction in the hybrid topology state. Default False. """ - softcore_LJ_v2 = True + softcore_LJ: Literal['gapsys', 'beutler'] """ - Whether to use the LJ softcore function as defined by - Gapsys et al. JCTC 2012 Default True. + Whether to use the LJ softcore function as defined by Gapsys et al. + JCTC 2012, or the one by Beutler et al. Chem. Phys. Lett. 1994. + Default 'gapsys'. """ softcore_alpha = 0.85 """Softcore alpha parameter. Default 0.85""" - interpolate_old_and_new_14s = False + turn_off_core_unique_exceptions = False """ Whether to turn off interactions for new exceptions (not just 1,4s) - at lambda 0 and old exceptions at lambda 1. If False they are present - in the nonbonded force. Default False. + at lambda 0 and old exceptions at lambda 1 between unique atoms and core + atoms. If False they are present in the nonbonded force. Default False. """ explicit_charge_correction = False """ @@ -105,8 +106,19 @@ class Config: class RelativeHybridTopologyProtocolSettings(Settings): - class Config: - arbitrary_types_allowed = True + protocol_repeats: int + """ + The number of completely independent repeats of the entire sampling + process. The mean of the repeats defines the final estimate of FE + difference, while the variance between repeats is used as the uncertainty. + """ + + @validator('protocol_repeats') + def must_be_positive(cls, v): + if v <= 0: + errmsg = f"protocol_repeats must be a positive value, got {v}." + raise ValueError(errmsg) + return v # Inherited things @@ -116,8 +128,6 @@ class Config: """Settings for thermodynamic parameters.""" # Things for creating the systems - system_settings: SystemSettings - """Simulation system settings including the long-range non-bonded method.""" solvation_settings: SolvationSettings """Settings for solvating the system.""" @@ -130,9 +140,9 @@ class Config: """ Alchemical protocol settings including soft core scaling. """ - alchemical_sampler_settings: AlchemicalSamplerSettings + simulation_settings: MultiStateSimulationSettings """ - Settings for sampling alchemical space, including the number of repeats. + Settings for alchemical sampler. """ # MD Engine things @@ -143,8 +153,7 @@ class Config: integrator_settings: IntegratorSettings """Settings for the integrator such as timestep and barostat settings.""" - # Simulation run settings - simulation_settings: SimulationSettings + output_settings: OutputSettings """ - Simulation control settings, including simulation lengths and record-keeping. + Simulation output control settings. """ diff --git a/openfe/protocols/openmm_utils/omm_settings.py b/openfe/protocols/openmm_utils/omm_settings.py index 9bcbcb41b..fdafe4110 100644 --- a/openfe/protocols/openmm_utils/omm_settings.py +++ b/openfe/protocols/openmm_utils/omm_settings.py @@ -12,7 +12,6 @@ from typing import Optional from openff.units import unit from openff.models.types import FloatQuantity -import os from gufe.settings import ( Settings, @@ -28,42 +27,6 @@ from pydantic import validator # type: ignore[assignment] -class SystemSettings(SettingsBaseModel): - """Settings describing the simulation system settings.""" - - class Config: - arbitrary_types_allowed = True - - nonbonded_method = 'PME' - """ - Method for treating nonbonded interactions, currently only PME and - NoCutoff are allowed. Default PME. - """ - nonbonded_cutoff: FloatQuantity['nanometer'] = 1.0 * unit.nanometer - """ - Cutoff value for short range nonbonded interactions. - Default 1.0 * unit.nanometer. - """ - - @validator('nonbonded_method') - def allowed_nonbonded(cls, v): - if v.lower() not in ['pme', 'nocutoff']: - errmsg = ("Only PME and NoCutoff are allowed nonbonded_methods") - raise ValueError(errmsg) - return v - - @validator('nonbonded_cutoff') - def is_positive_distance(cls, v): - # these are time units, not simulation steps - if not v.is_compatible_with(unit.nanometer): - raise ValueError("nonbonded_cutoff must be in distance units " - "(i.e. nanometers)") - if v < 0: - errmsg = "nonbonded_cutoff must be a positive value" - raise ValueError(errmsg) - return v - - class SolvationSettings(SettingsBaseModel): """Settings for solvating the system @@ -106,120 +69,6 @@ def is_positive_distance(cls, v): return v -class AlchemicalSamplerSettings(SettingsBaseModel): - """Settings for the Equilibrium Alchemical sampler, currently supporting - either MultistateSampler, SAMSSampler or ReplicaExchangeSampler. - - """ - - """ - TODO - ---- - * It'd be great if we could pass in the sampler object rather than using - strings to define which one we want. - * Make n_replicas optional such that: If `None` or greater than the number - of lambda windows set in :class:`AlchemicalSettings`, this will default - to the number of lambda windows. If less than the number of lambda - windows, the replica lambda states will be picked at equidistant - intervals along the lambda schedule. - """ - class Config: - arbitrary_types_allowed = True - - sampler_method = "repex" - """ - Alchemical sampling method, must be one of; - `repex` (Hamiltonian Replica Exchange), - `sams` (Self-Adjusted Mixture Sampling), - or `independent` (independently sampled lambda windows). - Default `repex`. - """ - online_analysis_interval: Optional[int] = 250 - """ - MCMC steps (i.e. ``IntegratorSettings.n_steps``) interval at which - to perform an analysis of the free energies. - - At each interval, real time analysis data (e.g. current free energy - estimate and timing data) will be written to a yaml file named - ``_real_time_analysis.yaml``. The - current error in the estimate will also be assed and if it drops - below ``AlchemicalSamplerSettings.online_analysis_target_error`` - the simulation will be terminated. - - If ``None``, no real time analysis will be performed and the yaml - file will not be written. - - Must be a multiple of ``SimulationSettings.checkpoint_interval`` - - Default `250`. - """ - online_analysis_target_error: FloatQuantity = 0.0 * unit.boltzmann_constant * unit.kelvin - """ - Target error for the online analysis measured in kT. Once the free energy - is at or below this value, the simulation will be considered complete. A - suggested value of 0.2 * `unit.boltzmann_constant` * `unit.kelvin` has - shown to be effective in both hydration and binding free energy benchmarks. - Default 0.0 * `unit.boltzmann_constant` * `unit.kelvin`, i.e. no early - termination will occur. - """ - online_analysis_minimum_iterations = 500 - """ - Number of iterations which must pass before online analysis is - carried out. Default 500. - """ - n_repeats: int = 3 - """ - Number of independent repeats to run. Default 3 - """ - flatness_criteria = 'logZ-flatness' - """ - SAMS only. Method for assessing when to switch to asymptomatically - optimal scheme. - One of ['logZ-flatness', 'minimum-visits', 'histogram-flatness']. - Default 'logZ-flatness'. - """ - gamma0 = 1.0 - """SAMS only. Initial weight adaptation rate. Default 1.0.""" - n_replicas = 11 - """Number of replicas to use. Default 11.""" - - @validator('flatness_criteria') - def supported_flatness(cls, v): - supported = [ - 'logz-flatness', 'minimum-visits', 'histogram-flatness' - ] - if v.lower() not in supported: - errmsg = ("Only the following flatness_criteria are " - f"supported: {supported}") - raise ValueError(errmsg) - return v - - @validator('sampler_method') - def supported_sampler(cls, v): - supported = ['repex', 'sams', 'independent'] - if v.lower() not in supported: - errmsg = ("Only the following sampler_method values are " - f"supported: {supported}") - raise ValueError(errmsg) - return v - - @validator('n_repeats', 'n_replicas') - def must_be_positive(cls, v): - if v <= 0: - errmsg = "n_repeats and n_replicas must be positive values" - raise ValueError(errmsg) - return v - - @validator('online_analysis_target_error', 'n_repeats', - 'online_analysis_minimum_iterations', 'gamma0', 'n_replicas') - def must_be_zero_or_positive(cls, v): - if v < 0: - errmsg = ("Online analysis target error, minimum iteration " - "and SAMS gamm0 must be 0 or positive values.") - raise ValueError(errmsg) - return v - - class OpenMMEngineSettings(SettingsBaseModel): """OpenMM MD engine settings""" @@ -244,13 +93,8 @@ class Config: timestep: FloatQuantity['femtosecond'] = 4 * unit.femtosecond """Size of the simulation timestep. Default 4 * unit.femtosecond.""" - collision_rate: FloatQuantity['1/picosecond'] = 1.0 / unit.picosecond + langevin_collision_rate: FloatQuantity['1/picosecond'] = 1.0 / unit.picosecond """Collision frequency. Default 1.0 / unit.pisecond.""" - n_steps = 250 * unit.timestep # todo: IntQuantity - """ - Number of integration timesteps between each time the MCMC move - is applied. Default 250 * unit.timestep. - """ reassign_velocities = False """ If True, velocities are reassigned from the Maxwell-Boltzmann @@ -268,20 +112,24 @@ class Config: Frequency at which volume scaling changes should be attempted. Default 25 * unit.timestep. """ + remove_com: bool = False + """ + Whether or not to remove the center of mass motion. Default False. + """ - @validator('collision_rate', 'n_restart_attempts') + @validator('langevin_collision_rate', 'n_restart_attempts') def must_be_positive_or_zero(cls, v): if v < 0: - errmsg = ("collision_rate, and n_restart_attempts must be " - "zero or positive values") + errmsg = ("langevin_collision_rate, and n_restart_attempts must be" + f" zero or positive values, got {v}.") raise ValueError(errmsg) return v - @validator('timestep', 'n_steps', 'constraint_tolerance') + @validator('timestep', 'constraint_tolerance') def must_be_positive(cls, v): if v <= 0: - errmsg = ("timestep, n_steps, constraint_tolerance " - "must be positive values") + errmsg = ("timestep, and constraint_tolerance " + f"must be positive values, got {v}.") raise ValueError(errmsg) return v @@ -293,40 +141,22 @@ def is_time(cls, v): "(i.e. picoseconds)") return v - @validator('collision_rate') + @validator('langevin_collision_rate') def must_be_inverse_time(cls, v): if not v.is_compatible_with(1 / unit.picosecond): - raise ValueError("collision_rate must be in inverse time " + raise ValueError("langevin collision_rate must be in inverse time " "(i.e. 1/picoseconds)") return v -class SimulationSettings(SettingsBaseModel): +class OutputSettings(SettingsBaseModel): """ - Settings for simulation control, including lengths, + Settings for simulation output settings, writing to disk, etc... """ class Config: arbitrary_types_allowed = True - minimization_steps = 5000 - """Number of minimization steps to perform. Default 5000.""" - equilibration_length: FloatQuantity['nanosecond'] - """ - Length of the equilibration phase in units of time. The total number of - steps from this equilibration length - (i.e. ``equilibration_length`` / :class:`IntegratorSettings.timestep`) - must be a multiple of the value defined for - :class:`IntegratorSettings.n_steps`. - """ - production_length: FloatQuantity['nanosecond'] - """ - Length of the production phase in units of time. The total number of - steps from this production length (i.e. - ``production_length`` / :class:`IntegratorSettings.timestep`) must be - a multiple of the value defined for :class:`IntegratorSettings.nsteps`. - """ - # reporter settings output_filename = 'simulation.nc' """Path to the trajectory storage file. Default 'simulation.nc'.""" @@ -341,14 +171,14 @@ class Config: Selection string for which part of the system to write coordinates for. Default 'not water'. """ - checkpoint_interval = 250 * unit.timestep # todo: Needs IntQuantity + checkpoint_interval: FloatQuantity['picosecond'] = 1 * unit.picosecond """ - Frequency to write the checkpoint file. Default 250 * unit.timestep. + Frequency to write the checkpoint file. Default 1 * unit.picosecond. """ - checkpoint_storage = 'checkpoint.nc' + checkpoint_storage_filename = 'checkpoint.chk' """ Separate filename for the checkpoint file. Note, this should - not be a full path, just a filename. Default 'checkpoint.nc'. + not be a full path, just a filename. Default 'checkpoint.chk'. """ forcefield_cache: Optional[str] = 'db.json' """ @@ -356,6 +186,39 @@ class Config: later reused. """ + @validator('checkpoint_interval') + def must_be_positive(cls, v): + if v <= 0: + errmsg = f"Checkpoint intervals must be positive, got {v}." + raise ValueError(errmsg) + return v + + +class SimulationSettings(SettingsBaseModel): + """ + Settings for simulation control, including lengths, etc... + """ + class Config: + arbitrary_types_allowed = True + + minimization_steps = 5000 + """Number of minimization steps to perform. Default 5000.""" + equilibration_length: FloatQuantity['nanosecond'] + """ + Length of the equilibration phase in units of time. The total number of + steps from this equilibration length + (i.e. ``equilibration_length`` / :class:`IntegratorSettings.timestep`) + must be a multiple of the value defined for + :class:`AlchemicalSamplerSettings.steps_per_iteration`. + """ + production_length: FloatQuantity['nanosecond'] + """ + Length of the production phase in units of time. The total number of + steps from this production length (i.e. + ``production_length`` / :class:`IntegratorSettings.timestep`) must be + a multiple of the value defined for :class:`IntegratorSettings.nsteps`. + """ + @validator('equilibration_length', 'production_length') def is_time(cls, v): # these are time units, not simulation steps @@ -364,19 +227,143 @@ def is_time(cls, v): return v @validator('minimization_steps', 'equilibration_length', - 'production_length', 'checkpoint_interval') + 'production_length') + def must_be_positive(cls, v): + if v <= 0: + errmsg = ("Minimization steps, and MD lengths must be positive, " + f"got {v}") + raise ValueError(errmsg) + return v + + +class MultiStateSimulationSettings(SimulationSettings): + """ + Settings for simulation control for multistate simulations, + including simulation length and details of the alchemical sampler. + """ + + """ + TODO + ---- + * It'd be great if we could pass in the sampler object rather than using + strings to define which one we want. + * Make n_replicas optional such that: If `None` or greater than the number + of lambda windows set in :class:`AlchemicalSettings`, this will default + to the number of lambda windows. If less than the number of lambda + windows, the replica lambda states will be picked at equidistant + intervals along the lambda schedule. + """ + + class Config: + arbitrary_types_allowed = True + + sampler_method = "repex" + """ + Alchemical sampling method, must be one of; + `repex` (Hamiltonian Replica Exchange), + `sams` (Self-Adjusted Mixture Sampling), + or `independent` (independently sampled lambda windows). + Default `repex`. + """ + time_per_iteration: FloatQuantity['picosecond'] = 1 * unit.picosecond + # todo: Add validators in the protocol + """ + Simulation time between each MCMC move attempt. Default 1 * unit.picosecond. + """ + real_time_analysis_interval: Optional[FloatQuantity['picosecond']] = 250 * unit.picosecond + # todo: Add validators in the protocol + """ + Time interval at which to perform an analysis of the free energies. + + At each interval, real time analysis data (e.g. current free energy + estimate and timing data) will be written to a yaml file named + ``_real_time_analysis.yaml``. The + current error in the estimate will also be assessed and if it drops + below ``MultiStateSimulationSettings.early_termination_target_error`` + the simulation will be terminated. + + If ``None``, no real time analysis will be performed and the yaml + file will not be written. + + Must be a multiple of ``OutputSettings.checkpoint_interval`` + + Default `250`. + + """ + early_termination_target_error: Optional[FloatQuantity['kcal/mol']] = 0.0 * unit.kilocalorie_per_mole + # todo: have default ``None`` or ``0.0 * unit.kilocalorie_per_mole`` + # (later would give an example of unit). + """ + Target error for the real time analysis measured in kcal/mol. Once the MBAR + error of the free energy is at or below this value, the simulation will be + considered complete. + A suggested value of 0.12 * `unit.kilocalorie_per_mole` has + shown to be effective in both hydration and binding free energy benchmarks. + Default ``None``, i.e. no early termination will occur. + """ + real_time_analysis_minimum_time: FloatQuantity['picosecond'] = 500 * unit.picosecond + # todo: Add validators in the protocol + """ + Simulation time which must pass before real time analysis is + carried out. + + Default 500 * unit.picosecond. + """ + + sams_flatness_criteria = 'logZ-flatness' + """ + SAMS only. Method for assessing when to switch to asymptomatically + optimal scheme. + One of ['logZ-flatness', 'minimum-visits', 'histogram-flatness']. + Default 'logZ-flatness'. + """ + sams_gamma0 = 1.0 + """SAMS only. Initial weight adaptation rate. Default 1.0.""" + n_replicas = 11 + """Number of replicas to use. Default 11.""" + + @validator('sams_flatness_criteria') + def supported_flatness(cls, v): + supported = [ + 'logz-flatness', 'minimum-visits', 'histogram-flatness' + ] + if v.lower() not in supported: + errmsg = ("Only the following sams_flatness_criteria are " + f"supported: {supported}") + raise ValueError(errmsg) + return v + + @validator('sampler_method') + def supported_sampler(cls, v): + supported = ['repex', 'sams', 'independent'] + if v.lower() not in supported: + errmsg = ("Only the following sampler_method values are " + f"supported: {supported}") + raise ValueError(errmsg) + return v + + @validator('n_replicas', 'time_per_iteration') def must_be_positive(cls, v): if v <= 0: - errmsg = ("Minimization steps, MD lengths, and checkpoint " - "intervals must be positive") + errmsg = "n_replicas and steps_per_iteration must be positive " \ + f"values, got {v}." + raise ValueError(errmsg) + return v + + @validator('early_termination_target_error', + 'real_time_analysis_minimum_time', 'sams_gamma0', + 'n_replicas') + def must_be_zero_or_positive(cls, v): + if v < 0: + errmsg = ("Early termination target error, minimum iteration and" + f" SAMS gamma0 must be 0 or positive values, got {v}.") raise ValueError(errmsg) return v -class SimulationSettingsMD(SimulationSettings): +class MDSimulationSettings(SimulationSettings): """ - Settings for simulation control for plain MD simulations, including - writing outputs + Settings for simulation control for plain MD simulations """ class Config: arbitrary_types_allowed = True @@ -385,10 +372,15 @@ class Config: """ Length of the equilibration phase in the NVT ensemble in units of time. The total number of steps from this equilibration length - (i.e. ``equilibration_length_nvt`` / :class:`IntegratorSettings.timestep`) - must be a multiple of the value defined for - :class:`IntegratorSettings.n_steps`. + (i.e. ``equilibration_length_nvt`` / :class:`IntegratorSettings.timestep`). """ + + +class MDOutputSettings(OutputSettings): + """ Settings for simulation output settings for plain MD simulations.""" + class Config: + arbitrary_types_allowed = True + # reporter settings production_trajectory_filename = 'simulation.xtc' """Path to the storage file for analysis. Default 'simulation.xtc'.""" @@ -397,7 +389,8 @@ class Config: Frequency to write the xtc file. Default 5000 * unit.timestep. """ preminimized_structure = 'system.pdb' - """Path to the pdb file of the full pre-minimized system. Default 'system.pdb'.""" + """Path to the pdb file of the full pre-minimized system. + Default 'system.pdb'.""" minimized_structure = 'minimized.pdb' """Path to the pdb file of the system after minimization. Only the specified atom subset is saved. Default 'minimized.pdb'.""" @@ -407,11 +400,6 @@ class Config: equil_NPT_structure = 'equil_NPT.pdb' """Path to the pdb file of the system after NPT equilibration. Only the specified atom subset is saved. Default 'equil_NPT.pdb'.""" - checkpoint_storage_filename = 'checkpoint.chk' - """ - Separate filename for the checkpoint file. Note, this should - not be a full path, just a filename. Default 'checkpoint.chk'. - """ log_output = 'simulation.log' """ Filename for writing the log of the MD simulation, including timesteps, diff --git a/openfe/protocols/openmm_utils/settings_validation.py b/openfe/protocols/openmm_utils/settings_validation.py index 976bc4963..2ebf6187e 100644 --- a/openfe/protocols/openmm_utils/settings_validation.py +++ b/openfe/protocols/openmm_utils/settings_validation.py @@ -4,11 +4,11 @@ Reusable utility methods to validate input settings to OpenMM-based alchemical Protocols. """ -import warnings -from typing import Dict, List, Tuple from openff.units import unit -from gufe import ( - Component, ChemicalSystem, SolventComponent, ProteinComponent +from typing import Optional +from .omm_settings import ( + IntegratorSettings, + MultiStateSimulationSettings, ) @@ -71,3 +71,110 @@ def get_simsteps(sim_length: unit.Quantity, raise ValueError(errmsg) return sim_steps + + +def convert_steps_per_iteration( + simulation_settings: MultiStateSimulationSettings, + integrator_settings: IntegratorSettings, +) -> int: + """Convert time per iteration to steps + + Parameters + ---------- + simulation_settings: MultiStateSimulationSettings + integrator_settings: IntegratorSettings + + Returns + ------- + steps_per_iteration : int + suitable for input to Integrator + """ + tpi_fs = round(simulation_settings.time_per_iteration.to(unit.attosecond).m) + ts_fs = round(integrator_settings.timestep.to(unit.attosecond).m) + steps_per_iteration, rem = divmod(tpi_fs, ts_fs) + + if rem: + raise ValueError(f"time_per_iteration ({simulation_settings.time_per_iteration}) " + f"not divisible by timestep ({integrator_settings.timestep})") + + return steps_per_iteration + + +def convert_real_time_analysis_iterations( + simulation_settings: MultiStateSimulationSettings, +) -> tuple[Optional[int], Optional[int]]: + """Convert time units in Settings to various other units + + Interally openmmtools uses various quantities with units of time, + steps, and iterations. + + Our Settings objects instead have things defined in time (fs or ps). + + This function generates suitable inputs for the openmmtools objects + + Parameters + ---------- + simulation_settings: MultiStateSimulationSettings + + Returns + ------- + real_time_analysis_iterations : Optional[int] + suitable for input to online_analysis_interval + real_time_analysis_minimum_iterations : Optional[int] + suitable for input to real_time_analysis_minimum_iterations + """ + if simulation_settings.real_time_analysis_interval is None: + # option to turn off real time analysis + return None, None + + tpi_fs = round(simulation_settings.time_per_iteration.to(unit.attosecond).m) + + # convert real_time_analysis time to interval + # rta_its must be number of MCMC iterations + # i.e. rta_fs / tpi_fs -> number of iterations + rta_fs = round(simulation_settings.real_time_analysis_interval.to(unit.attosecond).m) + + rta_its, rem = divmod(rta_fs, tpi_fs) + if rem: + raise ValueError(f"real_time_analysis_interval ({simulation_settings.real_time_analysis_interval}) " + f"is not divisible by time_per_iteration ({simulation_settings.time_per_iteration})") + + # convert RTA_minimum_time to iterations + rta_min_fs = round(simulation_settings.real_time_analysis_minimum_time.to(unit.attosecond).m) + rta_min_its, rem = divmod(rta_min_fs, tpi_fs) + if rem: + raise ValueError(f"real_time_analysis_minimum_time ({simulation_settings.real_time_analysis_minimum_time}) " + f"is not divisible by time_per_iteration ({simulation_settings.time_per_iteration})") + + return rta_its, rta_min_its + + +def convert_target_error_from_kcal_per_mole_to_kT( + temperature, + target_error, +) -> float: + """Convert kcal/mol target error to kT units + + If target_error is 0.0, returns 0.0 + + Parameters + ---------- + temperature : unit.Quantity + temperature in K + target_error : unit.Quantity + error in kcal/mol + + Returns + ------- + early_termination_target_error : float + in units of kT, suitable for input as "online_analysis_target_error" in a + Sampler + """ + if target_error: + kB = 0.001987204 * unit.kilocalorie_per_mole / unit.kelvin + kT = temperature * kB + early_termination_target_error = target_error / kT + else: + return 0.0 + + return early_termination_target_error.m diff --git a/openfe/protocols/openmm_utils/system_creation.py b/openfe/protocols/openmm_utils/system_creation.py index fa28ba8d0..e67a6862c 100644 --- a/openfe/protocols/openmm_utils/system_creation.py +++ b/openfe/protocols/openmm_utils/system_creation.py @@ -18,8 +18,7 @@ Component, ProteinComponent, SolventComponent, SmallMoleculeComponent ) from ..openmm_rfe.equil_rfe_settings import ( - SystemSettings, SimulationSettings, SolvationSettings, - IntegratorSettings, + SolvationSettings, IntegratorSettings, ) @@ -27,7 +26,6 @@ def get_system_generator( forcefield_settings: OpenMMSystemGeneratorFFSettings, thermo_settings: ThermoSettings, integrator_settings: IntegratorSettings, - system_settings: SystemSettings, cache: Optional[Path], has_solvent: bool, ) -> SystemGenerator: @@ -38,16 +36,15 @@ def get_system_generator( --------- forcefield_settings : OpenMMSystemGeneratorFFSettings Force field settings, including necessary information - for constraints, hydrogen mass, rigid waters, COM removal, + for constraints, hydrogen mass, rigid waters, non-ligand FF xmls, and the ligand FF name. + integrator_settings: IntegratorSettings + Integrator settings, including COM removal. thermo_settings : ThermoSettings Thermodynamic settings, including necessary settings for defining the ensemble conditions. integrator_settings : IntegratorSettings Integrator settings, including barostat control variables. - system_settings : SystemSettings - System settings including all necessary information for - the nonbonded methods. cache : Optional[pathlib.Path] Path to openff force field cache. has_solvent : bool @@ -76,7 +73,7 @@ def get_system_generator( forcefield_kwargs = { 'constraints': constraints, 'rigidWater': forcefield_settings.rigid_water, - 'removeCMMotion': forcefield_settings.remove_com, + 'removeCMMotion': integrator_settings.remove_com, 'hydrogenMass': forcefield_settings.hydrogen_mass * omm_unit.amu, } @@ -87,10 +84,10 @@ def get_system_generator( 'cutoffnonperiodic': app.CutoffNonPeriodic, 'cutoffperiodic': app.CutoffPeriodic, 'ewald': app.Ewald - }[system_settings.nonbonded_method.lower()] + }[forcefield_settings.nonbonded_method.lower()] nonbonded_cutoff = to_openmm( - system_settings.nonbonded_cutoff + forcefield_settings.nonbonded_cutoff, # type: ignore ) # create the periodic_kwarg entry diff --git a/openfe/setup/alchemical_network_planner/relative_alchemical_network_planner.py b/openfe/setup/alchemical_network_planner/relative_alchemical_network_planner.py index 4299de44e..cda2cb720 100644 --- a/openfe/setup/alchemical_network_planner/relative_alchemical_network_planner.py +++ b/openfe/setup/alchemical_network_planner/relative_alchemical_network_planner.py @@ -209,7 +209,7 @@ def _build_transformation( # Todo: Another dirty hack! - START protocol_settings = transformation_protocol.settings.unfrozen_copy() if "vacuum" in transformation_name: - protocol_settings.system_settings.nonbonded_method = "nocutoff" + protocol_settings.forcefield_settings.nonbonded_method = "nocutoff" transformation_protocol = transformation_protocol.__class__( settings=protocol_settings diff --git a/openfe/tests/data/openmm_afe/AHFEProtocol_json_results.gz b/openfe/tests/data/openmm_afe/AHFEProtocol_json_results.gz index 1cc67f05a..79640ee29 100644 Binary files a/openfe/tests/data/openmm_afe/AHFEProtocol_json_results.gz and b/openfe/tests/data/openmm_afe/AHFEProtocol_json_results.gz differ diff --git a/openfe/tests/data/openmm_md/MDProtocol_json_results.gz b/openfe/tests/data/openmm_md/MDProtocol_json_results.gz index 44d3bc24b..fda867320 100644 Binary files a/openfe/tests/data/openmm_md/MDProtocol_json_results.gz and b/openfe/tests/data/openmm_md/MDProtocol_json_results.gz differ diff --git a/openfe/tests/data/openmm_rfe/RHFEProtocol_json_results.gz b/openfe/tests/data/openmm_rfe/RHFEProtocol_json_results.gz index 09fc1c8c7..43b064681 100644 Binary files a/openfe/tests/data/openmm_rfe/RHFEProtocol_json_results.gz and b/openfe/tests/data/openmm_rfe/RHFEProtocol_json_results.gz differ diff --git a/openfe/tests/protocols/conftest.py b/openfe/tests/protocols/conftest.py index fcaf13fbc..67cc21d1c 100644 --- a/openfe/tests/protocols/conftest.py +++ b/openfe/tests/protocols/conftest.py @@ -195,7 +195,10 @@ def toluene_many_solv_system(benzene_modifications): @pytest.fixture def rfe_transformation_json() -> str: - """string of a RFE results similar to quickrun""" + """string of a RFE results similar to quickrun + + generated with gen-serialized-results.py + """ d = resources.files('openfe.tests.data.openmm_rfe') with gzip.open((d / 'RHFEProtocol_json_results.gz').as_posix(), 'r') as f: # type: ignore @@ -206,6 +209,8 @@ def rfe_transformation_json() -> str: def afe_solv_transformation_json() -> str: """ string of a Absolute Solvation result (CN in water) generated by quickrun + + generated with gen-serialized-results.py """ d = resources.files('openfe.tests.data.openmm_afe') fname = "AHFEProtocol_json_results.gz" @@ -218,6 +223,8 @@ def afe_solv_transformation_json() -> str: def md_json() -> str: """ string of a MD result (TYK ligand lig_ejm_31 in water) generated by quickrun + + generated with gen-serialized-results.py """ d = resources.files('openfe.tests.data.openmm_md') fname = "MDProtocol_json_results.gz" diff --git a/openfe/tests/protocols/test_openmm_afe_slow.py b/openfe/tests/protocols/test_openmm_afe_slow.py index bf2ce4f73..8dfad8fa8 100644 --- a/openfe/tests/protocols/test_openmm_afe_slow.py +++ b/openfe/tests/protocols/test_openmm_afe_slow.py @@ -46,16 +46,16 @@ def test_openmm_run_engine(platform, # Run a really short calculation to check everything is going well s = openmm_afe.AbsoluteSolvationProtocol.default_settings() s.alchemsampler_settings.n_repeats = 1 - s.solvent_simulation_settings.output_indices = "resname UNK" + s.solvent_output_settings.output_indices = "resname UNK" s.vacuum_simulation_settings.equilibration_length = 0.1 * unit.picosecond s.vacuum_simulation_settings.production_length = 0.1 * unit.picosecond s.solvent_simulation_settings.equilibration_length = 0.1 * unit.picosecond s.solvent_simulation_settings.production_length = 0.1 * unit.picosecond s.vacuum_engine_settings.compute_platform = platform s.solvent_engine_settings.compute_platform = platform - s.integrator_settings.n_steps = 5 * unit.timestep - s.vacuum_simulation_settings.checkpoint_interval = 5 * unit.timestep - s.solvent_simulation_settings.checkpoint_interval = 5 * unit.timestep + s.alchemsampler_settings.steps_per_iteration = 5 * unit.timestep + s.vacuum_output_settings.checkpoint_interval = 5 * unit.timestep + s.solvent_output_settings.checkpoint_interval = 5 * unit.timestep s.alchemsampler_settings.n_replicas = 20 s.lambda_settings.lambda_elec = \ [0.0, 0.25, 0.5, 0.75, 1.0, 1.0, 1.0, 1.0, 1.0, 1.0, diff --git a/openfe/tests/protocols/test_openmm_afe_solvation_protocol.py b/openfe/tests/protocols/test_openmm_afe_solvation_protocol.py index 1c821742f..a81c350b3 100644 --- a/openfe/tests/protocols/test_openmm_afe_solvation_protocol.py +++ b/openfe/tests/protocols/test_openmm_afe_solvation_protocol.py @@ -70,11 +70,17 @@ def test_validate_lambda_schedule_naked_charge(val, default_settings): default_settings.lambda_settings.lambda_elec = val['elec'] default_settings.lambda_settings.lambda_vdw = val['vdw'] default_settings.lambda_settings.lambda_restraints = val['restraints'] - default_settings.alchemsampler_settings.n_replicas = 2 + default_settings.vacuum_simulation_settings.n_replicas = 2 + default_settings.solvent_simulation_settings.n_replicas = 2 with pytest.raises(ValueError, match=errmsg): AbsoluteSolvationProtocol._validate_lambda_schedule( default_settings.lambda_settings, - default_settings.alchemsampler_settings, + default_settings.vacuum_simulation_settings, + ) + with pytest.raises(ValueError, match=errmsg): + AbsoluteSolvationProtocol._validate_lambda_schedule( + default_settings.lambda_settings, + default_settings.solvent_simulation_settings, ) @@ -86,13 +92,13 @@ def test_validate_lambda_schedule_nreplicas(val, default_settings): default_settings.lambda_settings.lambda_vdw = val['vdw'] default_settings.lambda_settings.lambda_restraints = val['restraints'] n_replicas = 3 - default_settings.alchemsampler_settings.n_replicas = n_replicas + default_settings.vacuum_simulation_settings.n_replicas = n_replicas errmsg = (f"Number of replicas {n_replicas} does not equal the" f" number of lambda windows {len(val['vdw'])}") with pytest.raises(ValueError, match=errmsg): AbsoluteSolvationProtocol._validate_lambda_schedule( default_settings.lambda_settings, - default_settings.alchemsampler_settings, + default_settings.vacuum_simulation_settings, ) @@ -104,7 +110,7 @@ def test_validate_lambda_schedule_nwindows(val, default_settings): default_settings.lambda_settings.lambda_vdw = val['vdw'] default_settings.lambda_settings.lambda_restraints = val['restraints'] n_replicas = 3 - default_settings.alchemsampler_settings.n_replicas = n_replicas + default_settings.vacuum_simulation_settings.n_replicas = n_replicas errmsg = ( "Components elec and vdw must have equal amount" f" of lambda windows. Got {len(val['elec'])} elec lambda" @@ -112,7 +118,7 @@ def test_validate_lambda_schedule_nwindows(val, default_settings): with pytest.raises(ValueError, match=errmsg): AbsoluteSolvationProtocol._validate_lambda_schedule( default_settings.lambda_settings, - default_settings.alchemsampler_settings, + default_settings.vacuum_simulation_settings, ) @@ -126,11 +132,11 @@ def test_validate_lambda_schedule_nonzero_restraints(val, default_settings): default_settings.lambda_settings.lambda_elec = val['elec'] default_settings.lambda_settings.lambda_vdw = val['vdw'] default_settings.lambda_settings.lambda_restraints = val['restraints'] - default_settings.alchemsampler_settings.n_replicas = 2 + default_settings.vacuum_simulation_settings.n_replicas = 2 with pytest.warns(UserWarning, match=wmsg): AbsoluteSolvationProtocol._validate_lambda_schedule( default_settings.lambda_settings, - default_settings.alchemsampler_settings, + default_settings.vacuum_simulation_settings, ) @@ -206,7 +212,7 @@ def test_validate_solvent_endstates_nosolvcomp_stateB( with pytest.raises( ValueError, match="No SolventComponent found in stateB" ): - comps = AbsoluteSolvationProtocol._validate_solvent_endstates(stateA, stateB) + AbsoluteSolvationProtocol._validate_solvent_endstates(stateA, stateB) def test_validate_alchem_comps_appearingB(benzene_modifications): @@ -262,7 +268,7 @@ def test_validate_alchem_nonsmc(benzene_modifications): def test_vac_bad_nonbonded(benzene_modifications): settings = openmm_afe.AbsoluteSolvationProtocol.default_settings() - settings.vacuum_system_settings.nonbonded_method = 'pme' + settings.vacuum_forcefield_settings.nonbonded_method = 'pme' protocol = openmm_afe.AbsoluteSolvationProtocol(settings=settings) stateA = ChemicalSystem({ @@ -284,8 +290,8 @@ def test_vac_bad_nonbonded(benzene_modifications): def test_dry_run_vac_benzene(benzene_modifications, method, tmpdir): s = openmm_afe.AbsoluteSolvationProtocol.default_settings() - s.alchemsampler_settings.n_repeats = 1 - s.alchemsampler_settings.sampler_method = method + s.protocol_repeats = 1 + s.vacuum_simulation_settings.sampler_method = method protocol = openmm_afe.AbsoluteSolvationProtocol( settings=s, @@ -327,7 +333,7 @@ def test_dry_run_vac_benzene(benzene_modifications, def test_confgen_fail_AFE(benzene_modifications, tmpdir): # check system parametrisation works even if confgen fails s = openmm_afe.AbsoluteSolvationProtocol.default_settings() - s.alchemsampler_settings.n_repeats = 1 + s.protocol_repeats = 1 protocol = openmm_afe.AbsoluteSolvationProtocol( settings=s, @@ -362,8 +368,8 @@ def test_confgen_fail_AFE(benzene_modifications, tmpdir): def test_dry_run_solv_benzene(benzene_modifications, tmpdir): s = openmm_afe.AbsoluteSolvationProtocol.default_settings() - s.alchemsampler_settings.n_repeats = 1 - s.solvent_simulation_settings.output_indices = "resname UNK" + s.protocol_repeats = 1 + s.solvent_output_settings.output_indices = "resname UNK" protocol = openmm_afe.AbsoluteSolvationProtocol( settings=s, @@ -407,8 +413,13 @@ def test_dry_run_solv_benzene(benzene_modifications, tmpdir): def test_dry_run_solv_benzene_tip4p(benzene_modifications, tmpdir): s = AbsoluteSolvationProtocol.default_settings() - s.alchemsampler_settings.n_repeats = 1 - s.forcefield_settings.forcefields = [ + s.protocol_repeats = 1 + s.vacuum_forcefield_settings.forcefields = [ + "amber/ff14SB.xml", # ff14SB protein force field + "amber/tip4pew_standard.xml", # FF we are testsing with the fun VS + "amber/phosaa10.xml", # Handles THE TPO + ] + s.solvent_forcefield_settings.forcefields = [ "amber/ff14SB.xml", # ff14SB protein force field "amber/tip4pew_standard.xml", # FF we are testsing with the fun VS "amber/phosaa10.xml", # Handles THE TPO @@ -453,7 +464,7 @@ def test_dry_run_solv_user_charges_benzene(benzene_modifications, tmpdir): alchemical system. """ s = openmm_afe.AbsoluteSolvationProtocol.default_settings() - s.alchemsampler_settings.n_repeats = 1 + s.protocol_repeats = 1 protocol = openmm_afe.AbsoluteSolvationProtocol( settings=s, @@ -534,8 +545,9 @@ def assign_fictitious_charges(offmol): def test_high_timestep(benzene_modifications, tmpdir): s = AbsoluteSolvationProtocol.default_settings() - s.alchemsampler_settings.n_repeats = 1 - s.forcefield_settings.hydrogen_mass = 1.0 + s.protocol_repeats = 1 + s.solvent_forcefield_settings.hydrogen_mass = 1.0 + s.vacuum_forcefield_settings.hydrogen_mass = 1.0 protocol = AbsoluteSolvationProtocol( settings=s, diff --git a/openfe/tests/protocols/test_openmm_equil_rfe_protocols.py b/openfe/tests/protocols/test_openmm_equil_rfe_protocols.py index 6258b6f35..869672efe 100644 --- a/openfe/tests/protocols/test_openmm_equil_rfe_protocols.py +++ b/openfe/tests/protocols/test_openmm_equil_rfe_protocols.py @@ -66,7 +66,7 @@ def test_append_topology(benzene_complex_system, toluene_complex_system): ) assert len(list(top2.atoms())) == 2625 + 3 # added methyl - assert len(list(top2.bonds())) == 2645 + 4 - 1 # add methyl bonds, minus hydrogen + assert len(list(top2.bonds())) == 2645 + 4 - 1 # add methyl bonds, minus hydrogen assert appended_resids[0] == len(list(top1.residues())) - 1 @@ -95,7 +95,7 @@ def test_append_topology_no_exclude(benzene_complex_system, ) assert len(list(top2.atoms())) == 2625 + 15 # added toluene - assert len(list(top2.bonds())) == 2645 + 15 # 15 bonds in toluene + assert len(list(top2.bonds())) == 2645 + 15 # 15 bonds in toluene assert appended_resids[0] == len(list(top1.residues())) @@ -167,7 +167,7 @@ def test_validate_alchemical_components_wrong_mappings(mapping): def test_validate_alchemical_components_missing_alchem_comp( benzene_to_toluene_mapping): - alchem_comps = {'stateA': [openfe.SolventComponent(),], 'stateB': []} + alchem_comps = {'stateA': [openfe.SolventComponent(), ], 'stateB': []} with pytest.raises(ValueError, match="Unmapped alchemical component"): _validate_alchemical_components( alchem_comps, {'ligand': benzene_to_toluene_mapping}, @@ -181,9 +181,9 @@ def test_dry_run_default_vacuum(benzene_vacuum_system, toluene_vacuum_system, benzene_to_toluene_mapping, method, tmpdir): vac_settings = openmm_rfe.RelativeHybridTopologyProtocol.default_settings() - vac_settings.system_settings.nonbonded_method = 'nocutoff' - vac_settings.alchemical_sampler_settings.sampler_method = method - vac_settings.alchemical_sampler_settings.n_repeats = 1 + vac_settings.forcefield_settings.nonbonded_method = 'nocutoff' + vac_settings.simulation_settings.sampler_method = method + vac_settings.protocol_repeats = 1 protocol = openmm_rfe.RelativeHybridTopologyProtocol( settings=vac_settings, @@ -226,7 +226,7 @@ def test_dry_run_default_vacuum(benzene_vacuum_system, toluene_vacuum_system, def test_dry_run_gaff_vacuum(benzene_vacuum_system, toluene_vacuum_system, benzene_to_toluene_mapping, tmpdir): vac_settings = openmm_rfe.RelativeHybridTopologyProtocol.default_settings() - vac_settings.system_settings.nonbonded_method = 'nocutoff' + vac_settings.forcefield_settings.nonbonded_method = 'nocutoff' vac_settings.forcefield_settings.small_molecule_forcefield = 'gaff-2.11' protocol = openmm_rfe.RelativeHybridTopologyProtocol( @@ -348,15 +348,15 @@ def test_dry_core_element_change(tmpdir): ) settings = openmm_rfe.RelativeHybridTopologyProtocol.default_settings() - settings.system_settings.nonbonded_method = 'nocutoff' + settings.forcefield_settings.nonbonded_method = 'nocutoff' protocol = openmm_rfe.RelativeHybridTopologyProtocol( settings=settings, ) dag = protocol.create( - stateA=openfe.ChemicalSystem({'ligand': benz,}), - stateB=openfe.ChemicalSystem({'ligand': pyr,}), + stateA=openfe.ChemicalSystem({'ligand': benz, }), + stateB=openfe.ChemicalSystem({'ligand': pyr, }), mapping={'whatamapping': mapping}, ) @@ -383,9 +383,9 @@ def test_dry_run_ligand(benzene_system, toluene_system, benzene_to_toluene_mapping, method, tmpdir): # this might be a bit time consuming settings = openmm_rfe.RelativeHybridTopologyProtocol.default_settings() - settings.alchemical_sampler_settings.sampler_method = method - settings.alchemical_sampler_settings.n_repeats = 1 - settings.simulation_settings.output_indices = 'resname UNK' + settings.simulation_settings.sampler_method = method + settings.protocol_repeats = 1 + settings.output_settings.output_indices = 'resname UNK' protocol = openmm_rfe.RelativeHybridTopologyProtocol( settings=settings, @@ -416,7 +416,7 @@ def test_confgen_mocked_fail(benzene_system, toluene_system, Check that even if conformer generation fails, we can still perform a sim """ settings = openmm_rfe.RelativeHybridTopologyProtocol.default_settings() - settings.alchemical_sampler_settings.n_repeats = 1 + settings.protocol_repeats = 1 protocol = openmm_rfe.RelativeHybridTopologyProtocol(settings=settings) @@ -442,11 +442,11 @@ def tip4p_hybrid_factory( settings = openmm_rfe.RelativeHybridTopologyProtocol.default_settings() settings.forcefield_settings.forcefields = [ "amber/ff14SB.xml", # ff14SB protein force field - "amber/tip4pew_standard.xml", # FF we are testsing with the fun VS + "amber/tip4pew_standard.xml", # FF we are testsing with the fun VS "amber/phosaa10.xml", # Handles THE TPO ] settings.solvation_settings.solvent_padding = 1.0 * unit.nanometer - settings.system_settings.nonbonded_cutoff = 0.9 * unit.nanometer + settings.forcefield_settings.nonbonded_cutoff = 0.9 * unit.nanometer settings.solvation_settings.solvent_model = 'tip4pew' settings.integrator_settings.reassign_velocities = True @@ -594,8 +594,8 @@ def test_dry_run_user_charges(benzene_modifications, tmpdir): hybrid topology. """ vac_settings = openmm_rfe.RelativeHybridTopologyProtocol.default_settings() - vac_settings.system_settings.nonbonded_method = 'nocutoff' - vac_settings.alchemical_sampler_settings.n_repeats = 1 + vac_settings.forcefield_settings.nonbonded_method = 'nocutoff' + vac_settings.protocol_repeats = 1 protocol = openmm_rfe.RelativeHybridTopologyProtocol( settings=vac_settings, @@ -641,8 +641,8 @@ def check_propchgs(smc, charge_array): # create DAG from protocol and take first (and only) work unit from within dag = protocol.create( - stateA=openfe.ChemicalSystem({'l': benzene_smc,}), - stateB=openfe.ChemicalSystem({'l': toluene_smc,}), + stateA=openfe.ChemicalSystem({'l': benzene_smc, }), + stateB=openfe.ChemicalSystem({'l': toluene_smc, }), mapping={'ligand': mapping}, ) dag_unit = list(dag.protocol_units)[0] @@ -723,11 +723,11 @@ def test_virtual_sites_no_reassign(benzene_system, toluene_system, settings = openmm_rfe.RelativeHybridTopologyProtocol.default_settings() settings.forcefield_settings.forcefields = [ "amber/ff14SB.xml", # ff14SB protein force field - "amber/tip4pew_standard.xml", # FF we are testsing with the fun VS + "amber/tip4pew_standard.xml", # FF we are testsing with the fun VS "amber/phosaa10.xml", # Handles THE TPO ] settings.solvation_settings.solvent_padding = 1.0 * unit.nanometer - settings.system_settings.nonbonded_cutoff = 0.9 * unit.nanometer + settings.forcefield_settings.nonbonded_cutoff = 0.9 * unit.nanometer settings.solvation_settings.solvent_model = 'tip4pew' settings.integrator_settings.reassign_velocities = False @@ -754,8 +754,8 @@ def test_dry_run_complex(benzene_complex_system, toluene_complex_system, # this will be very time consuming settings = openmm_rfe.RelativeHybridTopologyProtocol.default_settings() settings.alchemical_sampler_settings.sampler_method = method - settings.alchemical_sampler_settings.n_repeats = 1 - settings.simulation_settings.output_indices = 'protein or resname UNK' + settings.protocol_repeats = 1 + settings.output_settings.output_indices = 'protein or resname UNK' protocol = openmm_rfe.RelativeHybridTopologyProtocol( settings=settings, @@ -797,7 +797,7 @@ def test_hightimestep(benzene_vacuum_system, benzene_to_toluene_mapping, tmpdir): settings = openmm_rfe.RelativeHybridTopologyProtocol.default_settings() settings.forcefield_settings.hydrogen_mass = 1.0 - settings.system_settings.nonbonded_method = 'nocutoff' + settings.forcefield_settings.nonbonded_method = 'nocutoff' p = openmm_rfe.RelativeHybridTopologyProtocol( settings=settings, @@ -823,8 +823,8 @@ def test_n_replicas_not_n_windows(benzene_vacuum_system, # equals the numbers of replicas used - TODO: remove limitation settings = openmm_rfe.RelativeHybridTopologyProtocol.default_settings() # default lambda windows is 11 - settings.alchemical_sampler_settings.n_replicas = 13 - settings.system_settings.nonbonded_method = 'nocutoff' + settings.simulation_settings.n_replicas = 13 + settings.forcefield_settings.nonbonded_method = 'nocutoff' errmsg = ("Number of replicas 13 does not equal the number of " "lambda windows 11") @@ -967,7 +967,7 @@ def test_too_many_specified_mappings(benzene_system, toluene_system, stateA=benzene_system, stateB=toluene_system, mapping={'solvent': benzene_to_toluene_mapping, - 'ligand': benzene_to_toluene_mapping,} + 'ligand': benzene_to_toluene_mapping, } ) @@ -1022,7 +1022,7 @@ def test_element_change_warning(atom_mapping_basic_test_files): def test_ligand_overlap_warning(benzene_vacuum_system, toluene_vacuum_system, benzene_to_toluene_mapping, tmpdir): vac_settings = openmm_rfe.RelativeHybridTopologyProtocol.default_settings() - vac_settings.system_settings.nonbonded_method = 'nocutoff' + vac_settings.forcefield_settings.nonbonded_method = 'nocutoff' protocol = openmm_rfe.RelativeHybridTopologyProtocol( settings=vac_settings, @@ -1062,11 +1062,9 @@ def test_ligand_overlap_warning(benzene_vacuum_system, toluene_vacuum_system, @pytest.fixture def solvent_protocol_dag(benzene_system, toluene_system, benzene_to_toluene_mapping): settings = openmm_rfe.RelativeHybridTopologyProtocol.default_settings() - protocol = openmm_rfe.RelativeHybridTopologyProtocol( settings=settings, ) - return protocol.create( stateA=benzene_system, stateB=toluene_system, mapping={'ligand': benzene_to_toluene_mapping}, @@ -1082,7 +1080,6 @@ def test_unit_tagging(solvent_protocol_dag, tmpdir): for u in dag_units: ret = u.execute(context=gufe.Context(tmpdir, tmpdir)) results.append(ret) - repeats = set() for ret in results: assert isinstance(ret, gufe.ProtocolUnitResult) @@ -1336,7 +1333,7 @@ def tyk2_xml(tmp_path_factory): settings.forcefield_settings.small_molecule_forcefield = 'openff-2.0.0' settings.system_settings.nonbonded_method = 'nocutoff' settings.forcefield_settings.hydrogen_mass = 3.0 - settings.alchemical_sampler_settings.n_repeats = 1 + settings.protocol_repeats = 1 protocol = openmm_rfe.RelativeHybridTopologyProtocol(settings) @@ -1553,9 +1550,8 @@ def benzene_solvent_openmm_system(benzene_modifications): system_generator = system_creation.get_system_generator( forcefield_settings=settings.forcefield_settings, - thermo_settings=settings.thermo_settings, integrator_settings=settings.integrator_settings, - system_settings=settings.system_settings, + thermo_settings=settings.thermo_settings, cache=None, has_solvent=True, ) @@ -1595,9 +1591,8 @@ def benzene_tip4p_solvent_openmm_system(benzene_modifications): system_generator = system_creation.get_system_generator( forcefield_settings=settings.forcefield_settings, - thermo_settings=settings.thermo_settings, integrator_settings=settings.integrator_settings, - system_settings=settings.system_settings, + thermo_settings=settings.thermo_settings, cache=None, has_solvent=True, ) diff --git a/openfe/tests/protocols/test_openmm_plain_md_protocols.py b/openfe/tests/protocols/test_openmm_plain_md_protocols.py index b96bd94a4..17af7d379 100644 --- a/openfe/tests/protocols/test_openmm_plain_md_protocols.py +++ b/openfe/tests/protocols/test_openmm_plain_md_protocols.py @@ -51,7 +51,7 @@ def test_create_independent_repeat_ids(benzene_system): # this allows multiple DAGs in flight for one Transformation that don't clash on gather settings = PlainMDProtocol.default_settings() # Default protocol is 1 repeat, change to 3 repeats - settings.repeat_settings.n_repeats = 3 + settings.protocol_repeats = 3 protocol = PlainMDProtocol( settings=settings, ) @@ -79,7 +79,7 @@ def test_create_independent_repeat_ids(benzene_system): def test_dry_run_default_vacuum(benzene_vacuum_system, tmpdir): vac_settings = PlainMDProtocol.default_settings() - vac_settings.system_settings.nonbonded_method = 'nocutoff' + vac_settings.forcefield_settings.nonbonded_method = 'nocutoff' protocol = PlainMDProtocol( settings=vac_settings, @@ -104,7 +104,7 @@ def test_dry_run_default_vacuum(benzene_vacuum_system, tmpdir): def test_dry_run_logger_output(benzene_vacuum_system, tmpdir, caplog): vac_settings = PlainMDProtocol.default_settings() - vac_settings.system_settings.nonbonded_method = 'nocutoff' + vac_settings.forcefield_settings.nonbonded_method = 'nocutoff' vac_settings.simulation_settings.equilibration_length_nvt = 1 * unit.picosecond vac_settings.simulation_settings.equilibration_length = 1 * unit.picosecond vac_settings.simulation_settings.production_length = 1 * unit.picosecond @@ -135,13 +135,13 @@ def test_dry_run_logger_output(benzene_vacuum_system, tmpdir, caplog): def test_dry_run_ffcache_none_vacuum(benzene_vacuum_system, tmpdir): vac_settings = PlainMDProtocol.default_settings() - vac_settings.system_settings.nonbonded_method = 'nocutoff' - vac_settings.simulation_settings.forcefield_cache = None + vac_settings.forcefield_settings.nonbonded_method = 'nocutoff' + vac_settings.output_settings.forcefield_cache = None protocol = PlainMDProtocol( settings=vac_settings, ) - assert protocol.settings.simulation_settings.forcefield_cache is None + assert protocol.settings.output_settings.forcefield_cache is None # create DAG from protocol and take first (and only) work unit from within dag = protocol.create( @@ -157,7 +157,7 @@ def test_dry_run_ffcache_none_vacuum(benzene_vacuum_system, tmpdir): def test_dry_run_gaff_vacuum(benzene_vacuum_system, tmpdir): vac_settings = PlainMDProtocol.default_settings() - vac_settings.system_settings.nonbonded_method = 'nocutoff' + vac_settings.forcefield_settings.nonbonded_method = 'nocutoff' vac_settings.forcefield_settings.small_molecule_forcefield = 'gaff-2.11' protocol = PlainMDProtocol( @@ -278,7 +278,7 @@ def test_dry_run_ligand_tip4p(benzene_system, tmpdir): "amber/phosaa10.xml", # Handles THE TPO ] settings.solvation_settings.solvent_padding = 1.0 * unit.nanometer - settings.system_settings.nonbonded_cutoff = 0.9 * unit.nanometer + settings.forcefield_settings.nonbonded_cutoff = 0.9 * unit.nanometer settings.solvation_settings.solvent_model = 'tip4pew' settings.integrator_settings.reassign_velocities = True @@ -326,7 +326,7 @@ def test_dry_run_complex(benzene_complex_system, tmpdir): def test_hightimestep(benzene_vacuum_system, tmpdir): settings = PlainMDProtocol.default_settings() settings.forcefield_settings.hydrogen_mass = 1.0 - settings.system_settings.nonbonded_method = 'nocutoff' + settings.forcefield_settings.nonbonded_method = 'nocutoff' p = PlainMDProtocol( settings=settings, @@ -362,7 +362,7 @@ def test_vaccuum_PME_error(benzene_vacuum_system): @pytest.fixture def solvent_protocol_dag(benzene_system): settings = PlainMDProtocol.default_settings() - settings.repeat_settings.n_repeats = 3 + settings.protocol_repeats = 3 protocol = PlainMDProtocol( settings=settings, ) @@ -402,7 +402,7 @@ def test_gather(solvent_protocol_dag, tmpdir): keep_shared=True) settings = PlainMDProtocol.default_settings() - settings.repeat_settings.n_repeats = 3 + settings.protocol_repeats = 3 prot = PlainMDProtocol( settings=settings ) diff --git a/openfe/tests/protocols/test_openmm_rfe_slow.py b/openfe/tests/protocols/test_openmm_rfe_slow.py index 2af72929c..4e8a8b522 100644 --- a/openfe/tests/protocols/test_openmm_rfe_slow.py +++ b/openfe/tests/protocols/test_openmm_rfe_slow.py @@ -47,11 +47,11 @@ def test_openmm_run_engine(benzene_vacuum_system, platform, s = openfe.protocols.openmm_rfe.RelativeHybridTopologyProtocol.default_settings() s.simulation_settings.equilibration_length = 0.1 * unit.picosecond s.simulation_settings.production_length = 0.1 * unit.picosecond - s.integrator_settings.n_steps = 5 * unit.timestep + s.alchemical_sampler_settings.steps_per_iteration = 5 * unit.timestep s.system_settings.nonbonded_method = 'nocutoff' - s.alchemical_sampler_settings.n_repeats = 1 + s.n_repeats = 1 s.engine_settings.compute_platform = platform - s.simulation_settings.checkpoint_interval = 5 * unit.timestep + s.output_settings.checkpoint_interval = 5 * unit.timestep p = openmm_rfe.RelativeHybridTopologyProtocol(s) @@ -105,9 +105,9 @@ def test_run_eg5_sim(eg5_protein, eg5_ligands, eg5_cofactor, tmpdir): s = openfe.protocols.openmm_rfe.RelativeHybridTopologyProtocol.default_settings() s.simulation_settings.equilibration_length = 0.1 * unit.picosecond s.simulation_settings.production_length = 0.1 * unit.picosecond - s.integrator_settings.n_steps = 5 * unit.timestep - s.alchemical_sampler_settings.n_repeats = 1 - s.simulation_settings.checkpoint_interval = 5 * unit.timestep + s.alchemical_sampler_settings.steps_per_iteration = 5 * unit.timestep + s.n_repeats = 1 + s.output_settings.checkpoint_interval = 5 * unit.timestep p = openmm_rfe.RelativeHybridTopologyProtocol(s) diff --git a/openfe/tests/protocols/test_openmm_settings.py b/openfe/tests/protocols/test_openmm_settings.py index f73a2907b..5c95e0e2d 100644 --- a/openfe/tests/protocols/test_openmm_settings.py +++ b/openfe/tests/protocols/test_openmm_settings.py @@ -12,7 +12,7 @@ class TestOMMSettingsFromStrings: # checks that we can set Settings fields via strings def test_system_settings(self): - s = omm_settings.SystemSettings() + s = omm_settings.OpenMMSystemGeneratorFFSettings() s.nonbonded_cutoff = '1.1 nm' @@ -26,7 +26,7 @@ def test_solvation_settings(self): assert s.solvent_padding == 1.1 * unit.nanometer def test_alchemical_sampler_settings(self): - # todo: online_analysis_target_error is in kT, how to pass this as string? + # todo: early_termination_target_error is in kT, how to pass this as string? pass def test_integator_settings(self): @@ -36,9 +36,9 @@ def test_integator_settings(self): assert s.timestep == 3.0 * unit.femtosecond - s.collision_rate = '1.1 / ps' + s.langevin_collision_rate = '1.1 / ps' - assert s.collision_rate == 1.1 / unit.picosecond + assert s.langevin_collision_rate == 1.1 / unit.picosecond # todo: nsteps, barostat frequency require IntQuantity @@ -59,7 +59,7 @@ def test_simulation_settings(self): class TestEquilRFESettingsFromString: def test_alchemical_settings(self): - s = equil_rfe_settings.AlchemicalSettings() + s = equil_rfe_settings.AlchemicalSettings(softcore_LJ='gapsys') s.explicit_charge_correction_cutoff = '0.85 nm' diff --git a/openfe/tests/protocols/test_openmmutils.py b/openfe/tests/protocols/test_openmmutils.py index 9cba0c25b..4776c444b 100644 --- a/openfe/tests/protocols/test_openmmutils.py +++ b/openfe/tests/protocols/test_openmmutils.py @@ -15,10 +15,10 @@ import openfe from openfe.protocols.openmm_utils import ( settings_validation, system_validation, system_creation, - multistate_analysis + multistate_analysis, omm_settings ) from openfe.protocols.openmm_rfe.equil_rfe_settings import ( - SystemSettings, SolvationSettings, IntegratorSettings, + SolvationSettings, IntegratorSettings, ) @@ -165,15 +165,13 @@ def test_components_complex(T4_protein_component, benzene_modifications): @pytest.fixture(scope='module') def get_settings(): forcefield_settings = OpenMMSystemGeneratorFFSettings() + integrator_settings = IntegratorSettings() thermo_settings = ThermoSettings( temperature=298.15 * unit.kelvin, pressure=1 * unit.bar, ) - system_settings = SystemSettings() - integrator_settings = IntegratorSettings() - return (forcefield_settings, thermo_settings, system_settings, - integrator_settings) + return forcefield_settings, integrator_settings, thermo_settings class TestFEAnalysis: @@ -269,9 +267,9 @@ def test_analyze_unknown_method_warning_and_error(self, reporter): class TestSystemCreation: def test_system_generator_nosolv_nocache(self, get_settings): - ffsets, thermosets, systemsets, intsets = get_settings + ffsets, intsets, thermosets = get_settings generator = system_creation.get_system_generator( - ffsets, thermosets, intsets, systemsets, None, False + ffsets, thermosets, intsets, None, False ) assert generator.barostat is None assert generator.template_generator._cache is None @@ -293,13 +291,13 @@ def test_system_generator_nosolv_nocache(self, get_settings): assert generator.periodic_forcefield_kwargs == periodic_kwargs def test_system_generator_solv_cache(self, get_settings): - ffsets, thermosets, systemsets, intsets = get_settings + ffsets, intsets, thermosets = get_settings thermosets.temperature = 320 * unit.kelvin thermosets.pressure = 1.25 * unit.bar intsets.barostat_frequency = 200 * unit.timestep generator = system_creation.get_system_generator( - ffsets, thermosets, intsets, systemsets, Path('./db.json'), True + ffsets, thermosets, intsets, Path('./db.json'), True ) # Check barostat conditions @@ -323,9 +321,9 @@ def test_system_generator_solv_cache(self, get_settings): def test_get_omm_modeller_complex(self, T4_protein_component, benzene_modifications, get_settings): - ffsets, thermosets, systemsets, intsets = get_settings + ffsets, intsets, thermosets = get_settings generator = system_creation.get_system_generator( - ffsets, thermosets, intsets, systemsets, None, True + ffsets, thermosets, intsets, None, True ) smc = benzene_modifications['toluene'] @@ -349,9 +347,9 @@ def test_get_omm_modeller_complex(self, T4_protein_component, np.linspace(165, len(resids)-1, len(resids)-165)) def test_get_omm_modeller_ligand_no_neutralize(self, get_settings): - ffsets, thermosets, systemsets, intsets = get_settings + ffsets, intsets, thermosets = get_settings generator = system_creation.get_system_generator( - ffsets, thermosets, intsets, systemsets, None, True + ffsets, thermosets, intsets, None, True ) offmol = OFFMol.from_smiles('[O-]C=O') @@ -386,3 +384,108 @@ def test_get_omm_modeller_ligand_no_neutralize(self, get_settings): charge = ensure_quantity(charge, 'openff') assert pytest.approx(charge.m) == -1.0 + + +def test_convert_steps_per_iteration(): + sim = omm_settings.MultiStateSimulationSettings( + equilibration_length='10 ps', + production_length='10 ps', + time_per_iteration='1.0 ps', + ) + inty = omm_settings.IntegratorSettings( + timestep='4 fs' + ) + + spi = settings_validation.convert_steps_per_iteration(sim, inty) + + assert spi == 250 + + +def test_convert_steps_per_iteration_failure(): + sim = omm_settings.MultiStateSimulationSettings( + equilibration_length='10 ps', + production_length='10 ps', + time_per_iteration='1.0 ps', + ) + inty = omm_settings.IntegratorSettings( + timestep='3 fs' + ) + + with pytest.raises(ValueError, match="not divisible"): + settings_validation.convert_steps_per_iteration(sim, inty) + + +def test_convert_real_time_analysis_iterations(): + sim = omm_settings.MultiStateSimulationSettings( + equilibration_length='10 ps', + production_length='10 ps', + time_per_iteration='1.0 ps', + real_time_analysis_interval='250 ps', + real_time_analysis_minimum_time='500 ps', + ) + + rta_its, rta_min_its = settings_validation.convert_real_time_analysis_iterations(sim) + + assert rta_its == 250, 500 + + +def test_convert_real_time_analysis_iterations_interval_fail(): + # shouldn't like 250.5 ps / 1.0 ps + sim = omm_settings.MultiStateSimulationSettings( + equilibration_length='10 ps', + production_length='10 ps', + time_per_iteration='1.0 ps', + real_time_analysis_interval='250.5 ps', + real_time_analysis_minimum_time='500 ps', + ) + + with pytest.raises(ValueError, match='not divisible'): + settings_validation.convert_real_time_analysis_iterations(sim) + + +def test_convert_real_time_analysis_iterations_min_interval_fail(): + # shouldn't like 500.5 ps / 1 ps + sim = omm_settings.MultiStateSimulationSettings( + equilibration_length='10 ps', + production_length='10 ps', + time_per_iteration='1.0 ps', + real_time_analysis_interval='250 ps', + real_time_analysis_minimum_time='500.5 ps', + ) + + with pytest.raises(ValueError, match='not divisible'): + settings_validation.convert_real_time_analysis_iterations(sim) + + +def test_convert_real_time_analysis_iterations_None(): + sim = omm_settings.MultiStateSimulationSettings( + equilibration_length='10 ps', + production_length='10 ps', + time_per_iteration='1.0 ps', + real_time_analysis_interval=None, + real_time_analysis_minimum_time='500 ps', + ) + + rta_its, rta_min_its = settings_validation.convert_real_time_analysis_iterations(sim) + + assert rta_its is None + assert rta_min_its is None + + +def test_convert_target_error_from_kcal_per_mole_to_kT(): + kT = settings_validation.convert_target_error_from_kcal_per_mole_to_kT( + temperature=298.15 * unit.kelvin, + target_error=0.12 * unit.kilocalorie_per_mole, + ) + + assert kT == pytest.approx(0.20253681663365392) + + +def test_convert_target_error_from_kcal_per_mole_to_kT_zero(): + # special case, 0 input gives 0 output + kT = settings_validation.convert_target_error_from_kcal_per_mole_to_kT( + temperature=298.15 * unit.kelvin, + target_error=0.0 * unit.kilocalorie_per_mole, + ) + + assert kT == 0.0 diff --git a/openfe/tests/protocols/test_rfe_tokenization.py b/openfe/tests/protocols/test_rfe_tokenization.py index 928b5786d..044f3a56d 100644 --- a/openfe/tests/protocols/test_rfe_tokenization.py +++ b/openfe/tests/protocols/test_rfe_tokenization.py @@ -39,7 +39,7 @@ def instance(self): class TestRelativeHybridTopologyProtocol(GufeTokenizableTestsMixin): cls = openmm_rfe.RelativeHybridTopologyProtocol - key = "RelativeHybridTopologyProtocol-02ad30e0e0332ab6a8ce97bcdf3f287c" + key = "RelativeHybridTopologyProtocol-a7946099ab9756c22f3a79fe9082585a" repr = f"<{key}>" @pytest.fixture() diff --git a/openfe/tests/protocols/test_solvation_afe_tokenization.py b/openfe/tests/protocols/test_solvation_afe_tokenization.py index dfce2ae97..2169023c6 100644 --- a/openfe/tests/protocols/test_solvation_afe_tokenization.py +++ b/openfe/tests/protocols/test_solvation_afe_tokenization.py @@ -49,7 +49,7 @@ def protocol_result(afe_solv_transformation_json): class TestAbsoluteSolvationProtocol(GufeTokenizableTestsMixin): cls = openmm_afe.AbsoluteSolvationProtocol - key = "AbsoluteSolvationProtocol-045abf8f41dbb14db338c2c8583da308" + key = "AbsoluteSolvationProtocol-9a18332b06a721da1b0fcaf5cc86dd25" repr = f"<{key}>" @pytest.fixture() @@ -85,7 +85,7 @@ def test_key_stable(self): class TestAbsoluteSolvationProtocolResult(GufeTokenizableTestsMixin): cls = openmm_afe.AbsoluteSolvationProtocolResult - key = "AbsoluteSolvationProtocolResult-bd4bb7aadc07066d532596543ec193a2" + key = "AbsoluteSolvationProtocolResult-55ac1971317176d6e84549601d1eed5e" repr = f"<{key}>" @pytest.fixture() diff --git a/openfecli/tests/commands/test_gather.py b/openfecli/tests/commands/test_gather.py index 2eac0c460..7bb71d7ed 100644 --- a/openfecli/tests/commands/test_gather.py +++ b/openfecli/tests/commands/test_gather.py @@ -209,8 +209,8 @@ def test_missing_leg_allow_partial(results_dir): RBFE_RESULTS = pooch.create( pooch.os_cache('openfe'), - base_url="doi:10.6084/m9.figshare.24542059", - registry={"results.tar.gz": None}, + base_url="doi:10.6084/m9.figshare.25148945", + registry={"results.tar.gz": "bf27e728935b31360f95188f41807558156861f6d89b8a47854502a499481da3"}, ) diff --git a/openfecli/tests/data/bad_transformation.json b/openfecli/tests/data/bad_transformation.json index bcfc8107c..fe0bc82f0 100644 --- a/openfecli/tests/data/bad_transformation.json +++ b/openfecli/tests/data/bad_transformation.json @@ -1 +1 @@ -{"__module__": "gufe.transformations.transformation", "__qualname__": "Transformation", "mapping": {"__module__": "gufe.mapping.ligandatommapping", "__qualname__": "LigandAtomMapping", "annotations": "{}", "componentA": {"__module__": "gufe.components.smallmoleculecomponent", "__qualname__": "SmallMoleculeComponent", "atoms": [[6, 0, 0, true, 0, 0, {}], [6, 0, 0, true, 0, 0, {}], [6, 0, 0, true, 0, 0, {}], [6, 0, 0, true, 0, 0, {}], [6, 0, 0, true, 0, 0, {}], [6, 0, 0, true, 0, 0, {}], [1, 0, 0, false, 0, 0, {}], [1, 0, 0, false, 0, 0, {}], [1, 0, 0, false, 0, 0, {}], [1, 0, 0, false, 0, 0, {}], [1, 0, 0, false, 0, 0, {}], [1, 0, 0, false, 0, 0, {}]], "bonds": [[0, 1, 12, 0, {}], [0, 5, 12, 0, {}], [0, 6, 1, 0, {}], [1, 2, 12, 0, {}], [1, 7, 1, 0, {}], [2, 3, 12, 0, {}], [2, 8, 1, 0, {}], [3, 4, 12, 0, {}], [3, 9, 1, 0, {}], [4, 5, 12, 0, {}], [4, 10, 1, 0, {}], [5, 11, 1, 0, {}]], "conformer": ["\u0093NUMPY\u0001\u0000v\u0000{'descr': '\u00e8<@[\u00b1\u00bf\u00ec\u009e|!@\u00b0rh\u0091\u00ed|\u0014@\u00c3d\u00aa`T2<@\u009a\b\u001b\u009e^I @\u00caT\u00c1\u00a8\u00a4\u008e\u001a@?\u00c6\u00dc\u00b5\u0084\u00fc;@\u00da\u001b|a2\u00d5 @$(~\u008c\u00b9k\u0016@n\u00a3\u0001\u00bc\u0005B;@\u00c0\u00ec\u009e<,t\"@\u0084\u009e\u00cd\u00aa\u00cfU\u0016@\u00ee|?5^\u00fa9@\u0002+\u0087\u0016\u00d9N\u0015@\u0004V\u000e-\u00b2\u001d\u0013@\u0085\u00ebQ\u00b8\u001ee:@\u00b2\u009d\u00ef\u00a7\u00c6K\u0014@\u00cb\u00a1E\u00b6\u00f3\u00fd\u000b@\u00d7\u00a3p=\nW;@q=\n\u00d7\u00a3p\u0017@\u009e\u00ef\u00a7\u00c6K7\u0007@\u0083\u00c0\u00ca\u00a1E\u00d6;@\u00c9v\u00be\u009f\u001a\u00af\u001b@Zd;\u00dfO\u008d\f@\u00ecQ\u00b8\u001e\u0085k;@b\u0010X9\u00b4\u00c8\u001c@\u0006\u0081\u0095C\u008bl\u0013@sh\u0091\u00ed|\u007f:@j\u00bct\u0093\u0018\u0084\u0019@\u00c7K7\u0089A\u00e0\u0015@\u00ed\u009e<,\u00d4:9@\u00e8<@[\u00b1\u00bf\u00ec\u009e|!@\u00b0rh\u0091\u00ed|\u0014@\u00c3d\u00aa`T2<@\u009a\b\u001b\u009e^I @\u00caT\u00c1\u00a8\u00a4\u008e\u001a@?\u00c6\u00dc\u00b5\u0084\u00fc;@\u00da\u001b|a2\u00d5 @$(~\u008c\u00b9k\u0016@n\u00a3\u0001\u00bc\u0005B;@\u00c0\u00ec\u009e<,t\"@\u0084\u009e\u00cd\u00aa\u00cfU\u0016@\u00ee|?5^\u00fa9@\u0002+\u0087\u0016\u00d9N\u0015@\u0004V\u000e-\u00b2\u001d\u0013@\u0085\u00ebQ\u00b8\u001ee:@\u00b2\u009d\u00ef\u00a7\u00c6K\u0014@\u00cb\u00a1E\u00b6\u00f3\u00fd\u000b@\u00d7\u00a3p=\nW;@q=\n\u00d7\u00a3p\u0017@\u009e\u00ef\u00a7\u00c6K7\u0007@\u0083\u00c0\u00ca\u00a1E\u00d6;@\u00c9v\u00be\u009f\u001a\u00af\u001b@Zd;\u00dfO\u008d\f@\u00ecQ\u00b8\u001e\u0085k;@b\u0010X9\u00b4\u00c8\u001c@\u0006\u0081\u0095C\u008bl\u0013@sh\u0091\u00ed|\u007f:@j\u00bct\u0093\u0018\u0084\u0019@\u00c7K7\u0089A\u00e0\u0015@\u00ed\u009e<,\u00d4:9@\u00e8<@[\u00b1\u00bf\u00ec\u009e|!@\u00b0rh\u0091\u00ed|\u0014@\u00c3d\u00aa`T2<@\u009a\b\u001b\u009e^I @\u00caT\u00c1\u00a8\u00a4\u008e\u001a@?\u00c6\u00dc\u00b5\u0084\u00fc;@\u00da\u001b|a2\u00d5 @$(~\u008c\u00b9k\u0016@n\u00a3\u0001\u00bc\u0005B;@\u00c0\u00ec\u009e<,t\"@\u0084\u009e\u00cd\u00aa\u00cfU\u0016@\u00ee|?5^\u00fa9@\u0002+\u0087\u0016\u00d9N\u0015@\u0004V\u000e-\u00b2\u001d\u0013@\u0085\u00ebQ\u00b8\u001ee:@\u00b2\u009d\u00ef\u00a7\u00c6K\u0014@\u00cb\u00a1E\u00b6\u00f3\u00fd\u000b@\u00d7\u00a3p=\nW;@q=\n\u00d7\u00a3p\u0017@\u009e\u00ef\u00a7\u00c6K7\u0007@\u0083\u00c0\u00ca\u00a1E\u00d6;@\u00c9v\u00be\u009f\u001a\u00af\u001b@Zd;\u00dfO\u008d\f@\u00ecQ\u00b8\u001e\u0085k;@b\u0010X9\u00b4\u00c8\u001c@\u0006\u0081\u0095C\u008bl\u0013@sh\u0091\u00ed|\u007f:@j\u00bct\u0093\u0018\u0084\u0019@\u00c7K7\u0089A\u00e0\u0015@\u00ed\u009e<,\u00d4:9@\u00e8<@[\u00b1\u00bf\u00ec\u009e|!@\u00b0rh\u0091\u00ed|\u0014@\u00c3d\u00aa`T2<@\u009a\b\u001b\u009e^I @\u00caT\u00c1\u00a8\u00a4\u008e\u001a@?\u00c6\u00dc\u00b5\u0084\u00fc;@\u00da\u001b|a2\u00d5 @$(~\u008c\u00b9k\u0016@n\u00a3\u0001\u00bc\u0005B;@\u00c0\u00ec\u009e<,t\"@\u0084\u009e\u00cd\u00aa\u00cfU\u0016@\u00ee|?5^\u00fa9@\u0002+\u0087\u0016\u00d9N\u0015@\u0004V\u000e-\u00b2\u001d\u0013@\u0085\u00ebQ\u00b8\u001ee:@\u00b2\u009d\u00ef\u00a7\u00c6K\u0014@\u00cb\u00a1E\u00b6\u00f3\u00fd\u000b@\u00d7\u00a3p=\nW;@q=\n\u00d7\u00a3p\u0017@\u009e\u00ef\u00a7\u00c6K7\u0007@\u0083\u00c0\u00ca\u00a1E\u00d6;@\u00c9v\u00be\u009f\u001a\u00af\u001b@Zd;\u00dfO\u008d\f@\u00ecQ\u00b8\u001e\u0085k;@b\u0010X9\u00b4\u00c8\u001c@\u0006\u0081\u0095C\u008bl\u0013@sh\u0091\u00ed|\u007f:@j\u00bct\u0093\u0018\u0084\u0019@\u00c7K7\u0089A\u00e0\u0015@\u00ed\u009e<,\u00d4:9@\u00e8<@[\u00b1\u00bf\u00ec\u009e|!@\u00b0rh\u0091\u00ed|\u0014@\u00c3d\u00aa`T2<@\u009a\b\u001b\u009e^I @\u00caT\u00c1\u00a8\u00a4\u008e\u001a@?\u00c6\u00dc\u00b5\u0084\u00fc;@\u00da\u001b|a2\u00d5 @$(~\u008c\u00b9k\u0016@n\u00a3\u0001\u00bc\u0005B;@\u00c0\u00ec\u009e<,t\"@\u0084\u009e\u00cd\u00aa\u00cfU\u0016@\u00ee|?5^\u00fa9@\u0002+\u0087\u0016\u00d9N\u0015@\u0004V\u000e-\u00b2\u001d\u0013@\u0085\u00ebQ\u00b8\u001ee:@\u00b2\u009d\u00ef\u00a7\u00c6K\u0014@\u00cb\u00a1E\u00b6\u00f3\u00fd\u000b@\u00d7\u00a3p=\nW;@q=\n\u00d7\u00a3p\u0017@\u009e\u00ef\u00a7\u00c6K7\u0007@\u0083\u00c0\u00ca\u00a1E\u00d6;@\u00c9v\u00be\u009f\u001a\u00af\u001b@Zd;\u00dfO\u008d\f@\u00ecQ\u00b8\u001e\u0085k;@b\u0010X9\u00b4\u00c8\u001c@\u0006\u0081\u0095C\u008bl\u0013@sh\u0091\u00ed|\u007f:@j\u00bct\u0093\u0018\u0084\u0019@\u00c7K7\u0089A\u00e0\u0015@\u00ed\u009e<,\u00d4:9@\u00e8<@[\u00b1\u00bf\u00ec\u009e|!@\u00b0rh\u0091\u00ed|\u0014@\u00c3d\u00aa`T2<@\u009a\b\u001b\u009e^I @\u00caT\u00c1\u00a8\u00a4\u008e\u001a@?\u00c6\u00dc\u00b5\u0084\u00fc;@\u00da\u001b|a2\u00d5 @$(~\u008c\u00b9k\u0016@n\u00a3\u0001\u00bc\u0005B;@\u00c0\u00ec\u009e<,t\"@\u0084\u009e\u00cd\u00aa\u00cfU\u0016@\u00ee|?5^\u00fa9@\u0002+\u0087\u0016\u00d9N\u0015@\u0004V\u000e-\u00b2\u001d\u0013@\u0085\u00ebQ\u00b8\u001ee:@\u00b2\u009d\u00ef\u00a7\u00c6K\u0014@\u00cb\u00a1E\u00b6\u00f3\u00fd\u000b@\u00d7\u00a3p=\nW;@q=\n\u00d7\u00a3p\u0017@\u009e\u00ef\u00a7\u00c6K7\u0007@\u0083\u00c0\u00ca\u00a1E\u00d6;@\u00c9v\u00be\u009f\u001a\u00af\u001b@Zd;\u00dfO\u008d\f@\u00ecQ\u00b8\u001e\u0085k;@b\u0010X9\u00b4\u00c8\u001c@\u0006\u0081\u0095C\u008bl\u0013@sh\u0091\u00ed|\u007f:@j\u00bct\u0093\u0018\u0084\u0019@\u00c7K7\u0089A\u00e0\u0015@\u00ed\u009e<,\u00d4:9@\u00e8<@[\u00b1\u00bf\u00ec\u009e|!@\u00b0rh\u0091\u00ed|\u0014@\u00c3d\u00aa`T2<@\u009a\b\u001b\u009e^I @\u00caT\u00c1\u00a8\u00a4\u008e\u001a@?\u00c6\u00dc\u00b5\u0084\u00fc;@\u00da\u001b|a2\u00d5 @$(~\u008c\u00b9k\u0016@n\u00a3\u0001\u00bc\u0005B;@\u00c0\u00ec\u009e<,t\"@\u0084\u009e\u00cd\u00aa\u00cfU\u0016@\u00ee|?5^\u00fa9@\u0002+\u0087\u0016\u00d9N\u0015@\u0004V\u000e-\u00b2\u001d\u0013@\u0085\u00ebQ\u00b8\u001ee:@\u00b2\u009d\u00ef\u00a7\u00c6K\u0014@\u00cb\u00a1E\u00b6\u00f3\u00fd\u000b@\u00d7\u00a3p=\nW;@q=\n\u00d7\u00a3p\u0017@\u009e\u00ef\u00a7\u00c6K7\u0007@\u0083\u00c0\u00ca\u00a1E\u00d6;@\u00c9v\u00be\u009f\u001a\u00af\u001b@Zd;\u00dfO\u008d\f@\u00ecQ\u00b8\u001e\u0085k;@b\u0010X9\u00b4\u00c8\u001c@\u0006\u0081\u0095C\u008bl\u0013@sh\u0091\u00ed|\u007f:@j\u00bct\u0093\u0018\u0084\u0019@\u00c7K7\u0089A\u00e0\u0015@\u00ed\u009e<,\u00d4:9@\u00e8<@[\u00b1\u00bf\u00ec\u009e|!@\u00b0rh\u0091\u00ed|\u0014@\u00c3d\u00aa`T2<@\u009a\b\u001b\u009e^I @\u00caT\u00c1\u00a8\u00a4\u008e\u001a@?\u00c6\u00dc\u00b5\u0084\u00fc;@\u00da\u001b|a2\u00d5 @$(~\u008c\u00b9k\u0016@n\u00a3\u0001\u00bc\u0005B;@\u00c0\u00ec\u009e<,t\"@\u0084\u009e\u00cd\u00aa\u00cfU\u0016@\u00ee|?5^\u00fa9@\u0002+\u0087\u0016\u00d9N\u0015@\u0004V\u000e-\u00b2\u001d\u0013@\u0085\u00ebQ\u00b8\u001ee:@\u00b2\u009d\u00ef\u00a7\u00c6K\u0014@\u00cb\u00a1E\u00b6\u00f3\u00fd\u000b@\u00d7\u00a3p=\nW;@q=\n\u00d7\u00a3p\u0017@\u009e\u00ef\u00a7\u00c6K7\u0007@\u0083\u00c0\u00ca\u00a1E\u00d6;@\u00c9v\u00be\u009f\u001a\u00af\u001b@Zd;\u00dfO\u008d\f@\u00ecQ\u00b8\u001e\u0085k;@b\u0010X9\u00b4\u00c8\u001c@\u0006\u0081\u0095C\u008bl\u0013@sh\u0091\u00ed|\u007f:@j\u00bct\u0093\u0018\u0084\u0019@\u00c7K7\u0089A\u00e0\u0015@\u00ed\u009e<,\u00d4:9@