Structures

StructureContainer

class hiphive.StructureContainer(cs, fit_structure_list=None)[source]

This class serves as a container for structures as well as associated fit properties and fit matrices.

Parameters:
  • cs (ClusterSpace object) – cluster space that is the basis for the container
  • fit_structure_list (list of FitStructure objects (optional)) – structures to be added to the container
add_structure(atoms, user_tag=None, compute_fit_matrix=True)[source]

Add a structure to the container.

Note that custom information about the atoms object may not be stored inside, for example a SinglePointCalculator will not be kept.

Parameters:
  • atoms (ASE Atoms object) – the structure to be added; the Atoms object must contain supplementary per-atom arrays with displacements and forces
  • user_tag (string) – custom user tag to identify structure
  • compute_fit_matrix (boolean) – if True the fit matrix of the structure is computed
cluster_space

ClusterSpace object – copy of the cluster space the structure container is based on

data_shape

tuple – tuple of integers representing the shape of the fit data matrix

get_fit_data(structures=None)[source]

Return fit data for structures. The fit matrices and target forces for the structures are stacked into NumPy arrays.

Parameters:structures (list, tuple) – list of integers corresponding to structure indices. Defaults to None and in that case returns all fit data available.
Returns:stacked fit matrices, stacked target forces for the structures
Return type:NumPy array, NumPy array
get_structure_indices(user_tag)[source]

Get structure indices via user tag.

Parameters:user_tag (string) – user tag used for selecting structures
Returns:list of structures with specified user tag
Return type:list
print_cluster_space(include_orbits=False)[source]

Print information concerning the cluster space this structure container is based on.

Parameters:include_orbits (boolean) – if True also print the list of orbits associated with the cluster space
static read(fileobj, read_fit_matrices=True)[source]

Restore a StructureContainer object from file.

Parameters:
  • f (string or file object) – name of input file (string) or stream to load from (file object)
  • read_fit_matrices (bool) – Whether or not to load the fit matrices.
update_fit_data()[source]

Compute the fit matrices for structures for which they are not yet available.

write(f)[source]

Write a StructureContainer instance to a file.

Parameters:f (str or file object) – name of input file (str) or stream to write to (file object)

FitStructure

class hiphive.core.structure_container.FitStructure(atoms, user_tag=None, fit_matrix=None)[source]

This class holds a structure with displacements and forces as well as the fit matrix.

Parameters:
  • atoms (ASE Atoms object) – supercell structure
  • user_tag (string) – custom user tag
  • fit_matrix (NumPy array) – fit matrix; NumPy (N, M) array with N = 3 * len(atoms)
atoms

ASE Atoms object – supercell structure

displacements

NumPy array – atomic displacements

fit_matrix

NumPy array – the fit matrix

fit_ready

boolean – True if the structure is prepared for fitting, i.e. the fit matrix is available

forces

NumPy array – forces

static read(fileobj, read_fit_matrix=True)[source]
set_fit_matrix(fit_matrix)[source]

Set the fit matrix.

Parameters:fit_matrix (NumPy array) – fit matrix for this structure; NumPy (N, M) array with N = 3 * len(atoms)
user_tag
write(fileobj)[source]

Structure generation

Module for generating displaced structures

Todo

Add harmonic phonon eigenmode displacement method

hiphive.structure_generation.generate_mc_rattled_structures(atoms, n_configs, rattle_std, d_min, seed=42, **kwargs)[source]

Return list of Monte Carlo rattled configurations.

Rattling atom i is a Monte Carlo move and is accepted with a probability that is determined from the minimum interatomic distance d_ij. If min(d_ij) is smaller than d_min the move is only accepted with a low probability.

This process is repeated for each atom a number of times meaning the magnitude of the final displacements are not connected to rattle_std.

Warning

Repeatedly calling this function without providing different seeds will yield identical or correlated results to avoid this please specify a different seed for each call in such a case.

Notes

Note that this might not generate a symmetric distribution for the displacements kwargs will be forwarded to mc_rattle (see doc for this function for detailed explanation)

Parameters:
  • atoms (ASE Atoms object) – prototype ase atoms object
  • n_structures (int) – number of structures to generate
  • rattle_std (float) – rattle amplitude (standard deviation in normal distribution). Note this value is not connected to the final average displacement for the structures.
  • d_min (float) – interatomic distance used for computing the probability for each rattle move.
  • seed (int) – seed for setting up NumPy random state from which random numbers are generated
hiphive.structure_generation.generate_rattled_structures(atoms, n_structures, rattle_std, seed=42)[source]

Return list of rattled configurations.

Displacements are drawn from normal distributions in x,y,z directions for each atom independently.

Warning

Repeatedly calling this function without providing different seeds will yield identical or correlated results to avoid this please specify a different seed for each call in such a case.

Parameters:
  • atoms (ASE Atoms object) – prototype ase atoms object
  • n_structures (int) – number of structures to generate
  • rattle_std (float) – rattle amplitude (standard deviation in normal distribution)
  • seed (int) – seed for setting up NumPy random state from which random numbers are generated