Tools¶
Mapping structures¶
- icet.tools.map_structure_to_reference(structure, reference, inert_species=None, tol_positions=0.0001, suppress_warnings=False, assume_no_cell_relaxation=False)[source]¶
Maps a structure onto a reference structure. This is often desirable when, for example, a structure has been relaxed using DFT, and one wants to use it as a training structure in a cluster expansion.
The function returns a tuple comprising the ideal supercell most closely matching the input structure and a dictionary with supplementary information concerning the mapping. The latter includes for example the largest deviation of any position in the input structure from its reference position (
drmax
), the average deviation of the positions in the input structure from the reference positions (dravg
), and the strain tensor for the input structure relative to the reference structure (strain_tensor
).The returned Atoms object provides further supplemental information via custom per-atom arrays including the atomic displacements (
Displacement
,Displacement_Magnitude
),the distances to the three closest sites (Minimum_Distances
), as well as a mapping between the indices of the returned Atoms object and those of the input structure (IndexMapping
).Note
The input and reference structures should adhere to the same crystallographic setting. In other words, they should not be rotated or translated with respect to each other.
- Parameters:
structure (
Atoms
) – Input structure, typically a relaxed structure.reference (
Atoms
) – Reference structure, which can but need not be the primitive structure.inert_species (
Optional
[List
[str
]]) – List of chemical symbols (e.g.,['Au', 'Pd']
) that are never substituted for a vacancy. The number of inert sites is used to rescale the volume of the input structure to match the reference structure.tol_positions (
float
) – Tolerance factor applied when scanning for overlapping positions in Ångstrom (forwarded toase.build.make_supercell()
).suppress_warnings (
bool
) – ifTrue
print no warnings of large strain or relaxation distances.assume_no_cell_relaxation (
bool
) –If
False
volume and cell metric of the input structure are rescaled to match the reference structure. This can be unnecessary (and counterproductive) for some structures, e.g., with many vacancies.Note
When setting this parameter to
False
the reference cell metric must be obtainable via an integer transformation matrix from the reference cell metric. In other words the input structure should not involve relaxations of the volume or the cell metric.
- Return type:
Tuple
[Atoms
,dict
]
Example
The following code snippet illustrates the general usage. It first creates a primitive FCC cell, which is latter used as reference structure. To emulate a relaxed structure obtained from, e.g., a density functional theory calculation, the code then creates a 4x4x4 conventional FCC supercell, which is populated with two different atom types, has distorted cell vectors, and random displacements to the atoms. Finally, the present function is used to map the structure back the ideal lattice:
>>> from ase.build import bulk >>> reference = bulk('Au', a=4.09) >>> structure = bulk('Au', cubic=True, a=4.09).repeat(4) >>> structure.set_chemical_symbols(10 * ['Ag'] + (len(structure) - 10) * ['Au']) >>> structure.set_cell(structure.cell * 1.02, scale_atoms=True) >>> structure.rattle(0.1) >>> mapped_structure, info = map_structure_to_reference(structure, reference)
Structure enumeration¶
- icet.tools.enumerate_structures(structure, sizes, chemical_symbols, concentration_restrictions=None, niggli_reduce=None, symprec=1e-05, position_tolerance=None)[source]¶
Yields a sequence of enumerated structures. The function generates all inequivalent structures that are permissible given a certain lattice. Using the
chemical_symbols
andconcentration_restrictions
keyword arguments it is possible to specify which chemical_symbols are to be included on which site and in which concentration range.The function is sensitive to the boundary conditions of the input structure. An enumeration of, for example, a surface can thus be performed by setting
structure.pbc = [True, True, False]
.The algorithm implemented here was developed by Gus L. W. Hart and Rodney W. Forcade in Phys. Rev. B 77, 224115 (2008) [HarFor08] and Phys. Rev. B 80, 014120 (2009) [HarFor09].
- Parameters:
structure (
Atoms
) – Primitive structure from which derivative superstructures should be generated.sizes (
Union
[List
[int
],range
]) – Number of sites (included in enumeration).chemical_symbols (
list
) – Chemical species with which to decorate the structure, e.g.,['Au', 'Ag']
; see below for more examples.concentration_restrictions (
Optional
[dict
]) – Allowed concentration range for one or more element inchemical_symbols
, e.g.,{'Au': (0, 0.2)}
will only enumerate structures in which the Au content is between 0 and 20 %. Here, concentration is always defined as the number of atoms of the specified kind divided by the number of all atoms.niggli_reduction – If
True
perform a Niggli reduction with spglib for each structure. The default isTrue
ifstructure
is periodic in all directions,False
otherwise.symprec (
float
) – Tolerance imposed when analyzing the symmetry using spglib.position_tolerance (
Optional
[float
]) – Tolerance applied when comparing positions in Cartesian coordinates; by default this value is set equal tosymprec
.
- Return type:
Examples
The following code snippet illustrates how to enumerate structures with up to 6 atoms in the unit cell for a binary alloy without any constraints:
>>> from ase.build import bulk >>> prim = bulk('Ag') >>> for structure in enumerate_structures(structure=prim, ... sizes=range(1, 5), ... chemical_symbols=['Ag', 'Au']): ... pass # Do something with the structure
To limit the concentration range to 10 to 40% Au the code should be modified as follows:
>>> conc_restr = {'Au': (0.1, 0.4)} >>> for structure in enumerate_structures(structure=prim, ... sizes=range(1, 5), ... chemical_symbols=['Ag', 'Au'], ... concentration_restrictions=conc_restr): ... pass # Do something with the structure
Often one would like to consider mixing on only one sublattice. This can be achieved as illustrated for a Ga(1-x)Al(x)As alloy as follows:
>>> prim = bulk('GaAs', crystalstructure='zincblende', a=5.65) >>> for structure in enumerate_structures(structure=prim, ... sizes=range(1, 9), ... chemical_symbols=[['Ga', 'Al'], ['As']]): ... pass # Do something with the structure
- icet.tools.enumerate_supercells(structure, sizes, niggli_reduce=None, symprec=1e-05, position_tolerance=None)[source]¶
Yields a sequence of enumerated supercells. The function generates all inequivalent supercells that are permissible given a certain lattice. Any supercell can be reduced to one of the supercells generated.
The function is sensitive to the boundary conditions of the input structure. An enumeration of, for example, a surface can thus be performed by setting
structure.pbc = [True, True, False]
.The algorithm is based on Gus L. W. Hart and Rodney W. Forcade in Phys. Rev. B 77, 224115 (2008) [HarFor08] and Phys. Rev. B 80, 014120 (2009) [HarFor09].
- Parameters:
structure (
Atoms
) – Primitive structure from which supercells should be generated.sizes (
Union
[List
[int
],range
]) – Number of sites (included in enumeration).niggli_reduction – If
True
perform a Niggli reduction with spglib for each supercell. The default isTrue
ifstructure
is periodic in all directions,False
otherwise.symprec (
float
) – Tolerance imposed when analyzing the symmetry using spglib.position_tolerance (
Optional
[float
]) – Tolerance applied when comparing positions in Cartesian coordinates. By default this value is set equal tosymprec
.
- Return type:
Examples
The following code snippet illustrates how to enumerate supercells with up to 6 atoms in the unit cell:
>>> from ase.build import bulk >>> prim = bulk('Ag') >>> for supercell in enumerate_supercells(structure=prim, sizes=range(1, 7)): ... pass # Do something with the supercell
Generation of training structures¶
- icet.tools.training_set_generation.structure_selection_annealing(cluster_space, monte_carlo_structures, n_structures_to_add, n_steps, base_structures=None, cooling_start=5, cooling_stop=0.001, cooling_function='exponential', initial_indices=None)[source]¶
Given a cluster space, a base pool of structures, and a new pool of structures, this function uses a Monte Carlo inspired annealing method to find a good structure pool for training.
- Return type:
Tuple
[List
[int
],List
[float
]]- Returns:
A tuple comprising the indices of the optimal structures in the
monte_carlo_structures
pool and a list of accepted metric values.- Parameters:
cluster_space (
ClusterSpace
) – A cluster space defining the lattice to be occupied.monte_carlo_structures (
List
[Atoms
]) – A list of candidate training structures.n_structures_to_add (
int
) – How many of the structures in themonte_carlo_structures
pool that should be kept for training.n_steps (
int
) – Number of steps in the annealing algorithm.base_structures (
Optional
[List
[Atoms
]]) – A list of structures that is already in your training pool; can beNone
if you do not have any structures yet.cooling_start (
float
) – Initial value of thecooling_function
.cooling_stop (
float
) – Last value of thecooling_function
.cooling_function (
Union
[str
,Callable
]) – Artificial number that rescales the difference between the metric value between two iterations. Available options are'linear'
and'exponential'
.initial_indices (
Optional
[List
[int
]]) – Picks out the starting structure from themonte_carlo_structures
pool. Can be used if you want to continue from an old run for example.
Example
The following snippet demonstrates the use of this function for generating an optimized structure pool. Here, we first set up a pool of candidate structures by randomly occupying a FCC supercell with Au and Pd:
>>> from ase.build import bulk >>> from icet.tools.structure_generation import occupy_structure_randomly >>> prim = bulk('Au', a=4.0) >>> cs = ClusterSpace(prim, [6.0], [['Au', 'Pd']]) >>> structure_pool = [] >>> for _ in range(500): >>> # Create random supercell. >>> supercell = np.random.randint(1, 4, size=3) >>> structure = prim.repeat(supercell) >>> >>> # Randomize concentrations in the supercell >>> n_atoms = len(structure) >>> n_Au = np.random.randint(0, n_atoms) >>> n_Pd = n_atoms - n_Au >>> concentration = {'Au': n_Au / n_atoms, 'Pd': n_Pd / n_atoms} >>> >>> # Occupy the structure randomly and store it. >>> occupy_structure_randomly(structure, cs, concentration) >>> structure_pool.append(structure) >>> start_inds = [f for f in range(10)]
Now we can use the
structure_selection_annealing()
function to find an optimized structure pool:>>> inds, cond = structure_selection_annealing(cs, >>> structure_pool, >>> n_structures_to_add=10, >>> n_steps=100) >>> training_structures = [structure_pool[ind] for ind in inds] >>> print(training_structures)
Generation of special structures¶
- icet.tools.structure_generation.generate_sqs(cluster_space, max_size, target_concentrations, include_smaller_cells=True, pbc=None, T_start=5.0, T_stop=0.001, n_steps=None, optimality_weight=1.0, random_seed=None, tol=1e-05)[source]¶
Given a
cluster_space
, generate a special quasirandom structure (SQS), i.e., a structure that for a given supercell size provides the best possible approximation to a random alloy [ZunWeiFer90].In the present case, this means that the generated structure will have a cluster vector that as closely as possible matches the cluster vector of an infintely large randomly occupied supercell. Internally the function uses a simulated annealing algorithm and the difference between two cluster vectors is calculated with the measure suggested by A. van de Walle et al. in Calphad 42, 13-18 (2013) [WalTiwJon13] (for more information, see
mchammer.calculators.TargetVectorCalculator
).- Parameters:
cluster_space (
ClusterSpace
) – Cluster space defining the lattice to be occupied.max_size (
int
) – Maximum supercell size.target_concentrations (
dict
) – Concentration of each species in the target structure, per sublattice (for example{'Au': 0.5, 'Pd': 0.5}
for a single sublattice Au-Pd structure, or{'A': {'Au': 0.5, 'Pd': 0.5}, 'B': {'H': 0.25, 'X': 0.75}}
for a system with two sublattices. The symbols defining sublattices (‘A’, ‘B’ etc) can be found by printing thecluster_space
.include_smaller_cells (
bool
) – IfTrue
, search among all supercell sizes includingmax_size
, else search only among those exactly matchingmax_size
pbc (
Union
[Tuple
[bool
,bool
,bool
],Tuple
[int
,int
,int
],None
]) – Periodic boundary conditions for each direction, e.g.,(True, True, False)
. The axes are defined by the cell ofcluster_space.primitive_structure
. Default is periodic boundary in all directions.T_start (
float
) – Artificial temperature at which the simulated annealing starts.T_stop (
float
) – Artifical temperature at which the simulated annealing stops.n_steps (
Optional
[int
]) – Total number of Monte Carlo steps in the simulation.optimality_weight (
float
) – Controls weighting \(L\) of perfect correlations, seemchammer.calculators.TargetVectorCalculator
.random_seed (
Optional
[int
]) – Seed for the random number generator used in the Monte Carlo simulation.tol (
float
) – Numerical tolerance.
- Return type:
- icet.tools.structure_generation.generate_sqs_by_enumeration(cluster_space, max_size, target_concentrations, include_smaller_cells=True, pbc=None, optimality_weight=1.0, tol=1e-05)[source]¶
Given a
cluster_space
, generate a special quasirandom structure (SQS), i.e., a structure that for a given supercell size provides the best possible approximation to a random alloy [ZunWeiFer90].In the present case, this means that the generated structure will have a cluster vector that as closely as possible matches the cluster vector of an infintely large randomly occupied supercell. Internally the function uses a simulated annealing algorithm and the difference between two cluster vectors is calculated with the measure suggested by A. van de Walle et al. in Calphad 42, 13-18 (2013) [WalTiwJon13] (for more information, see
mchammer.calculators.TargetVectorCalculator
).This functions generates SQS cells by exhaustive enumeration, which means that the generated SQS cell is guaranteed to be optimal with regard to the specified measure and cell size.
- Parameters:
cluster_space (
ClusterSpace
) – Cluster space defining the lattice to be occupied.max_size (
int
) – Maximum supercell size.target_concentrations (
dict
) – Concentration of each species in the target structure, per sublattice (for example{'Au': 0.5, 'Pd': 0.5}
for a single sublattice Au-Pd structure, or{'A': {'Au': 0.5, 'Pd': 0.5}, 'B': {'H': 0.25, 'X': 0.75}}
for a system with two sublattices. The symbols defining sublattices (‘A’, ‘B’ etc) can be found by printing thecluster_space
.include_smaller_cells (
bool
) – ifTrue
search among all supercell sizes includingmax_size
, else search only among those exactly matchingmax_size
.pbc (
Union
[Tuple
[bool
,bool
,bool
],Tuple
[int
,int
,int
],None
]) – Periodic boundary conditions for each direction, e.g.,(True, True, False)
. The axes are defined by the cell ofcluster_space.primitive_structure
. Default is periodic boundary in all directions.optimality_weight (
float
) – Controls weighting \(L\) of perfect correlations, seemchammer.calculators.TargetVectorCalculator
.tol (
float
) – Numerical tolerance.
- Return type:
- icet.tools.structure_generation.generate_sqs_from_supercells(cluster_space, supercells, target_concentrations, T_start=5.0, T_stop=0.001, n_steps=None, optimality_weight=1.0, random_seed=None, random_start=True, tol=1e-05)[source]¶
Given a
cluster_space`
and one or moresupercells
, generate a special quasirandom structure (SQS), i.e., a structure that for the provided supercells size provides the best possible approximation to a random alloy [ZunWeiFer90].In the present case, this means that the generated structure will have a cluster vector that as closely as possible matches the cluster vector of an infintely large randomly occupied supercell. Internally the function uses a simulated annealing algorithm and the difference between two cluster vectors is calculated with the measure suggested by A. van de Walle et al. in Calphad 42, 13-18 (2013) [WalTiwJon13] (for more information, see
mchammer.calculators.TargetVectorCalculator
).- Parameters:
cluster_space (
ClusterSpace
) – Cluster space defining the lattice to be occupied.supercells (
List
[Atoms
]) – List of one or more supercells among which an optimal structure will be searched for.target_concentrations (
dict
) – Concentration of each species in the target structure, per sublattice (for example{'Au': 0.5, 'Pd': 0.5}
for a single sublattice Au-Pd structure, or{'A': {'Au': 0.5, 'Pd': 0.5}, 'B': {'H': 0.25, 'X': 0.75}}
for a system with two sublattices. The symbols defining sublattices (‘A’, ‘B’ etc) can be found by printing thecluster_space
.T_start (
float
) – Artificial temperature at which the simulated annealing starts.T_stop (
float
) – Artifical temperature at which the simulated annealing stops.n_steps (
Optional
[int
]) – Total number of Monte Carlo steps in the simulation.optimality_weight (
float
) – Controls weighting \(L\) of perfect correlations, seemchammer.calculators.TargetVectorCalculator
.random_seed (
Optional
[int
]) – Seed for the random number generator used in the Monte Carlo simulation and used for initializing the occupation of the supercells if random_start isTrue
.random_start (
bool
) – Randomly occupy starting structure, can be disabled if the user prefers to pass an initial structure.tol (
float
) – Numerical tolerance.
- Return type:
- icet.tools.structure_generation.generate_target_structure(cluster_space, max_size, target_concentrations, target_cluster_vector, include_smaller_cells=True, pbc=None, T_start=5.0, T_stop=0.001, n_steps=None, optimality_weight=1.0, random_seed=None, tol=1e-05)[source]¶
Given a
cluster_space
and atarget_cluster_vector
, generate a structure that as closely as possible matches that cluster vector. The search is performed among all inequivalent supercells shapes up to a certain size.Internally the function uses a simulated annealing algorithm and the difference between two cluster vectors is calculated with the measure suggested by A. van de Walle et al. in Calphad 42, 13-18 (2013) [WalTiwJon13] (for more information, see
mchammer.calculators.TargetVectorCalculator
).- Parameters:
cluster_space (
ClusterSpace
) – Cluster space defining the lattice to be occupied.max_size (
int
) – Maximum supercell size.target_concentrations (
dict
) – Concentration of each species in the target structure, per sublattice (for example{'Au': 0.5, 'Pd': 0.5}
for a single sublattice Au-Pd structure, or{'A': {'Au': 0.5, 'Pd': 0.5}, 'B': {'H': 0.25, 'X': 0.75}}
for a system with two sublattices. The symbols defining sublattices (‘A’, ‘B’ etc) can be found by printing thecluster_space
.target_cluster_vector (
List
[float
]) – Cluster vector that the generated structure should match as closely as possible.include_smaller_cells (
bool
) – IfTrue
, search among all supercell sizes includingmax_size
, else search only among those exactly matchingmax_size
pbc (
Union
[Tuple
[bool
,bool
,bool
],Tuple
[int
,int
,int
],None
]) – Periodic boundary conditions for each direction, e.g.,(True, True, False)
. The axes are defined by the cell ofcluster_space.primitive_structure`
. Default is periodic boundary in all directions.T_start (
float
) – Artificial temperature at which the simulated annealing starts.T_stop (
float
) – Artifical temperature at which the simulated annealing stops.n_steps (
Optional
[int
]) – Total number of Monte Carlo steps in the simulation.optimality_weight (
float
) – Controls weighting \(L\) of perfect correlations, seemchammer.calculators.TargetVectorCalculator
.random_seed (
Optional
[int
]) – Seed for the random number generator used in the Monte Carlo simulation.tol (
float
) – Numerical tolerance.
- Return type:
- icet.tools.structure_generation.generate_target_structure_from_supercells(cluster_space, supercells, target_concentrations, target_cluster_vector, T_start=5.0, T_stop=0.001, n_steps=None, optimality_weight=1.0, random_seed=None, random_start=True, tol=1e-05)[source]¶
Given a
cluster_space
and atarget_cluster_vector
and one or moresupercells
, generate a structure that as closely as possible matches that cluster vector.Internally the function uses a simulated annealing algorithm and the difference between two cluster vectors is calculated with the measure suggested by A. van de Walle et al. in Calphad 42, 13-18 (2013) [WalTiwJon13] (for more information, see
mchammer.calculators.TargetVectorCalculator
).- Parameters:
cluster_space (
ClusterSpace
) – A cluster space defining the lattice to be occupied.supercells (
List
[Atoms
]) – List of one or more supercells among which an optimal structure will be searched for.target_concentrations (
dict
) – Concentration of each species in the target structure, per sublattice (for example{'Au': 0.5, 'Pd': 0.5}
for a single sublattice Au-Pd structure, or{'A': {'Au': 0.5, 'Pd': 0.5}, 'B': {'H': 0.25, 'X': 0.75}}
for a system with two sublattices. The symbols defining sublattices (‘A’, ‘B’ etc) can be found by printing thecluster_space
.target_cluster_vector (
List
[float
]) – Cluster vector that the generated structure should match as closely as possible.T_start (
float
) – Artificial temperature at which the simulated annealing starts.T_stop (
float
) – Artifical temperature at which the simulated annealing stops.n_steps (
Optional
[int
]) – Total number of Monte Carlo steps in the simulation.optimality_weight (
float
) – Controls weighting \(L\) of perfect correlations, seemchammer.calculators.TargetVectorCalculator
.random_seed (
Optional
[int
]) – Seed for the random number generator used in the Monte Carlo simulation and used for initializing the occupation of the supercells if random_start isTrue
.random_start (
bool
) – Randomly occupy starting structure, can be disabled if the user prefers to pass an initial structure.tol (
float
) – Numerical tolerance.
- Return type:
- icet.tools.structure_generation.occupy_structure_randomly(structure, cluster_space, target_concentrations, random_seed=None)[source]¶
Occupy a structure with quasirandom order but fulfilling
target_concentrations
.- Parameters:
structure (
Atoms
) – ASE Atoms object that will be occupied randomly.cluster_space (
ClusterSpace
) – Cluster space (needed as it carries information about sublattices).target_concentrations (
dict
) – Concentration of each species in the target structure, per sublattice (for example{'Au': 0.5, 'Pd': 0.5}
for a single sublattice Au-Pd structure, or{'A': {'Au': 0.5, 'Pd': 0.5}, 'B': {'H': 0.25, 'X': 0.75}}
for a system with two sublattices. The symbols defining sublattices (‘A’, ‘B’ etc) can be found by printing thecluster_space
.random_seed (
Optional
[int
]) – Seed for the random number generator.
- Return type:
None
Ground state finder¶
- class icet.tools.ground_state_finder.GroundStateFinder(cluster_expansion, structure, solver_name=None, verbose=True)[source]¶
This class provides functionality for determining the ground states using a binary cluster expansion. This is efficiently achieved through the use of mixed integer programming (MIP) as developed by Larsen et al. in Phys. Rev. Lett. 120, 256101 (2018).
This class relies on the Python-MIP package. Python-MIP can be used together with Gurobi, which is not open source but issues academic licenses free of charge. Pleaase note that Gurobi needs to be installed separately. The
GroundStateFinder
works also without Gurobi, but if performance is critical, Gurobi is highly recommended.Warning
In order to be able to use Gurobi with python-mip one must ensure that GUROBI_HOME should point to the installation directory (
<installdir>
):export GUROBI_HOME=<installdir>
Note
The current implementation only works for binary systems.
- Parameters:
cluster_expansion (
ClusterExpansion
) – Cluster expansion for which to find ground states.structure (
Atoms
) – Atomic configuration.solver_name (
Optional
[str
]) –'gurobi'
,'grb'
or'cbc'
. Searches for available solvers if no value is provided.verbose (
bool
) – IfTrue
print solver messages to stdout.
Example
The following snippet illustrates how to determine the ground state for a Au-Ag alloy. Here, the parameters of the cluster expansion are set to emulate a simple Ising model in order to obtain an example that can be run without modification. In practice, one should of course use a proper cluster expansion:
>>> from ase.build import bulk >>> from icet import ClusterExpansion, ClusterSpace >>> # prepare cluster expansion >>> # the setup emulates a second nearest-neighbor (NN) Ising model >>> # (zerolet and singlet parameters are zero; only first and second neighbor >>> # pairs are included) >>> prim = bulk('Au') >>> chemical_symbols = ['Ag', 'Au'] >>> cs = ClusterSpace(prim, cutoffs=[4.3], chemical_symbols=chemical_symbols) >>> ce = ClusterExpansion(cs, [0, 0, 0.1, -0.02]) >>> # prepare initial configuration >>> structure = prim.repeat(3) >>> # set up the ground state finder and calculate the ground state energy >>> gsf = GroundStateFinder(ce, structure) >>> ground_state = gsf.get_ground_state({'Ag': 5}) >>> print('Ground state energy:', ce.predict(ground_state))
- get_ground_state(species_count=None, max_seconds=inf, threads=0)[source]¶
Finds the ground state for a given structure and species count. If
species_count
is not provided when initializing the instance of this class the first species in the list of chemical symbols for the active sublattice will be used.- Parameters:
species_count (
Optional
[Dict
[str
,int
]]) – Dictionary with count for one of the species on each active sublattice. If no count is provided for a sublattice, the concentration is allowed to vary.max_seconds (
float
) – Maximum runtime in seconds.threads (
int
) – Number of threads to be used when solving the problem, given that a positive integer has been provided. If set to \(0\) the solver default configuration is used while \(-1\) corresponds to all available processing cores.
- Return type:
- property model: Model¶
Python-MIP model.
- property optimization_status: OptimizationStatus¶
Optimization status.
Convex hull construction¶
- class icet.tools.ConvexHull(concentrations, energies)[source]¶
This class provides functionality for extracting the convex hull of the (free) energy of mixing. It is based on the convex hull calculator in SciPy.
- Parameters:
concentrations (
Union
[List
[float
],List
[List
[float
]]]) – Concentrations for each structure listed as[[c1, c2], [c1, c2], ...]
; for binaries, in which case there is only one independent concentration, the format[c1, c2, c3, ...]
works as well.energies (
List
[float
]) – Energy (or energy of mixing) for each structure.
- concentrations¶
Concentrations of the structures on the convex hull.
- Type:
np.ndarray
- energies¶
Energies of the structures on the convex hull.
- Type:
np.ndarray
- dimensions¶
Number of independent concentrations needed to specify a point in concentration space (1 for binaries, 2 for ternaries and so on).
- Type:
int
- structures¶
Indices of structures that constitute the convex hull (indices are defined by the order of their concentrations and energies are fed when initializing the
ConvexHull
object).- Type:
List[int]
Examples
A
ConvexHull
object is easily initialized by providing lists of concentrations and energies:>>> data = {'concentration': [0, 0.2, 0.2, 0.3, 0.4, 0.5, 0.8, 1.0], ... 'mixing_energy': [0.1, -0.2, -0.1, -0.2, 0.2, -0.4, -0.2, -0.1]} >>> hull = ConvexHull(data['concentration'], data['mixing_energy'])
Now one can for example access the points along the convex hull directly:
>>> for c, e in zip(hull.concentrations, hull.energies): ... print(c, e) 0.0 0.1 0.2 -0.2 0.5 -0.4 1.0 -0.1
or plot the convex hull along with the original data using e.g., matplotlib:
>>> import matplotlib.pyplot as plt >>> plt.scatter(data['concentration'], data['mixing_energy'], color='darkred') >>> plt.plot(hull.concentrations, hull.energies) >>> plt.show(block=False)
It is also possible to extract structures at or close to the convex hull:
>>> low_energy_structures = hull.extract_low_energy_structures( ... data['concentration'], data['mixing_energy'], ... energy_tolerance=0.005)
A complete example can be found in the basic tutorial.
- extract_low_energy_structures(concentrations, energies, energy_tolerance)[source]¶
Returns the indices of energies that lie within a certain tolerance of the convex hull.
- Parameters:
concentrations (
Union
[List
[float
],List
[List
[float
]]]) –Concentrations of candidate structures.
If there is one independent concentration, a list of floats is sufficient. Otherwise, the concentrations must be provided as a list of lists, such as
[[0.1, 0.2], [0.3, 0.1], ...]
.energies (
List
[float
]) – Energies of candidate structures.energy_tolerance (
float
) – Include structures with an energy that is at most this far from the convex hull.
- Return type:
List
[int
]
- get_energy_at_convex_hull(target_concentrations)[source]¶
Returns the energy of the convex hull at specified concentrations. If any concentration is outside the allowed range, NaN is returned.
- Parameters:
target_concentrations (
Union
[List
[float
],List
[List
[float
]]]) –Concentrations at target points.
If there is one independent concentration, a list of floats is sufficient. Otherwise, the concentrations ought to be provided as a list of lists, such as
[[0.1, 0.2], [0.3, 0.1], ...]
.- Return type:
Fitting with constraints¶
- class icet.tools.constraints.Constraints(n_params)[source]¶
Class for handling linear constraints with right-hand-side equal to zero.
- Parameters:
n_params (
int
) – Number of parameters in model.
Example
The following example demonstrates fitting of a cluster expansion under the constraint that parameter 2 and parameter 4 should be equal:
>>> import numpy as np >>> from icet.tools import Constraints >>> from trainstation import Optimizer >>> # Set up random sensing matrix and target "energies" >>> n_params = 10 >>> n_energies = 20 >>> A = np.random.random((n_energies, n_params)) >>> y = np.random.random(n_energies) >>> # Define constraints >>> c = Constraints(n_params=n_params) >>> M = np.zeros((1, n_params)) >>> M[0, [2, 4]] = 1 >>> c.add_constraint(M) >>> # Do the actual fit and finally extract parameters >>> A_constrained = c.transform(A) >>> opt = Optimizer((A_constrained, y), fit_method='ridge') >>> opt.train() >>> parameters = c.inverse_transform(opt.parameters)
- add_constraint(M)[source]¶
Add a constraint matrix and resolve for the constraint space.
- Parameters:
M (
ndarray
) – Constraint matrix with each constraint as a row. Can (but need not be) cluster vectors.- Return type:
None
- icet.tools.constraints.get_mixing_energy_constraints(cluster_space)[source]¶
A cluster expansion of the mixing energy should ideally predict zero energy for concentrations 0 and 1. This function constructs a
Constraints
object that enforces that condition during training.- Parameters:
cluster_space – Cluster space corresponding to cluster expansion for which constraints should be imposed.
- Return type:
Example
This example demonstrates how to constrain the mixing energy to zero at the pure phases in a toy example with random cluster vectors and random target energies:
>>> import numpy as np >>> from ase.build import bulk >>> from icet import ClusterSpace >>> from icet.tools import get_mixing_energy_constraints >>> from trainstation import Optimizer >>> # Set up cluster space along with random sensing matrix and target "energies" >>> prim = bulk('Au') >>> cs = ClusterSpace(prim, cutoffs=[6.0, 5.0], chemical_symbols=['Au', 'Ag']) >>> n_params = len(cs) >>> n_energies = 20 >>> A = np.random.random((n_energies, n_params)) >>> y = np.random.random(n_energies) >>> # Define constraints >>> c = get_mixing_energy_constraints(cs) >>> # Do the actual fit and finally extract parameters >>> A_constrained = c.transform(A) >>> opt = Optimizer((A_constrained, y), fit_method='ridge') >>> opt.train() >>> parameters = c.inverse_transform(opt.parameters)
Warning
Constraining the energy of one structure is always done at the expense of the fit quality of the others. Always expect that your cross-validation scores will increase somewhat when using this function.
Constituent strain¶
- class icet.tools.ConstituentStrain(supercell, primitive_structure, chemical_symbols, concentration_symbol, strain_energy_function, k_to_parameter_function=None, damping=1.0, tol=1e-06)[source]¶
Class for handling constituent strain in cluster expansions (see Laks et al., Phys. Rev. B 46, 12587 (1992) [LakFerFro92]). This makes it possible to use cluster expansions to describe systems with strain due to, for example, coherent phase separation. For an extensive example on how to use this module, please see this example.
- Parameters:
supercell (
Atoms
) – Defines supercell that will be used when calculating constituent strain.primitive_structure (
Atoms
) – Primitive structure the supercell is based on.chemical_symbols (
List
[str
]) – List with chemical symbols involved, such as['Ag', 'Cu']
.concentration_symbol (
str
) – Chemical symbol used to define concentration, such as'Ag'
.strain_energy_function (
Callable
[[float
,List
[float
]],float
]) – A function that takes two arguments, a list of parameters and concentration (e.g.,[0.5, 0.5, 0.5]
and0.3
), and returns the corresponding strain energy. The parameters are in turn determined byk_to_parameter_function
(see below). Ifk_to_parameter_function
is None, the parameters list will be the k-point. For more information, see this example.k_to_parameter_function (
Optional
[Callable
[[List
[float
]],List
[float
]]]) – A function that takes a k-point as a list of three floats and returns a parameter vector that will be fed intostrain_energy_function
(see above). IfNone
, the k-point itself will be the parameter vector tostrain_energy_function
. The purpose of this function is to be able to precompute any factor in the strain energy that depends on the k-point but not the concentration. For more information, see this example.damping (
float
) – Damping factor \(\eta\) used to suppress impact of large-magnitude k-points by multiplying strain with \(\exp(-(\eta \mathbf{k})^2)\) (unit Ångstrom).tol (
float
) – Numerical tolerance when comparing k-points (units of inverse Ångstrom).
- accept_change()[source]¶
Update structure factor for each kpoint to the value in
structure_factor_after
. This makes it possible to efficiently calculate changes in constituent strain with theget_constituent_strain_change()
function; this function should be called if the last occupations used to callget_constituent_strain_change()
should be the starting point for the next call ofget_constituent_strain_change()
. This is taken care of automatically by the Monte Carlo simulations in mchammer.- Return type:
None
- get_concentration(occupations)[source]¶
Calculate current concentration.
- Parameters:
occupations (
ndarray
) – Current occupations.- Return type:
float
- get_constituent_strain(occupations)[source]¶
Calculate total constituent strain.
- Parameters:
occupations (
List
[int
]) – Current occupations.- Return type:
float
- get_constituent_strain_change(occupations, atom_index)[source]¶
Calculate change in constituent strain upon change of the occupation of one site.
Warning
This function is dependent on the internal state of the
ConstituentStrain
object and should typically only be used internally by mchammer. Specifically, the structure factor is saved internally to speed up computation. The first time this function is called,occupations
must be the same array as was used to initialize theConstituentStrain
object, or the same as was last used whenget_constituent_strain()
was called. After the present function has been called, the same occupations vector need to be used the next time as well, unlessaccept_change()
has been called, in which caseoccupations
should incorporate the changes implied by the previous call to the function.- Parameters:
occupations (
ndarray
) – Occupations before change.atom_index (
int
) – Index of site the occupation of which is to be changed.
- Return type:
float
Other structure tools¶
- icet.tools.get_primitive_structure(structure, no_idealize=True, to_primitive=True, symprec=1e-05)[source]¶
Returns the primitive structure using spglib.
- icet.tools.get_wyckoff_sites(structure, map_occupations=None, symprec=1e-05)[source]¶
Returns the Wyckoff symbols of the input structure. The Wyckoff sites are of general interest for symmetry analysis but can be especially useful when setting up, e.g., a
SiteOccupancyObserver
. The Wyckoff labels can be conveniently attached as an array to the structure object as demonstrated in the examples section below.By default the occupation of the sites is part of the symmetry analysis. If a chemically disordered structure is provided this will usually reduce the symmetry substantially. If one is interested in the symmetry of the underlying structure one can control how occupations are handled. To this end, one can provide the
map_occupations
keyword argument. The latter must be a list, each entry of which is a list of species that should be treated as indistinguishable. As a shortcut, if all species should be treated as indistinguishable one can provide an empty list. Examples that illustrate the usage of the keyword are given below.- Parameters:
structure (
Atoms
) – Input structure. Note that the occupation of the sites is included in the symmetry analysis.map_occupations (
Optional
[List
[List
[str
]]]) – Each sublist in this list specifies a group of chemical species that shall be treated as indistinguishable for the purpose of the symmetry analysis.symprec (
float
) – Tolerance imposed when analyzing the symmetry using spglib.
- Return type:
List
[str
]
Examples
Wyckoff sites of a hexagonal-close packed structure:
>>> from ase.build import bulk >>> structure = bulk('Ti') >>> wyckoff_sites = get_wyckoff_sites(structure) >>> print(wyckoff_sites) ['2d', '2d']
The Wyckoff labels can also be attached as an array to the structure, in which case the information is also included when storing the Atoms object:
>>> from ase.io import write >>> structure.new_array('wyckoff_sites', wyckoff_sites, str) >>> write('structure.xyz', structure)
The function can also be applied to supercells:
>>> structure = bulk('GaAs', crystalstructure='zincblende', a=3.0).repeat(2) >>> wyckoff_sites = get_wyckoff_sites(structure) >>> print(wyckoff_sites) ['4a', '4c', '4a', '4c', '4a', '4c', '4a', '4c', '4a', '4c', '4a', '4c', '4a', '4c', '4a', '4c']
Now assume that one is given a supercell of a (Ga,Al)As alloy. Applying the function directly yields much lower symmetry since the symmetry of the original structure is broken:
>>> structure.set_chemical_symbols( ... ['Ga', 'As', 'Al', 'As', 'Ga', 'As', 'Al', 'As', ... 'Ga', 'As', 'Ga', 'As', 'Al', 'As', 'Ga', 'As']) >>> print(get_wyckoff_sites(structure)) ['8g', '8i', '4e', '8i', '8g', '8i', '2c', '8i', '2d', '8i', '8g', '8i', '4e', '8i', '8g', '8i']
Since Ga and Al occupy the same sublattice, they should, however, be treated as indistinguishable for the purpose of the symmetry analysis, which can be achieved via the
map_occupations
keyword:>>> print(get_wyckoff_sites(structure, map_occupations=[['Ga', 'Al'], ['As']])) ['4a', '4c', '4a', '4c', '4a', '4c', '4a', '4c', '4a', '4c', '4a', '4c', '4a', '4c', '4a', '4c']
If occupations are to ignored entirely, one can simply provide an empty list. In the present case, this turns the zincblende lattice into a diamond lattice, on which case there is only one Wyckoff site:
>>> print(get_wyckoff_sites(structure, map_occupations=[])) ['8a', '8a', '8a', '8a', '8a', '8a', '8a', '8a', '8a', '8a', '8a', '8a', '8a', '8a', '8a', '8a']