aiida.orm.nodes.data.array package#
Module with Node sub classes for array based data structures.
Submodules#
AiiDA ORM data class storing (numpy) arrays
- class aiida.orm.nodes.data.array.array.ArrayData(arrays: ndarray | dict[str, ndarray] | None = None, **kwargs)[source]#
Bases:
Data
Store a set of arrays on disk (rather than on the database) in an efficient way using numpy.save() (therefore, this class requires numpy to be installed).
Each array is stored within the Node folder as a different .npy file.
- Note:
Before storing, no caching is done: if you perform a
get_array()
call, the array will be re-read from disk. If instead the ArrayData node has already been stored, the array is cached in memory after the first read, and the cached array is used thereafter. If too much RAM memory is used, you can clear the cache with theclear_internal_cache()
method.
- __abstractmethods__ = frozenset({})#
- __init__(arrays: ndarray | dict[str, ndarray] | None = None, **kwargs)[source]#
Construct a new instance and set one or multiple numpy arrays.
- Parameters:
arrays – An optional single numpy array, or dictionary of numpy arrays to store.
- __module__ = 'aiida.orm.nodes.data.array.array'#
- __parameters__ = ()#
- _abc_impl = <_abc._abc_data object>#
- _arraynames_from_files() list[str] [source]#
Return a list of all arrays stored in the node, listing the files (and not relying on the properties).
- _arraynames_from_properties() list[str] [source]#
Return a list of all arrays stored in the node, listing the attributes starting with the correct prefix.
- _get_array_entries() dict[str, Any] [source]#
Return a dictionary with the different array entries.
The idea is that this dictionary contains the array name as a key and the value is the numpy array transformed into a list. This is so that it can be transformed into a json object.
- _prepare_json(main_file_name='', comments=True) tuple[bytes, dict] [source]#
Dump the content of the arrays stored in this node into JSON format.
- Parameters:
comments – if True, includes comments (if it makes sense for the given format)
- _validate() bool [source]#
Check if the list of .npy files stored inside the node and the list of properties match. Just a name check, no check on the size since this would require to reload all arrays and this may take time and memory.
- array_prefix = 'array|'#
- clear_internal_cache() None [source]#
Clear the internal memory cache where the arrays are stored after being read from disk (used in order to reduce at minimum the readings from disk). This function is useful if you want to keep the node in memory, but you do not want to waste memory to cache the arrays in RAM.
- default_array_name = 'default'#
- delete_array(name: str) None [source]#
Delete an array from the node. Can only be called before storing.
- Parameters:
name – The name of the array to delete from the node.
- get_array(name: str | None = None) ndarray [source]#
Return an array stored in the node
- Parameters:
name – The name of the array to return. The name can be omitted in case the node contains only a single array, which will be returned in that case. If
name
isNone
and the node contains multiple arrays or no arrays at all aValueError
is raised.- Raises:
ValueError – If
name
isNone
and the node contains more than one arrays or no arrays at all.
- get_arraynames() list[str] [source]#
Return a list of all arrays stored in the node, listing the files (and not relying on the properties).
New in version 0.7: Renamed from arraynames
- get_iterarrays() Iterator[tuple[str, ndarray]] [source]#
Iterator that returns tuples (name, array) for each array stored in the node.
New in version 1.0: Renamed from iterarrays
- aiida.orm.nodes.data.array.array.clean_array(array: ndarray) list [source]#
Replacing np.nan and np.inf/-np.inf for Nones.
The function will also sanitize the array removing
np.nan
andnp.inf
forNone
of this way the resulting JSON is always valid. Bothnp.nan
andnp.inf
/-np.inf
are set to None to be in accordance with the ECMA-262 standard.- Parameters:
array – input array to be cleaned
- Returns:
cleaned list to be serialized
- Return type:
This module defines the classes related to band structures or dispersions in a Brillouin zone, and how to operate on them.
- class aiida.orm.nodes.data.array.bands.BandsData(arrays: ndarray | dict[str, ndarray] | None = None, **kwargs)[source]#
Bases:
KpointsData
Class to handle bands data
- __abstractmethods__ = frozenset({})#
- __module__ = 'aiida.orm.nodes.data.array.bands'#
- __parameters__ = ()#
- _abc_impl = <_abc._abc_data object>#
- _get_bandplot_data(cartesian, prettify_format=None, join_symbol=None, get_segments=False, y_origin=0.0)[source]#
Get data to plot a band structure
- Parameters:
cartesian – if True, distances (for the x-axis) are computed in cartesian coordinates, otherwise they are computed in reciprocal coordinates. cartesian=True will fail if no cell has been set.
prettify_format – by default, strings are not prettified. If you want to prettify them, pass a valid prettify_format string (see valid options in the docstring of :py:func:prettify_labels).
join_symbols – by default, strings are not joined. If you pass a string, this is used to join strings that are much closer than a given threshold. The most typical string is the pipe symbol:
|
.get_segments – if True, also computes the band split into segments
y_origin – if present, shift bands so to set the value specified at
y=0
- Returns:
a plot_info dictiorary, whose keys are
x
(array of distances for the x axis of the plot);y
(array of bands),labels
(list of tuples in the format (float x value of the label, label string),band_type_idx
(array containing an index for each band: if there is only one spin, then it’s an array of zeros, of length equal to the number of bands at each point; if there are two spins, then it’s an array of zeros or ones depending on the type of spin; the length is always equalt to the total number of bands per kpoint).
- _matplotlib_get_dict(main_file_name='', comments=True, title='', legend=None, legend2=None, y_max_lim=None, y_min_lim=None, y_origin=0.0, prettify_format=None, **kwargs)[source]#
Prepare the data to send to the python-matplotlib plotting script.
- Parameters:
comments – if True, print comments (if it makes sense for the given format)
plot_info – a dictionary
setnumber_offset – an offset to be applied to all set numbers (i.e. s0 is replaced by s[offset], s1 by s[offset+1], etc.)
color_number – the color number for lines, symbols, error bars and filling (should be less than the parameter MAX_NUM_AGR_COLORS defined below)
title – the title
legend – the legend (applied only to the first of the set)
legend2 – the legend for second-type spins (applied only to the first of the set)
y_max_lim – the maximum on the y axis (if None, put the maximum of the bands)
y_min_lim – the minimum on the y axis (if None, put the minimum of the bands)
y_origin – the new origin of the y axis -> all bands are replaced by bands-y_origin
prettify_format – if None, use the default prettify format. Otherwise specify a string with the prettifier to use.
kwargs – additional customization variables; only a subset is accepted, see internal variable ‘valid_additional_keywords
- _prepare_agr(main_file_name='', comments=True, setnumber_offset=0, color_number=1, color_number2=2, legend='', title='', y_max_lim=None, y_min_lim=None, y_origin=0.0, prettify_format=None)[source]#
Prepare an xmgrace agr file.
- Parameters:
comments – if True, print comments (if it makes sense for the given format)
plot_info – a dictionary
setnumber_offset – an offset to be applied to all set numbers (i.e. s0 is replaced by s[offset], s1 by s[offset+1], etc.)
color_number – the color number for lines, symbols, error bars and filling (should be less than the parameter MAX_NUM_AGR_COLORS defined below)
color_number2 – the color number for lines, symbols, error bars and filling for the second-type spins (should be less than the parameter MAX_NUM_AGR_COLORS defined below)
legend – the legend (applied only to the first set)
title – the title
y_max_lim – the maximum on the y axis (if None, put the maximum of the bands); applied after shifting the origin by
y_origin
y_min_lim – the minimum on the y axis (if None, put the minimum of the bands); applied after shifting the origin by
y_origin
y_origin – the new origin of the y axis -> all bands are replaced by bands-y_origin
prettify_format – if None, use the default prettify format. Otherwise specify a string with the prettifier to use.
- _prepare_agr_batch(main_file_name='', comments=True, prettify_format=None)[source]#
Prepare two files, data and batch, to be plot with xmgrace as: xmgrace -batch file.dat
- Parameters:
main_file_name – if the user asks to write the main content on a file, this contains the filename. This should be used to infer a good filename for the additional files. In this case, we remove the extension, and add ‘_data.dat’
comments – if True, print comments (if it makes sense for the given format)
prettify_format – if None, use the default prettify format. Otherwise specify a string with the prettifier to use.
- _prepare_dat_blocks(main_file_name='', comments=True)[source]#
Format suitable for gnuplot using blocks. Columns with x and y (path and band energy). Several blocks, separated by two empty lines, one per energy band.
- Parameters:
comments – if True, print comments (if it makes sense for the given format)
- _prepare_dat_multicolumn(main_file_name='', comments=True)[source]#
Write an N x M matrix. First column is the distance between kpoints, The other columns are the bands. Header contains number of kpoints and the number of bands (commented).
- Parameters:
comments – if True, print comments (if it makes sense for the given format)
- _prepare_gnuplot(main_file_name=None, title='', comments=True, prettify_format=None, y_max_lim=None, y_min_lim=None, y_origin=0.0)[source]#
Prepare an gnuplot script to plot the bands, with the .dat file returned as an independent file.
- Parameters:
main_file_name – if the user asks to write the main content on a file, this contains the filename. This should be used to infer a good filename for the additional files. In this case, we remove the extension, and add ‘_data.dat’
title – if specified, add a title to the plot
comments – if True, print comments (if it makes sense for the given format)
prettify_format – if None, use the default prettify format. Otherwise specify a string with the prettifier to use.
- _prepare_json(main_file_name='', comments=True)[source]#
Prepare a json file in a format compatible with the AiiDA band visualizer
- Parameters:
comments – if True, print comments (if it makes sense for the given format)
- _prepare_mpl_pdf(main_file_name='', *args, **kwargs)[source]#
Prepare a python script using matplotlib to plot the bands, with the JSON returned as an independent file.
For the possible parameters, see documentation of
_matplotlib_get_dict()
- _prepare_mpl_png(main_file_name='', *args, **kwargs)[source]#
Prepare a python script using matplotlib to plot the bands, with the JSON returned as an independent file.
For the possible parameters, see documentation of
_matplotlib_get_dict()
- _prepare_mpl_singlefile(*args, **kwargs)[source]#
Prepare a python script using matplotlib to plot the bands
For the possible parameters, see documentation of
_matplotlib_get_dict()
- _prepare_mpl_withjson(main_file_name='', *args, **kwargs)[source]#
Prepare a python script using matplotlib to plot the bands, with the JSON returned as an independent file.
For the possible parameters, see documentation of
_matplotlib_get_dict()
- _validate_bands_occupations(bands, occupations=None, labels=None)[source]#
Validate the list of bands and of occupations before storage. Kpoints must be set in advance. Bands and occupations must be convertible into arrays of Nkpoints x Nbands floats or Nspins x Nkpoints x Nbands; Nkpoints must correspond to the number of kpoints.
- property array_labels#
Get the labels associated with the band arrays
- get_bands(also_occupations=False, also_labels=False)[source]#
Returns an array (nkpoints x num_bands or nspins x nkpoints x num_bands) of energies. :param also_occupations: if True, returns also the occupations array. Default = False
- set_bands(bands, units=None, occupations=None, labels=None)[source]#
Set an array of band energies of dimension (nkpoints x nbands). Kpoints must be set in advance. Can contain floats or None. :param bands: a list of nkpoints lists of nbands bands, or a 2D array of shape (nkpoints x nbands), with band energies for each kpoint :param units: optional, energy units :param occupations: optional, a 2D list or array of floats of same shape as bands, with the occupation associated to each band
- set_kpointsdata(kpointsdata)[source]#
Load the kpoints from a kpoint object. :param kpointsdata: an instance of KpointsData class
- show_mpl(**kwargs)[source]#
Call a show() command for the band structure using matplotlib. This uses internally the ‘mpl_singlefile’ format, with empty main_file_name.
Other kwargs are passed to self._exportcontent.
- property units#
Units in which the data in bands were stored. A string
- aiida.orm.nodes.data.array.bands._extract_formula(akinds, asites, args)[source]#
Extract formula from the structure object.
- Parameters:
akinds – list of kinds, e.g. [{‘mass’: 55.845, ‘name’: ‘Fe’, ‘symbols’: [‘Fe’], ‘weights’: [1.0]}, {‘mass’: 15.9994, ‘name’: ‘O’, ‘symbols’: [‘O’], ‘weights’: [1.0]}]
asites – list of structure sites e.g. [{‘position’: [0.0, 0.0, 0.0], ‘kind_name’: ‘Fe’}, {‘position’: [2.0, 2.0, 2.0], ‘kind_name’: ‘O’}]
args (dict) – a namespace with parsed command line parameters, here only ‘element’ and ‘element_only’ are used
- Returns:
a string with formula if the formula is found
- aiida.orm.nodes.data.array.bands.find_bandgap(bandsdata, number_electrons=None, fermi_energy=None)[source]#
Tries to guess whether the bandsdata represent an insulator. This method is meant to be used only for electronic bands (not phonons) By default, it will try to use the occupations to guess the number of electrons and find the Fermi Energy, otherwise, it can be provided explicitely. Also, there is an implicit assumption that the kpoints grid is “sufficiently” dense, so that the bandsdata are not missing the intersection between valence and conduction band if present. Use this function with care!
- Parameters:
number_electrons – (optional, float) number of electrons in the unit cell
fermi_energy – (optional, float) value of the fermi energy.
- Note:
By default, the algorithm uses the occupations array to guess the number of electrons and the occupied bands. This is to be used with care, because the occupations could be smeared so at a non-zero temperature, with the unwanted effect that the conduction bands might be occupied in an insulator. Prefer to pass the number_of_electrons explicitly
- Note:
Only one between number_electrons and fermi_energy can be specified at the same time.
- Returns:
(is_insulator, gap), where is_insulator is a boolean, and gap a float. The gap is None in case of a metal, zero when the homo is equal to the lumo (e.g. in semi-metals).
- aiida.orm.nodes.data.array.bands.get_bands_and_parents_structure(args, backend=None)[source]#
Search for bands and return bands and the closest structure that is a parent of the instance.
- Returns:
- A list of sublists, each latter containing (in order):
pk as string, formula as string, creation date, bandsdata-label
- aiida.orm.nodes.data.array.bands.prepare_header_comment(uuid, plot_info, comment_char='#')[source]#
Prepare the header.
Module of the KpointsData class, defining the AiiDA data type for storing lists and meshes of k-points (i.e., points in the reciprocal space of a periodic crystal structure).
- class aiida.orm.nodes.data.array.kpoints.KpointsData(arrays: ndarray | dict[str, ndarray] | None = None, **kwargs)[source]#
Bases:
ArrayData
Class to handle array of kpoints in the Brillouin zone. Provide methods to generate either user-defined k-points or path of k-points along symmetry lines. Internally, all k-points are defined in terms of crystal (fractional) coordinates. Cell and lattice vector coordinates are in Angstroms, reciprocal lattice vectors in Angstrom^-1 . :note: The methods setting and using the Bravais lattice info assume the PRIMITIVE unit cell is provided in input to the set_cell or set_cell_from_structure methods.
- __abstractmethods__ = frozenset({})#
- __module__ = 'aiida.orm.nodes.data.array.kpoints'#
- __parameters__ = ()#
- _abc_impl = <_abc._abc_data object>#
- _change_reference(kpoints, to_cartesian=True)[source]#
Change reference system, from cartesian to crystal coordinates (units of b1,b2,b3) or viceversa. :param kpoints: a list of (3) point coordinates :return kpoints: a list of (3) point coordinates in the new reference
- property _dimension#
Dimensionality of the structure, found from its pbc (i.e. 1 if it’s a 1D structure, 2 if its 2D, 3 if it’s 3D …). :return dimensionality: 0, 1, 2 or 3 :note: will return 3 if pbc has not been set beforehand
- _set_cell(value)[source]#
Validate if ‘value’ is a allowed crystal unit cell :param value: something compatible with a 3x3 tuple of floats
- _set_labels(value)[source]#
set label names. Must pass in input a list like:
[[0,'X'],[34,'L'],... ]
- _validate_kpoints_weights(kpoints, weights)[source]#
Validate the list of kpoints and of weights before storage. Kpoints and weights must be convertible respectively to an array of N x dimension and N floats
- property cell#
The crystal unit cell. Rows are the crystal vectors in Angstroms. :return: a 3x3 numpy.array
- get_description()[source]#
Returns a string with infos retrieved from kpoints node’s properties. :param node: :return: retstr
- get_kpoints(also_weights=False, cartesian=False)[source]#
Return the list of kpoints
- Parameters:
also_weights – if True, returns also the list of weights. Default = False
cartesian – if True, returns points in cartesian coordinates, otherwise, returns in crystal coordinates. Default = False.
- get_kpoints_mesh(print_list=False)[source]#
Get the mesh of kpoints.
- Parameters:
print_list – default=False. If True, prints the mesh of kpoints as a list
- Raises:
AttributeError – if no mesh has been set
- Return mesh,offset:
(if print_list=False) a list of 3 integers and a list of three floats 0<x<1, representing the mesh and the offset of kpoints
- Return kpoints:
(if print_list = True) an explicit list of kpoints coordinates, similar to what returned by get_kpoints()
- property labels#
Labels associated with the list of kpoints. List of tuples with kpoint index and kpoint name:
[(0,'G'),(13,'M'),...]
- property pbc#
The periodic boundary conditions along the vectors a1,a2,a3.
- Returns:
a tuple of three booleans, each one tells if there are periodic boundary conditions for the i-th real-space direction (i=1,2,3)
- property reciprocal_cell#
Compute reciprocal cell from the internally set cell.
- Returns:
reciprocal cell in units of 1/Angstrom with cell vectors stored as rows. Use e.g. reciprocal_cell[0] to access the first reciprocal cell vector.
- set_cell(cell, pbc=None)[source]#
Set a cell to be used for symmetry analysis. To set a cell from an AiiDA structure, use “set_cell_from_structure”.
- Parameters:
cell – 3x3 matrix of cell vectors. Orientation: each row represent a lattice vector. Units are Angstroms.
pbc – list of 3 booleans, True if in the nth crystal direction the structure is periodic. Default = [True,True,True]
- set_cell_from_structure(structuredata)[source]#
Set a cell to be used for symmetry analysis from an AiiDA structure. Inherits both the cell and the pbc’s. To set manually a cell, use “set_cell”
- Parameters:
structuredata – an instance of StructureData
- set_kpoints(kpoints, cartesian=False, labels=None, weights=None, fill_values=0)[source]#
Set the list of kpoints. If a mesh has already been stored, raise a ModificationNotAllowed
- Parameters:
kpoints –
a list of kpoints, each kpoint being a list of one, two or three coordinates, depending on self.pbc: if structure is 1D (only one True in self.pbc) one allows singletons or scalars for each k-point, if it’s 2D it can be a length-2 list, and in all cases it can be a length-3 list. Examples:
[[0.,0.,0.],[0.1,0.1,0.1],…] for 1D, 2D or 3D
[[0.,0.],[0.1,0.1,],…] for 1D or 2D
[[0.],[0.1],…] for 1D
[0., 0.1, …] for 1D (list of scalars)
For 0D (all pbc are False), the list can be any of the above or empty - then only Gamma point is set. The value of k for the non-periodic dimension(s) is set by fill_values
cartesian – if True, the coordinates given in input are treated as in cartesian units. If False, the coordinates are crystal, i.e. in units of b1,b2,b3. Default = False
labels – optional, the list of labels to be set for some of the kpoints. See labels for more info
weights – optional, a list of floats with the weight associated to the kpoint list
fill_values – scalar to be set to all non-periodic dimensions (indicated by False in self.pbc), or list of values for each of the non-periodic dimensions.
- set_kpoints_mesh(mesh, offset=None)[source]#
Set KpointsData to represent a uniformily spaced mesh of kpoints in the Brillouin zone. This excludes the possibility of set/get kpoints
- Parameters:
mesh – a list of three integers, representing the size of the kpoint mesh along b1,b2,b3.
offset – (optional) a list of three floats between 0 and 1. [0.,0.,0.] is Gamma centered mesh [0.5,0.5,0.5] is half shifted [1.,1.,1.] by periodicity should be equivalent to [0.,0.,0.] Default = [0.,0.,0.].
- set_kpoints_mesh_from_density(distance, offset=None, force_parity=False)[source]#
Set a kpoints mesh using a kpoints density, expressed as the maximum distance between adjacent points along a reciprocal axis
- Parameters:
distance – distance (in 1/Angstrom) between adjacent kpoints, i.e. the number of kpoints along each reciprocal axis i is \(|b_i|/distance\) where \(|b_i|\) is the norm of the reciprocal cell vector.
offset – (optional) a list of three floats between 0 and 1. [0.,0.,0.] is Gamma centered mesh [0.5,0.5,0.5] is half shifted Default = [0.,0.,0.].
force_parity – (optional) if True, force each integer in the mesh to be even (except for the non-periodic directions).
- Note:
a cell should be defined first.
- Note:
the number of kpoints along non-periodic axes is always 1.
Data plugin to represet arrays of projected wavefunction components.
- class aiida.orm.nodes.data.array.projection.ProjectionData(arrays: ndarray | dict[str, ndarray] | None = None, **kwargs)[source]#
Bases:
OrbitalData
,ArrayData
A class to handle arrays of projected wavefunction data. That is projections of a orbitals, usually an atomic-hydrogen orbital, onto a given bloch wavefunction, the bloch wavefunction being indexed by s, n, and k. E.g. the elements are the projections described as < orbital | Bloch wavefunction (s,n,k) >
- __abstractmethods__ = frozenset({})#
- __module__ = 'aiida.orm.nodes.data.array.projection'#
- __parameters__ = ()#
- _abc_impl = <_abc._abc_data object>#
- _check_projections_bands(projection_array)[source]#
Checks to make sure that a reference bandsdata is already set, and that projection_array is of the same shape of the bands data
- Parameters:
projwfc_arrays – nk x nb x nwfc array, to be checked against bands
- Raise:
AttributeError if energy is not already set
- Raise:
AttributeError if input_array is not of same shape as dos_energy
- _find_orbitals_and_indices(**kwargs)[source]#
Finds all the orbitals and their indicies associated with kwargs essential for retrieving the other indexed array parameters
- Parameters:
kwargs – kwargs that can call orbitals as in get_orbitals()
- Returns:
retrieve_indexes, list of indicicies of orbitals corresponding to the kwargs
- Returns:
all_orbitals, list of orbitals to which the indexes correspond
- get_pdos(**kwargs)[source]#
Retrieves all the pdos arrays corresponding to the input kwargs
- Parameters:
kwargs – inputs describing the orbitals associated with the pdos arrays
- Returns:
a list of tuples containing the orbital, energy array and pdos array associated with all orbitals that correspond to kwargs
- get_projections(**kwargs)[source]#
Retrieves all the pdos arrays corresponding to the input kwargs
- Parameters:
kwargs – inputs describing the orbitals associated with the pdos arrays
- Returns:
a list of tuples containing the orbital, and projection arrays associated with all orbitals that correspond to kwargs
- get_reference_bandsdata()[source]#
Returns the reference BandsData, using the set uuid via set_reference_bandsdata
- Returns:
a BandsData instance
- Raises:
AttributeError – if the bandsdata has not been set yet
exceptions.NotExistent – if the bandsdata uuid did not retrieve bandsdata
- set_orbitals(**kwargs)[source]#
This method is inherited from OrbitalData, but is blocked here. If used will raise a NotImplementedError
- set_projectiondata(list_of_orbitals, list_of_projections=None, list_of_energy=None, list_of_pdos=None, tags=None, bands_check=True)[source]#
Stores the projwfc_array using the projwfc_label, after validating both.
- Parameters:
list_of_orbitals – list of orbitals, of class orbital data. They should be the ones up on which the projection array corresponds with.
list_of_projections – list of arrays of projections of a atomic wavefunctions onto bloch wavefunctions. Since the projection is for every bloch wavefunction which can be specified by its spin (if used), band, and kpoint the dimensions must be nspin x nbands x nkpoints for the projwfc array. Or nbands x nkpoints if spin is not used.
energy_axis – list of energy axis for the list_of_pdos
list_of_pdos – a list of projected density of states for the atomic wavefunctions, units in states/eV
tags – A list of tags, not supported currently.
bands_check – if false, skips checks of whether the bands has been already set, and whether the sizes match. For use in parsers, where the BandsData has not yet been stored and therefore get_reference_bandsdata cannot be called
- set_reference_bandsdata(value)[source]#
Sets a reference bandsdata, creates a uuid link between this data object and a bandsdata object, must be set before any projection arrays
- Parameters:
value – a BandsData instance, a uuid or a pk
- Raise:
exceptions.NotExistent if there was no BandsData associated with uuid or pk
AiiDA class to deal with crystal structure trajectories.
- class aiida.orm.nodes.data.array.trajectory.TrajectoryData(structurelist=None, **kwargs)[source]#
Bases:
ArrayData
Stores a trajectory (a sequence of crystal structures with timestamps, and possibly with velocities).
- __abstractmethods__ = frozenset({})#
- __init__(structurelist=None, **kwargs)[source]#
Construct a new instance and set one or multiple numpy arrays.
- Parameters:
arrays – An optional single numpy array, or dictionary of numpy arrays to store.
- __module__ = 'aiida.orm.nodes.data.array.trajectory'#
- __parameters__ = ()#
- _abc_impl = <_abc._abc_data object>#
- _internal_validate(stepids, cells, symbols, positions, times, velocities)[source]#
Internal function to validate the type and shape of the arrays. See the documentation of py:meth:.set_trajectory for a description of the valid shape and type of the parameters.
- _parse_xyz_pos(inputstring)[source]#
Load positions from a XYZ file.
Note
The steps and symbols must be set manually before calling this import function as a consistency measure. Even though the symbols and steps could be extracted from the XYZ file, the data present in the XYZ file may or may not be correct and the same logic would have to be present in the XYZ-velocities function. It was therefore decided not to implement it at all but require it to be set explicitly.
Usage:
from aiida.orm.nodes.data.array.trajectory import TrajectoryData t = TrajectoryData() # get sites and number of timesteps t.set_array('steps', arange(ntimesteps)) t.set_array('symbols', array([site.kind for site in s.sites])) t.importfile('some-calc/AIIDA-PROJECT-pos-1.xyz', 'xyz_pos')
- _parse_xyz_vel(inputstring)[source]#
Load velocities from a XYZ file.
Note
The steps and symbols must be set manually before calling this import function as a consistency measure. See also comment for
_parse_xyz_pos()
- _prepare_cif(trajectory_index=None, main_file_name='')[source]#
Write the given trajectory to a string of format CIF.
- _prepare_xsf(index=None, main_file_name='')[source]#
Write the given trajectory to a string of format XSF (for XCrySDen).
- _validate()[source]#
Verify that the required arrays are present and that their type and dimension are correct.
- get_cells()[source]#
Return the array of cells, if it has already been set.
- Raises:
KeyError – if the trajectory has not been set yet.
- get_cif(index=None, **kwargs)[source]#
Creates
aiida.orm.nodes.data.cif.CifData
New in version 1.0: Renamed from _get_cif
- get_index_from_stepid(stepid)[source]#
Given a value for the stepid (i.e., a value among those of the
steps
array), return the array index of that stepid, that can be used in other methods such asget_step_data()
orget_step_structure()
.New in version 0.7: Renamed from get_step_index
Note
Note that this function returns the first index found (i.e. if multiple steps are present with the same value, only the index of the first one is returned).
- Raises:
ValueError – if no step with the given value is found.
- get_positions()[source]#
Return the array of positions, if it has already been set.
- Raises:
KeyError – if the trajectory has not been set yet.
- get_step_data(index)[source]#
Return a tuple with all information concerning the stepid with given index (0 is the first step, 1 the second step and so on). If you know only the step value, use the
get_index_from_stepid()
method to get the corresponding index.If no velocities, cells, or times were specified, None is returned as the corresponding element.
- Returns:
A tuple in the format
(stepid, time, cell, symbols, positions, velocities)
, wherestepid
is an integer,time
is a float,cell
is a \(3 imes 3\) matrix,symbols
is an array of lengthn
, positions is a \(n imes 3\) array, and velocities is eitherNone
or a \(n imes 3\) array- Parameters:
index – The index of the step that you want to retrieve, from 0 to
self.numsteps - 1
.- Raises:
IndexError – if you require an index beyond the limits.
KeyError – if you did not store the trajectory yet.
- get_step_structure(index, custom_kinds=None)[source]#
Return an AiiDA
aiida.orm.nodes.data.structure.StructureData
node (not stored yet!) with the coordinates of the given step, identified by its index. If you know only the step value, use theget_index_from_stepid()
method to get the corresponding index.Note
The periodic boundary conditions are always set to True.
New in version 0.7: Renamed from step_to_structure
- Parameters:
index – The index of the step that you want to retrieve, from 0 to
self.numsteps- 1
.custom_kinds – (Optional) If passed must be a list of
aiida.orm.nodes.data.structure.Kind
objects. There must be one kind object for each different string in thesymbols
array, withkind.name
set to this string. If this parameter is omitted, the automatic kind generation of AiiDAaiida.orm.nodes.data.structure.StructureData
nodes is used, meaning that the strings in thesymbols
array must be valid chemical symbols.
- Returns:
- get_stepids()[source]#
Return the array of steps, if it has already been set.
New in version 0.7: Renamed from get_steps
- Raises:
KeyError – if the trajectory has not been set yet.
- get_structure(store=False, **kwargs)[source]#
Creates
aiida.orm.nodes.data.structure.StructureData
.New in version 1.0: Renamed from _get_aiida_structure
- Parameters:
store – If True, intermediate calculation gets stored in the AiiDA database for record. Default False.
index – The index of the step that you want to retrieve, from 0 to
self.numsteps- 1
.custom_kinds – (Optional) If passed must be a list of
aiida.orm.nodes.data.structure.Kind
objects. There must be one kind object for each different string in thesymbols
array, withkind.name
set to this string. If this parameter is omitted, the automatic kind generation of AiiDAaiida.orm.nodes.data.structure.StructureData
nodes is used, meaning that the strings in thesymbols
array must be valid chemical symbols.custom_cell – (Optional) The cell matrix of the structure. If omitted, the cell will be read from the trajectory, if present, otherwise the default cell of
aiida.orm.nodes.data.structure.StructureData
will be used.
- Returns:
- get_times()[source]#
Return the array of times (in ps), if it has already been set.
- Raises:
KeyError – if the trajectory has not been set yet.
- get_velocities()[source]#
Return the array of velocities, if it has already been set.
Note
This function (differently from all other
get_*
functions, will not raise an exception if the velocities are not set, but rather returnNone
(both if no trajectory was not set yet, and if it the trajectory was set but no velocities were specified).
- property numsites#
Return the number of stored sites, or zero if nothing has been stored yet.
- property numsteps#
Return the number of stored steps, or zero if nothing has been stored yet.
- set_structurelist(structurelist)[source]#
Create trajectory from the list of
aiida.orm.nodes.data.structure.StructureData
instances.- Parameters:
structurelist – a list of
aiida.orm.nodes.data.structure.StructureData
instances.- Raises:
ValueError – if symbol lists of supplied structures are different
- set_trajectory(symbols, positions, stepids=None, cells=None, times=None, velocities=None)[source]#
Store the whole trajectory, after checking that types and dimensions are correct.
Parameters
stepids
,cells
andvelocities
are optional variables. If nothing is passed forcells
orvelocities
nothing will be stored. However, if no input is given forstepids
a consecutive sequence [0,1,2,…,len(positions)-1] will be assumed.- Parameters:
symbols – string list with dimension
n
, wheren
is the number of atoms (i.e., sites) in the structure. The same list is used for each step. Normally, the string should be a valid chemical symbol, but actually any unique string works and can be used as the name of the atomic kind (see also theget_step_structure()
method).positions – float array with dimension \(s \times n \times 3\), where
s
is the length of thestepids
array andn
is the length of thesymbols
array. Units are angstrom. In particular,positions[i,j,k]
is thek
-th component of thej
-th atom (or site) in the structure at the time step with indexi
(identified by step numberstep[i]
and with timestamptimes[i]
).stepids – integer array with dimension
s
, wheres
is the number of steps. Typically represents an internal counter within the code. For instance, if you want to store a trajectory with one step every 10, starting from step 65, the array will be[65,75,85,...]
. No checks are done on duplicate elements or on the ordering, but anyway this array should be sorted in ascending order, without duplicate elements. (If not specified, stepids will be set tonumpy.arange(s)
by default) It is internally stored as an array named ‘steps’.cells – if specified float array with dimension \(s \times 3 \times 3\), where
s
is the length of thestepids
array. Units are angstrom. In particular,cells[i,j,k]
is thek
-th component of thej
-th cell vector at the time step with indexi
(identified by step numberstepid[i]
and with timestamptimes[i]
).times – if specified, float array with dimension
s
, wheres
is the length of thestepids
array. Contains the timestamp of each step in picoseconds (ps).velocities – if specified, must be a float array with the same dimensions of the
positions
array. The array contains the velocities in the atoms.
- show_mpl_pos(**kwargs)[source]#
Shows the positions as a function of time, separate for XYZ coordinates
- Parameters:
stepsize (int) – The stepsize for the trajectory, set higher than 1 to reduce number of points
mintime (int) – Time to start from
maxtime (int) – Maximum time
elements (list) – A list of atomic symbols that should be displayed. If not specified, all atoms are displayed.
indices (list) – A list of indices of that atoms that can be displayed. If not specified, all atoms of the correct species are displayed.
dont_block (bool) – If True, interpreter is not blocked when figure is displayed.
- aiida.orm.nodes.data.array.trajectory.plot_positions_XYZ(times, positions, indices_to_show, color_list, label, positions_unit='A', times_unit='ps', dont_block=False, mintime=None, maxtime=None, n_labels=10)[source]#
Plot with matplotlib the positions of the coordinates of the atoms over time for a trajectory
- Parameters:
times – array of times
positions – array of positions
indices_to_show – list of indices of to show (0, 1, 2 for X, Y, Z)
color_list – list of valid color specifications for matplotlib
label – label for this plot to put in the title
positions_unit – label for the units of positions (for the x label)
times_unit – label for the units of times (for the y label)
dont_block – passed to plt.show() as
block=not dont_block
mintime – if specified, cut the time axis at the specified min value
maxtime – if specified, cut the time axis at the specified max value
n_labels – how many labels (t, coord) to put
This module defines the classes related to Xy data. That is data that contains collections of y-arrays bound to a single x-array, and the methods to operate on them.
- class aiida.orm.nodes.data.array.xy.XyData(x_array: 'ndarray' | None = None, y_arrays: 'ndarray' | list['ndarray'] | None = None, *, x_name: str | None = None, x_units: str | None = None, y_names: str | list[str] | None = None, y_units: str | list[str] | None = None, **kwargs)[source]#
Bases:
ArrayData
A subclass designed to handle arrays that have an “XY” relationship to each other. That is there is one array, the X array, and there are several Y arrays, which can be considered functions of X.
- __abstractmethods__ = frozenset({})#
- __init__(x_array: 'ndarray' | None = None, y_arrays: 'ndarray' | list['ndarray'] | None = None, *, x_name: str | None = None, x_units: str | None = None, y_names: str | list[str] | None = None, y_units: str | list[str] | None = None, **kwargs)[source]#
Construct a new instance, optionally setting the x and y arrays.
Note
If the
x_array
is specified, all other keywords need to be specified as well.- Parameters:
x_array – The x array.
y_arrays – The y arrays.
x_name – The name of the x array.
x_units – The unit of the x array.
y_names – The names of the y arrays.
y_units – The units of the y arrays.
- __module__ = 'aiida.orm.nodes.data.array.xy'#
- __parameters__ = ()#
- _abc_impl = <_abc._abc_data object>#
- static _arrayandname_validator(array: ndarray, name: str, units: str) None [source]#
Validates that the array is an numpy.ndarray and that the name is of type str. Raises TypeError or ValueError if this not the case.
- get_x() tuple[str, ndarray, str] [source]#
Tries to retrieve the x array and x name raises a NotExistent exception if no x array has been set yet. :return x_name: the name set for the x_array :return x_array: the x array set earlier :return x_units: the x units set earlier
- get_y() list[tuple[str, ndarray, str]] [source]#
Tries to retrieve the y arrays and the y names, raises a NotExistent exception if they have not been set yet, or cannot be retrieved :return y_names: list of strings naming the y_arrays :return y_arrays: list of y_arrays :return y_units: list of strings giving the units for the y_arrays
- set_x(x_array: ndarray, x_name: str, x_units: str) None [source]#
Sets the array and the name for the x values.
- Parameters:
x_array – A numpy.ndarray, containing only floats
x_name – a string for the x array name
x_units – the units of x
- set_y(y_arrays: 'ndarray' | Sequence['ndarray'], y_names: str | Sequence[str], y_units: str | Sequence[str]) None [source]#
Set array(s) for the y part of the dataset. Also checks if the x_array has already been set, and that, the shape of the y_arrays agree with the x_array. :param y_arrays: A list of y_arrays, numpy.ndarray :param y_names: A list of strings giving the names of the y_arrays :param y_units: A list of strings giving the units of the y_arrays
- aiida.orm.nodes.data.array.xy.check_convert_single_to_tuple(item: Any | Sequence[Any]) Sequence[Any] [source]#
Checks if the item is a list or tuple, and converts it to a list if it is not already a list or tuple
- Parameters:
item – an object which may or may not be a list or tuple
- Returns:
item_list: the input item unchanged if list or tuple and [item] otherwise