Velocity Autocorrelation Function
Velocity Autocorrelation Function —
transport_analysis.velocityautocorr
This module offers a class to efficiently calculate a velocity autocorrelation function (VACF). Averaging over all atoms in the atom group
ag
, regardless of type, it will calculate\[C(j \Delta t) = {1 \over N - j} \sum_{i=0}^{N-1-j} v(i \Delta t)v((i+j)\Delta t)\]where \(N\) is the number of time frames and \(\Delta t\) are discrete time intervals between data points.
Basic usage
This example uses the files
PRM_NCBOX
andTRJ_NCBOX
from the MDAnalysis test suite. To get started, execute# imports >>> import MDAnalysis as mda >>> from transport_analysis.velocityautocorr import VelocityAutocorr # test data for this example >>> from MDAnalysis.tests.datafiles import PRM_NCBOX, TRJ_NCBOXWe will calculate the VACF of an atom group of all water atoms in residues 1-5. To select these atoms:
>>> u = mda.Universe(PRM_NCBOX, TRJ_NCBOX) >>> ag = u.select_atoms("resname WAT and resid 1-5")We can run the calculation using any variable of choice such as
wat_vacf
and access our results usingwat_vacf.results.timeseries
:>>> wat_vacf = VelocityAutocorr(ag).run() >>> wat_vacf.results.timeseries array([275.62075467, -18.42008255, -23.94383428, 41.41415381, -2.3164344 , -35.66393559, -22.66874897, -3.97575003, 6.57888933, -5.29065096])Notice that this example data is insufficient to provide a well-defined VACF. When working with real data, ensure that the frames are captured frequently enough to obtain a VACF suitable for your needs.
- class transport_analysis.velocityautocorr.VelocityAutocorr(atomgroup: AtomGroup, dim_type: str = 'xyz', fft: bool = True, **kwargs)[source]
Class to calculate a velocity autocorrelation function (VACF).
- Parameters
atomgroup (AtomGroup) – An MDAnalysis
AtomGroup
. Note thatUpdatingAtomGroup
instances are not accepted.dim_type ({‘xyz’, ‘xy’, ‘yz’, ‘xz’, ‘x’, ‘y’, ‘z’}) – Desired dimensions to be included in the VACF. Defaults to ‘xyz’.
fft (bool) – If
True
, uses a fast FFT based algorithm for computation of the VACF. Otherwise, use the simple “windowed” algorithm. The tidynamics package is required for fft=True. Defaults toTrue
.- Variables
atomgroup (
AtomGroup
) – The atoms to which this analysis is applieddim_fac (int) – Dimensionality \(d\) of the VACF.
results.timeseries (
numpy.ndarray
) – The averaged VACF over all the particles with respect to lag-time. Obtained after callingVelocityAutocorr.run()
results.vacf_by_particle (
numpy.ndarray
) – The VACF of each individual particle with respect to lag-time.start (Optional[int]) – The first frame of the trajectory used to compute the analysis.
stop (Optional[int]) – The frame to stop at for the analysis, non-inclusive.
step (Optional[int]) – Number of frames to skip between each analyzed frame.
n_frames (int) – Number of frames analysed in the trajectory.
n_particles (int) – Number of particles VACF was calculated over.
- plot_running_integral(start=0, stop=0, step=1, initial=0, xlabel='Time (ps)', ylabel='Running Integral of the VACF (Å^2 / ps)')[source]
Returns a plot of the running integral of the velocity autocorrelation function (VACF) via
Matplotlib
. In this case, the running integral is the integral of the VACF divided by the dimensionality. Analysis must be run prior to plotting.
- Parameters
start (Optional[int]) – The first frame of
self.results.timeseries
used for the plot.stop (Optional[int]) – The frame of
self.results.timeseries
to stop at for the plot, non-inclusive.step (Optional[int]) – Number of frames to skip between each plotted frame.
initial (Optional[float]) – Inserted value at the beginning of the integrated result array. Defaults to 0.
xlabel (Optional[str]) – The x-axis label text. Defaults to “Time (ps)”.
ylabel (Optional[str]) – The y-axis label text. Defaults to “Running Integral of the VACF (Å^2 / ps)”.
- Returns
A
matplotlib.lines.Line2D
instance with the desired VACF plotting information.- Return type
matplotlib.lines.Line2D
- plot_vacf(start=0, stop=0, step=1, xlabel='Time (ps)', ylabel='Velocity Autocorrelation Function (Å^2 / ps^2)')[source]
Returns a velocity autocorrelation function (VACF) plot via
Matplotlib
. Analysis must be run prior to plotting.
- Parameters
start (Optional[int]) – The first frame of
self.results.timeseries
used for the plot.stop (Optional[int]) – The frame of
self.results.timeseries
to stop at for the plot, non-inclusive.step (Optional[int]) – Number of frames to skip between each plotted frame.
xlabel (Optional[str]) – The x-axis label text. Defaults to “Time (ps)”.
ylabel (Optional[str]) – The y-axis label text. Defaults to “Velocity Autocorrelation Function (Å^2 / ps^2)”.
- Returns
A
matplotlib.lines.Line2D
instance with the desired VACF plotting information.- Return type
matplotlib.lines.Line2D
- self_diffusivity_gk(start=0, stop=0, step=1)[source]
Returns a self-diffusivity value using
scipy.integrate.trapezoid
. Analysis must be run prior to computing self-diffusivity.
- Parameters
start (Optional[int]) – The first frame of
self.results.timeseries
used for the calculation.stop (Optional[int]) – The frame of
self.results.timeseries
to stop at for the calculation, non-inclusive.step (Optional[int]) – Number of frames to skip between each frame used for the calculation.
- Returns
The calculated self-diffusivity value for the analysis.
- Return type
numpy.float64
- self_diffusivity_gk_odd(start=0, stop=0, step=1)[source]
Returns a self-diffusivity value using
scipy.integrate.simpson
. Recommended for use with an odd number of evenly spaced data points. Analysis must be run prior to computing self-diffusivity.
- Parameters
start (Optional[int]) – The first frame of
self.results.timeseries
used for the calculation.stop (Optional[int]) – The frame of
self.results.timeseries
to stop at for the calculation, non-inclusive.step (Optional[int]) – Number of frames to skip between each frame used for the calculation.
- Returns
The calculated self-diffusivity value for the analysis.
- Return type
numpy.float64