12.2.1. Fast distance array computation — MDAnalysis.lib.distances
¶
Fast C-routines to calculate distance arrays from coordinate arrays. Many of the functions also exist in parallel versions, that typically provide higher performance than the serial code. The boolean attribute MDAnalysis.lib.distances.USED_OPENMP can be checked to see if OpenMP was used in the compilation of MDAnalysis.
12.2.1.1. Selection of acceleration (“backend”)¶
All functions take the optional keyword backend, which determines the type of acceleration. Currently, the following choices are implemented (backend is case-insensitive):
backend | module | description |
---|---|---|
“serial” | c_distances |
serial implementation in C/Cython |
“OpenMP” | c_distances_openmp |
parallel implementation in C/Cython with OpenMP |
New in version 0.13.0.
12.2.1.2. Functions¶
-
MDAnalysis.lib.distances.
distance_array
(reference, configuration[, box[, result[, backend]]])[source]¶ Calculate all distances between a reference set and another configuration.
If there are i positions in reference, and j positions in configuration, will calculate a i x j array of distances If an box is supplied then a minimum image convention is used when calculating distances.
If a 2D numpy array of dtype
numpy.float64
with the shape(len(reference), len(configuration))
is provided in result then this preallocated array is filled. This can speed up calculations.Parameters: - reference (numpy.ndarray) – Reference coordinate array of shape
(n, 3)
(dtype
is arbitrary, will be converted todtype=numpy.float32
internally) - configuration (numpy.ndarray) – Configuration coordinate array of shape
(m, 3)
(dtype
is arbitrary, will be converted todtype=numpy.float32
internally) - box (numpy.ndarray or None) – Dimensions of the cell; if provided, the minimum image convention is
applied. The dimensions must be provided in the same format as returned
by by
MDAnalysis.coordinates.base.Timestep.dimensions
:[lx, ly, lz, alpha, beta, gamma]
. - result (numpy.ndarray(dtype=numpy.float64), optional) – Preallocated result array which must have the
shape
(len(ref), len(conf))
anddtype=numpy.float64
. Avoids creating the array which saves time when the function is called repeatedly. [None
] - backend (str) – Select the type of acceleration; “serial” is always available. Other possibilities are “OpenMP” (OpenMP).
Returns: d –
(len(reference),len(configuration))
numpy array with the distancesd[i,j]
between reference coordinates i and configuration coordinates j.Return type: Note
This method is slower than it could be because internally we need to make copies of the ref and conf arrays.
Changed in version 0.13.0: Added backend keyword.
Changed in version 0.19.0: Internal dtype conversion of input coordinates to
numpy.float32
.- reference (numpy.ndarray) – Reference coordinate array of shape
-
MDAnalysis.lib.distances.
self_distance_array
(reference[, box[, result[, backend]]])[source]¶ Calculate all distances within a configuration reference.
If a box is supplied then a minimum image convention is used before calculating distances.
If a 1D numpy array of dtype
numpy.float64
with the shape(N*(N-1)/2)
is provided in result then this preallocated array is filled. This can speed up calculations.Parameters: - reference (numpy.ndarray) – Reference coordinate array with
N=len(ref)
coordinates (dtype
is arbitrary, will be converted todtype=numpy.float32
internally) - box (numpy.ndarray or None) – Dimensions of the cell; if provided, the minimum image convention is
applied. The dimensions must be provided in the same format as returned
by
MDAnalysis.coordinates.base.Timestep.dimensions
:[lx, ly, lz, alpha, beta, gamma]
. - result (numpy.ndarray(dtype=numpy.float64), optional) – Preallocated result array which must have the shape
(N*(N-1)/2,)
and dtypenumpy.float64
. Avoids creating the array which saves time when the function is called repeatedly. [None
] - backend (str) – Select the type of acceleration; “serial” is always available. Other possibilities are “OpenMP” (OpenMP).
Returns: d –
N*(N-1)/2
numpy 1D array with the distances dist[i,j] between ref coordinates i and j at position d[k]. Loop through d:for i in range(N): for j in range(i+1, N): k += 1 dist[i,j] = d[k]
Return type: Note
This method is slower than it could be because internally we need to make copies of the coordinate array.
Changed in version 0.13.0: Added backend keyword.
Changed in version 0.19.0: Internal dtype conversion of input coordinates to
numpy.float32
.- reference (numpy.ndarray) – Reference coordinate array with
-
MDAnalysis.lib.distances.
calc_bonds
(atom1, atom2[, box[, result[, backend]]])[source]¶ Calculate all distances between a pair of atoms. atom1 and atom2 are both arrays of coordinates, where atom1[i] and atom2[i] represent a bond.
In comparison to distance_array and self_distance_array which calculate distances between all combinations of coordinates, calc_bonds can be used to calculate distance between pairs of objects, similar to:
numpy.linalg.norm(a - b) for a, b in zip(coords1, coords2)
The optional argument box applies minimum image convention if supplied. box can be either orthogonal or triclinic
If a 1D numpy array of dtype
numpy.float64
withlen(atom1)
elements is provided in result then this preallocated array is filled. This can speed up calculations.bondlengths = calc_bonds(coords1, coords2 [, box [,result=bondlengths]])
Parameters: - coords1 (array) – An array of coordinates for one half of the bond (
dtype
is arbitrary, will be converted todtype=numpy.float32
internally) - coords2 (array) – An array of coordinates for the other half of bond (
dtype
is arbitrary, will be converted todtype=numpy.float32
internally) - box (array) – The unitcell dimesions for this system.
The dimensions must be provided in the same format as returned
by
MDAnalysis.coordinates.base.Timestep.dimensions
:[lx, ly, lz, alpha, beta, gamma]
. - result (array, optional) – Preallocated result array which must be same length as coord
arrays and
dtype=numpy.float64
. Avoids creating the array which saves time when the function is called repeatedly. [None] - backend (str) – Select the type of acceleration; “serial” is always available. Other possibilities are “OpenMP” (OpenMP).
Returns: bondlengths – The length between each pair in coords1 and coords2
Return type: array
New in version 0.8.
Changed in version 0.13.0: Added backend keyword.
Changed in version 0.19.0: Internal dtype conversion of input coordinates to
numpy.float32
.- coords1 (array) – An array of coordinates for one half of the bond (
-
MDAnalysis.lib.distances.
calc_angles
(atom1, atom2, atom3[, box[, result[, backend]]])[source]¶ Calculates the angle formed between three atoms, over a list of coordinates. All atom inputs are lists of coordinates of equal length, with atom2 representing the apex of the angle.
If a 1D numpy array of dtype
numpy.float64
withlen(atom1)
elements is provided in result then this preallocated array is filled. This can speed up calculations.The optional argument
box
ensures that periodic boundaries are taken into account when constructing the connecting vectors between atoms, ie that the vector between atoms 1 & 2 goes between coordinates in the same image.angles = calc_angles(coords1, coords2, coords3, [[box=None],result=angles])
Parameters: - coords1 (numpy.ndarray) – Coordinate array of one side of angles (
dtype
is arbitrary, will be converted todtype=numpy.float32
internally) - coords2 (numpy.ndarray) – Coordinate array of apex of angles (
dtype
is arbitrary, will be converted todtype=numpy.float32
internally) - coords3 (numpy.ndarray) – Coordinate array of other side of angles (
dtype
is arbitrary, will be converted todtype=numpy.float32
internally) - box (numpy.ndarray, optional) – The unitcell dimesions for this system.
The dimensions must be provided in the same format as returned
by
MDAnalysis.coordinates.base.Timestep.dimensions
:[lx, ly, lz, alpha, beta, gamma]
. - result (numpy.ndarray, optional) – Preallocated result array which must be same length as coord
arrays and
dtype=numpy.float64
. Avoids creating the array which saves time when the function is called repeatedly. [None] - backend (str) – Select the type of acceleration; “serial” is always available. Other possibilities are “OpenMP” (OpenMP).
Returns: angles – An array of angles in radians.
Return type: New in version 0.8.
Changed in version 0.9.0: Added optional box argument to account for periodic boundaries in calculation
Changed in version 0.13.0: Added backend keyword.
Changed in version 0.19.0: Internal dtype conversion of input coordinates to
numpy.float32
.- coords1 (numpy.ndarray) – Coordinate array of one side of angles (
-
MDAnalysis.lib.distances.
calc_dihedrals
(atom1, atom2, atom3, atom4[, box[, result[, backend]]])[source]¶ Calculate the dihedral angle formed by four atoms, over a list of coordinates.
Dihedral angle around axis connecting atoms 1 and 2 (i.e. the angle between the planes spanned by atoms (0,1,2) and (1,2,3)):
3 | 1-----2 / 0
If a 1D numpy array of dtype
numpy.float64
withlen(atom1)
elements is provided in result then this preallocated array is filled. This can speed up calculations.The optional argument
box
ensures that periodic boundaries are taken into account when constructing the connecting vectors between atoms, ie that the vector between atoms 1 & 2 goes between coordinates in the same image:angles = calc_dihedrals(coords1, coords2, coords3, coords4 [,box=box, result=angles])
Parameters: - coords1 (array) – Coordinate array of 1st atom in dihedrals (
dtype
is arbitrary, will be converted todtype=numpy.float32
internally) - coords2 (array) – Coordinate array of 2nd atom in dihedrals (
dtype
is arbitrary, will be converted todtype=numpy.float32
internally) - coords3 (array) – Coordinate array of 3rd atom in dihedrals (
dtype
is arbitrary, will be converted todtype=numpy.float32
internally) - coords4 (array) – Coordinate array of 4th atom in dihedrals (
dtype
is arbitrary, will be converted todtype=numpy.float32
internally) - box (array) – The unitcell dimesions for this system.
The dimensions must be provided in the same format as returned
by
MDAnalysis.coordinates.base.Timestep.dimensions
:[lx, ly, lz, alpha, beta, gamma]
. - result (array, optional) – Preallocated result array which must be same length as coord
arrays and
dtype=numpy.float64
. Avoids creating the array which saves time when the function is called repeatedly. [None] - backend (str) – Select the type of acceleration; “serial” is always available. Other possibilities are “OpenMP” (OpenMP).
Returns: angles – A numpy.array of angles in radians.
Return type: array
New in version 0.8.
Changed in version 0.9.0: Added optional box argument to account for periodic boundaries in calculation
Changed in version 0.11.0: Renamed from calc_torsions to calc_dihedrals
Changed in version 0.13.0: Added backend keyword.
Changed in version 0.19.0: Internal dtype conversion of input coordinates to
numpy.float32
.- coords1 (array) – Coordinate array of 1st atom in dihedrals (
-
MDAnalysis.lib.distances.
apply_PBC
(coordinates, box[, backend])[source]¶ Moves a set of coordinates to all be within the primary unit cell
newcoords = apply_PBC(coords, box)
Parameters: - incoords (numpy.ndarray) – Coordinate array of shape
(n, 3)
(dtype
is arbitrary, will be converted todtype=numpy.float32
internally) - box (array) – The unitcell dimesions for this system; can be either orthogonal or
triclinic information. The dimensions must be provided in the same
format as returned by
MDAnalysis.coordinates.base.Timestep.dimensions
:[lx, ly, lz, alpha, beta, gamma]
. - backend (str) – Select the type of acceleration;
"serial"
is always available. Other possibilities are"OpenMP"
(OpenMP).
Returns: newcoords – Coordinates that are now all within the primary unit cell, as defined by box.
Return type: numpy.ndarray(dtype=numpy.float32)
New in version 0.8.
Changed in version 0.13.0: Added backend keyword.
Changed in version 0.19.0: Internal dtype conversion of input coordinates to
numpy.float32
.- incoords (numpy.ndarray) – Coordinate array of shape
-
MDAnalysis.lib.distances.
capped_distance
(reference, configuration, max_cutoff[, min_cutoff[, box[, method]]])[source]¶ Calculates the pairs and distances within a specified distance
If a box is supplied, then a minimum image convention is used to evaluate the distances.
An automatic guessing of optimized method to calculate the distances is included in the function. An optional keyword for the method is also provided. Users can override the method with this functionality. Currently pkdtree and bruteforce are implemented.
Parameters: - reference (array) – reference coordinates array with shape
reference.shape = (3,)
orreference.shape = (len(reference), 3)
- configuration (array) – Configuration coordinate array with shape
reference.shape = (3,)
orreference.shape = (len(reference), 3)
- max_cutoff (float) – Maximum cutoff distance between the reference and configuration
- min_cutoff ((optional) float) – Minimum cutoff distance between reference and configuration [None]
- box ((optional) array or None) – The dimensions, if provided, must be provided in the same
The unitcell dimesions for this system format as returned
by
MDAnalysis.coordinates.base.Timestep.dimensions
:[lx,ly, lz, alpha, beta, gamma]
. Minimum image convention is applied if the box is provided [None] - method ((optional) 'bruteforce' or 'pkdtree' or 'None') – Keyword to override the automatic guessing of method built-in in the function [None]
Returns: pairs (array) – Pair of indices, one from each reference and configuration such that distance between them is within the
max_cutoff
andmin_cutoff
pairs[i,j] contains the indices i from reference coordinates, and j from configurationdistances (array) – Distances corresponding to each pair of indices. d[k] corresponding to the pairs[i,j] gives the distance between i-th and j-th coordinate in reference and configuration respectively
pairs, distances = capped_distances(reference, coordinates, max_cutoff) for indx, [a,b] in enumerate(pairs): coord1 = reference[a] coord2 = configuration[b] distance = distances[indx]
Note
Currently only supports brute force and Periodic KDtree
See also
:func:’MDAnalysis.lib.distances.distance_array’
See also
:func:’MDAnalysis.lib.pkdtree.PeriodicKDTree’
- reference (array) – reference coordinates array with shape
-
MDAnalysis.lib.distances.
self_capped_distance
(reference, max_cutoff[, min_cutoff[, box[, method]]])[source]¶ Finds all the pairs and respective distances within a specified cutoff for a configuration reference
If a box is supplied, then a minimum image convention is used to evaluate the distances.
An automatic guessing of optimized method to calculate the distances is included in the function. An optional keyword for the method is also provided. Users can override the method with this functionality. Currently pkdtree and bruteforce are implemented.
Parameters: - reference (array) – reference coordinates array with shape
reference.shape = (3,)
orreference.shape = (len(reference), 3)
- max_cutoff (float) – Maximum cutoff distance to check the neighbors with itself
- min_cutoff ((optional) float) – Minimum cutoff distance [None]
- box ((optional) array or None) – The dimensions, if provided, must be provided in the same
The unitcell dimesions for this system format as returned
by
MDAnalysis.coordinates.base.Timestep.dimensions
:[lx,ly, lz, alpha, beta, gamma]
. Minimum image convention is applied if the box is provided [None] - method ((optional) 'bruteforce' or 'pkdtree' or 'None') – Keyword to override the automatic guessing of method built-in in the function [None]
Returns: pairs (array) – Pair of indices such that distance between them is within the
max_cutoff
andmin_cutoff
distances (array) – Distances corresponding to each pair of indices. d[k] corresponding to the pairs[i,j] gives the distance between i-th and j-th coordinate in reference
pairs, distances = self_capped_distances(reference, max_cutoff) for indx, [a,b] in enumerate(pairs): coord1, coords2 = reference[a], reference[b] distance = distances[indx]
Note
Currently only supports brute force and Periodic KDtree
See also
:func:’MDAnalysis.lib.distances.self_distance_array’
See also
:func:’MDAnalysis.lib.pkdtree.PeriodicKDTree’
- reference (array) – reference coordinates array with shape
-
MDAnalysis.lib.distances.
transform_RtoS
(coordinates, box[, backend])[source]¶ Transform an array of coordinates from real space to S space (aka lambda space)
S space represents fractional space within the unit cell for this system
Reciprocal operation to
transform_StoR()
Parameters: - inputcoords (array) – A (3,) coordinate or (n x 3) array of coordinates (
dtype
is arbitrary, will be converted todtype=numpy.float32
internally) - box (array) – The unitcell dimesions for this system.
The dimensions must be provided in the same format as returned
by
MDAnalysis.coordinates.base.Timestep.dimensions
:[lx, ly, lz, alpha, beta, gamma]
. - backend (str, optional) – Select the type of acceleration; “serial” is always available. Other possibilities are “OpenMP” (OpenMP).
Returns: outcoords – A n x 3 array of fractional coordiantes.
Return type: array
Changed in version 0.13.0: Added backend keyword.
Changed in version 0.19.0: Internal dtype conversion of input coordinates to
numpy.float32
.- inputcoords (array) – A (3,) coordinate or (n x 3) array of coordinates (
-
MDAnalysis.lib.distances.
transform_StoR
(coordinates, box[, backend])[source]¶ Transform an array of coordinates from S space into real space.
S space represents fractional space within the unit cell for this system
Reciprocal operation to
transform_RtoS()
Parameters: - inputcoords (array) – A (3,) coordinate or (n x 3) array of coordinates (
dtype
is arbitrary, will be converted todtype=numpy.float32
internally) - box (array) – The unitcell dimesions for this system.
The dimensions must be provided in the same format as returned
by
MDAnalysis.coordinates.base.Timestep.dimensions
:[lx, ly, lz, alpha, beta, gamma]
. - backend (str, optional) – Select the type of acceleration; “serial” is always available. Other possibilities are “OpenMP” (OpenMP).
Returns: outcoords – A n x 3 array of fractional coordiantes.
Return type: array
Changed in version 0.13.0: Added backend keyword.
Changed in version 0.19.0: Internal dtype conversion of input coordinates to
numpy.float32
.- inputcoords (array) – A (3,) coordinate or (n x 3) array of coordinates (
-
MDAnalysis.lib._augment.
augment_coordinates
(coordinates, box, radius)¶ Calculates the relevant images of particles which are within a distance ‘r’ from the box walls
The algorithm works by generating explicit periodic images of interior atoms residing close to any of the six box walls. The steps involved in generating images involves evaluation of reciprocal vectors for the given box vectors followed by calculation of projection distance of atom along the reciprocal vectors. If the distance is less than a specified cutoff distance, relevant periodic images are generated using box translation vectors i.e.
l[a] + m[b] + n[c]
, where[l, m, n]
are the neighbouring cell indices relative to the central cell, and[a, b, c]
are the box vectors. For instance, an atom close toXY
plane containing origin will generate a periodic image outside the central cell and close to the opposite XY plane of the box i.e. at0[a] + 0[b] + 1[c]
. Similarly, if the particle is close to more than one box walls, images along the diagonals are also generated| x x +---------------+ | +---------------+ | | | | | | | | | | | | | | | | o | | x | o | +---------------+ | +---------------+ |
Parameters: - coordinates (array) – Input coordinate array to generate duplicate images in the vicinity of the central cell. All the coordinates must be within the primary unit cell. (dtype = numpy.float32)
- box (array) – Box dimension of shape (6, ). The dimensions must be
provided in the same format as returned
by
MDAnalysis.coordinates.base.Timestep.dimensions
:[lx, ly, lz, alpha, beta, gamma]
(dtype = numpy.float32) - r (float) – thickness of cutoff region for duplicate image generation
Returns: - output (array) – coordinates of duplicate(augmented) particles (dtype = numpy.float32)
- indices (array) – original indices of the augmented coordinates (dtype = numpy.int64)
A map which translates the indices of augmented particles
to their original particle index such that
indices[augmentedindex] = originalindex
Note
Output doesnot return coordinates from the initial array. To merge the particles with their respective images, following operation needs to be superseded after generating the images:
images, mapping = augment_coordinates(coordinates, box, max_cutoff) all_coords = np.concatenate([coordinates, images])
New in version 0.19.0.
-
MDAnalysis.lib._augment.
undo_augment
(indices, translation, nreal)¶ Translate augmented indices back to original indices
Parameters: - results (numpy.ndarray) – indices of coordinates, including “augmented” indices (dtype = numpy.int64)
- translation (numpy.ndarray) – Map to link the augmented indices to the original particle indices
such that
translation[augmentedindex] = originalindex
(dtype = numpy.int64) - nreal (int) – number of real coordinates, i.e. values in results equal or larger than this need to be translated to their real counterpart
Returns: results – modified input results with all the augmented indices translated to their corresponding initial original indices (dtype = numpy.int64)
Return type: Note
Modifies the results array in place
See also
‘MDAnalysis.lib._augment.augment_coordinates’
New in version 0.19.0.
-
MDAnalysis.lib.distances.
apply_PBC
(incoords, box, backend='serial')[source] Moves a set of coordinates to all be within the primary unit cell
newcoords = apply_PBC(coords, box)
Parameters: - incoords (numpy.ndarray) – Coordinate array of shape
(n, 3)
(dtype
is arbitrary, will be converted todtype=numpy.float32
internally) - box (array) – The unitcell dimesions for this system; can be either orthogonal or
triclinic information. The dimensions must be provided in the same
format as returned by
MDAnalysis.coordinates.base.Timestep.dimensions
:[lx, ly, lz, alpha, beta, gamma]
. - backend (str) – Select the type of acceleration;
"serial"
is always available. Other possibilities are"OpenMP"
(OpenMP).
Returns: newcoords – Coordinates that are now all within the primary unit cell, as defined by box.
Return type: numpy.ndarray(dtype=numpy.float32)
New in version 0.8.
Changed in version 0.13.0: Added backend keyword.
Changed in version 0.19.0: Internal dtype conversion of input coordinates to
numpy.float32
.- incoords (numpy.ndarray) – Coordinate array of shape
-
MDAnalysis.lib.distances.
calc_angle
(a, b, c, box=None)[source]¶ Angle (in degrees) between a, b and c, where b is apex of angle
Parameters: - b, c (a,) – single coordinate vectors
- box (numpy.ndarray) – simulation box if given periodic boundary conditions will be applied to the vectors between atoms
New in version 0.18.1.
-
MDAnalysis.lib.distances.
calc_angles
(coords1, coords2, coords3, box=None, result=None, backend='serial')[source] Calculates the angle formed between three atoms, over a list of coordinates. All atom inputs are lists of coordinates of equal length, with atom2 representing the apex of the angle.
If a 1D numpy array of dtype
numpy.float64
withlen(atom1)
elements is provided in result then this preallocated array is filled. This can speed up calculations.The optional argument
box
ensures that periodic boundaries are taken into account when constructing the connecting vectors between atoms, ie that the vector between atoms 1 & 2 goes between coordinates in the same image.angles = calc_angles(coords1, coords2, coords3, [[box=None],result=angles])
Parameters: - coords1 (numpy.ndarray) – Coordinate array of one side of angles (
dtype
is arbitrary, will be converted todtype=numpy.float32
internally) - coords2 (numpy.ndarray) – Coordinate array of apex of angles (
dtype
is arbitrary, will be converted todtype=numpy.float32
internally) - coords3 (numpy.ndarray) – Coordinate array of other side of angles (
dtype
is arbitrary, will be converted todtype=numpy.float32
internally) - box (numpy.ndarray, optional) – The unitcell dimesions for this system.
The dimensions must be provided in the same format as returned
by
MDAnalysis.coordinates.base.Timestep.dimensions
:[lx, ly, lz, alpha, beta, gamma]
. - result (numpy.ndarray, optional) – Preallocated result array which must be same length as coord
arrays and
dtype=numpy.float64
. Avoids creating the array which saves time when the function is called repeatedly. [None] - backend (str) – Select the type of acceleration; “serial” is always available. Other possibilities are “OpenMP” (OpenMP).
Returns: angles – An array of angles in radians.
Return type: New in version 0.8.
Changed in version 0.9.0: Added optional box argument to account for periodic boundaries in calculation
Changed in version 0.13.0: Added backend keyword.
Changed in version 0.19.0: Internal dtype conversion of input coordinates to
numpy.float32
.- coords1 (numpy.ndarray) – Coordinate array of one side of angles (
-
MDAnalysis.lib.distances.
calc_bonds
(coords1, coords2, box=None, result=None, backend='serial')[source] Calculate all distances between a pair of atoms. atom1 and atom2 are both arrays of coordinates, where atom1[i] and atom2[i] represent a bond.
In comparison to distance_array and self_distance_array which calculate distances between all combinations of coordinates, calc_bonds can be used to calculate distance between pairs of objects, similar to:
numpy.linalg.norm(a - b) for a, b in zip(coords1, coords2)
The optional argument box applies minimum image convention if supplied. box can be either orthogonal or triclinic
If a 1D numpy array of dtype
numpy.float64
withlen(atom1)
elements is provided in result then this preallocated array is filled. This can speed up calculations.bondlengths = calc_bonds(coords1, coords2 [, box [,result=bondlengths]])
Parameters: - coords1 (array) – An array of coordinates for one half of the bond (
dtype
is arbitrary, will be converted todtype=numpy.float32
internally) - coords2 (array) – An array of coordinates for the other half of bond (
dtype
is arbitrary, will be converted todtype=numpy.float32
internally) - box (array) – The unitcell dimesions for this system.
The dimensions must be provided in the same format as returned
by
MDAnalysis.coordinates.base.Timestep.dimensions
:[lx, ly, lz, alpha, beta, gamma]
. - result (array, optional) – Preallocated result array which must be same length as coord
arrays and
dtype=numpy.float64
. Avoids creating the array which saves time when the function is called repeatedly. [None] - backend (str) – Select the type of acceleration; “serial” is always available. Other possibilities are “OpenMP” (OpenMP).
Returns: bondlengths – The length between each pair in coords1 and coords2
Return type: array
New in version 0.8.
Changed in version 0.13.0: Added backend keyword.
Changed in version 0.19.0: Internal dtype conversion of input coordinates to
numpy.float32
.- coords1 (array) – An array of coordinates for one half of the bond (
-
MDAnalysis.lib.distances.
calc_dihedral
(a, b, c, d, box=None)[source]¶ Dihedral angle (in degrees) between planes (a, b, c) and (b, c, d)
Parameters: - b, c, d (a,) – single coordinate vectors
- box (numpy.ndarray, optional) – simulation box, if given periodic boundary conditions will be applied
New in version 0.18.1.
-
MDAnalysis.lib.distances.
calc_dihedrals
(coords1, coords2, coords3, coords4, box=None, result=None, backend='serial')[source] Calculate the dihedral angle formed by four atoms, over a list of coordinates.
Dihedral angle around axis connecting atoms 1 and 2 (i.e. the angle between the planes spanned by atoms (0,1,2) and (1,2,3)):
3 | 1-----2 / 0
If a 1D numpy array of dtype
numpy.float64
withlen(atom1)
elements is provided in result then this preallocated array is filled. This can speed up calculations.The optional argument
box
ensures that periodic boundaries are taken into account when constructing the connecting vectors between atoms, ie that the vector between atoms 1 & 2 goes between coordinates in the same image:angles = calc_dihedrals(coords1, coords2, coords3, coords4 [,box=box, result=angles])
Parameters: - coords1 (array) – Coordinate array of 1st atom in dihedrals (
dtype
is arbitrary, will be converted todtype=numpy.float32
internally) - coords2 (array) – Coordinate array of 2nd atom in dihedrals (
dtype
is arbitrary, will be converted todtype=numpy.float32
internally) - coords3 (array) – Coordinate array of 3rd atom in dihedrals (
dtype
is arbitrary, will be converted todtype=numpy.float32
internally) - coords4 (array) – Coordinate array of 4th atom in dihedrals (
dtype
is arbitrary, will be converted todtype=numpy.float32
internally) - box (array) – The unitcell dimesions for this system.
The dimensions must be provided in the same format as returned
by
MDAnalysis.coordinates.base.Timestep.dimensions
:[lx, ly, lz, alpha, beta, gamma]
. - result (array, optional) – Preallocated result array which must be same length as coord
arrays and
dtype=numpy.float64
. Avoids creating the array which saves time when the function is called repeatedly. [None] - backend (str) – Select the type of acceleration; “serial” is always available. Other possibilities are “OpenMP” (OpenMP).
Returns: angles – A numpy.array of angles in radians.
Return type: array
New in version 0.8.
Changed in version 0.9.0: Added optional box argument to account for periodic boundaries in calculation
Changed in version 0.11.0: Renamed from calc_torsions to calc_dihedrals
Changed in version 0.13.0: Added backend keyword.
Changed in version 0.19.0: Internal dtype conversion of input coordinates to
numpy.float32
.- coords1 (array) – Coordinate array of 1st atom in dihedrals (
-
MDAnalysis.lib.distances.
calc_distance
(a, b, box=None)[source]¶ Distance between a and b
Parameters: - b (a,) – single coordinate vectors
- box (numpy.ndarray, optional) – simulation box, if given periodic boundary conditions will be applied
New in version 0.18.1.
-
MDAnalysis.lib.distances.
capped_distance
(reference, configuration, max_cutoff, min_cutoff=None, box=None, method=None)[source] Calculates the pairs and distances within a specified distance
If a box is supplied, then a minimum image convention is used to evaluate the distances.
An automatic guessing of optimized method to calculate the distances is included in the function. An optional keyword for the method is also provided. Users can override the method with this functionality. Currently pkdtree and bruteforce are implemented.
Parameters: - reference (array) – reference coordinates array with shape
reference.shape = (3,)
orreference.shape = (len(reference), 3)
- configuration (array) – Configuration coordinate array with shape
reference.shape = (3,)
orreference.shape = (len(reference), 3)
- max_cutoff (float) – Maximum cutoff distance between the reference and configuration
- min_cutoff ((optional) float) – Minimum cutoff distance between reference and configuration [None]
- box ((optional) array or None) – The dimensions, if provided, must be provided in the same
The unitcell dimesions for this system format as returned
by
MDAnalysis.coordinates.base.Timestep.dimensions
:[lx,ly, lz, alpha, beta, gamma]
. Minimum image convention is applied if the box is provided [None] - method ((optional) 'bruteforce' or 'pkdtree' or 'None') – Keyword to override the automatic guessing of method built-in in the function [None]
Returns: pairs (array) – Pair of indices, one from each reference and configuration such that distance between them is within the
max_cutoff
andmin_cutoff
pairs[i,j] contains the indices i from reference coordinates, and j from configurationdistances (array) – Distances corresponding to each pair of indices. d[k] corresponding to the pairs[i,j] gives the distance between i-th and j-th coordinate in reference and configuration respectively
pairs, distances = capped_distances(reference, coordinates, max_cutoff) for indx, [a,b] in enumerate(pairs): coord1 = reference[a] coord2 = configuration[b] distance = distances[indx]
Note
Currently only supports brute force and Periodic KDtree
See also
:func:’MDAnalysis.lib.distances.distance_array’
See also
:func:’MDAnalysis.lib.pkdtree.PeriodicKDTree’
- reference (array) – reference coordinates array with shape
-
MDAnalysis.lib.distances.
distance_array
(reference, configuration, box=None, result=None, backend='serial')[source] Calculate all distances between a reference set and another configuration.
If there are i positions in reference, and j positions in configuration, will calculate a i x j array of distances If an box is supplied then a minimum image convention is used when calculating distances.
If a 2D numpy array of dtype
numpy.float64
with the shape(len(reference), len(configuration))
is provided in result then this preallocated array is filled. This can speed up calculations.Parameters: - reference (numpy.ndarray) – Reference coordinate array of shape
(n, 3)
(dtype
is arbitrary, will be converted todtype=numpy.float32
internally) - configuration (numpy.ndarray) – Configuration coordinate array of shape
(m, 3)
(dtype
is arbitrary, will be converted todtype=numpy.float32
internally) - box (numpy.ndarray or None) – Dimensions of the cell; if provided, the minimum image convention is
applied. The dimensions must be provided in the same format as returned
by by
MDAnalysis.coordinates.base.Timestep.dimensions
:[lx, ly, lz, alpha, beta, gamma]
. - result (numpy.ndarray(dtype=numpy.float64), optional) – Preallocated result array which must have the
shape
(len(ref), len(conf))
anddtype=numpy.float64
. Avoids creating the array which saves time when the function is called repeatedly. [None
] - backend (str) – Select the type of acceleration; “serial” is always available. Other possibilities are “OpenMP” (OpenMP).
Returns: d –
(len(reference),len(configuration))
numpy array with the distancesd[i,j]
between reference coordinates i and configuration coordinates j.Return type: Note
This method is slower than it could be because internally we need to make copies of the ref and conf arrays.
Changed in version 0.13.0: Added backend keyword.
Changed in version 0.19.0: Internal dtype conversion of input coordinates to
numpy.float32
.- reference (numpy.ndarray) – Reference coordinate array of shape
-
MDAnalysis.lib.distances.
self_capped_distance
(reference, max_cutoff, min_cutoff=None, box=None, method=None)[source] Finds all the pairs and respective distances within a specified cutoff for a configuration reference
If a box is supplied, then a minimum image convention is used to evaluate the distances.
An automatic guessing of optimized method to calculate the distances is included in the function. An optional keyword for the method is also provided. Users can override the method with this functionality. Currently pkdtree and bruteforce are implemented.
Parameters: - reference (array) – reference coordinates array with shape
reference.shape = (3,)
orreference.shape = (len(reference), 3)
- max_cutoff (float) – Maximum cutoff distance to check the neighbors with itself
- min_cutoff ((optional) float) – Minimum cutoff distance [None]
- box ((optional) array or None) – The dimensions, if provided, must be provided in the same
The unitcell dimesions for this system format as returned
by
MDAnalysis.coordinates.base.Timestep.dimensions
:[lx,ly, lz, alpha, beta, gamma]
. Minimum image convention is applied if the box is provided [None] - method ((optional) 'bruteforce' or 'pkdtree' or 'None') – Keyword to override the automatic guessing of method built-in in the function [None]
Returns: pairs (array) – Pair of indices such that distance between them is within the
max_cutoff
andmin_cutoff
distances (array) – Distances corresponding to each pair of indices. d[k] corresponding to the pairs[i,j] gives the distance between i-th and j-th coordinate in reference
pairs, distances = self_capped_distances(reference, max_cutoff) for indx, [a,b] in enumerate(pairs): coord1, coords2 = reference[a], reference[b] distance = distances[indx]
Note
Currently only supports brute force and Periodic KDtree
See also
:func:’MDAnalysis.lib.distances.self_distance_array’
See also
:func:’MDAnalysis.lib.pkdtree.PeriodicKDTree’
- reference (array) – reference coordinates array with shape
-
MDAnalysis.lib.distances.
self_distance_array
(reference, box=None, result=None, backend='serial')[source] Calculate all distances within a configuration reference.
If a box is supplied then a minimum image convention is used before calculating distances.
If a 1D numpy array of dtype
numpy.float64
with the shape(N*(N-1)/2)
is provided in result then this preallocated array is filled. This can speed up calculations.Parameters: - reference (numpy.ndarray) – Reference coordinate array with
N=len(ref)
coordinates (dtype
is arbitrary, will be converted todtype=numpy.float32
internally) - box (numpy.ndarray or None) – Dimensions of the cell; if provided, the minimum image convention is
applied. The dimensions must be provided in the same format as returned
by
MDAnalysis.coordinates.base.Timestep.dimensions
:[lx, ly, lz, alpha, beta, gamma]
. - result (numpy.ndarray(dtype=numpy.float64), optional) – Preallocated result array which must have the shape
(N*(N-1)/2,)
and dtypenumpy.float64
. Avoids creating the array which saves time when the function is called repeatedly. [None
] - backend (str) – Select the type of acceleration; “serial” is always available. Other possibilities are “OpenMP” (OpenMP).
Returns: d –
N*(N-1)/2
numpy 1D array with the distances dist[i,j] between ref coordinates i and j at position d[k]. Loop through d:for i in range(N): for j in range(i+1, N): k += 1 dist[i,j] = d[k]
Return type: Note
This method is slower than it could be because internally we need to make copies of the coordinate array.
Changed in version 0.13.0: Added backend keyword.
Changed in version 0.19.0: Internal dtype conversion of input coordinates to
numpy.float32
.- reference (numpy.ndarray) – Reference coordinate array with
-
MDAnalysis.lib.distances.
transform_RtoS
(inputcoords, box, backend='serial')[source] Transform an array of coordinates from real space to S space (aka lambda space)
S space represents fractional space within the unit cell for this system
Reciprocal operation to
transform_StoR()
Parameters: - inputcoords (array) – A (3,) coordinate or (n x 3) array of coordinates (
dtype
is arbitrary, will be converted todtype=numpy.float32
internally) - box (array) – The unitcell dimesions for this system.
The dimensions must be provided in the same format as returned
by
MDAnalysis.coordinates.base.Timestep.dimensions
:[lx, ly, lz, alpha, beta, gamma]
. - backend (str, optional) – Select the type of acceleration; “serial” is always available. Other possibilities are “OpenMP” (OpenMP).
Returns: outcoords – A n x 3 array of fractional coordiantes.
Return type: array
Changed in version 0.13.0: Added backend keyword.
Changed in version 0.19.0: Internal dtype conversion of input coordinates to
numpy.float32
.- inputcoords (array) – A (3,) coordinate or (n x 3) array of coordinates (
-
MDAnalysis.lib.distances.
transform_StoR
(inputcoords, box, backend='serial')[source] Transform an array of coordinates from S space into real space.
S space represents fractional space within the unit cell for this system
Reciprocal operation to
transform_RtoS()
Parameters: - inputcoords (array) – A (3,) coordinate or (n x 3) array of coordinates (
dtype
is arbitrary, will be converted todtype=numpy.float32
internally) - box (array) – The unitcell dimesions for this system.
The dimensions must be provided in the same format as returned
by
MDAnalysis.coordinates.base.Timestep.dimensions
:[lx, ly, lz, alpha, beta, gamma]
. - backend (str, optional) – Select the type of acceleration; “serial” is always available. Other possibilities are “OpenMP” (OpenMP).
Returns: outcoords – A n x 3 array of fractional coordiantes.
Return type: array
Changed in version 0.13.0: Added backend keyword.
Changed in version 0.19.0: Internal dtype conversion of input coordinates to
numpy.float32
.- inputcoords (array) – A (3,) coordinate or (n x 3) array of coordinates (
12.2.2. Low-level modules for MDAnalysis.lib.distances
¶
MDAnalysis.lib._distances
contains low level access to the
serial MDAnalysis Cython functions in distances. These have
little to no error checking done on inputs so should be used with
caution. Similarly, MDAnalysis.lib._distances_openmp
contains
the OpenMP-enable functions.