Restraints

Restraints are geometry-dependent potentials that are employed to bias simulated systems towards selected structural features. These are used whenever it is desirable to encourage rather than enforce a particular geometry; rigid constraints, by contrast, eliminate the corresponding degrees of freedom entirely. By introducing supplementary energy terms that penalise deviations from specified reference values, such as inter-atomic distances or bond angles, experimental observations or prior theoretical insights may be honoured, enhanced convergence can be promoted, and rare events may be explored within practicable timescales; all whilst permitting natural fluctuations. Typically, restraints are applied only to a limited selection of atoms, leaving the remainder of the model to evolve according to the underlying force-field. In this way, a balance is struck between preserving the natural dynamical character of the system and focusing computational effort on the configurational regions of greatest significance to the study at hand.

Restraints are defined and applied at the tasks layer. Upon execution, each task—such as SP or Opt—automatically iterates over the supplied restraints, evaluating them in sequence. All concrete restraint classes inherit from the abstract base chemsh.base.restraint.Restraint. Accordingly, an instance of the appropriate concrete class is required for every restraint to be used, and a list of these instances should be supplied to the task via the restraints keyword argument.

The chemsh package offers harmonic restraint potentials for the six following geometric properties:

  • Bond: bond distance restraint between two atoms.

  • Angle: bond angle restraint involving three atoms.

  • Torsion: torsion angle restraint defined by four atoms.

  • BondDifference2: difference of two bond distances.

  • BondDifference3: difference of three bond distances.

  • BondDifference4: difference of four bond distances.

A minimum working example is provided below demonstrating the application of three harmonic restraints to a single-point task:

import chemsh
from chemsh import Fragment, NWChem, SP
from chemsh.base.restraint import Bond, Angle, Torsion

fragment = Fragment("some_molecule_file.xyz")

theory = NWChem(frag=fragment, method="hf", basis="3-21g")

restraints = [
    Bond((0, 8), 2.4, 2.0),
    Angle((8, 0, 2), 115.0, 2.0),
    Torsion((0, 2, 3, 1), -50.0, 2.0)]

task = SP(frag=fragment, theory=theory,
          restraints=restraints, gradients=True)

task.run()

Abstract base

class chemsh.base.restraint.Restraint[source]

Bases: ABC

Abstract base class for geometric restraints.

This abstract base class defines the interface method “accumulate()” which all restraint classes must implement. It is through this interface that restraints will be evaluated. The interface method is tasked with 1) evaluating the restraint energy and, optionally, its forces, and 2) accumulate the results into the supplied array(s).

Methods

accumulate(frag, energy[, gradient])

Accumulate the restraint energy and, optionally, its forces.

abstract accumulate(frag: Fragment, energy: ndarray[Any, dtype[float]], gradient: ndarray[Any, dtype[float]] | None = None) None[source]

Accumulate the restraint energy and, optionally, its forces.

Parameters:
  • frag (Fragment) – The molecular fragment whose Cartesian coordinates are to be restrained.

  • energy (NDArray[float]) – A length‑1 array into which the restraint energy is to be added.

  • gradient (Optional[NDArray[float]]) – An optional array, of shape frag.coords, into which the associated gradients are to be added. If no array is provided, as is the default, then the gradient values will not be computed.

Harmonic restraints

Bond

class chemsh.base.restraint.Bond(indices: tuple[int, int], r0: float, k: float)[source]

Bases: Restraint

A harmonic bond-distance restraint between two atoms.

The potential

\[E = \tfrac12 k \left(r - r_{0}\right)^{2}\]

is applied, where r denotes the present inter‑atomic distance, \(r_0\) the target distance, and \(k\) the spring constant.

r0

Equilibrium distance \(r_{0}\) in Bohr.

Type:

float

k

Harmonic spring constant \(k\) in units of \(\text{Hartree}/\text{Bohr}^{2}\).

Type:

float

Attributes:
indices

Atomic indices of the two target atoms.

Methods

accumulate(frag, energy[, gradient])

Evaluate and accumulate the bond-distance restraint.

__init__(indices: tuple[int, int], r0: float, k: float)[source]

Initialise harmonic bond distance restraint instance.

Parameters:
  • indices (Indices2) – A tuple specifying the atomic indices of the two atoms to be restrained.

  • r0 (float) – Equilibrium distance \(r_{0}\) (positive, Bohr).

  • k (float) – Harmonic spring constant \(k\) in units of \(\text{Hartree}/\text{Bohr}^{2}\).

accumulate(frag: Fragment, energy: ndarray[Any, dtype[float]], gradient: ndarray[Any, dtype[float]] | None = None) None[source]

Evaluate and accumulate the bond-distance restraint.

Parameters:
  • frag (Fragment) – The molecular fragment whose Cartesian coordinates are to be restrained.

  • energy (NDArray[float]) – A length‑1 array into which the restraint energy is to be added.

  • gradient (Optional[NDArray[float]]) – An optional array, of shape frag.coords, into which the associated gradients are to be added. If no array is provided, as is the default, then the gradient values will not be computed.

property indices: tuple[int, int]

Atomic indices of the two target atoms.

Angle

class chemsh.base.restraint.Angle(indices: tuple[int, int, int], theta0: float, k: float)[source]

Bases: Restraint

A harmonic bond‑angle restraint involving three atoms.

The potential

\[E = \tfrac12 k (\theta - \theta_{0})^{2}\]

penalises deviations of the bond angle \(\theta\) from the target angle \(\theta_{0}\).

k

Harmonic spring constant \(k\) in units of \(\text{Hartree}/\text{rad}^{2}\).

Type:

float

Attributes:
indices

Atomic indices of the three target atoms.

theta0

Equilibrium angle, \(\theta_{0}\).

Methods

accumulate(frag, energy[, gradient])

Evaluate and accumulate the bond-angle restraint.

__init__(indices: tuple[int, int, int], theta0: float, k: float)[source]

Initialise harmonic bond angle restraint instance.

Parameters:
  • indices (Indices3) – A tuple specifying the atomic indices of the three atoms to be restrained. The bond angle is defined as the angle between the vectors i-j & j-k.

  • theta0 (float) – Target angle \(\theta_{0}\) in degrees, restricted to the interval [0°, 180°].

  • k (float) – Harmonic spring constant \(k\) in units of \(\text{Hartree}/\text{rad}^{2}\).

accumulate(frag: Fragment, energy: ndarray[Any, dtype[float]], gradient: ndarray[Any, dtype[float]] | None = None) None[source]

Evaluate and accumulate the bond-angle restraint.

Parameters:
  • frag (Fragment) – The molecular fragment whose Cartesian coordinates are to be restrained.

  • energy (NDArray[float]) – A length‑1 array into which the restraint energy is to be added.

  • gradient (Optional[NDArray[float]]) – An optional array, of shape frag.coords, into which the associated gradients are to be added. If no array is provided, as is the default, then the gradient values will not be computed.

property indices: tuple[int, int, int]

Atomic indices of the three target atoms.

property theta0: float

Equilibrium angle, \(\theta_{0}\).

Torsion

class chemsh.base.restraint.Torsion(indices: tuple[int, int, int, int], phi0: float, k: float)[source]

Bases: Restraint

A harmonic torsion‑angle restraint involving four atoms.

The potential

\[E = \tfrac12 k (\phi - \phi_{0})^{2}\]

is used, where \(\phi\) is the proper dihedral angle and \(\phi_{0}\) the target value.

k

Harmonic spring constant \(k\) in units of \(\text{Hartree}/\text{rad}^{2}\).

Type:

float

Attributes:
indices

Atomic indices of the four target atoms.

phi0

Equilibrium torsion angle, \(\phi_{0}\).

Methods

accumulate(frag, energy[, gradient])

Evaluate and accumulate the torsion-angle restraint.

torsion_angle(r)

Proper dihedral angle of atoms A–B–C–D.

torsion_angle_and_derivative(r)

Proper dihedral angle and its cartesian gradient.

__init__(indices: tuple[int, int, int, int], phi0: float, k: float)[source]

Initialise harmonic torsion-angle restraint instance.

Parameters:
  • indices (Indices4) – A tuple specifying the atomic indices of the four atoms to be restrained. The torsion angle is defined as the angle between the plains i-j-k & j-k-l.

  • phi0 (float) – Target dihedral \(\phi_{0}\) in degrees, within the interval [-180°, +180°].

  • k (float) – Harmonic spring constant \(k\) in units of \(\text{Hartree}/\text{rad}^{2}\).

accumulate(frag: Fragment, energy: ndarray[Any, dtype[float]], gradient: ndarray[Any, dtype[float]] | None = None) None[source]

Evaluate and accumulate the torsion-angle restraint.

Parameters:
  • frag (Fragment) – The molecular fragment whose Cartesian coordinates are to be restrained.

  • energy (NDArray[float]) – A length‑1 array into which the restraint energy is to be added.

  • gradient (Optional[NDArray[float]]) – An optional array, of shape frag.coords, into which the associated gradients are to be added. If no array is provided, as is the default, then the gradient values will not be computed.

property indices: tuple[int, int, int, int]

Atomic indices of the four target atoms.

property phi0: float

Equilibrium torsion angle, \(\phi_{0}\).

static torsion_angle(r: ndarray[Any, dtype[float]])[source]

Proper dihedral angle of atoms A–B–C–D.

Parameters:

r (NDArray[float]) – A 4×3 array specifying the cartesian coordinates of atoms A, B, C, and D.

Returns:

phi – Dihedral angle in radians, in the range [-π, π].

Return type:

float

static torsion_angle_and_derivative(r: ndarray[Any, dtype[float]])[source]

Proper dihedral angle and its cartesian gradient.

Parameters:

r (NDArray[float]) – A 4×3 array specifying the cartesian coordinates of atoms A, B, C, and D.

Returns:

  • phi (float) – Dihedral angle in radians, in the range [-π, π].

  • dphi_dr (NDArray[float]) – A 4×3 array storing the analytical gradients \(\tfrac{\partial\phi}{r}\), row-wise for A, B, C, and D.

Distance-difference restraints

class chemsh.base.restraint.BondDifference2(indices: tuple[int, int, int, int], r0: float, k: float)[source]

Bases: Restraint

A harmonic restraint on the difference of two bond distances.

The potential of the two bond difference restraint is evaluated as follows:

\[E = \tfrac12 k\left(\left(r_{12}-r_{34}\right)-r_0\right)^2\]

Where \(r_{ij}\) is the distance between the two points \(\bf{r}_i\) and \(\bf{r}_i\), \(k\) is the spring constant, and \(r_0\) is the equilibrium distance difference.

r0

Equilibrium distance difference, in units of Bohr, computed as \(r_{12}-r_{34}\).

Type:

float

k

Harmonic spring constant \(k\) in units of \(\text{Hartree}/\text{Bohr}^{2}\).

Type:

float

Attributes:
indices

Atomic indices of the four target atoms.

Methods

accumulate(frag, energy[, gradient])

Evaluate & accumulate two bond distance difference restraint.

__init__(indices: tuple[int, int, int, int], r0: float, k: float)[source]

Initialise BondDifference2 restraint instance.

Parameters:
  • indices (Indices4) – A tuple specifying the atomic indices of the four target atoms ordered such that the first two indices designate the atoms of the first bond, and the last two the second bond.

  • r0 (float) – Equilibrium distance difference, in units of Bohr, computed as \(r_{12}-r_{34}\).

  • k (float) – Harmonic spring constant \(k\) in units of :math:` ext{Hartree}/ ext{Bohr}^{2}`.

accumulate(frag: Fragment, energy: ndarray[Any, dtype[float]], gradient: ndarray[Any, dtype[float]] | None = None) None[source]

Evaluate & accumulate two bond distance difference restraint.

Parameters:
  • frag (Fragment) – The molecular fragment whose Cartesian coordinates are to be restrained.

  • energy (NDArray[float]) – A length‑1 array into which the restraint energy is to be added.

  • gradient (Optional[NDArray[float]]) – An optional array, of shape frag.coords, into which the associated gradients are to be added. If no array is provided, as is the default, then the gradient values will not be computed.

property indices: tuple[int, int, int, int]

Atomic indices of the four target atoms.

class chemsh.base.restraint.BondDifference3(indices: tuple[int, int, int, int, int, int], r0: float, k: float)[source]

Bases: Restraint

A harmonic restraint on the difference of three bond distances.

The potential of the three bond difference restraint is evaluated as follows:

\[E = \tfrac12 k\left(\left(r_{12}-r_{34}-r_{56}\right)-r_0\right)^2\]

Where \(r_{ij}\) is the distance between the two points \(\bf{r}_i\) and \(\bf{r}_i\), \(k\) is the spring constant, and \(r_0\) is the equilibrium distance difference.

r0

Equilibrium distance difference, in units of Bohr, computed as \(r_{12}-r_{34}-r_{56}\).

Type:

float

k

Harmonic spring constant \(k\) in units of \(\text{Hartree}/\text{Bohr}^{2}\).

Type:

float

Attributes:
indices

Atomic indices of the six target atoms.

Methods

accumulate(frag, energy[, gradient])

Evaluate & accumulate three bond distance difference restraint.

__init__(indices: tuple[int, int, int, int, int, int], r0: float, k: float)[source]

Initialise BondDifference3 restraint instance.

Parameters:
  • indices (Indices6) – A tuple specifying the atomic indices of the six target atoms ordered such that the first two indices designate the atoms of the first bond, the next two the second bond, and the final pair those in the third bond.

  • r0 (float) – Equilibrium distance difference, in units of Bohr, computed as \(r_{12}-r_{34}-r_{56}\).

  • k (float) – Harmonic spring constant \(k\) in units of :math:` ext{Hartree}/ ext{Bohr}^{2}`.

accumulate(frag: Fragment, energy: ndarray[Any, dtype[float]], gradient: ndarray[Any, dtype[float]] | None = None) None[source]

Evaluate & accumulate three bond distance difference restraint.

Parameters:
  • frag (Fragment) – The molecular fragment whose Cartesian coordinates are to be restrained.

  • energy (NDArray[float]) – A length‑1 array into which the restraint energy is to be added.

  • gradient (Optional[NDArray[float]]) – An optional array, of shape frag.coords, into which the associated gradients are to be added. If no array is provided, as is the default, then the gradient values will not be computed.

property indices: tuple[int, int, int, int, int, int]

Atomic indices of the six target atoms.

class chemsh.base.restraint.BondDifference4(indices: tuple[int, int, int, int, int, int, int, int], r0: float, k: float)[source]

Bases: Restraint

A harmonic restraint on the difference of four bond distances.

The potential of the four bond difference restraint is evaluated as follows:

\[E = \tfrac12 k\left(\left(r_{12}+r_{34}-r_{56}-r_{78}\right)-r_0\right)^2\]

Where \(r_{ij}\) is the distance between the two points \(\bf{r}_i\) and \(\bf{r}_i\), \(k\) is the spring constant, and \(r_0\) is the equilibrium distance difference.

r0

Equilibrium distance difference, in units of Bohr, computed as \(r_{12}+r_{34}-r_{56}-r_{78}\).

Type:

float

k

Harmonic spring constant \(k\) in units of \(\text{Hartree}/\text{Bohr}^{2}\).

Type:

float

Attributes:
indices

Atomic indices of the eight target atoms.

Methods

accumulate(frag, energy[, gradient])

Evaluate & accumulate four bond distance difference restraint.

__init__(indices: tuple[int, int, int, int, int, int, int, int], r0: float, k: float)[source]

Initialise BondDifference4 restraint instance.

Parameters:
  • indices (Indices8) – A tuple specifying the atomic indices of the six target atoms ordered such that the first two indices designate the atoms of the first bond, the next two the second bond, and so on.

  • r0 (float) – Equilibrium distance difference, in units of Bohr, computed as \(r_{12}+r_{34}-r_{56}-r_{78}\).

  • k (float) – Harmonic spring constant \(k\) in units of :math:` ext{Hartree}/ ext{Bohr}^{2}`.

accumulate(frag: Fragment, energy: ndarray[Any, dtype[float]], gradient: ndarray[Any, dtype[float]] | None = None) None[source]

Evaluate & accumulate four bond distance difference restraint.

Parameters:
  • frag (Fragment) – The molecular fragment whose Cartesian coordinates are to be restrained.

  • energy (NDArray[float]) – A length‑1 array into which the restraint energy is to be added.

  • gradient (Optional[NDArray[float]]) – An optional array, of shape frag.coords, into which the associated gradients are to be added. If no array is provided, as is the default, then the gradient values will not be computed.

property indices: tuple[int, int, int, int, int, int, int, int]

Atomic indices of the eight target atoms.

Type aliases