IrregularTimeMeanSquareDisplacement

class IrregularTimeMeanSquareDisplacement(trajectories, time_measurement=None, position_measurement=None, maximum_time=None, atom_selection=None, anisotropy=None, resolution=None, prune_data=None, verbose=None)

Calculate the mean square displacement of a group of atomic positions. It is assumed that the times are unevenly spaced, and thus this algorithm is suited to being used with hypertime or kinetic Monte Carlo trajectories.

Parameters:
  • trajectories (MDTrajectory | sequence of MDTrajectory) – The trajectories to be analyzed

  • time_measurement (str) – The key used for the simulation time in measurements. If None is given the normal molecular dynamics time will be used.

  • position_measurement (str) – The key used for the atom positions in measurements. If None is given the trajectory atomic positions will be used.

  • maximum_time (str) – The maximum time to measure the mean square displacement over. If None is given the average length of the input trajectories is used.

  • atom_selection (None | PeriodicTableElement | str | sequence of PeriodicTableElement or str) – List of either elements or tags used to select atoms from the trajectory. This can only be used when using the configurations in the trajectory.

  • anisotropy (int | sequence of int | None) – The Cartesian directions (x=0, y=1, z=2) to calculate the anisotropic MSD.
    Default: None.

  • resolution (int) – The number of divisions in the mean square displacement histogram.

  • prune_data (int) – Use every nth timestep for the analysis. A value of 1 uses all data.

  • verbose (PhysicalQuantity of type length) – Print the progress of the analysis as it runs.

dataCount()
Returns:

The number of data points used in each time bin.

Return type:

Sequence of int

dataIndices()
Returns:

The indices of the times that are used in calculating the mean square displacement. Times in which no displacements are found are removed from the final histogram.

Return type:

ndarray

displacementStandardDeviation()
Returns:

The standard deviation of the displacement histogram.

Return type:

PhysicalQuantity of type square length

displacementStandardError()
Returns:

The standard error of the displacement histogram.

Return type:

PhysicalQuantity of type square length

displacementVariance()
Returns:

The variance of the displacement histogram.

Return type:

PhysicalQuantity of type quartic length

static fromData(times, positions, maximum_time=None, atom_selection=None, anisotropy=None, resolution=None, prune_data=None, verbose=None)

Create a IrregularTimeMeanSquareDisplacement object from time and position data.

Parameters:
  • times (PhysicalQuantity of type time) – The times of each state.

  • positions (PhysicalQuantity of type distance) – The atomic positions of each state.

  • maximum_time (str) – The maximum time to measure the mean square displacement over. IfNone is given the average length of the input trajectories is used.

  • atom_selection (None | PeriodicTableElement | str | sequence of PeriodicTableElement or str) – List of either elements or tags used to select atoms from the trajectory. This can only be used when using the configurations in the trajectory.

  • anisotropy (int | sequence of int | None) – The Cartesian directions (x=0, y=1, z=2) to calculate the anisotropic MSD.
    Default: None.

  • resolution (int) – The number of divisions in the mean square displacement histogram.

  • prune_data (int) – Use every nth timestep for the analysis. A value of 1 uses all data.

  • verbose (bool) – Print the progress of the analysis as it runs.

Returns:

An IrregularTimeMeanSquareDisplacement object.

Return type:

IrregularTimeMeanSquareDisplacement

meanSquareDisplacement()
Returns:

The mean square displacement histogram.

Return type:

PhysicalQuantity of type square length

times()
Returns:

The times at which the mean square displacement histogram was calculated.

Return type:

PhysicalQuantity of type time

uniqueString()

Return a unique string representing the state of the object.

Usage Examples

In this example the mean square displacement (MSD) is calculated for a trajectory from a hyperdynamics simulation. Such a trajectory might come from the diffusion example in CollectiveVariableHyperdynamics documentation. Here the Hypertime measurement is used for the time data and the tracked_atoms measurement for the position data. The data is also pruned so that only every second data point is used.

msd = IrregularTimeMeanSquareDisplacement(
    trajectories=trajectory,
    time_measurement='Hypertime',
    position_measurement='tracked_atoms',
    prune_data=2,
    verbose=True,
)

Once the MSD is calculated it can be plotted and a line fitted to calculate the diffusion coefficient. This is 1/6th the gradient of the MSD.

# Create the plot model
model = Plot.PlotModel(ps, Angstrom**2)
model.xAxis().setLabel('Time')
model.yAxis().setLabel('Mean Square Displacement')
model.legend().setVisible(True)

# Add the MSD data and a linear fit to the plot
line = Plot.Line(msd.times(), msd.meanSquareDisplacement())
line.setLabel('MSD Data')
line.setColor('mediumpurple')
fit = Plot.PolynomialFit(order=1)
line.addItem(fit)
model.addItem(line)

# Set the range of the fit so that it can be easily adjusted
max_time = msd.times()[-1]
fit.setBounds(0.1 * max_time, 0.9 * max_time)

# Show the plot
model.setLimits()
Plot.show(model)

Notes

The IrregularTimeMeanSquareDisplacement class implements the calculation of the mean square displacement (MSD) of atoms in cases where the time steps between atomic positions are not evenly spaced. These kinds of trajectories can be generated by simulations such as hyperdynamics or KineticMonteCarlo. In hyperdynamics the elapsed time between frames depends on the applied bias potential. Similarly in kinetic Monte Carlo the elapsed time between images depends on the rates to other states from the current state.

In order to calculated the MSD the IrregularTimeMeanSquareDisplacement selects a period of time to calculate the MSD over, and then divides that up into a number of segments. A histogram of atomic square displacements is then constructed by binning each square displacement into the nearest time segment. With a regularly spaced trajectory this can be done exactly by using the time step size as the time segment length. With irregular time steps each step must be rounded into the nearest time segment. From the histogram the average square displacement is calculated at each time segment. Additionally the variance, standard deviation and standard error of the mean are also calculated. These give a measure of the statistical uncertainty in the MSD calculation. These do not account for systematic errors, due largely to insufficient sampling of the of the system.

To calculate the MSD one or more MDTrajectory objects must be given. If multiple trajectories are given then the MSD is averaged over all the trajectories. The atomic time and position data can be taken either directly from the trajectory or from measurements contained within the trajectory. The keyword arguments time_measurement and position_measurement can be used to specify the measurements for the time and position respectively. If not given then the data is taken from the trajectory directly. Regardless of the source, the time and position data must have the same number of entries. If the position data is taken from the trajectory then the atoms used can be selected either by tags or element type. This is done with the keyword argument atom_selection. Calculating the MSD in specific cartesian directions is also possible using the anisotropy keyword argument. If this is not given then the radial MSD is calculated.

A number of options control how the histogram is constructed. The total histgram time can be selected using the maximum_time keyword argument. If not given then the average length of the trajectories is used. Setting this explicitly can be useful in cases where the MSD is calculated from multiple trajectories with substantially different lengths. The number of segments is likewise given with the resolution keyword argument. Formally calculating the MSD scales as \(O(N^2)\), where \(N\) is the number of data points. This can take a significant amount of time if a very large number of data points are used. The prune_data keyword allows giving the spacing between data points. Giving a value of 1 will use all data points, while a value of 2 will only use every second data point. Finally, additional information tracking the progress of the calculation can be printed using the verbose keyword argument. This is useful in estimating the time required to calculate the MSD when using large trajectories or many steps.