restraints

Create OpenMM restraint forces.

Functions:

• create_flat_bottom_restraint –

  Creates a flat bottom position restraint.

• create_position_restraints –

  Creates position restraints for all ligand atoms.

• create_boresch_restraint –

  Creates a 'Boresch' style restraint force, useful for aligning a receptor and

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: Topology,
    mask: str,
    config: FlatBottomRestraint,
) -> CustomExternalForce

Creates position restraints for all ligand atoms.

Parameters:

• topology (Topology) –

  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: mdtop.Topology,
    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_idxs = topology.select(mask)
    assert len(selection_idxs) > 0, "no atoms were found to restrain"

    assert topology.xyz is not None, "topology must have coordinates to restrain"
    xyz = topology.xyz.value_in_unit(openmm.unit.angstrom)

    coords = {i: openmm.Vec3(*xyz[i]) * 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' style restraint force, useful for aligning a receptor and ligand.

It applies to three receptor atoms (r1, r2, r3), and three ligand atoms (l1, l2, l3).

Namely, the following will be restrained: * r3 -- l1 distance. * (θ_a) r2 -- r3 -- l1 angle. * (θ_b) r3 -- l1 -- l2 angle. * (φ_a) r1 -- r2 -- r3 -- l1 dih. * (φ_b) r2 -- r3 -- l1 -- l2 dih. * (φ_c) r3 -- l1 -- l2 -- l3 dih.

Parameters:

• config (BoreschRestraint) –

  The restraint configuration.

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

  The indices of the three receptor atoms to restrain.

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

  The indices of the three 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' style restraint force, useful for aligning a receptor and
    ligand.

    It applies to three receptor atoms (r1, r2, r3), and three ligand atoms
    (l1, l2, l3).

    Namely, the following will be restrained:
        * ``r3`` -- ``l1`` distance.
        * (θ_a) ``r2`` -- ``r3`` -- ``l1`` angle.
        * (θ_b) ``r3`` -- ``l1`` -- ``l2`` angle.
        * (φ_a) ``r1`` -- ``r2`` -- ``r3`` -- ``l1`` dih.
        * (φ_b) ``r2`` -- ``r3`` -- ``l1`` -- ``l2`` dih.
        * (φ_c) ``r3`` -- ``l1`` -- ``l2`` -- ``l3`` dih.

    Args:
        config: The restraint configuration.
        receptor_atoms: The indices of the three receptor atoms to restrain.
        ligand_atoms: The indices of the three 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