dynsight.analysis.spatialaverage

dynsight.analysis.spatialaverage(universe, descriptor_array, selection, r_cut, trajslice=None, n_jobs=1)[source]

Compute spatially averaged descriptor values over neighboring particles.

This function takes a molecular dynamics (MD) simulation stored in an MDAnalysis.Universe object and a NumPy array of descriptor values (such as a physical property for each particle in each frame of the simulation). For each particle in the system, the function calculates the average of its descriptor values with the descriptor values of its neighboring particles within a specified r_cut radius. The calculation is parallelized across multiple processes for efficiency.

Caution

This function utilizes multiprocessing and must be called from within a main() function. To avoid runtime errors, ensure that the script includes the following guard:

if __name__ == "__main__":
    main()

Failure to follow this structure may result in unexpected behavior or crashes, especially on Windows and MacOS.

Important

  • Supports both scalar descriptors (2D) and vector descriptors (3D).

  • Utilizes multiprocessing to speed up the computation.

Parameters:

  • universe (Universe) – The MDAnalysis Universe object containing the molecular dynamics simulation data, including atom positions and trajectory.

  • descriptor_array (NDArray[np.float64]) – NumPy array containing the descriptor values. The array should have dimensions corresponding to either (n_atoms, n_frames) for scalar descriptors, or (n_atoms, n_frames, n_features) for vector descriptors.

  • selection (str) – An atom selection string compatible with MDAnalysis. This defines the subset of atoms for which the spatial averaging will be computed.

  • r_cut (float) – The distance r_cut (in the same units as the trajectory) that defines the neighborhood radius within which particles are considered as neighbors.

  • traj_cut – The number of frames to exclude from the end of the trajectory.

  • n_jobs (int) – The number of processes to use for parallel computation. Warning: Adjust this based on the available cores.

  • trajslice (slice | None)

Returns:

A NumPy array of the same shape as the input descriptor array, containing the spatially averaged descriptor values. If the input array is 2D (n_atoms, n_frames), the output will be a 2D array of the same shape with spatially averaged values. Otherwise, if the input is 3D (n_atoms, n_frames, n_features), the output will also be a 3D array of the same shape with averaged vector values.

Raises:

ValueError – If the input descriptor array does not have 2 or 3 dimensions, an error is raised.

Return type:

NDArray[np.float64]

Example

from dynsight.analysis import spatialaverage
import numpy as np

u = MDAnalysis.Universe('topology.gro', 'trajectory.xtc')
descriptor = np.load('descriptor_array.npy')

averaged_values = spatialaverage(
    universe=u,
    descriptor_array=descriptor,
    selection='name CA',
    r_cut=5.0,
    n_jobs=8,
)

This example computes the spatial averages of the descriptor values for atoms selected as CA atoms, within a r_cut of 5.0 units, using 8 processes in parallel. The result is stored in averaged_values, a NumPy array. All supported input file formats by MDAnalysis can be used to set up the Universe.