pyiron.atomistics.structure.atoms module¶
-
class
pyiron.atomistics.structure.atoms.
Atoms
(symbols=None, positions=None, numbers=None, tags=None, momenta=None, masses=None, magmoms=None, charges=None, scaled_positions=None, cell=None, pbc=None, celldisp=None, constraint=None, calculator=None, info=None, indices=None, elements=None, dimension=None, species=None, **qwargs)[source]¶ Bases:
ase.atoms.Atoms
The Atoms class represents all the information required to describe a structure at the atomic scale. This class is derived from the ASE atoms class.
- Parameters
elements (list/numpy.ndarray) – List of strings containing the elements or a list of atomistics.structure.periodic_table.ChemicalElement instances
numbers (list/numpy.ndarray) – List of atomic numbers of elements
symbols (list/numpy.ndarray) – List of chemical symbols
positions (list/numpy.ndarray) – List of positions
scaled_positions (list/numpy.ndarray) – List of scaled positions (relative coordinates)
pbc (list/numpy.ndarray/boolean) – Tells if periodic boundary conditions should be applied on the three axes
cell (list/numpy.ndarray instance) – A 3x3 array representing the lattice vectors of the structure
Note: Only one of elements/symbols or numbers should be assigned during initialization
-
indices
¶ A list of size N which gives the species index of the structure which has N atoms
- Type
numpy.ndarray
-
add_high_symmetry_path
(path)[source]¶ Adds a new path to the dictionary of pathes for band structure calculations.
- Parameters
path (dict) – dictionary of lists of tuples with start and end point. E.G. {“my_path”: [(‘Gamma’, ‘X’), (‘X’, ‘Y’)]}
-
add_high_symmetry_points
(new_points)[source]¶ Adds new points to the dict of existing high symmetry points.
- Parameters
new_points (dict) – Points to add
-
add_tag
(*args, **qwargs)[source]¶ Add tags to the atoms object.
Examples
For selective dynamics:
>>> self.add_tag(selective_dynamics=[False, False, False])
-
property
analyse
¶
-
analyse_ovito_centro_symmetry
(num_neighbors=12)[source]¶ Analyse centrosymmetry parameter
- Parameters
num_neighbors (int) – number of neighbors
- Returns
list of centrosymmetry parameter
- Return type
list
-
analyse_ovito_cna_adaptive
(mode='total')[source]¶ Use common neighbor analysis
- Parameters
mode ("total"/"numeric"/"str") –
Controls the style and level of detail of the output. - total : return number of atoms belonging to each structure - numeric : return a per atom list of numbers- 0 for unknown,
1 fcc, 2 hcp, 3 bcc and 4 icosa
str : return a per atom string of sructures
ovito_compatibility (bool) – use ovito compatiblity mode
- Returns
(depends on mode)
-
analyse_pyscal_centro_symmetry
(num_neighbors=12)[source]¶ Analyse centrosymmetry parameter
- Parameters
num_neighbors (int) – number of neighbors
- Returns
list of centrosymmetry parameter
- Return type
list
-
analyse_pyscal_cna_adaptive
(mode='total', ovito_compatibility=False)[source]¶ Use common neighbor analysis
- Parameters
mode ("total"/"numeric"/"str") –
Controls the style and level of detail of the output. - total : return number of atoms belonging to each structure - numeric : return a per atom list of numbers- 0 for unknown,
1 fcc, 2 hcp, 3 bcc and 4 icosa
str : return a per atom string of sructures
ovito_compatibility (bool) – use ovito compatiblity mode
- Returns
(depends on mode)
-
analyse_pyscal_diamond_structure
(mode='total', ovito_compatibility=False)[source]¶ Analyse diamond structure
- Parameters
mode ("total"/"numeric"/"str") – Controls the style and level
detail of the output. (of) –
total : return number of atoms belonging to each structure
- numericreturn a per atom list of numbers- 0 for unknown,
1 fcc, 2 hcp, 3 bcc and 4 icosa
str : return a per atom string of sructures
ovito_compatibility (bool) – use ovito compatiblity mode
- Returns
(depends on mode)
-
analyse_pyscal_steinhardt_parameter
(neighbor_method='cutoff', cutoff=0, n_clusters=2, q=(4, 6), averaged=False, clustering=True)[source]¶ Calculate Steinhardts parameters
- Parameters
neighbor_method (str) – can be [‘cutoff’, ‘voronoi’]
cutoff (float) – can be 0 for adaptive cutoff or any other value
n_clusters (int) – number of clusters for K means clustering
q (list) – can be from 2-12, the required q values to be calculated
averaged (bool) – If True, calculates the averaged versions of the parameter
clustering (bool) – If True, cluster based on the q values
- Returns
calculated q parameters
- Return type
list
-
apply_strain
(epsilon, return_box=False)[source]¶ Apply a given strain on the structure
- Parameters
epsilon (float/list/ndarray) – epsilon matrix. If a single number is set, the same strain is applied in each direction. If a 3-dim vector is set, it will be multiplied by a unit matrix.
return_box (bool) – whether to return a box. If set to True, only the returned box will have the desired strain and the original box will stay unchanged.
-
center_coordinates_in_unit_cell
(origin=0, eps=0.0001)[source]¶ Wrap atomic coordinates within the supercell as given by a1, a2., a3
- Parameters
origin (float) – 0 to confine between 0 and 1, -0.5 to confine between -0.5 and 0.5
eps (float) – Tolerance to detect atoms at cell edges
- Returns
Wrapped structure
- Return type
-
cluster_analysis
(id_list, neighbors=None, radius=None, return_cluster_sizes=False)[source]¶ - Parameters
id_list –
neighbors –
radius –
return_cluster_sizes –
Returns:
-
convert_element
(el, pse=None)[source]¶ Convert a string or an atom instance into a ChemicalElement instance
- Parameters
el (str/atomistics.structure.atom.Atom) – String or atom instance from which the element should be generated
pse (atomistics.structure.periodictable.PeriodicTable) – PeriodicTable instance from which the element is generated (optional)
- Returns
The required chemical element
- Return type
atomistics.structure.periodictable.ChemicalElement
-
static
convert_formula
(elements)[source]¶ Convert a given chemical formula into a list of elements
- Parameters
elements (str) – A string of the required chamical formula (eg. H2O)
- Returns
A list of elements corresponding to the formula
- Return type
list
-
create_line_mode_structure
(with_time_reversal=True, recipe='hpkot', threshold=1e-07, symprec=1e-05, angle_tolerance=- 1.0)[source]¶ Uses ‘seekpath’ to create a new structure with high symmetry points and path for band structure calculations.
- Parameters
with_time_reversal (bool) – if False, and the group has no inversion symmetry, additional lines are returned as described in the HPKOT paper.
recipe (str) – choose the reference publication that defines the special points and paths. Currently, only ‘hpkot’ is implemented.
threshold (float) – the threshold to use to verify if we are in and edge case (e.g., a tetragonal cell, but a==c). For instance, in the tI lattice, if abs(a-c) < threshold, a EdgeCaseWarning is issued. Note that depending on the bravais lattice, the meaning of the threshold is different (angle, length, …)
symprec (float) – the symmetry precision used internally by SPGLIB
angle_tolerance (float) – the angle_tolerance used internally by SPGLIB
- Returns
new structure
- Return type
-
property
elements
¶ A size N list of atomistics.structure.periodic_table.ChemicalElement instances according to the ordering of the atoms in the instance
- Type
numpy.ndarray
-
extend
(other)[source]¶ Extend atoms object by appending atoms from other. (Extending the ASE function)
- Parameters
other (pyiron.atomistics.structure.atoms.Atoms/ase.atoms.Atoms) – Structure to append
- Returns
The extended structure
- Return type
-
find_mic
(v, vectors=True)[source]¶ Find vectors following minimum image convention (mic). In principle this function does the same as ase.geometry.find_mic
- Parameters
v (list/numpy.ndarray) – 3d vector or a list/array of 3d vectors
vectors (bool) – Whether to return vectors (distances are returned if False)
Returns: numpy.ndarray of the same shape as input with mic
-
find_neighbors_by_vector
(vector, return_deviation=False, num_neighbors=96)[source]¶ - Parameters
vector (list/np.ndarray) – vector by which positions are translated (and neighbors are searched)
return_deviation (bool) – whether to return distance between the expect positions and real positions
- Returns
list of id’s for the specified translation
- Return type
np.ndarray
Example
a_0 = 2.832 structure = pr.create_structure(‘Fe’, ‘bcc’, a_0) id_list = structure.find_neighbors_by_vector([0, 0, a_0]) # In this example, you get a list of neighbor atom id’s at z+=a_0 for each atom. # This is particularly powerful for SSA when the magnetic structure has to be translated # in each direction.
-
from_hdf
(hdf, group_name='structure')[source]¶ Retrieve the object from a HDF5 file
- Parameters
hdf (pyiron_base.generic.hdfio.FileHDFio) – HDF path to which the object is to be saved
group_name (str) – Group name from which the Atoms object is retreived.
- Returns
The retrieved atoms class
- Return type
pyiron_atomistic.structure.atoms.Atoms
-
get_atomic_numbers
()[source]¶ Returns the atomic numbers of all the atoms in the structure
- Returns
A list of atomic numbers
- Return type
numpy.ndarray
-
get_bonds
(radius=inf, max_shells=None, prec=0.1, num_neighbors=20)[source]¶ - Parameters
radius –
max_shells –
prec – minimum distance between any two clusters (if smaller considered to be single cluster)
num_neighbors –
Returns:
-
get_chemical_elements
()[source]¶ Returns the list of chemical element instances
- Returns
A list of chemical element instances
- Return type
numpy.ndarray
-
get_chemical_formula
()[source]¶ Returns the chemical formula of structure
- Returns
The chemical formula as a string
- Return type
str
-
get_chemical_indices
()[source]¶ Returns the list of chemical indices as ordered in self.species
- Returns
A list of chemical indices
- Return type
numpy.ndarray
-
get_chemical_symbols
()[source]¶ Returns the chemical symbols for all the atoms in the structure
- Returns
A list of chemical symbols
- Return type
numpy.ndarray
-
get_density
()[source]¶ Returns the density in g/cm^3
- Returns
Density of the structure
- Return type
float
-
get_distance
(a0, a1, mic=True, vector=False)[source]¶ Return distance between two atoms. Use mic=True to use the Minimum Image Convention. vector=True gives the distance vector (from a0 to a1).
- Parameters
a0 (int/numpy.ndarray/list) – position or atom ID
a1 (int/numpy.ndarray/list) – position or atom ID
mic (bool) – minimum image convention (True if periodic boundary conditions should be considered)
vector (bool) – True, if instead of distnce the vector connecting the two positions should be returned
- Returns
distance or vectors in length unit
- Return type
float
-
get_distances_array
(p1=None, p2=None, mic=True, vectors=False)[source]¶ Return distance matrix of every position in p1 with every position in p2. If p2 is not set, it is assumed that distances between all positions in p1 are desired. p2 will be set to p1 in this case. If both p1 and p2 are not set, the distances between all atoms in the box are returned.
- Parameters
p1 (numpy.ndarray/list) – Nx3 array of positions
p2 (numpy.ndarray/list) – Nx3 array of positions
mic (bool) – minimum image convention
vectors (bool) – return vectors instead of distances
- Returns
NxN if vector=False and NxNx3 if vector=True
- Return type
numpy.ndarray
-
get_equivalent_points
(points, use_magmoms=False, use_elements=True, symprec=1e-05, angle_tolerance=- 1.0)[source]¶ - Parameters
points (list/ndarray) – 3d vector
use_magmoms (bool) – cf. get_symmetry()
use_elements (bool) – cf. get_symmetry()
symprec (float) – cf. get_symmetry()
angle_tolerance (float) – cf. get_symmetry()
- Returns
array of equivalent points with respect to box symmetries
- Return type
(ndarray)
-
get_equivalent_voronoi_vertices
(return_box=False, minimum_dist=0.1, symprec=1e-05, angle_tolerance=- 1.0)[source]¶ This function gives the positions of spatially equivalent Voronoi vertices in lists, which most likely represent interstitial points or vacancies (along with other high symmetry points) Each list item contains an array of positions which are spacially equivalent. This function does not work if there are Hs atoms in the box
- Parameters
return_box – True, if the box containing atoms on the positions of Voronoi vertices should be returned (which are represented by Hs atoms)
minimum_dist – Minimum distance between two Voronoi vertices to be considered as one
Returns: List of numpy array positions of spacially equivalent Voronoi vertices
-
get_extended_positions
(width, return_indices=False, norm_order=2)[source]¶ Get all atoms in the boundary around the supercell which have a distance to the supercell boundary of less than dist
- Parameters
width (float) – width of the buffer layer on every periodic box side within which all atoms across periodic boundaries are chosen.
return_indices (bool) – Whether or not return the original indices of the appended atoms.
norm_order (float) – Order of Lp-norm.
- Returns
- Positions of all atoms in the extended box, indices of atoms in
their original option (if return_indices=True)
- Return type
numpy.ndarray
-
get_high_symmetry_path
()[source]¶ Path used for band structure calculations
- Returns
dict of pathes with start and end points.
- Return type
dict
-
get_high_symmetry_points
()[source]¶ dictionary of high-symmetry points defined for this specific structure.
- Returns
high_symmetry_points
- Return type
dict
-
get_initial_magnetic_moments
()[source]¶ Get array of initial magnetic moments.
- Returns
numpy.array()
-
get_ir_reciprocal_mesh
(mesh, is_shift=array([0, 0, 0], dtype=int32), is_time_reversal=True, symprec=1e-05)[source]¶ - Parameters
mesh –
is_shift –
is_time_reversal –
symprec –
Returns:
-
get_majority_species
(return_count=False)[source]¶ This function returns the majority species and their number in the box
- Returns
number of atoms of the majority species, chemical symbol and chemical index
-
get_masses
()[source]¶ Gets the atomic masses of all atoms in the structure
- Returns
Array of masses
- Return type
numpy.ndarray
-
get_neighborhood
(positions, num_neighbors=12, t_vec=True, cutoff_radius=inf, width_buffer=1.2, norm_order=2)[source]¶ - Parameters
position – Position in a box whose neighborhood information is analysed
num_neighbors (int) – Number of nearest neighbors
t_vec (bool) – True: compute distance vectors (pbc are taken into account)
cutoff_radius (float) – Upper bound of the distance to which the search is to be done
width_buffer (float) – Width of the layer to be added to account for pbc.
norm_order (int) – Norm to use for the neighborhood search and shell recognition. The definition follows the conventional Lp norm (cf. https://en.wikipedia.org/wiki/Lp_space). This is an feature and for anything other than norm_order=2, there is no guarantee that this works flawlessly.
- Returns
- Neighbors instances with the neighbor
indices, distances and vectors
- Return type
pyiron.atomistics.structure.atoms.Tree
-
get_neighbors
(num_neighbors=12, t_vec=True, tolerance=2, id_list=None, cutoff_radius=inf, width_buffer=1.2, allow_ragged=False, norm_order=2)[source]¶ - Parameters
num_neighbors (int) – number of neighbors
t_vec (bool) – True: compute distance vectors (pbc are automatically taken into account)
tolerance (int) – tolerance (round decimal points) used for computing neighbor shells
id_list (list) – list of atoms the neighbors are to be looked for
cutoff_radius (float) – Upper bound of the distance to which the search must be done
width_buffer (float) – width of the layer to be added to account for pbc.
allow_ragged (bool) – Whether to allow ragged list of arrays or rectangular numpy.ndarray filled with np.inf for values outside cutoff_radius
norm_order (int) – Norm to use for the neighborhood search and shell recognition. The definition follows the conventional Lp norm (cf. https://en.wikipedia.org/wiki/Lp_space). This is an feature and for anything other than norm_order=2, there is no guarantee that this works flawlessly.
- Returns
- Neighbors instances with the neighbor
indices, distances and vectors
- Return type
pyiron.atomistics.structure.atoms.Neighbors
-
get_neighbors_by_distance
(cutoff_radius=5, num_neighbors=None, t_vec=True, tolerance=2, id_list=None, width_buffer=1.2, allow_ragged=True, norm_order=2)[source]¶ - Parameters
cutoff_radius (float) – Upper bound of the distance to which the search must be done
num_neighbors (int/None) – maximum number of neighbors found; if None this is estimated based on the density.
t_vec (bool) – True: compute distance vectors (pbc are automatically taken into account)
tolerance (int) – tolerance (round decimal points) used for computing neighbor shells
id_list (list) – list of atoms the neighbors are to be looked for
width_buffer (float) – width of the layer to be added to account for pbc.
allow_ragged (bool) – Whether to allow ragged list of arrays or rectangular numpy.ndarray filled with np.inf for values outside cutoff_radius
norm_order (int) – Norm to use for the neighborhood search and shell recognition. The definition follows the conventional Lp norm (cf. https://en.wikipedia.org/wiki/Lp_space). This is an feature and for anything other than norm_order=2, there is no guarantee that this works flawlessly.
- Returns
- Neighbors instances with the neighbor
indices, distances and vectors
- Return type
pyiron.atomistics.structure.atoms.Neighbors
-
get_number_species_atoms
()[source]¶ Returns a dictionary with the species in the structure and the corresponding count in the structure
- Returns
An ordered dictionary with the species and the corresponding count
- Return type
collections.OrderedDict
-
get_numbers_of_neighbors_in_sphere
(cutoff_radius=10, num_neighbors=None, id_list=None, width_buffer=1.2)[source]¶ Function to compute the maximum number of neighbors in a sphere around each atom. :param cutoff_radius: Upper bound of the distance to which the search must be done :type cutoff_radius: float :param num_neighbors: maximum number of neighbors found :type num_neighbors: int/None :param id_list: list of atoms the neighbors are to be looked for :type id_list: list :param width_buffer: width of the layer to be added to account for pbc. :type width_buffer: float
- Returns
- for each atom the number of neighbors found in the sphere of radius
cutoff_radius (<= num_neighbors if specified)
- Return type
(np.ndarray)
-
get_parent_basis
()[source]¶ Returns the basis with all user defined/special elements as the it’s parent
- Returns
Structure without any user defined elements
- Return type
-
get_parent_symbols
()[source]¶ Returns the chemical symbols for all the atoms in the structure even for user defined elements
- Returns
A list of chemical symbols
- Return type
numpy.ndarray
-
get_primitive_cell
(symprec=1e-05, angle_tolerance=- 1.0)[source]¶ - Parameters
symprec –
angle_tolerance –
Returns:
-
get_shell_matrix
(id_list=None, chemical_pair=None, num_neighbors=100, tolerance=2, cluster_by_distances=False, cluster_by_vecs=False)[source]¶ Shell matrices for pairwise interaction. Note: The matrices are always symmetric, meaning if you use them as bilinear operators, you have to divide the results by 2.
- Parameters
chemical_pair (list) – pair of chemical symbols (e.g. [‘Fe’, ‘Ni’])
- Returns
list of sparse matrices for different shells
Example
from pyiron import Project structure = Project(‘.’).create_structure(‘Fe’, ‘bcc’, 2.83).repeat(2) J = -0.1 # Ising parameter magmoms = 2*np.random.random((len(structure)), 3)-1 # Random magnetic moments between -1 and 1 neigh = structure.get_neighbors(num_neighbors=8) # Iron first shell shell_matrices = neigh.get_shell_matrix() print(‘Energy =’, 0.5*J*magmoms.dot(shell_matrices[0].dot(matmoms)))
-
get_spacegroup
(symprec=1e-05, angle_tolerance=- 1.0)[source]¶ - Parameters
symprec –
angle_tolerance –
Returns:
-
get_species_symbols
()[source]¶ Returns the symbols of the present species
- Returns
List of the symbols of the species
- Return type
numpy.ndarray
-
get_spherical_coordinates
(x=None)[source]¶ - Parameters
x (list/ndarray) – coordinates to transform. If not set, the positions in structure will be returned.
- Returns
array in spherical coordinates
-
get_symmetry
(use_magmoms=False, use_elements=True, symprec=1e-05, angle_tolerance=- 1.0)[source]¶ - Parameters
use_magmoms –
use_elements – True or False. If False, chemical elements will be ignored
symprec –
angle_tolerance –
Returns:
-
get_symmetry_dataset
(symprec=1e-05, angle_tolerance=- 1.0)[source]¶ - Parameters
symprec –
angle_tolerance –
Returns:
Returns the keys of the stored tags of the structure
- Returns
Keys of the stored tags
- Return type
dict_keys
-
get_volume
(per_atom=False)[source]¶ - Parameters
per_atom (bool) – True if volume per atom is to be returned
- Returns
Volume in A**3
- Return type
volume (float)
-
get_wrapped_coordinates
(positions, epsilon=1e-08)[source]¶ Return coordinates in wrapped in the periodic cell
- Parameters
positions (list/numpy.ndarray) – Positions
epsilon (float) – displacement to add to avoid wrapping of atoms at borders
- Returns
Wrapped positions
- Return type
numpy.ndarray
-
group_points_by_symmetry
(points)[source]¶ This function classifies the points into groups according to the box symmetry given by spglib.
- Parameters
points – (np.array/list) nx3 array which contains positions
Returns: list of arrays containing geometrically equivalent positions
It is possible that the original points are not found in the returned list, as the positions outsie the box will be projected back to the box.
-
numbers_to_elements
(numbers)[source]¶ Convert atomic numbers in element objects (needed for compatibility with ASE)
- Parameters
numbers (list) – List of Element Numbers (as Integers; default in ASE)
- Returns
A list of elements as needed for pyiron
- Return type
list
-
plot3d
(mode='NGLview', show_cell=True, show_axes=True, camera='orthographic', spacefill=True, particle_size=1.0, select_atoms=None, background='white', color_scheme=None, colors=None, scalar_field=None, scalar_start=None, scalar_end=None, scalar_cmap=None, vector_field=None, vector_color=None, magnetic_moments=False, view_plane=array([0, 0, 1]), distance_from_camera=1.0, opacity=1.0)[source]¶ Plot3d relies on NGLView or plotly to visualize atomic structures. Here, we construct a string in the “protein database”
The final widget is returned. If it is assigned to a variable, the visualization is suppressed until that variable is evaluated, and in the meantime more NGL operations can be applied to it to modify the visualization.
- Parameters
mode (str) – NGLView, plotly or ase
show_cell (bool) – Whether or not to show the frame. (Default is True.)
show_axes (bool) – Whether or not to show xyz axes. (Default is True.)
camera (str) – ‘perspective’ or ‘orthographic’. (Default is ‘perspective’.)
spacefill (bool) – Whether to use a space-filling or ball-and-stick representation. (Default is True, use space-filling atoms.)
particle_size (float) – Size of the particles. (Default is 1.)
select_atoms (numpy.ndarray) – Indices of atoms to show, either as integers or a boolean array mask. (Default is None, show all atoms.)
background (str) – Background color. (Default is ‘white’.)
color_scheme (str) – NGLView color scheme to use. (Default is None, color by element.)
colors (numpy.ndarray) – A per-atom array of HTML color names or hex color codes to use for atomic colors. (Default is None, use coloring scheme.)
scalar_field (numpy.ndarray) – Color each atom according to the array value (Default is None, use coloring scheme.)
scalar_start (float) – The scalar value to be mapped onto the low end of the color map (lower values are clipped). (Default is None, use the minimum value in scalar_field.)
scalar_end (float) – The scalar value to be mapped onto the high end of the color map (higher values are clipped). (Default is None, use the maximum value in scalar_field.)
scalar_cmap (matplotlib.cm) – The colormap to use. (Default is None, giving a blue-red divergent map.)
vector_field (numpy.ndarray) – Add vectors (3 values) originating at each atom. (Default is None, no vectors.)
vector_color (numpy.ndarray) – Colors for the vectors (only available with vector_field). (Default is None, vectors are colored by their direction.)
magnetic_moments (bool) – Plot magnetic moments as ‘scalar_field’ or ‘vector_field’.
view_plane (numpy.ndarray) – A Nx3-array (N = 1,2,3); the first 3d-component of the array specifies which plane of the system to view (for example, [1, 0, 0], [1, 1, 0] or the [1, 1, 1] planes), the second 3d-component (if specified, otherwise [1, 0, 0]) gives the horizontal direction, and the third component (if specified) is the vertical component, which is ignored and calculated internally. The orthonormality of the orientation is internally ensured, and therefore is not required in the function call. (Default is np.array([0, 0, 1]), which is view normal to the x-y plane.)
distance_from_camera (float) – Distance of the camera from the structure. Higher = farther away. (Default is 14, which also seems to be the NGLView default value.)
NGLView color schemes (Possible) – ” “, “picking”, “random”, “uniform”, “atomindex”, “residueindex”, “chainindex”, “modelindex”, “sstruc”, “element”, “resname”, “bfactor”, “hydrophobicity”, “value”, “volume”, “occupancy”
- Returns
The NGLView widget itself, which can be operated on further or viewed as-is.
- Return type
(nglview.NGLWidget)
Warning
Many features only work with space-filling atoms (e.g. coloring by a scalar field).
The colour interpretation of some hex codes is weird, e.g. ‘green’.
-
refine_cell
(symprec=1e-05, angle_tolerance=- 1.0)[source]¶ - Parameters
symprec –
angle_tolerance –
Returns:
-
repeat
(rep)[source]¶ Create new repeated atoms object.
The rep argument should be a sequence of three positive integers like (2,3,1) or a single integer (r) equivalent to (r,r,r).
-
repeat_points
(points, rep, centered=False)[source]¶ Return points with repetition given according to periodic boundary conditions
- Parameters
points (np.ndarray/list) – xyz vector or list/array of xyz vectors
rep (int/list/np.ndarray) – Repetition in each direction. If int is given, the same value is used for every direction
centered (bool) – Whether the original points should be in the center of repeated points.
- Returns
(np.ndarray) repeated points
-
rotate
(a=0.0, v=None, center=(0, 0, 0), rotate_cell=False, index_list=None)[source]¶ Rotate atoms based on a vector and an angle, or two vectors. This function is completely adopted from ASE code (https://wiki.fysik.dtu.dk/ase/_modules/ase/atoms.html#Atoms.rotate)
- Parameters
a (float/list) – Angle that the atoms is rotated around the vecor ‘v’. If an angle is not specified, the length of ‘v’ is used as the angle (default). The angle can also be a vector and then ‘v’ is rotated into ‘a’.
v (list/numpy.ndarray/string) – Vector to rotate the atoms around. Vectors can be given as strings: ‘x’, ‘-x’, ‘y’, … .
center (tuple/list/numpy.ndarray/str) – The center is kept fixed under the rotation. Use ‘COM’ to fix the center of mass, ‘COP’ to fix the center of positions or ‘COU’ to fix the center of cell.
= False (rotate_cell) – If true the cell is also rotated.
index_list (list/numpy.ndarray) – Indices of atoms to be rotated
Examples:
Rotate 90 degrees around the z-axis, so that the x-axis is rotated into the y-axis:
>>> atoms = Atoms() >>> atoms.rotate(90, 'z') >>> atoms.rotate(90, (0, 0, 1)) >>> atoms.rotate(-90, '-z') >>> atoms.rotate('x', 'y') >>> atoms.rotate((1, 0, 0), (0, 1, 0))
-
select_index
(el)[source]¶ Returns the indices of a given element in the structure
- Parameters
el (str/atomistics.structures.periodic_table.ChemicalElement/list) – Element for which the indices should be returned
- Returns
An array of indices of the atoms of the given element
- Return type
numpy.ndarray
-
select_parent_index
(el)[source]¶ Returns the indices of a given element in the structure ignoring user defined elements
- Parameters
el (str/atomistics.structures.periodic_table.ChemicalElement) – Element for which the indices should be returned
- Returns
An array of indices of the atoms of the given element
- Return type
numpy.ndarray
-
set_constraint
(constraint=None)[source]¶ Apply one or more constrains.
The constraint argument must be one constraint object or a list of constraint objects.
-
set_initial_magnetic_moments
(magmoms=None)[source]¶ Set array of initial magnetic moments.
- Parameters
magmoms (numpy.ndarray/list) – List of magneric moments
-
set_species
(value)[source]¶ Setting the species list
- Parameters
value (list) – A list atomistics.structure.periodic_table.ChemicalElement instances
-
property
species
¶ A list of atomistics.structure.periodic_table.ChemicalElement instances
- Type
list
-
symmetrize_vectors
(vectors, force_update=False, use_magmoms=False, use_elements=True, symprec=1e-05, angle_tolerance=- 1.0)[source]¶ Symmetrization of natom x 3 vectors according to box symmetries
- Parameters
vectors (ndarray/list) – natom x 3 array to symmetrize
force_update (bool) – whether to update the symmetry info
use_magmoms (bool) – cf. get_symmetry
use_elements (bool) – cf. get_symmetry
symprec (float) – cf. get_symmetry
angle_tolerance (float) – cf. get_symmetry
- Returns
(np.ndarray) symmetrized vectors
-
to_hdf
(hdf, group_name='structure')[source]¶ Save the object in a HDF5 file
- Parameters
hdf (pyiron_base.generic.hdfio.FileHDFio) – HDF path to which the object is to be saved
group_name (str) – Group name with which the object should be stored. This same name should be used to retrieve the object
-
property
visualize
¶
-
pyiron.atomistics.structure.atoms.
ase_to_pyiron
(ase_obj)[source]¶ Convert an ase.atoms.Atoms structure object to its equivalent pyiron structure
- Parameters
ase_obj (ase.atoms.Atoms) – The ase atoms instance to convert
- Returns
The equivalent pyiron structure
- Return type
-
pyiron.atomistics.structure.atoms.
default
(data, dflt)[source]¶ Helper function for setting default values.
- Parameters
data –
dflt –
Returns:
-
pyiron.atomistics.structure.atoms.
ovito_to_pyiron
(ovito_obj)[source]¶ - Parameters
ovito_obj –
Returns:
-
pyiron.atomistics.structure.atoms.
pyiron_to_pymatgen
(pyiron_obj)[source]¶ Convert pyiron atoms object to pymatgen atoms object
- Parameters
pyiron_obj – pyiron atoms object
- Returns
pymatgen atoms object
-
pyiron.atomistics.structure.atoms.
pymatgen_to_pyiron
(pymatgen_obj)[source]¶ Convert pymatgen atoms object to pyiron atoms object (pymatgen->ASE->pyiron)
- Parameters
pymatgen_obj – pymatgen atoms object
- Returns
pyiron atoms object