API Documentation
This section contains an overview of the waterdynamics API.
- class waterdynamics.AngularDistribution(universe, select, bins=40, nproc=1, axis='z')[source]
Angular distribution function analysis
The angular distribution function (AD) is defined as the distribution probability of the cosine of the \(\theta\) angle formed by the OH vector, HH vector or dipolar vector of water molecules and a vector \(\hat n\) parallel to chosen axis (z is the default value). The cosine is define as \(\cos \theta = \hat u \cdot \hat n\), where \(\hat u\) is OH, HH or dipole vector. The AD is also know as Angular Probability (AP).
- Parameters:
universe (Universe) – Universe object
select (str) – Selection string to evaluate its angular distribution [‘byres name OH2’]
bins (int (optional)) – Number of bins to create the histogram by means of
numpy.histogram()
axis ({‘x’, ‘y’, ‘z’} (optional)) – Axis to create angle with the vector (HH, OH or dipole) and calculate cosine theta [‘z’].
- class waterdynamics.WaterOrientationalRelaxation(universe, select, t0, tf, dtmax, nproc=1, order=2)[source]
Water orientation relaxation analysis
Function to evaluate the Water Orientational Relaxation proposed by Yu-ling Yeh and Chung-Yuan Mou [1]. WaterOrientationalRelaxation indicates “how fast” water molecules are rotating or changing direction. This is a time correlation function given by:
\[C_{\hat u}(\tau)=\langle \mathit{P}_2[\mathbf{\hat{u}}(t_0)\cdot\mathbf{\hat{u}}(t_0+\tau)]\rangle\]where \(P_2=(3x^2-1)/2\) is the second-order Legendre polynomial and \(\hat{u}\) is a unit vector along HH, OH or dipole vector. Another option is to select the first-order Legendre polynomial, \(P_1=x\).
- Parameters:
universe (Universe) – Universe object
selection (str) – Selection string for water [‘byres name OH2’].
t0 (int) – frame where analysis begins
tf (int) – frame where analysis ends
dtmax (int) – Maximum dt size, dtmax < tf or it will crash.
order (1 or 2 (default)) – first- or second-order Legendre polynomial
- class waterdynamics.MeanSquareDisplacement(universe, select, t0, tf, dtmax, nproc=1)[source]
Mean square displacement analysis
Function to evaluate the Mean Square Displacement (MSD). The MSD gives the average distance that particles travels. The MSD is given by:
\[\langle\Delta r(t)^2\rangle = 2nDt\]where \(r(t)\) is the position of particle in time \(t\), \(\Delta r(t)\) is the displacement after time lag \(t\), \(n\) is the dimensionality, in this case \(n=3\), \(D\) is the diffusion coefficient and \(t\) is the time.
- Parameters:
universe (Universe) – Universe object
select (str) – Selection string for water [‘byres name OH2’].
t0 (int) – frame where analysis begins
tf (int) – frame where analysis ends
dtmax (int) – Maximum dt size, dtmax < tf or it will crash.
- class waterdynamics.SurvivalProbability(universe, select, verbose=False)[source]
Survival Probability (SP) gives the probability for a group of particles to remain in a certain region. The SP is given by:
\[P(\tau) = \langle \frac{ N(t, t + \tau )} { N(t) }\rangle\]where \(\tau\) is the timestep, \(N(t)\) the number of particles at time \(t\), and \(N(t, t+\tau)\) is the number of particles at every frame from \(t\) to \(t + \tau\). The angular brackets represent an average over all time origins, \(t\). See
MDAnalysis.lib.correlations.autocorrelation()
for technical details.- Parameters:
universe (Universe) – Universe object
select (str) – Selection string; any selection is allowed. With this selection you define the region/zone where to analyze, e.g.: “resname SOL and around 5 (resid 10)”.
verbose (Boolean, optional) – When True, prints progress and comments to the console.
Notes
Currently
SurvivalProbability
is the only one inwaterdynamics
to support an exclusive behaviour (i.e. similar to the current behaviour ofAnalysisBase
to the stop keyword passed toSurvivalProbability.run()
. Unlike otherwaterdynamics
final frame definitions which are inclusive.- run(tau_max=20, start=None, stop=None, step=None, residues=False, intermittency=0, verbose=False)[source]
Computes and returns the Survival Probability (SP) timeseries
- Parameters:
start (int, optional) – Zero-based index of the first frame to be analysed, Default: None (first frame).
stop (int, optional) – Zero-based index of the last frame to be analysed (exclusive), Default: None (last frame).
step (int, optional) – Jump every step-th frame. This is compatible but independant of the taus used, and it is good to consider using the step equal to tau_max to remove the overlap. Note that step and tau_max work consistently with intermittency. Default: None (use every frame).
tau_max (int, optional) – Survival probability is calculated for the range 1 <= tau <= tau_max.
residues (Boolean, optional) – If true, the analysis will be carried out on the residues (.resids) rather than on atom (.ids). A single atom is sufficient to classify the residue as within the distance.
intermittency (int, optional) – The maximum number of consecutive frames for which an atom can leave but be counted as present if it returns at the next frame. An intermittency of 0 is equivalent to a continuous survival probability, which does not allow for the leaving and returning of atoms. For example, for intermittency=2, any given atom may leave a region of interest for up to two consecutive frames yet be treated as being present at all frames. The default is continuous (0).
verbose (Boolean, optional) – Print the progress to the console.
- Returns:
tau_timeseries (list) – tau from 1 to tau_max. Saved in the field tau_timeseries.
sp_timeseries (list) – survival probability for each value of tau. Saved in the field sp_timeseries.
sp_timeseries_data (list) – raw datapoints from which the average is taken (sp_timeseries). Time dependancy and distribution can be extracted.