kinisi.parser#
Parser functions, including implementation for pymatgen
compatible VASP files and MDAnalysis
compatible trajectories.
This parser borrows heavily from the 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_dt=None, max_dt=None, n_steps=100, spacing='linear', sampling='multi-origin', memory_limit=8.0, progress=True)#
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_dt (
Optional
[float
]) – Minimum timestep to be evaluated, in the simulation units. Optional, defaults to the produce oftime_step
andstep_skip
.max_dt (
Optional
[float
]) – Maximum timestep to be evaluated, in the simulation units. Optional, defaults to the length of the simulation.n_steps (
int
) – Number of steps to be used in the timestep function. Optional, defaults to100
unless this is fewer than the total number of steps in the trajectory when it defaults to this number.spacing (
str
) – The spacing of the steps that define the timestep, can be either'linear'
or'logarithmic'
. If'logarithmic'
the number of steps will be less than or equal to that in then_steps
argument, such that all values are unique. Optional, defaults tolinear
.sampling (
str
) – The ways that the time-windows are sampled. The options are'single-origin'
or'multi-origin'
with the former resulting in only one observation per atom per time-window and the latter giving the maximum number of trajectories. Optional, defaults to'multi-origin'
.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#
- Returns:
Volume of system, in cubic angstrom.
- static get_disp(coords, latt)#
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)#
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.
- get_timesteps(n_steps, spacing)#
Calculate the smoothed timesteps to be used.
- Parameters:
n_steps (
int
) – Number of time steps.step_spacing –
- Return type:
ndarray
- Returns:
Smoothed timesteps.
- get_disps(timesteps, drift_corrected, progress=True)#
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
,List
[ndarray
],ndarray
]- Returns:
Tuple containing: time step intervals and raw displacement.
- class kinisi.parser.ASEParser(atoms, specie, time_step, step_skip, sub_sample_traj=1, min_dt=None, max_dt=None, n_steps=100, spacing='linear', sampling='multi-origin', memory_limit=8.0, progress=True)#
Bases:
Parser
A parser for ASE Atoms objects
- Parameters:
atoms (
List
[ase.atoms.Atoms]) – Atoms ordered in sequence of run.specie (
str
) – symbol 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 frequency of the trajectory (time_step is multiplied by this number to get the real time between output from the simulation file).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 to the produce oftime_step
andstep_skip
.max_dt (
float
) – Maximum timestep to be evaluated, in the simulation units. Optional, defaults to the length of the simulation.n_steps (
int
) – Number of steps to be used in the timestep function. Optional, defaults to100
unless this is fewer than the total number of steps in the trajectory when it defaults to this number.spacing (
str
) – The spacing of the steps that define the timestep, can be either'linear'
or'logarithmic'
. If'logarithmic'
the number of steps will be less than or equal to that in then_steps
argument, such that all values are unique. Optional, defaults tolinear
.sampling (
str
) – The ways that the time-windows are sampled. The options are'single-origin'
or'multi-origin'
with the former resulting in only one observation per atom per time-window and the latter giving the maximum number of origins without sampling overlapping trajectories. Optional, defaults to'multi-origin'
.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(atoms, sub_sample_traj=1, progress=True)#
Obtain the initial structure and displacement from a
list
of :py:class`pymatgen.core.structure.Structure`.- Parameters:
structures – 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
[ase.atoms.Atoms,List
[ndarray
],List
[ndarray
]]- Returns:
Tuple containing: initial structure, fractional coordinates for all atoms, and lattice descriptions.
- static get_indices(structure, specie)#
Determine framework and non-framework indices for a
pymatgen
compatible file.- Parameters:
structure (ase.atoms.Atoms) – 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.
- class kinisi.parser.PymatgenParser(structures, specie, time_step, step_skip, sub_sample_traj=1, min_dt=None, max_dt=None, n_steps=100, spacing='linear', sampling='multi-origin', memory_limit=8.0, progress=True)#
Bases:
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).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 to the produce oftime_step
andstep_skip
.max_dt (
float
) – Maximum timestep to be evaluated, in the simulation units. Optional, defaults to the length of the simulation.n_steps (
int
) – Number of steps to be used in the timestep function. Optional, defaults to100
unless this is fewer than the total number of steps in the trajectory when it defaults to this number.spacing (
str
) – The spacing of the steps that define the timestep, can be either'linear'
or'logarithmic'
. If'logarithmic'
the number of steps will be less than or equal to that in then_steps
argument, such that all values are unique. Optional, defaults tolinear
.sampling (
str
) – The ways that the time-windows are sampled. The options are'single-origin'
or'multi-origin'
with the former resulting in only one observation per atom per time-window and the latter giving the maximum number of origins without sampling overlapping trajectories. Optional, defaults to'multi-origin'
.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)#
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)#
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,List
[pymatgen.core.periodic_table.Element],List
[pymatgen.core.periodic_table.Specie],List
[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, sub_sample_atoms=1, sub_sample_traj=1, min_dt=None, max_dt=None, n_steps=100, spacing='linear', sampling='multi-origin', memory_limit=8.0, progress=True)#
Bases:
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).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 to the produce oftime_step
andstep_skip
.max_dt (
float
) – Maximum timestep to be evaluated, in the simulation units. Optional, defaults to the length of the simulation.n_steps (
int
) – Number of steps to be used in the timestep function. Optional, defaults to100
unless this is fewer than the total number of steps in the trajectory when it defaults to this number.spacing (
str
) – The spacing of the steps that define the timestep, can be either'linear'
or'logarithmic'
. If'logarithmic'
the number of steps will be less than or equal to that in then_steps
argument, such that all values are unique. Optional, defaults tolinear
.sampling (
str
) – The ways that the time-windows are sampled. The options are'single-origin'
or'multi-origin'
with the former resulting in only one observation per atom per time-window and the latter giving the maximum number of origins without sampling overlapping trajectories. Optional, defaults to'multi-origin'
.memory_limit (
float
) – Upper limit in the amount of computer memory that the displacements can occupy in gigabytes (GB). Optional, defaults to8.
.n_steps – Number of steps to be used in the timestep function. Optional, defaults to all of them.
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)#
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)#
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.