Skip to content

restraints #

Create OpenMM restraint forces.

Modules:

  • femto

    A comprehensive toolkit for predicting free energies

Functions:

create_flat_bottom_restraint #

create_flat_bottom_restraint(
    config: FlatBottomRestraint, coords: dict[int, Quantity]
) -> CustomExternalForce

Creates a flat bottom position restraint.

Parameters:

  • config (FlatBottomRestraint) –

    The restraint configuration.

  • coords (dict[int, Quantity]) –

    A dictionary of indices of atoms to restrain and the corresponding coordinates to restrain them to.

Returns:

  • CustomExternalForce

    The restraint force.

Source code in femto/md/restraints.py
def create_flat_bottom_restraint(
    config: femto.md.config.FlatBottomRestraint,
    coords: dict[int, openmm.unit.Quantity],
) -> openmm.CustomExternalForce:
    """Creates a flat bottom position restraint.

    Args:
        config: The restraint configuration.
        coords: A dictionary of indices of atoms to restrain and the corresponding
            coordinates to restrain them to.

    Returns:
        The restraint force.
    """

    force = openmm.CustomExternalForce(_FLAT_BOTTOM_ENERGY_FN)

    force.addPerParticleParameter("x0")
    force.addPerParticleParameter("y0")
    force.addPerParticleParameter("z0")
    force.addPerParticleParameter("k")
    force.addPerParticleParameter("radius")

    k = config.k.value_in_unit_system(openmm.unit.md_unit_system)
    radius = config.radius.value_in_unit_system(openmm.unit.md_unit_system)

    for idx, coord in coords.items():
        force.addParticle(
            idx,
            [
                coord[0].value_in_unit_system(openmm.unit.md_unit_system),
                coord[1].value_in_unit_system(openmm.unit.md_unit_system),
                coord[2].value_in_unit_system(openmm.unit.md_unit_system),
                k,
                radius,
            ],
        )

    force.setForceGroup(femto.md.constants.OpenMMForceGroup.POSITION_RESTRAINT)
    force.setName(femto.md.constants.OpenMMForceName.POSITION_RESTRAINT)

    return force

create_position_restraints #

create_position_restraints(
    topology: Structure,
    mask: str,
    config: FlatBottomRestraint,
) -> CustomExternalForce

Creates position restraints for all ligand atoms.

Parameters:

  • topology (Structure) –

    The topology of the system being simulation.

  • mask (str) –

    The mask that defines which atoms to restrain.

  • config (FlatBottomRestraint) –

    The restraint configuration.

Returns:

  • CustomExternalForce

    The created restraint force.

Source code in femto/md/restraints.py
def create_position_restraints(
    topology: parmed.Structure,
    mask: str,
    config: femto.md.config.FlatBottomRestraint,
) -> openmm.CustomExternalForce:
    """Creates position restraints for all ligand atoms.

    Args:
        topology: The topology of the system being simulation.
        mask: The mask that defines which atoms to restrain.
        config: The restraint configuration.

    Returns:
        The created restraint force.
    """

    selection = parmed.amber.AmberMask(topology, mask).Selection()
    selection_idxs = tuple(i for i, matches in enumerate(selection) if matches)

    assert len(selection_idxs) > 0, "no atoms were found to restrain"

    coords = {
        i: openmm.Vec3(topology.atoms[i].xx, topology.atoms[i].xy, topology.atoms[i].xz)
        * openmm.unit.angstrom
        for i in selection_idxs
    }

    if not isinstance(config, femto.md.config.FlatBottomRestraint):
        raise NotImplementedError("only flat bottom restraints are currently supported")

    return femto.md.restraints.create_flat_bottom_restraint(config, coords)

create_boresch_restraint #

create_boresch_restraint(
    config: BoreschRestraint,
    receptor_atoms: tuple[int, int, int],
    ligand_atoms: tuple[int, int, int],
    coords: Quantity,
    ctx_parameter: str | None = None,
) -> CustomCompoundBondForce

Creates a Boresch restraint force useful in aligning a receptor and ligand.

Parameters:

  • config (BoreschRestraint) –

    The restraint configuration.

  • receptor_atoms (tuple[int, int, int]) –

    The indices of the receptor atoms to restrain.

  • ligand_atoms (tuple[int, int, int]) –

    The indices of the ligand atoms to restrain.

  • coords (Quantity) –

    The coordinates of the full system.

  • ctx_parameter (str | None, default: None ) –

    An optional context parameter to use to scale the strength of the restraint.

Returns:

  • CustomCompoundBondForce

    The restraint force.

Source code in femto/md/restraints.py
def create_boresch_restraint(
    config: femto.md.config.BoreschRestraint,
    receptor_atoms: tuple[int, int, int],
    ligand_atoms: tuple[int, int, int],
    coords: openmm.unit.Quantity,
    ctx_parameter: str | None = None,
) -> openmm.CustomCompoundBondForce:
    """Creates a Boresch restraint force useful in aligning a receptor and ligand.

    Args:
        config: The restraint configuration.
        receptor_atoms: The indices of the receptor atoms to restrain.
        ligand_atoms: The indices of the ligand atoms to restrain.
        coords: The coordinates of the *full* system.
        ctx_parameter: An optional context parameter to use to scale the strength of
            the restraint.

    Returns:
        The restraint force.
    """
    n_particles = 6  # 3 receptor + 3 ligand

    energy_fn = _BORESCH_ENERGY_FN

    if ctx_parameter is not None:
        energy_fn = f"{ctx_parameter} * {energy_fn}"

    force = openmm.CustomCompoundBondForce(n_particles, energy_fn)

    if ctx_parameter is not None:
        force.addGlobalParameter(ctx_parameter, 1.0)

    geometry = _compute_boresch_geometry(receptor_atoms, ligand_atoms, coords)

    parameters = []

    for key, value in [
        ("k_dist_a", config.k_distance),
        ("k_theta_a", config.k_angle_a),
        ("k_theta_b", config.k_angle_b),
        ("k_phi_a", config.k_dihedral_a),
        ("k_phi_b", config.k_dihedral_b),
        ("k_phi_c", config.k_dihedral_c),
        ("dist_0", geometry.dist_0),
        ("theta_a_0", geometry.theta_a_0),
        ("theta_b_0", geometry.theta_b_0),
        ("phi_a_0", geometry.phi_a_0),
        ("phi_b_0", geometry.phi_b_0),
        ("phi_c_0", geometry.phi_c_0),
    ]:
        force.addPerBondParameter(key)
        parameters.append(value.value_in_unit_system(openmm.unit.md_unit_system))

    force.addBond(receptor_atoms + ligand_atoms, parameters)
    force.setUsesPeriodicBoundaryConditions(False)
    force.setName(femto.md.constants.OpenMMForceName.ALIGNMENT_RESTRAINT)
    force.setForceGroup(femto.md.constants.OpenMMForceGroup.ALIGNMENT_RESTRAINT)

    return force