openmmtools.alchemy.AbsoluteAlchemicalFactory

class openmmtools.alchemy.AbsoluteAlchemicalFactory(consistent_exceptions=False, switch_width=Quantity(value=1.0, unit=angstrom), alchemical_pme_treatment='exact', alchemical_rf_treatment='switched', disable_alchemical_dispersion_correction=False, split_alchemical_forces=True)[source]

Factory of alchemically modified OpenMM Systems.

The context parameters created are: - softcore_alpha: factor controlling softcore lengthscale for Lennard-Jones - softcore_beta: factor controlling softcore lengthscale for Coulomb - softcore_a: softcore Lennard-Jones parameter from Eq. 13 of Ref [1] - softcore_b: softcore Lennard-Jones parameter from Eq. 13 of Ref [1] - softcore_c: softcore Lennard-Jones parameter from Eq. 13 of Ref [1] - softcore_d: softcore electrostatics parameter - softcore_e: softcore electrostatics parameter - softcore_f: softcore electrostatics parameter

Parameters:
consistent_exceptionsbool, optional, default = False

If True, the same functional form of the System’s nonbonded method will be use to determine the electrostatics contribution to the potential energy of 1,4 exceptions instead of the classical q1*q2/(4*epsilon*epsilon0*pi*r).

switch_widthfloat, optional, default = 1.0 * angstroms

Default switch width for electrostatics in periodic cutoff systems used in alchemical interactions only.

alchemical_pme_treatmentstr, optional, default = ‘exact’

Controls how alchemical region electrostatics are treated when PME is used. Options are [‘direct-space’, ‘coulomb’, ‘exact’]. - ‘direct-space’ only models the direct space contribution - ‘coulomb’ includes switched Coulomb interaction - ‘exact’ includes also the reciprocal space contribution, but it’s

only possible to annihilate the charges and the softcore parameters controlling the electrostatics are deactivated. Also, with this method, modifying the global variable lambda_electrostatics is not sufficient to control the charges. The recommended way to change them is through the AlchemicalState class.

alchemical_rf_treatmentstr, optional, default = ‘switched’

Controls how alchemical region electrostatics are treated when RF is used Options are [‘switched’, ‘shifted’] ‘switched’ sets c_rf = 0 for all reaction-field interactions and ensures continuity with a switch ‘shifted’ retains c_rf != 0 but can give erroneous results for hydration free energies

disable_alchemical_dispersion_correctionbool, optional, default=False

If True, the long-range dispersion correction will not be included for the alchemical region to avoid the need to recompute the correction (a CPU operation that takes ~ 0.5 s) every time ‘lambda_sterics’ is changed. If using nonequilibrium protocols, it is recommended that this be set to True since this can lead to enormous (100x) slowdowns if the correction must be recomputed every time step.

split_alchemical_forcesbool, optional, default=True

If True, forces that are altered to different alchemical variables will be split in different force groups. All non-alchemical forces will maintain their original force group. If more than 32 force groups are required, an error is thrown.

References

[1] Pham TT and Shirts MR. Identifying low variance pathways for free energy calculations of molecular transformations in solution phase. JCP 135:034114, 2011. http://dx.doi.org/10.1063/1.3607597

Examples

Create an alchemical factory to alchemically modify OpenMM System objects.

>>> factory = AbsoluteAlchemicalFactory(consistent_exceptions=False)

Create an alchemically modified version of p-xylene in T4 lysozyme L99A in GBSA.

>>> # Create a reference system.
>>> from openmmtools import testsystems
>>> reference_system = testsystems.LysozymeImplicit().system
>>> # Alchemically soften the pxylene atoms
>>> pxylene_atoms = range(2603,2621) # p-xylene
>>> alchemical_region = AlchemicalRegion(alchemical_atoms=pxylene_atoms)
>>> alchemical_system = factory.create_alchemical_system(reference_system, alchemical_region)

Alchemically modify one water in a water box.

>>> reference_system = testsystems.WaterBox().system
>>> alchemical_region = AlchemicalRegion(alchemical_atoms=[0, 1, 2])
>>> alchemical_system = factory.create_alchemical_system(reference_system, alchemical_region)

Alchemically modify some angles and torsions in alanine dipeptide and annihilate both sterics and electrostatics.

>>> reference_system = testsystems.AlanineDipeptideVacuum().system
>>> alchemical_region = AlchemicalRegion(alchemical_atoms=[0], alchemical_torsions=[0,1,2],
...                                      alchemical_angles=[0,1,2], annihilate_sterics=True,
...                                      annihilate_electrostatics=True)
>>> alchemical_system = factory.create_alchemical_system(reference_system, alchemical_region)

Alchemically modify a bond, angles, and torsions in toluene by automatically selecting bonds involving alchemical atoms.

>>> toluene_implicit = testsystems.TolueneImplicit()
>>> alchemical_region = AlchemicalRegion(alchemical_atoms=[0,1], alchemical_torsions=True,
...                                      alchemical_angles=True, annihilate_sterics=True)
>>> alchemical_system = factory.create_alchemical_system(reference_system, alchemical_region)

Once the alchemical system is created, you can modify its Hamiltonian through AlchemicalState

>>> alchemical_state = AlchemicalState.from_system(alchemical_system)
>>> alchemical_state.lambda_sterics
1.0
>>> alchemical_state.lambda_electrostatics = 0.5
>>> alchemical_state.apply_to_system(alchemical_system)

You can also modify its Hamiltonian directly into a context

>>> import openmm
>>> from openmm import unit
>>> integrator = openmm.VerletIntegrator(1.0*unit.femtosecond)
>>> context = openmm.Context(alchemical_system, integrator)
>>> alchemical_state.set_alchemical_parameters(0.0)  # Set all lambda to 0
>>> alchemical_state.apply_to_context(context)

Neglecting the long-range dispersion correction for the alchemical region (for nonequilibrium switching, for example) requires instantiating a factory with the appropriate options:

>>> new_factory = AbsoluteAlchemicalFactory(consistent_exceptions=False, disable_alchemical_dispersion_correction=True)
>>> reference_system = testsystems.WaterBox().system
>>> alchemical_region = AlchemicalRegion(alchemical_atoms=[0, 1, 2])
>>> alchemical_system = new_factory.create_alchemical_system(reference_system, alchemical_region)

Methods

create_alchemical_system(reference_system, ...)

Create an alchemically modified version of the reference system.

get_energy_components(alchemical_system, ...)

Compute potential energy of the alchemical system by Force.

__init__(consistent_exceptions=False, switch_width=Quantity(value=1.0, unit=angstrom), alchemical_pme_treatment='exact', alchemical_rf_treatment='switched', disable_alchemical_dispersion_correction=False, split_alchemical_forces=True)[source]

Methods

__init__([consistent_exceptions, ...])

create_alchemical_system(reference_system, ...)

Create an alchemically modified version of the reference system.

get_energy_components(alchemical_system, ...)

Compute potential energy of the alchemical system by Force.