MDMC.trajectory_analysis.observables package
Submodules
MDMC.trajectory_analysis.observables.concurrency_tools module
Tools for concurrency in observable calculation.
- MDMC.trajectory_analysis.observables.concurrency_tools.core_batch(generator: Generator) Generator[List, None, None] [source]
Turn a generator into a new generator that yields in batches according to the number of available cores, OMP_NUM_THREADS.
Example: >>> core_batch((i for i in range (10)) on 1 core produces [0], [1], [2], [3], [4], [5], [6], [7], [8], [9] on 4 cores produces [0, 1, 2, 3], [4, 5, 6, 7], [8, 9].
- Parameters:
generator (Generator) – The generator to batch.
- Returns:
A generator which iterates in batches of size OMP_NUM_THREADS.
- Return type:
Generator
- MDMC.trajectory_analysis.observables.concurrency_tools.create_executor() ThreadPoolExecutor [source]
Creates a ThreadPoolExecutor with the relevant number of workers (according to the number of cores specified).
- Returns:
A thread pool executor with max_workers=`OMP_NUM_THREADS`
- Return type:
ThreadPoolExecutor
MDMC.trajectory_analysis.observables.fqt module
Module for Intermediate Scattering Function class
- class MDMC.trajectory_analysis.observables.fqt.AbstractFQt[source]
Bases:
SQwMixins
,Observable
An abstract class for total, coherent, and incoherent intermediate scattering functions. Equations used for calculating this are based on Kneller et al. Comput. Phys. Commun. 91 (1995) 191-214.
Note that properties for MD frames & Q are found in the SQwMixins class.
- property FQt
- apply_resolution(resolution: Resolution) FQt [source]
Apply instrument resolution to an FQt object.
- Parameters:
resolution (Resolution) – The Resolution object to apply to FQt.
- Return type:
The FQt object with resolution applied.
- calculate_SQw(energy: list[float], resolution: Resolution = None) ndarray [source]
Calculates S(Q, w) from F(Q, t), accounting for instrument resolution.
In order to obtain
len(energy)
values in energy, we reflect the intermediate scattering function in time to give it dimensions of(len(self.Q), 2 * (len(self.t)) - 2)
. This uses the fact it is even in time, and the number of time points is chosen to be 1 greater than the number of energy points [Rapaport, The Art of Molecular Dynamics Simulation (2nd Edition), 2004, page 142].The numpy implementation of the FFT gives frequencies arranged so that the first
len(energy)
points in the energy dimension correspond to positive frequencies, and the remaining points have negative frequency.- Parameters:
energy (list of floats) – the list of energy (E) points at which S(Q, w) will be calculated.
resolution (Resolution (default None)) – The instrument resolution object which will be applied to FQt.
- Returns:
The S(Q, w) calculated from F(Q, t)
- Return type:
- calculate_from_MD(MD_input: CompactTrajectory, verbose: int = 0, **settings: dict) None [source]
Calculates the intermediate scattering function from a trajectory.
independent_variables
can either be set previously or defined within**settings
.- Parameters:
MD_input (CompactTrajectory) – a single
CompactTrajectory
object.verbose (int, optional) – The level of verbosity: Verbose level 0 gives no information. Verbose level 1 gives final time for the whole method. Verbose level 2 gives final time and also a progress bar. Verbose level 3 gives final time, a progress bar, and time per step.
**settings –
n_Q_vectors
(int)The maximum number of
Q_vectors
for anyQ
value. The greater the number ofQ_vectors
, the more accurate the calculation, but the longer it will take.dimensions
(list, tuple, numpy.ndarray)A 3 element tuple or
array
of float specifying the dimensions of theUniverse
in units ofAng
- calculate_rho_config(config: ndarray, single_Q_vectors: list) ndarray [source]
Calculate density over an entire configuration.
- Parameters:
config (np.ndarray) – An array of atom positions.
single_Q_vectors (list) – A list of Q-vectors for a single value of Q.
- Returns:
An array of rho values for each timestep summed over the atoms in the system, corresponding to each Q
- Return type:
np.ndarray
- property dependent_variables: dict
this is FQt, the intermediate scattering function (in
arb
)- Returns:
The dependent variables
- Return type:
- Type:
Get the dependent variables
- property dependent_variables_structure: dict[str, list]
The order in which the ‘FQt’ dependent variable is indexed in terms of ‘Q’ and ‘t’. Explicitly: we have that self.FQt[Q_index, t_index] is the data point for given indices of self.Q and self.t It also means that: np.shape(self.FQt)=(np.size(self.Q), np.size(self.t))
The purpose of this method is to ensure consistency between different readers/methods which create
FQt
objects.
- property errors: dict
Get or set the errors on the dependent variables, the intermediate scattering function (in
arb
)- Returns:
The errors on the
dependent_variables
- Return type:
- property independent_variables: dict
these are the frequency Q (in
Ang^-1
) and time t (infs
)- Returns:
The independent variables
- Return type:
- Type:
Get or set the independent variables
- property t: ndarray
Get or set the times of the intermediate scattering function in units of
fs
- Returns:
1D array of times in
fs
- Return type:
numpy.array
- property uniformity_requirements: dict[str, dict[str, bool]]
Captures the current limitations on the time ‘t’ and reciprocal lattice points ‘Q’ within the intermediate scattering function
Observables
. If using FFT, then ‘t’ must be uniform and start at zero, otherwise it has no restrictions. ‘Q’ must be uniform but does not need to start at zero.
- class MDMC.trajectory_analysis.observables.fqt.FQt[source]
Bases:
AbstractFQt
A class for containing, calculating and reading the intermediate scattering function for the total dynamic structure factor
- name = 'IntermediateScatteringFunction'
- MDMC.trajectory_analysis.observables.fqt.calc_incoherent_scatt_length(element)[source]
Takes the incoherent scattering cross section of element and calculates the incoherent scattering length to be returned.
- Parameters:
element (string) – a string representing the chemical symbol of the element, e.g ‘He’ for Helium.
- Returns:
Incoherent scattering length of chemical symbol passed in.
- Return type:
- MDMC.trajectory_analysis.observables.fqt.calculate_rho(positions: ndarray, Q_vectors: list) Generator[ndarray, None, None] [source]
Calculates
t
dependent number density in reciprocal space for all Q vectorsAs rho is the sum of the contributions for all of the specified Q vectors, these Q vectors should have the same Q value.
- Parameters:
positions (numpy.ndarray) – An
array
of atomic positions for which the reciprocal space number density should be calculatedQ_vectors (list) – A list of one or more Q vectors with the same Q value
- Returns:
A generator for the reciprocal space number density for each Q-vector.
- Return type:
Generator[np.ndarray]
- MDMC.trajectory_analysis.observables.fqt.get_point_group(dimensions: ndarray) str [source]
Gets the Hermann-Mauguin point group for the universe. Currently, MDMC can only create universes with mutually orthogonal sides, so this method can only produce cubic, tetragonal, or orthorhombic universe point groups.
For tetragonal universes, an additional identifier (x), (y), or (z) is added to identify which of the sides is unique.
- Parameters:
dimensions (numpy.ndarray) – An array of length 3, containing the dimensions of the universe.
- Returns:
The point group symbol for the universe.
- Return type:
MDMC.trajectory_analysis.observables.fqt_coh module
Module for coherent FQt class
- class MDMC.trajectory_analysis.observables.fqt_coh.FQtCoherent[source]
Bases:
AbstractFQt
A class for containing, calculating and reading the intermediate scattering function for the coherent dynamic structure factor
- name = 'FQt_coh'
MDMC.trajectory_analysis.observables.fqt_incoh module
Module for incoherent FQt class
- class MDMC.trajectory_analysis.observables.fqt_incoh.FQtIncoherent[source]
Bases:
AbstractFQt
A class for containing, calculating and reading the intermediate scattering function for the incoherent dynamic structure factor
- name = 'IncoherentIntermediateScatteringFunction'
MDMC.trajectory_analysis.observables.obs module
Module defining a class for storing, calculating and reading in observables from molecular dynamics trajectories.
- class MDMC.trajectory_analysis.observables.obs.Observable[source]
Bases:
ABC
Abstract class that defines methods common to all observable data containers
Observable data can either be from a file or calculated from MD and stored in the data property, along with the associated uncertainty. The bool property
from_MD
states the source of the information.- reader
The file reader for reading experimental data
- Type:
- abstract calculate_from_MD(MD_input: CompactTrajectory | list[CompactTrajectory], verbose: int = 0, **parameters: dict) None [source]
Calculates the observable using input from an MD simulation
- Parameters:
MD_input (Object) – Some input from an MD simulation, commonly a
CompactTrajectory
verbose (int) – Enables verbose printing of the calculation
**parameters – Additional parameters required for calculation specific
Observable
objects
- property data: dict
Get the independent, dependent and error data
- Returns:
The independent, dependent and error data
- Return type:
- abstract property dependent_variables: dict
The dependent variables
- Returns:
The dependent variables
- Return type:
- abstract property dependent_variables_structure: dict
The structure of the dependent variables with respect to the independent variables. Specifically, the order in which the dependent variables are indexed with regards to the independent variables. Example: if dep_var1[indep_var1_index, indep_var2_index, …] = data point for values of the independent_variables with the stated indices then the relevant entry in the returned dict should be: {‘dependent_variable1’: [independent_variable1, independent_variable2, …]} Note that this would also correspond to numpy.shape of the dependent variable being: np.shape(dependent_variable1)=(np.size(independent_variable1), np.size(independent_variable2), …)
The purpose of this method is to ensure that all ``Observable``s of a particular type are created with ‘dependent_variables’ that are consistent regardless of how they were created (e.g. by different ``Reader``s).
- Returns:
The np.shape of the dependent variables
- Return type:
- abstract property errors: dict
The errors on the dependent variables
- Returns:
The errors on the
dependent_variables
- Return type:
- abstract property independent_variables: dict
The independent variables
- Returns:
The independent variables
- Return type:
- abstract maximum_frames() int [source]
The maximum number of
CompactTrajectory
frames that can be used to calculate thedependent_variables
- Returns:
The maximum number of frames
- Return type:
- abstract minimum_frames(dt: float = None) int [source]
The minimum number of
CompactTrajectory
frames needed to calculate thedependent_variables
- property name: str
Get or set the module name that used for factory instantiation
- Returns:
The name of the module in which the
Observable
is located- Return type:
- property origin: str
Get or set the origin of the observable
- Returns:
The origin of the
Observable
, either'experiment'
or'MD'
- Return type:
- read_from_file(reader: str, file_name: str) None [source]
Reads in experimental data from a file using a specified reader
- abstract property uniformity_requirements: dict[str, dict[str, bool]]
Represents the current limitations on
independent_variables
of theObservable
. It captures if theindependent_variables
are required to be uniform or to start at zero The keys of the returned dictionary should be the variables that have such a restriction, with the associated values being a dictionary with booleans if the variables are ‘uniform’ or ‘zeroed’. Variables without any requirements do not need to be included, but can be included. If there are no uniformity requirements it is okay to return ‘None’.
MDMC.trajectory_analysis.observables.obs_factory module
Factory class for generating observables
- class MDMC.trajectory_analysis.observables.obs_factory.ObservableFactory[source]
Bases:
object
Provides a factory for creating an
Observable
. Any module within the observables submodule can be created with a string of the class name, as long as it is a subclass ofObservable
.- classmethod create_observable(name: str) Observable [source]
Creates an
Observable
object from a module nameThe
Observable
object must be registered with theObservableFactory
- Parameters:
name (str) – The name of the module with which the
Observable
is registered- Returns:
An
Observable
object- Return type:
- classmethod get_observable(name: str) Type[Observable] [source]
Gets an
Observable
class from a registry name- Parameters:
name (str) – The name of the module with which the
Observable
is registered- Returns:
A subclass of
Observable
- Return type:
cls
- classmethod register(names: str | Iterable) Callable [source]
A class level decorator for registering Observable classes
The names of the modules with which the Observable is registered should be the parameter passed to the decorator.
- Parameters:
names (str) – The names of the modules with which the Observable is registered
Example
To register the
SQw
class withObservableFactory
:@ObservableFactory.register('SQw') class SQw(Observable):
- registry: dict[str, Observable] = {'CoherentDynamicStructureFactor': <class 'MDMC.trajectory_analysis.observables.sqw.SQwCoherent'>, 'CoherentIntermediateScatteringFunction': <class 'MDMC.trajectory_analysis.observables.fqt_coh.FQtCoherent'>, 'DynamicStructureFactor': <class 'MDMC.trajectory_analysis.observables.sqw.SQw'>, 'FQt': <class 'MDMC.trajectory_analysis.observables.fqt.FQt'>, 'FQtCoh': <class 'MDMC.trajectory_analysis.observables.fqt_coh.FQtCoherent'>, 'FQtCoherent': <class 'MDMC.trajectory_analysis.observables.fqt_coh.FQtCoherent'>, 'FQtIncoherentFQtIncoh': <class 'MDMC.trajectory_analysis.observables.fqt_incoh.FQtIncoherent'>, 'FQt_coh': <class 'MDMC.trajectory_analysis.observables.fqt_coh.FQtCoherent'>, 'FQt_incoh': <class 'MDMC.trajectory_analysis.observables.fqt_incoh.FQtIncoherent'>, 'IncoherentDynamicStructureFactor': <class 'MDMC.trajectory_analysis.observables.sqw.SQwIncoherent'>, 'IncoherentIntermediateScatteringFunction': <class 'MDMC.trajectory_analysis.observables.fqt_incoh.FQtIncoherent'>, 'IntermediateScatteringFunction': <class 'MDMC.trajectory_analysis.observables.fqt.FQt'>, 'PDF': <class 'MDMC.trajectory_analysis.observables.pdf.PairDistributionFunction'>, 'PairDistributionFunction': <class 'MDMC.trajectory_analysis.observables.pdf.PairDistributionFunction'>, 'SQw': <class 'MDMC.trajectory_analysis.observables.sqw.SQw'>, 'SQwCoh': <class 'MDMC.trajectory_analysis.observables.sqw.SQwCoherent'>, 'SQwCoherent': <class 'MDMC.trajectory_analysis.observables.sqw.SQwCoherent'>, 'SQwIncoherentSQwIncoh': <class 'MDMC.trajectory_analysis.observables.sqw.SQwIncoherent'>, 'SQw_coh': <class 'MDMC.trajectory_analysis.observables.sqw.SQwCoherent'>, 'SQw_incoh': <class 'MDMC.trajectory_analysis.observables.sqw.SQwIncoherent'>}
MDMC.trajectory_analysis.observables.pdf module
Module for calculating the total pair distribution function (PDF)
- class MDMC.trajectory_analysis.observables.pdf.PairDistributionFunction[source]
Bases:
Observable
A class for containing, calculating and reading a pair distribution function (PDF).
We derive our definitions for this from the following publication: “A comparison of various commonly used correlation functions for describing total scattering” Keen, D. A. (2001). J. Appl. Cryst. 34, 172-177. DOI: https://doi.org/10.1107/S0021889800019993
We employ the following mathematical form for the total pair distribution function (
PDF
):\[G(r) = \sum_{i,j}^{N_{elements}} c_ic_jb_ib_j(g_{ij}(r) - 1)\]where \(c_i\) is the number concentration of element \(i\), \(b_i\) is the (coherent) scattering length of element \(i\). (This corresponds to equation 8 in the above publication)
The partial pair distribution, \(g_{ij}\), is:
\[g_{ij}(r) = \frac{h_{ij}(r)}{4 \pi r^2 \rho_{j} \Delta{r}}\]where \(h_{ij}`\) is the histogram of distances of \(j\) element atoms around atoms of element \(i\), with bins of size \(\Delta{r}\), and \(\rho_{j}\) is the number density of atoms of element \(j\). As \(g_{ij}(0) = 0\), it is evident that \(G(0) = -\sum_{i,j}^{N_{elements}} c_ic_jb_ib_j\). (This corresponds to equation 10 in the above publication)
The total PDF is contained in
PDF
and the partial pair PDFs (if calculated or imported) are contained inpartial_pdfs
.- property PDF
- property PDF_err
- calculate_from_MD(MD_input: CompactTrajectory, verbose: int = 0, **settings: dict)[source]
Calculate the pair distribution function, \(G(r)`\) from a
CompactTrajectory
The partial pair distribution for a pair i-j, \(g_{ij}\), is:
\[g_{ij}(r) = \frac{h_{ij}(r)}{4 \pi r^2 \rho_{j} \Delta{r}}\]where \(h_{ij}`\) is the histogram of distances of \(j\) element atoms around atoms of element \(i\), with bins of size \(\Delta{r}\), and \(\rho_{j}\) is the number density of atoms of element \(j\). As \(g_{ij}(0) = 0\), it is evident that \(G(0) = -\sum_{i,j}^{N_{elements}} c_ic_jb_ib_j\).
This corresponds to the equation (8) in the following paper: “A comparison of various commonly used correlation functions for
describing total scattering”
Keen, D. A. (2001). J. Appl. Cryst. 34, 172-177. DOI: https://doi.org/10.1107/S0021889800019993
The total pair distribution function (
pdf.PDF
) has the form:\[G(r) = \sum_{i,j}^{N_{elements}} c_ic_jb_ib_j(g_{ij}(r) - 1)\]where \(c_i\) is the proportion of element \(i\) in the material, \(b_i\) is the (coherent) scattering length of element \(i\)
This corresponds to the equation (10) in the above paper.
Independent variables can either be set previously or defined within settings.
A number of frames can be specified, from which the
PDF
and its error are calculated. If the number of frames is too large relative to the run length, the samples will be correlated, which will result in an underestimate of the error.- Parameters:
MD_input (CompactTrajectory) – A single
CompactTrajectory
object.verbose (int) – Verbose print settings. Not currently implemented for PDF.
**settings –
- n_framesint
The number of frames from which the pdf and its error are calculated. These frames are selected uniformly, and the step is taken to be n_frames / total number of frames rounded to the nearest positive integer. If this is not passed, 1% of the total number of frames are used (rounded up to the nearest positive integer).
- use_averagebool
Optional parameter. If set to True then the mean value for PDF is calculated across selected frames from the trajectory. Also, the errors are set to the standard deviation calculated over the multiple frames. If set to False (default), only the last frame of the trajectory is used and n_frames will be ignored.
- subsetlist of tuples
The subset of element pairs from which the PDF is calculated. This can be used to calculate the partial PDFs of a multicomponent system. If this is not passed, all combinations of elements are used i.e. the PDF is the total PDF.
- b_cohdict
A dictionary containing coherent scattering length values for one or more elements. This can be used to calculate the PDF of a system where one or more elements has a coherent scattering length different from the coherent scattering lengths from periodictable.
- r_minfloat
The minimum
r
(atomic separation) in Angstrom for which the PDF will be calculated. If this,r_max
, andr_step
are passed then these will create a range for the independent variabler
, which will overwrite anyr
which has previously been defined. This cannot be passed ifr
is passed.- r_maxfloat
The maximum
r
(atomic separation) in Angstrom for which the PDF will be calculated. If this,r_min
, andr_step
are passed then these will create a range for the independent variabler
, which will overwrite anyr
which has previously been defined. This cannot be passed ifr
is passed.- r_stepfloat
The step size of
r
(atomic separation) for which the PDF will be calculated. If this,r_min
, andr_max
are passed then these will create a range for the independent variabler
, which will overwrite anyr
which has previously been defined. This cannot be passed ifr
is passed.- rnumpy.ndarray
The uniform
r
values in Angstrom for which the PDF will be calculated. This cannot be passed ifr_min
,r_max
, andr_step
are passed.- dimensionsarray-like
A 3 element array-like (list, tuple) with the dimensions of the
Universe
in Angstrom.
Examples
To calculate the O-O partial PDF from a simulation of water, use the subset keyword:
pdf.calculate_from_MD(trajectory, subset=[(O, O)])
To calculate the sum of the H-O and O-O partial PDFs:
pdf.calculate_from_MD(trajectory, subset=[(O, O), (H, O)])
To calculate the total PDF for sodium chloride with 37Cl:
pdf.calculate_from_MD(trajectory, b_coh={‘Cl’:3.08})
To calculate the total PDF for r values of [1., 2., 3., 4.]:
pdf.calculate_from_MD(trajectory, r=[1., 2., 3., 4.])
To calculate the total PDF over an average of 5 frames:
pdf.calculate_from_MD(trajectory, use_average=True, n_frames=5)
- property dependent_variables: dict
these are PDF, the pair distribution function (in
barn
)- Returns:
The dependent variables
- Return type:
- Type:
Get the dependent variables
- property dependent_variables_structure: dict[str, list]
The shape of the ‘PDF’ dependent variable in terms of ‘r’’: np.shape(self.PDF)=(np.size(self.r))
- Returns:
The shape of the PDF dependent variable
- Return type:
- property errors: dict
Get or set the errors on the dependent variables, the pair distribution function (in
arb
)- Returns:
The errors on the
dependent_variables
- Return type:
- property independent_variables: dict
this is the atomic separation distance r (in
Ang
)- Returns:
The independent variables
- Return type:
- Type:
Get or set the independent variable
- maximum_frames() None [source]
There is no hard limit on the number of frames that can be used, so return None
- Return type:
None
- minimum_frames(dt: float = None) int [source]
The minimum number of
CompactTrajectory
frames needed to calculate thedependent_variables
is 1
- name = 'PairDistributionFunction'
MDMC.trajectory_analysis.observables.sqw module
Module for AbstractSQw and total SQw class
- class MDMC.trajectory_analysis.observables.sqw.AbstractSQw[source]
Bases:
SQwMixins
,Observable
An abstract class for total, coherent and incoherent dynamic structure factors. The equations used for calculating this are based on Kneller et al. Comput. Phys. Commun. 91 (1995) 191-214.
Note that properties for MD frames & Q are found in the SQwMixins class.
- property E
- property SQw
- property SQw_err
- static calculate_E(nE: int, dt: float) ndarray [source]
Calculates an array of
nE
uniformly spaced energy values from the time separation of theCompactTrajectory
frames,dt
. The frequencie are determined by the Fast Fourier Transform, as implemented by numpy, for2 * nE
points in time which we then crop to only includenE
positive frequencies. As we are dealing with frequency rather than angular frequency here, the relation to between energy is given by:\[E = h \nu\]- Parameters:
- Returns:
An
array
of float specifying the energy in units ofmeV
- Return type:
- calculate_dt() float [source]
Calculates the time separation of frames required by the experimental dataset, assuming uniform spacing. Note that this may be different from the time separation that the user has given as an input, as it only depends on the current values for
self.E
. The relationship between time and energy comes from the numpy implementation of the FFT for2 * nE
points where:\[\begin{split}\nu_{max} &=& \frac{n_E - 1}{2 n_E \Delta t} \\\\ \therefore \Delta t &=& \frac{h (n_E - 1)}{2 n_E E_{max}}\end{split}\]- Returns:
The time separation required by the current values of
self.E
- Return type:
- calculate_from_MD(MD_input: CompactTrajectory, verbose: int = 0, **settings: dict)[source]
Calculate the dynamic structure factor, S(Q, w) from a
CompactTrajectory
.If the
CompactTrajectory
has more frames than theself.maximum_frames()
that can be used to recreate the grid of energy points, it can slice theCompactTrajectory
into sub-trajectories of lengthself.maximum_frames()
, with the slicing specified through the settingsuse_average
andcont_slicing
.The
independent_variable
Q
can either be set previously or defined within**settings
.- Parameters:
MD_input (CompactTrajectory) – An MDMC
CompactTrajectory
from which to calculateSQw
verbose (int, optional) – The level of verbosity: Verbose level 0 gives no information. Verbose level 1 gives final time for the whole method. Verbose level 2 gives final time and also a progress bar. Verbose level 3 gives final time, a progress bar, and time per step.
**settings –
n_Q_vectors
(int)The maximum number of
Q_vectors
for anyQ
value. The greater the number ofQ_vectors
, the more accurate the calculation, but the longer it will take.dimensions
(list, tuple, numpy.ndarray)A 3 element tuple or
array
of float specifying the dimensions of theUniverse
in units ofAng
- ``energy_resolution` (dict)
Optionally specify energy resolution and function in units of ueV (micro eV), in the format of the one-line dict {‘function’: value}, where function is the resolution function and value is the desired FWHM. e.g. to pass a Gaussian resolution of 80ueV we use {‘gaussian’: 80}. Currently accepted functions are ‘gaussian’ and ‘lorentzian’ Can also be ‘lazily’ given as float, in which case it is assumed to be Gaussian.
Q_values
(array)1D array of Q float (in
Ang^-1
). (optional)use_average
(bool)Optional parameter if a list of more than one
Trajectory
is used. If set to True (default) then the mean value for S(Q, w) is calculated. Also, the errors are set to the standard deviation calculated over the list ofCompactTrajectory
objects.cont_slicing
(bool)Flag to decide between two possible behaviours when the number of
MD_steps
is larger than the minimum required to calculate the observables. IfFalse
(default) then theCompactTrajectory
is sliced into non-overlapping sub-CompactTrajectory
blocks for each of which the observable is calculated. IfTrue
, then theCompactTrajectory
is sliced into as many non-identical sub-CompactTrajectory
blocks as possible (with overlap allowed).
- calculate_resolution_functions(dt: float) dict [source]
Generates a resolution function in momentum and time that can be used in the calculation of SQw. Note that this uses the
SQw
values of theObservable
it is called from, and so should only be called for an observable which has been created from relevant resolution data, i.e. a vanadium sample.Note that if this resolution function is used on data outside its original range, then it will use nearest neighbour extrapolation. Additionally, the input will be reflected in the time/energy domain as symmetry about 0 is assumed. If for whatever reason this is not appropriate for the data in question, this function should not be used.
- Parameters:
dt (float) – The time spacing to use when performing the inverse Fourier transform in units of fs. Ideally this should be the same as the frame separation expected when applying this function.
- Returns:
A dictionary with the key ‘SQw’ corresponding to function which accepts arrays of time and momentum (respectively) and returns a 2D array of values for the instrument resolution.
- Return type:
- property dependent_variables: dict
this is SQw, the dynamic structure factor (in
arb
)- Returns:
The dependent variables
- Return type:
- Type:
Get the dependent variables
- property dependent_variables_structure: dict[str, list]
The order in which the ‘SQw’ dependent variable is indexed in terms of ‘Q’ and ‘E’. Explicitly: we have that self.SQw[Q_index, E_index] is the data point for given indices of self.Q and self.E It also means that: np.shape(self.SQw)=(np.size(self.Q), np.size(self.E))
The purpose of this method is to ensure consistency between different readers/methods which create
SQw
objects.
- property errors: dict
Get or set the errors on the dependent variables, the dynamic structure factor (in
arb
)- Returns:
The errors on the
dependent_variables
- Return type:
- property independent_variables: dict
these are the frequency Q (in
Ang^-1
) and energy E (in``meV``)- Returns:
The independent variables
- Return type:
- Type:
Get or set the independent variables
- property recreated_Q
Get the indices of the recreated Q_values.
- property uniformity_requirements: dict[str, dict[str, bool]]
Captures the current limitations on the energy ‘E’ and reciprocal lattice points ‘Q’ within the dynamic structure factor
Observables
. If using FFT, then ‘E’ must be uniform and start at zero, otherwise it has no restrictions. ‘Q’ must be uniform but does not need to start at zero.
- validate_energy(time_step: float = None)[source]
Asserts that the user set frame separation
dt
leads to energy separation that matches that of the experiment. If not, it changes the time step and trajectory step to fix this. The time step value is prioritised here.- Parameters:
time_step (float, optional) – User specified length of time for each update of the atoms trajectories in the simulation, default is None.
- Returns:
contains a boolean, and two floats (if traj_step and time_step have values) or two NoneType (if traj_step and time_step were passed in as None).
- Return type:
- Raises:
- property w
- class MDMC.trajectory_analysis.observables.sqw.SQw[source]
Bases:
AbstractSQw
A class for the total dynamic structure factor
Calculation is done in the respective FQt object, and this is just a reference to get the correct FQt object.
- name = 'SQw'
- class MDMC.trajectory_analysis.observables.sqw.SQwCoherent[source]
Bases:
AbstractSQw
A class for the coherent dynamic structure factor
- name = 'SQw_coh'
- class MDMC.trajectory_analysis.observables.sqw.SQwIncoherent[source]
Bases:
AbstractSQw
A class for the incoherent dynamic structure factor
- name = 'SQw_incoh'
- class MDMC.trajectory_analysis.observables.sqw.SQwMixins[source]
Bases:
object
A mixin class for properties used by both SQw and FQt objects
- property Q
- maximum_frames() int | None [source]
The maximum number of
CompactTrajectory
frames that can be used to calculate thedependent_variables
depends onself.use_FFT
.If True, it is the number of energy steps + 1, in order to allow for a reflection in time which only counts the end points once.
Otherwise, there is no limit and all frames will contribute to the calculation.
- Returns:
The maximum number of frames
- Return type:
- minimum_frames(dt: float = None) int [source]
The minimum number of
CompactTrajectory
frames needed to calculate thedependent_variables
depends onself.use_FFT
.If self.use_FFT == True, it is the number of energy steps + 1, in order to allow for a reflection in time which only counts the end points once.
If self.use_FFT == False, there is not a hard minimum on number of frames. However, to distinguish our smallest differences in energy \(F(Q,t)\) needs to cover at least a time period \(T_{min}\) such that:
\[T_{min} \sim \frac{h}{\Delta E_{min}}\]Due to the aforementioned reflection in the time domain, to cover a period of \(T_{min}\) we only need \(N\) frames:
\[N = \frac{T_{min}}{2 \Delta t} + 1 = \frac{h}{2 \Delta t \Delta E_{min}} + 1\]
Module contents
Modules calculating observables from molecular dynamics trajectories
Observables (/ indicates alias)
DynamicStructureFactor / SQw CoherentDynamicStructureFactor / SQwCoherent / SQwCoh / SQw_coh IncoherentDynamicStructureFactor / SQwIncoherentSQwIncoh / SQw_incoh IntermediateScatteringFunction / FQt IncoherentIntermediateScatteringFunction / FQtIncoherentFQtIncoh / FQt_incoh PDF / PairDistributionFunction CoherentIntermediateScatteringFunction / FQtCoherent / FQtCoh / FQt_coh
Examples
Observables can be instantiated using the names above. For instance an SQw observable can be instantiated using either of the aliases:
from MDMC.trajectory_analysis import observables
sqw = observables.SQw() # This line...
sqw = observables.DynamicStructureFactor() # ...is equivalent to this line
- MDMC.trajectory_analysis.observables.CoherentDynamicStructureFactor
alias of
SQwCoherent
- MDMC.trajectory_analysis.observables.CoherentIntermediateScatteringFunction
alias of
FQtCoherent
- class MDMC.trajectory_analysis.observables.FQt[source]
Bases:
AbstractFQt
A class for containing, calculating and reading the intermediate scattering function for the total dynamic structure factor
- name = 'IntermediateScatteringFunction'
- MDMC.trajectory_analysis.observables.FQtCoh
alias of
FQtCoherent
- class MDMC.trajectory_analysis.observables.FQtCoherent[source]
Bases:
AbstractFQt
A class for containing, calculating and reading the intermediate scattering function for the coherent dynamic structure factor
- name = 'FQt_coh'
- MDMC.trajectory_analysis.observables.FQtIncoherentFQtIncoh
alias of
FQtIncoherent
- MDMC.trajectory_analysis.observables.FQt_coh
alias of
FQtCoherent
- MDMC.trajectory_analysis.observables.FQt_incoh
alias of
FQtIncoherent
- MDMC.trajectory_analysis.observables.IncoherentDynamicStructureFactor
alias of
SQwIncoherent
- MDMC.trajectory_analysis.observables.IncoherentIntermediateScatteringFunction
alias of
FQtIncoherent
- MDMC.trajectory_analysis.observables.PDF
alias of
PairDistributionFunction
- class MDMC.trajectory_analysis.observables.PairDistributionFunction[source]
Bases:
Observable
A class for containing, calculating and reading a pair distribution function (PDF).
We derive our definitions for this from the following publication: “A comparison of various commonly used correlation functions for describing total scattering” Keen, D. A. (2001). J. Appl. Cryst. 34, 172-177. DOI: https://doi.org/10.1107/S0021889800019993
We employ the following mathematical form for the total pair distribution function (
PDF
):\[G(r) = \sum_{i,j}^{N_{elements}} c_ic_jb_ib_j(g_{ij}(r) - 1)\]where \(c_i\) is the number concentration of element \(i\), \(b_i\) is the (coherent) scattering length of element \(i\). (This corresponds to equation 8 in the above publication)
The partial pair distribution, \(g_{ij}\), is:
\[g_{ij}(r) = \frac{h_{ij}(r)}{4 \pi r^2 \rho_{j} \Delta{r}}\]where \(h_{ij}`\) is the histogram of distances of \(j\) element atoms around atoms of element \(i\), with bins of size \(\Delta{r}\), and \(\rho_{j}\) is the number density of atoms of element \(j\). As \(g_{ij}(0) = 0\), it is evident that \(G(0) = -\sum_{i,j}^{N_{elements}} c_ic_jb_ib_j\). (This corresponds to equation 10 in the above publication)
The total PDF is contained in
PDF
and the partial pair PDFs (if calculated or imported) are contained inpartial_pdfs
.- property PDF
- property PDF_err
- calculate_from_MD(MD_input: CompactTrajectory, verbose: int = 0, **settings: dict)[source]
Calculate the pair distribution function, \(G(r)`\) from a
CompactTrajectory
The partial pair distribution for a pair i-j, \(g_{ij}\), is:
\[g_{ij}(r) = \frac{h_{ij}(r)}{4 \pi r^2 \rho_{j} \Delta{r}}\]where \(h_{ij}`\) is the histogram of distances of \(j\) element atoms around atoms of element \(i\), with bins of size \(\Delta{r}\), and \(\rho_{j}\) is the number density of atoms of element \(j\). As \(g_{ij}(0) = 0\), it is evident that \(G(0) = -\sum_{i,j}^{N_{elements}} c_ic_jb_ib_j\).
This corresponds to the equation (8) in the following paper: “A comparison of various commonly used correlation functions for
describing total scattering”
Keen, D. A. (2001). J. Appl. Cryst. 34, 172-177. DOI: https://doi.org/10.1107/S0021889800019993
The total pair distribution function (
pdf.PDF
) has the form:\[G(r) = \sum_{i,j}^{N_{elements}} c_ic_jb_ib_j(g_{ij}(r) - 1)\]where \(c_i\) is the proportion of element \(i\) in the material, \(b_i\) is the (coherent) scattering length of element \(i\)
This corresponds to the equation (10) in the above paper.
Independent variables can either be set previously or defined within settings.
A number of frames can be specified, from which the
PDF
and its error are calculated. If the number of frames is too large relative to the run length, the samples will be correlated, which will result in an underestimate of the error.- Parameters:
MD_input (CompactTrajectory) – A single
CompactTrajectory
object.verbose (int) – Verbose print settings. Not currently implemented for PDF.
**settings –
- n_framesint
The number of frames from which the pdf and its error are calculated. These frames are selected uniformly, and the step is taken to be n_frames / total number of frames rounded to the nearest positive integer. If this is not passed, 1% of the total number of frames are used (rounded up to the nearest positive integer).
- use_averagebool
Optional parameter. If set to True then the mean value for PDF is calculated across selected frames from the trajectory. Also, the errors are set to the standard deviation calculated over the multiple frames. If set to False (default), only the last frame of the trajectory is used and n_frames will be ignored.
- subsetlist of tuples
The subset of element pairs from which the PDF is calculated. This can be used to calculate the partial PDFs of a multicomponent system. If this is not passed, all combinations of elements are used i.e. the PDF is the total PDF.
- b_cohdict
A dictionary containing coherent scattering length values for one or more elements. This can be used to calculate the PDF of a system where one or more elements has a coherent scattering length different from the coherent scattering lengths from periodictable.
- r_minfloat
The minimum
r
(atomic separation) in Angstrom for which the PDF will be calculated. If this,r_max
, andr_step
are passed then these will create a range for the independent variabler
, which will overwrite anyr
which has previously been defined. This cannot be passed ifr
is passed.- r_maxfloat
The maximum
r
(atomic separation) in Angstrom for which the PDF will be calculated. If this,r_min
, andr_step
are passed then these will create a range for the independent variabler
, which will overwrite anyr
which has previously been defined. This cannot be passed ifr
is passed.- r_stepfloat
The step size of
r
(atomic separation) for which the PDF will be calculated. If this,r_min
, andr_max
are passed then these will create a range for the independent variabler
, which will overwrite anyr
which has previously been defined. This cannot be passed ifr
is passed.- rnumpy.ndarray
The uniform
r
values in Angstrom for which the PDF will be calculated. This cannot be passed ifr_min
,r_max
, andr_step
are passed.- dimensionsarray-like
A 3 element array-like (list, tuple) with the dimensions of the
Universe
in Angstrom.
Examples
To calculate the O-O partial PDF from a simulation of water, use the subset keyword:
pdf.calculate_from_MD(trajectory, subset=[(O, O)])
To calculate the sum of the H-O and O-O partial PDFs:
pdf.calculate_from_MD(trajectory, subset=[(O, O), (H, O)])
To calculate the total PDF for sodium chloride with 37Cl:
pdf.calculate_from_MD(trajectory, b_coh={‘Cl’:3.08})
To calculate the total PDF for r values of [1., 2., 3., 4.]:
pdf.calculate_from_MD(trajectory, r=[1., 2., 3., 4.])
To calculate the total PDF over an average of 5 frames:
pdf.calculate_from_MD(trajectory, use_average=True, n_frames=5)
- property dependent_variables: dict
these are PDF, the pair distribution function (in
barn
)- Returns:
The dependent variables
- Return type:
- Type:
Get the dependent variables
- property dependent_variables_structure: dict[str, list]
The shape of the ‘PDF’ dependent variable in terms of ‘r’’: np.shape(self.PDF)=(np.size(self.r))
- Returns:
The shape of the PDF dependent variable
- Return type:
- property errors: dict
Get or set the errors on the dependent variables, the pair distribution function (in
arb
)- Returns:
The errors on the
dependent_variables
- Return type:
- property independent_variables: dict
this is the atomic separation distance r (in
Ang
)- Returns:
The independent variables
- Return type:
- Type:
Get or set the independent variable
- maximum_frames() None [source]
There is no hard limit on the number of frames that can be used, so return None
- Return type:
None
- minimum_frames(dt: float = None) int [source]
The minimum number of
CompactTrajectory
frames needed to calculate thedependent_variables
is 1
- name = 'PairDistributionFunction'
- class MDMC.trajectory_analysis.observables.SQw[source]
Bases:
AbstractSQw
A class for the total dynamic structure factor
Calculation is done in the respective FQt object, and this is just a reference to get the correct FQt object.
- name = 'SQw'
- MDMC.trajectory_analysis.observables.SQwCoh
alias of
SQwCoherent
- class MDMC.trajectory_analysis.observables.SQwCoherent[source]
Bases:
AbstractSQw
A class for the coherent dynamic structure factor
- name = 'SQw_coh'
- MDMC.trajectory_analysis.observables.SQwIncoherentSQwIncoh
alias of
SQwIncoherent
- MDMC.trajectory_analysis.observables.SQw_coh
alias of
SQwCoherent
- MDMC.trajectory_analysis.observables.SQw_incoh
alias of
SQwIncoherent