11.2.5. Mathematical helper functions — MDAnalysis.lib.mdamath

Helper functions for common mathematical operations

MDAnalysis.lib.mdamath.normal(vec1, vec2)[source]

Returns the unit vector normal to two vectors.

\[\hat{\mathbf{n}} = \frac{\mathbf{v}_1 \times \mathbf{v}_2}{|\mathbf{v}_1 \times \mathbf{v}_2|}\]

If the two vectors are collinear, the vector \(\mathbf{0}\) is returned.

Changed in version 0.11.0: Moved into lib.mdamath

MDAnalysis.lib.mdamath.norm(v)[source]

Calculate the norm of a vector v.

\[v = \sqrt{\mathbf{v}\cdot\mathbf{v}}\]

This version is faster then numpy.linalg.norm because it only works for a single vector and therefore can skip a lot of the additional fuss linalg.norm does.

Parameters:v (array_like) – 1D array of shape (N) for a vector of length N
Returns:norm of the vector
Return type:float
MDAnalysis.lib.mdamath.angle(a, b)[source]

Returns the angle between two vectors in radians

Changed in version 0.11.0: Moved into lib.mdamath

MDAnalysis.lib.mdamath.dihedral(ab, bc, cd)[source]

Returns the dihedral angle in radians between vectors connecting A,B,C,D.

The dihedral measures the rotation around bc:

  ab
A---->B
       \ bc
       _\'
         C---->D
           cd

The dihedral angle is restricted to the range -π <= x <= π.

New in version 0.8.

Changed in version 0.11.0: Moved into lib.mdamath

MDAnalysis.lib.mdamath.stp(vec1, vec2, vec3)[source]

Takes the scalar triple product of three vectors.

Returns the volume V of the parallel epiped spanned by the three vectors

\[V = \mathbf{v}_3 \cdot (\mathbf{v}_1 \times \mathbf{v}_2)\]

Changed in version 0.11.0: Moved into lib.mdamath

MDAnalysis.lib.mdamath.triclinic_box(x, y, z)[source]

Convert the three triclinic box vectors to [A,B,C,alpha,beta,gamma].

Angles are in degrees.

  • alpha = angle(y,z)
  • beta = angle(x,z)
  • gamma = angle(x,y)

Note

Definition of angles: http://en.wikipedia.org/wiki/Lattice_constant

MDAnalysis.lib.mdamath.triclinic_vectors(dimensions)[source]

Convert [A,B,C,alpha,beta,gamma] to a triclinic box representation.

Original code by Tsjerk Wassenaar posted on the Gromacs mailinglist.

If dimensions indicates a non-periodic system (i.e. all lengths 0) then null-vectors are returned.

Parameters:dimensions ([A, B, C, alpha, beta, gamma]) – list of box lengths and angles (in degrees) such as [A,B,C,alpha,beta,gamma]
Returns:numpy 3x3 array B, with B[0] as first box vector, B[1]``as second vector, ``B[2] as third box vector.
Return type:numpy.array

Notes

The first vector is always pointing along the X-axis, i.e., parallel to (1, 0, 0).

Changed in version 0.7.6: Null-vectors are returned for non-periodic (or missing) unit cell.

MDAnalysis.lib.mdamath.box_volume(dimensions)[source]

Return the volume of the unitcell described by dimensions.

The volume is computed as det(x1,x2,x2) where the xi are the triclinic box vectors from triclinic_vectors().

Parameters:dimensions ([A, B, C, alpha, beta, gamma]) – list of box lengths and angles (in degrees) such as [A,B,C,alpha,beta,gamma]
Returns:volume
Return type:float
MDAnalysis.lib.mdamath.make_whole(atomgroup, reference_atom=None)[source]

Move all atoms in a single molecule so that bonds don’t split over images

Atoms are modified in place.

This function is most useful when atoms have been packed into the primary unit cell, causing breaks mid molecule, with the molecule then appearing on either side of the unit cell. This is problematic for operations such as calculating the center of mass of the molecule.

+-----------+     +-----------+
|           |     |           |
| 6       3 |     |         3 | 6
| !       ! |     |         ! | !
|-5-8   1-2-| ->  |       1-2-|-5-8
| !       ! |     |         ! | !
| 7       4 |     |         4 | 7
|           |     |           |
+-----------+     +-----------+
Parameters:
  • atomgroup (AtomGroup) – The MDAnalysis.core.groups.AtomGroup to work with. The positions of this are modified in place. All these atoms must belong in the same molecule or fragment.
  • reference_atom (Atom) – The atom around which all other atoms will be moved. Defaults to atom 0 in the atomgroup.
Raises:
  • NoDataError – There are no bonds present. (See guess_bonds())
  • ValueError – The algorithm fails to work. This is usually caused by the atomgroup not being a single fragment. (ie the molecule can’t be traversed by following bonds)

Note

Only orthogonal boxes are currently supported.

Example

Make fragments whole:

from MDAnalysis.util.mdamath import make_whole

# This algorithm requires bonds, these can be guessed!
u = mda.Universe(......, guess_bonds=True)

# MDAnalysis can split molecules into their fragments
# based on bonding information.
# Note that this function will only handle a single fragment
# at a time, necessitating a loop.
for frag in u.fragments:
  make_whole(frag)

Alternatively, to keep a single atom in place as the anchor:

# This will mean that atomgroup[10] will NOT get moved,
# and all other atoms will move (if necessary).
make_whole(atomgroup, reference_atom=atomgroup[10])

New in version 0.11.0.

New in version 0.11.0.

MDAnalysis.lib.mdamath.angle(a, b)[source]

Returns the angle between two vectors in radians

Changed in version 0.11.0: Moved into lib.mdamath

MDAnalysis.lib.mdamath.box_volume(dimensions)[source]

Return the volume of the unitcell described by dimensions.

The volume is computed as det(x1,x2,x2) where the xi are the triclinic box vectors from triclinic_vectors().

Parameters:dimensions ([A, B, C, alpha, beta, gamma]) – list of box lengths and angles (in degrees) such as [A,B,C,alpha,beta,gamma]
Returns:volume
Return type:float
MDAnalysis.lib.mdamath.dihedral(ab, bc, cd)[source]

Returns the dihedral angle in radians between vectors connecting A,B,C,D.

The dihedral measures the rotation around bc:

  ab
A---->B
       \ bc
       _\'
         C---->D
           cd

The dihedral angle is restricted to the range -π <= x <= π.

New in version 0.8.

Changed in version 0.11.0: Moved into lib.mdamath

MDAnalysis.lib.mdamath.make_whole(atomgroup, reference_atom=None)[source]

Move all atoms in a single molecule so that bonds don’t split over images

Atoms are modified in place.

This function is most useful when atoms have been packed into the primary unit cell, causing breaks mid molecule, with the molecule then appearing on either side of the unit cell. This is problematic for operations such as calculating the center of mass of the molecule.

+-----------+     +-----------+
|           |     |           |
| 6       3 |     |         3 | 6
| !       ! |     |         ! | !
|-5-8   1-2-| ->  |       1-2-|-5-8
| !       ! |     |         ! | !
| 7       4 |     |         4 | 7
|           |     |           |
+-----------+     +-----------+
Parameters:
  • atomgroup (AtomGroup) – The MDAnalysis.core.groups.AtomGroup to work with. The positions of this are modified in place. All these atoms must belong in the same molecule or fragment.
  • reference_atom (Atom) – The atom around which all other atoms will be moved. Defaults to atom 0 in the atomgroup.
Raises:
  • NoDataError – There are no bonds present. (See guess_bonds())
  • ValueError – The algorithm fails to work. This is usually caused by the atomgroup not being a single fragment. (ie the molecule can’t be traversed by following bonds)

Note

Only orthogonal boxes are currently supported.

Example

Make fragments whole:

from MDAnalysis.util.mdamath import make_whole

# This algorithm requires bonds, these can be guessed!
u = mda.Universe(......, guess_bonds=True)

# MDAnalysis can split molecules into their fragments
# based on bonding information.
# Note that this function will only handle a single fragment
# at a time, necessitating a loop.
for frag in u.fragments:
  make_whole(frag)

Alternatively, to keep a single atom in place as the anchor:

# This will mean that atomgroup[10] will NOT get moved,
# and all other atoms will move (if necessary).
make_whole(atomgroup, reference_atom=atomgroup[10])

New in version 0.11.0.

MDAnalysis.lib.mdamath.norm(v)[source]

Calculate the norm of a vector v.

\[v = \sqrt{\mathbf{v}\cdot\mathbf{v}}\]

This version is faster then numpy.linalg.norm because it only works for a single vector and therefore can skip a lot of the additional fuss linalg.norm does.

Parameters:v (array_like) – 1D array of shape (N) for a vector of length N
Returns:norm of the vector
Return type:float
MDAnalysis.lib.mdamath.normal(vec1, vec2)[source]

Returns the unit vector normal to two vectors.

\[\hat{\mathbf{n}} = \frac{\mathbf{v}_1 \times \mathbf{v}_2}{|\mathbf{v}_1 \times \mathbf{v}_2|}\]

If the two vectors are collinear, the vector \(\mathbf{0}\) is returned.

Changed in version 0.11.0: Moved into lib.mdamath

MDAnalysis.lib.mdamath.one_to_many_pointers(Ni, Nj, i2j)[source]

Based on a many to one mapping of i to j, create the reverse mapping

Parameters:
  • Ni (int) – number of i components
  • Nj (int) – number of j components
  • i2j (numpy.ndarray) – relates i to parent js
Returns:

  • ordered (list) – an ordered list of i indices [size (i,)]
  • ptrs (list) – the start and end index for each j [size (Nj, 2)]

Example

# Residx - the resid of each Atom
ordered, ptrs = one_to_many_pointers(Natoms, Nres, Residx)

# Returns an array of the atom indices that are in resid 7
atoms = ordered[ptrs[7,0]:ptrs[7,1]]
MDAnalysis.lib.mdamath.stp(vec1, vec2, vec3)[source]

Takes the scalar triple product of three vectors.

Returns the volume V of the parallel epiped spanned by the three vectors

\[V = \mathbf{v}_3 \cdot (\mathbf{v}_1 \times \mathbf{v}_2)\]

Changed in version 0.11.0: Moved into lib.mdamath

MDAnalysis.lib.mdamath.triclinic_box(x, y, z)[source]

Convert the three triclinic box vectors to [A,B,C,alpha,beta,gamma].

Angles are in degrees.

  • alpha = angle(y,z)
  • beta = angle(x,z)
  • gamma = angle(x,y)

Note

Definition of angles: http://en.wikipedia.org/wiki/Lattice_constant

MDAnalysis.lib.mdamath.triclinic_vectors(dimensions)[source]

Convert [A,B,C,alpha,beta,gamma] to a triclinic box representation.

Original code by Tsjerk Wassenaar posted on the Gromacs mailinglist.

If dimensions indicates a non-periodic system (i.e. all lengths 0) then null-vectors are returned.

Parameters:dimensions ([A, B, C, alpha, beta, gamma]) – list of box lengths and angles (in degrees) such as [A,B,C,alpha,beta,gamma]
Returns:numpy 3x3 array B, with B[0] as first box vector, B[1]``as second vector, ``B[2] as third box vector.
Return type:numpy.array

Notes

The first vector is always pointing along the X-axis, i.e., parallel to (1, 0, 0).

Changed in version 0.7.6: Null-vectors are returned for non-periodic (or missing) unit cell.