kinisi.parser
kinisi.parser#
Parser functions, including implementation for pymatgen
compatible VASP files and MDAnalysis
compatible trajectories.
This parser borrows heavily from the :py:class`pymatgen.analysis.diffusion_analyzer.DiffusionAnalyzer` class, originally authored by Will Richards (wrichard@mit.edu) and Shyue Ping Ong. We include this statement to note that we make no claim to authorship of that code and make no attack on the original authors.
In fact, we love pymatgen!
- class kinisi.parser.Parser(disp, indices, framework_indices, time_step, step_skip, min_obs=30, min_dt=1, ndelta_t=75, memory_limit=8.0, progress=True)[source]#
Bases:
object
The base class for parsing.
- Parameters
disp (
ndarray
) – Displacements of atoms with the shape [site, time step, axis].indices (
ndarray
) – Indices for the atoms in the trajectory used in the diffusion calculation.framework_indices (
ndarray
) – Indices for the atoms in the trajectory that should not be used in the diffusion calculation.time_step (
float
) – Time step, in picoseconds, between steps in trajectory.step_skip (
int
) – Sampling freqency of the trajectory (time_step is multiplied by this number to get the real time between output from the simulation file).min_obs (
int
) – Minimum number of observations of an atom before including it in the MSD vs dt calculation. E.g. If a structure has 10 diffusing atoms, andmin_obs=30
, the MSD vs dt will be calculated up todt = total_run_time / 3
, so that each diffusing atom is measured at least 3 uncorrelated times. Optional, defaults to30
.min_dt (
int
) – Minimum timestep to be evaluated, in the simulation units. Optional, defaults to100
.ndelta_t (
int
) – The number ofdelta_t
values to calculate the MSD over. Optional, defaults to75
.memory_limit (
float
) – Upper limit in the amount of computer memory that the displacements can occupy in gigabytes (GB). Optional, defaults to8.
.progress (
bool
) – Print progress bars to screen. Optional, defaults toTrue
.
- property volume: float#
Volume of system, in cubic angstrom.
- Type
return
- Return type
float
- static get_disp(coords, latt)[source]#
Calculate displacements.
- Parameters
coords (
List
[ndarray
]) – Fractional coordinates for all atoms.latt (
List
[ndarray
]) – Lattice descriptions.
- Return type
ndarray
- Returns
Numpy array of with shape [site, time step, axis] describing displacements.
- static correct_drift(framework_indices, disp)[source]#
Perform drift correction, such that the displacement is calculated normalised to any framework drift.
- Parameters
framework_indices (
ndarray
) – Indices for the atoms in the trajectory that should not be used in the diffusion calculation.disp (
ndarray
) – Numpy array of with shape [site, time step, axis] that describes the displacements.
- Return type
ndarray
- Returns
Displacements corrected to account for drift of a framework.
- smoothed_timesteps(nsteps, min_obs, indices)[source]#
Calculate the smoothed timesteps to be used.
- Parameters
nsteps (
int
) – Number of time steps.min_obs (
int
) – Minimum number of observations to have before including in the MSD vs dt calculation. E.g. If a structure has 10 diffusing atoms, andmin_obs=30
, the MSD vs dt will be calculated up todt = total_run_time / 3
, so that each diffusing atom is measured at least 3 uncorrelated times.indices (
ndarray
) – Indices for the atoms in the trajectory used in the calculation of the diffusion.
- Return type
ndarray
- Returns
Smoothed timesteps.
- get_disps(timesteps, drift_corrected, progress=True)[source]#
Calculate the displacement at each timestep.
- Parameters
timesteps (
ndarray
) – Smoothed timesteps.drift_corrected (
ndarray
) – Drift of framework corrected disp.progress (
bool
) – Print progress bars to screen. Defaults toTrue
.
- Return type
Tuple
[ndarray
,ndarray
]- Returns
Tuple containing: time step intervals and raw displacement.
- class kinisi.parser.PymatgenParser(structures, specie, time_step, step_skip, min_obs=30, sub_sample_traj=1, min_dt=100, ndelta_t=75, memory_limit=8.0, progress=True)[source]#
Bases:
kinisi.parser.Parser
A parser for pymatgen structures.
- Parameters
structures (
List
[pymatgen.core.structure.Structure]) – Structures ordered in sequence of run.specie (
Union
[pymatgen.core.periodic_table.Element, pymatgen.core.periodic_table.Specie]) – Specie to calculate diffusivity for as a String, e.g.'Li'
.time_step (
float
) – Time step, in picoseconds, between steps in trajectory.step_skip (
int
) – Sampling freqency of the trajectory (time_step is multiplied by this number to get the real time between output from the simulation file).min_obs (
int
) – Minimum number of observations of an atom before including it in the MSD vs dt calculation. E.g. If a structure has 10 diffusing atoms, andmin_obs=30
, the MSD vs dt will be calculated up todt = total_run_time / 3
, so that each diffusing atom is measured at least 3 uncorrelated times. Optional, defaults to30
.sub_sample_traj (
int
) – Multiple of thetime_step
to sub sample at. Optional, defaults to1
where all timesteps are used.min_dt (
float
) – Minimum timestep to be evaluated, in the simulation units. Optional, defaults to100
.ndelta_t (
int
) – The number ofdelta_t
values to calculate the MSD over. Optional, defaults to75
.memory_limit (
float
) – Upper limit in the amount of computer memory that the displacements can occupy in gigabytes (GB). Optional, defaults to8.
.progress (
bool
) – Print progress bars to screen. Optional, defaults toTrue
.
- static get_structure_coords_latt(structures, sub_sample_traj=1, progress=True)[source]#
Obtain the initial structure and displacement from a
list
of :py:class`pymatgen.core.structure.Structure`.- Parameters
structures (
List
[pymatgen.core.structure.Structure]) – Structures ordered in sequence of run.sub_sample_traj (
int
) – Multiple of thetime_step
to sub sample at. Optional, default is1
.progress (
bool
) – Print progress bars to screen. Optional, defaults toTrue
.
- Return type
Tuple
[pymatgen.core.structure.Structure,List
[ndarray
],List
[ndarray
]]- Returns
Tuple containing: initial structure, fractional coordinates for all atoms, and lattice descriptions.
- static get_indices(structure, specie)[source]#
Determine framework and non-framework indices for a
pymatgen
compatible file.- Parameters
structure (pymatgen.core.structure.Structure) – Initial structure.
specie (
Union
[pymatgen.core.periodic_table.Element, pymatgen.core.periodic_table.Specie, pymatgen.core.periodic_table.Species]) – Specie to calculate diffusivity for as a String, e.g.'Li'
.
- Return type
Tuple
[ndarray
,ndarray
]- Returns
Tuple containing: indices for the atoms in the trajectory used in the calculation of the diffusion and indices of framework atoms.
- class kinisi.parser.MDAnalysisParser(universe, specie, time_step, step_skip, min_obs=30, sub_sample_atoms=1, sub_sample_traj=1, min_dt=100, ndelta_t=75, memory_limit=8.0, progress=True)[source]#
Bases:
kinisi.parser.Parser
A parser that consumes an MDAnalysis.Universe object.
- Parameters
universe (MDAnalysis.core.universe.Universe) – The MDAnalysis object of interest.
specie (
str
) – Specie to calculate diffusivity for as a String, e.g.'Li'
.time_step (
float
) – Time step, in picoseconds, between steps in trajectory.step_skip (
int
) – Sampling freqency of the trajectory (time_step is multiplied by this number to get the real time between output from the simulation file).min_obs (
int
) – Minimum number of observations of an atom before including it in the MSD vs dt calculation. E.g. If a structure has 10 diffusing atoms, andmin_obs=30
, the MSD vs dt will be calculated up todt = total_run_time / 3
, so that each diffusing atom is measured at least 3 uncorrelated times. Optional, defaults to30
.sub_sample_atoms (
int
) – The sampling rate to sample the atoms in the system. Optional, defaults to1
where all atoms are used.sub_sample_traj (
int
) – Multiple of thetime_step
to sub sample at. Optional, defaults to1
where all timesteps are used.min_dt (
float
) – Minimum timestep to be evaluated, in the simulation units. Optional, defaults to100
.ndelta_t (
int
) – The number ofdelta_t
values to calculate the MSD over. Optional, defaults to75
.memory_limit (
float
) – Upper limit in the amount of computer memory that the displacements can occupy in gigabytes (GB). Optional, defaults to8.
.progress (
bool
) – Print progress bars to screen. Optional, defaults toTrue
.
- static get_structure_coords_latt(universe, sub_sample_atoms=1, sub_sample_traj=1, progress=True)[source]#
Obtain the initial structure and displacement from a
MDAnalysis.universe.Universe
file.- Parameters
universe (MDAnalysis.core.universe.Universe) – Universe for analysis.
sub_sample_atoms (
int
) – Frequency to sub sample the number of atoms. Optional, default is1
.sub_sample_traj (
int
) – Multiple of thetime_step
to sub sample at. Optional, default is1
.progress (
bool
) – Print progress bars to screen. Optional, defaults toTrue
.
- Return type
Tuple
[MDAnalysis.core.groups.AtomGroup,List
[ndarray
],List
[ndarray
],float
]- Returns
Tuple containing: initial structure, fractional coordinates for all atoms, lattice descriptions, and the cell volume
- static get_indices(structure, specie)[source]#
Determine framework and non-framework indices for an
MDAnalysis
compatible file.- Parameters
structure (MDAnalysis.universe.Universe) – Initial structure.
specie (
str
) – Specie to calculate diffusivity for as a String, e.g.'Li'
.
- Return type
Tuple
[ndarray
,ndarray
]- Returns
Tuple containing: indices for the atoms in the trajectory used in the calculation of the diffusion and indices of framework atoms.