pymatgen.transformations.advanced_transformations module

class ChargeBalanceTransformation(charge_balance_sp)[source]

Bases: pymatgen.transformations.transformation_abc.AbstractTransformation

This is a transformation that disorders a structure to make it charge balanced, given an oxidation state-decorated structure.

Parameters:charge_balance_sp – specie to add or remove. Currently only removal is supported
apply_transformation(structure)[source]
inverse
is_one_to_many
class DopingTransformation(dopant, ionic_radius_tol=inf, min_length=10, alio_tol=0, codopant=False, max_structures_per_enum=100, allowed_doping_species=None, **kwargs)[source]

Bases: pymatgen.transformations.transformation_abc.AbstractTransformation

A transformation that performs doping of a structure.

Parameters:
  • dopant (Specie-like) – E.g., Al3+. Must have oxidation state.
  • ionic_radius_tol (float) – E.g., Fractional allowable ionic radii mismatch for dopant to fit into a site. Default of inf means that any dopant with the right oxidation state is allowed.
  • min_Length (float) – Min. lattice parameter between periodic images of dopant. Defaults to 10A for now.
  • alio_tol (int) – If this is not 0, attempt will be made to dope sites with oxidation_states +- alio_tol of the dopant. E.g., 1 means that the ions like Ca2+ and Ti4+ are considered as potential doping sites for Al3+.
  • codopant (bool) – If True, doping will be carried out with a codopant to maintain charge neutrality. Otherwise, vacancies will be used.
  • max_structures_per_enum (float) – Maximum number of structures to return per enumeration. Note that there can be more than one candidate doping site, and each site enumeration will return at max max_structures_per_enum structures. Defaults to 100.
  • allowed_doping_species (list) – Species that are allowed to be doping sites. This is an inclusionary list. If specified, any sites which are not
  • **kwargs – Same keyword args as EnumerateStructureTransformation, i.e., min_cell_size, etc.
apply_transformation(structure, return_ranked_list=False)[source]
Parameters:structure (Structure) – Input structure to dope
Returns:Structure, “energy”: float}]
Return type:[{“structure”
inverse
is_one_to_many
class EnumerateStructureTransformation(min_cell_size=1, max_cell_size=1, symm_prec=0.1, refine_structure=False, enum_precision_parameter=0.001, check_ordered_symmetry=True, max_disordered_sites=None, sort_criteria='ewald')[source]

Bases: pymatgen.transformations.transformation_abc.AbstractTransformation

Order a disordered structure using enumlib. For complete orderings, this generally produces fewer structures that the OrderDisorderedStructure transformation, and at a much faster speed.

Parameters:
  • min_cell_size – The minimum cell size wanted. Must be an int. Defaults to 1.
  • max_cell_size – The maximum cell size wanted. Must be an int. Defaults to 1.
  • symm_prec – Tolerance to use for symmetry.
  • refine_structure – This parameter has the same meaning as in enumlib_caller. If you are starting from a structure that has been relaxed via some electronic structure code, it is usually much better to start with symmetry determination and then obtain a refined structure. The refined structure have cell parameters and atomic positions shifted to the expected symmetry positions, which makes it much less sensitive precision issues in enumlib. If you are already starting from an experimental cif, refinment should have already been done and it is not necessary. Defaults to False.
  • enum_precision_parameter (float) – Finite precision parameter for enumlib. Default of 0.001 is usually ok, but you might need to tweak it for certain cells.
  • check_ordered_symmetry (bool) – Whether to check the symmetry of the ordered sites. If the symmetry of the ordered sites is lower, the lowest symmetry ordered sites is included in the enumeration. This is important if the ordered sites break symmetry in a way that is important getting possible structures. But sometimes including ordered sites slows down enumeration to the point that it cannot be completed. Switch to False in those cases. Defaults to True.
  • max_disordered_sites (int) – An alternate parameter to max_cell size. Will sequentially try larger and larger cell sizes until (i) getting a result or (ii) the number of disordered sites in the cell exceeds max_disordered_sites. Must set max_cell_size to None when using this parameter.
  • sort_criteria (str) – Sort by Ewald energy (“ewald”, must have oxidation states and slow) or by number of sites (“nsites”, much faster).
apply_transformation(structure, return_ranked_list=False)[source]

Return either a single ordered structure or a sequence of all ordered structures.

Parameters:
  • structure – Structure to order.
  • return_ranked_list (bool) – Whether or not multiple structures are returned. If return_ranked_list is a number, that number of structures is returned.
Returns:

Depending on returned_ranked list, either a transformed structure or a list of dictionaries, where each dictionary is of the form {“structure” = …. , “other_arguments”}

The list of ordered structures is ranked by ewald energy / atom, if the input structure is an oxidation state decorated structure. Otherwise, it is ranked by number of sites, with smallest number of sites first.

inverse
is_one_to_many
class MagOrderParameterConstraint(order_parameter, species_constraints=None, site_constraint_name=None, site_constraints=None)[source]

Bases: monty.json.MSONable

This class can be used to supply MagOrderingTransformation to just a specific subset of species or sites that satisfy the provided constraints. This can be useful for setting an order parameters for, for example, ferrimagnetic structures which might order on certain motifs, with the global order parameter dependent on how many sites satisfy that motif.

Parameters:(float) (order_parameter) – any number from 0.0 to 1.0,

typically 0.5 (antiferromagnetic) or 1.0 (ferromagnetic) :param species_constraint (list): str or list of strings of Specie symbols that the constraint should apply to :param site_constraint_name (str): name of the site property that the constraint should apply to, e.g. “coordination_no” :param site_constraints (list): list of values of the site property that the constraints should apply to

satisfies_constraint(site)[source]

Checks if a periodic site satisfies the constraint.

class MagOrderingTransformation(mag_species_spin, order_parameter=0.5, energy_model=<pymatgen.analysis.energy_models.SymmetryModel object>, **kwargs)[source]

Bases: pymatgen.transformations.transformation_abc.AbstractTransformation

This transformation takes a structure and returns a list of collinear magnetic orderings. For disordered structures, make an ordered approximation first.

Parameters:mag_species_spin – A mapping of elements/species to their

spin magnitudes, e.g. {“Fe3+”: 5, “Mn3+”: 4} :param order_parameter (float or list): if float, a specifies a global order parameter and can take values from 0.0 to 1.0 (e.g. 0.5 for antiferromagnetic or 1.0 for ferromagnetic), if list has to be a list of :class: pymatgen.transformations.advanced_transformations.MagOrderParameterConstraint to specify more complicated orderings, see documentation for MagOrderParameterConstraint more details on usage :param energy_model: Energy model to rank the returned structures, see :mod: pymatgen.analysis.energy_models for more information (note that this is not necessarily a physical energy). By default, returned structures use SymmetryModel() which ranks structures from most symmetric to least. :param kwargs: Additional kwargs that are passed to EnumerateStructureTransformation such as min_cell_size etc.

apply_transformation(structure, return_ranked_list=False)[source]

Apply MagOrderTransformation to an input structure. :param structure: Any ordered structure. :param return_ranked_list: As in other Transformations. :return:

static determine_min_cell(disordered_structure)[source]

Determine the smallest supercell that is able to enumerate the provided structure with the given order parameter

inverse
is_one_to_many
class MultipleSubstitutionTransformation(sp_to_replace, r_fraction, substitution_dict, charge_balance_species=None, order=True)[source]

Bases: object

Performs multiple substitutions on a structure. For example, can do a fractional replacement of Ge in LiGePS with a list of species, creating one structure for each substitution. Ordering is done using a dummy element so only one ordering must be done per substitution oxidation state. Charge balancing of the structure is optionally performed.

Note

There are no checks to make sure that removal fractions are possible and rounding may occur. Currently charge balancing only works for removal of species.

Performs multiple fractional substitutions on a transmuter.

Parameters:
  • sp_to_replace – species to be replaced
  • r_fraction – fraction of that specie to replace
  • substitution_dict – dictionary of the format {2: [“Mg”, “Ti”, “V”, “As”, “Cr”, “Ta”, “N”, “Nb”], 3: [“Ru”, “Fe”, “Co”, “Ce”, “As”, “Cr”, “Ta”, “N”, “Nb”], 4: [“Ru”, “V”, “Cr”, “Ta”, “N”, “Nb”], 5: [“Ru”, “W”, “Mn”] } The number is the charge used for each of the list of elements (an element can be present in multiple lists)
  • charge_balance_species – If specified, will balance the charge on the structure using that specie.
apply_transformation(structure, return_ranked_list=False)[source]
inverse
is_one_to_many
class SlabTransformation(miller_index, min_slab_size, min_vacuum_size, lll_reduce=False, center_slab=False, primitive=True, max_normal_search=None, shift=0, tol=0.1)[source]

Bases: pymatgen.transformations.transformation_abc.AbstractTransformation

A transformation that creates a slab from a structure.

Parameters:
  • miller_index (3-tuple or list) – miller index of slab
  • min_slab_size (float) – minimum slab size in angstroms
  • min_vacuum_size (float) – minimum size of vacuum
  • lll_reduce (bool) – whether to apply LLL reduction
  • center_slab (bool) – whether to center the slab
  • primitive (bool) – whether to reduce slabs to most primitive cell
  • max_normal_search (int) – maximum index to include in linear combinations of indices to find c lattice vector orthogonal to slab surface
  • shift (float) – shift to get termination
  • tol (float) – tolerance for primitive cell finding
apply_transformation(structure)[source]
inverse
is_one_to_many
class SubstitutionPredictorTransformation(threshold=0.01, **kwargs)[source]

Bases: pymatgen.transformations.transformation_abc.AbstractTransformation

This transformation takes a structure and uses the structure prediction module to find likely site substitutions.

Parameters:
  • threshold – Threshold for substitution.
  • **kwargs – Args for SubstitutionProbability class lambda_table, alpha
apply_transformation(structure, return_ranked_list=False)[source]
inverse
is_one_to_many
class SuperTransformation(transformations, nstructures_per_trans=1)[source]

Bases: pymatgen.transformations.transformation_abc.AbstractTransformation

This is a transformation that is inherently one-to-many. It is constructed from a list of transformations and returns one structure for each transformation. The primary use for this class is extending a transmuter object.

Parameters:
  • transformations ([transformations]) – List of transformations to apply to a structure. One transformation is applied to each output structure.
  • nstructures_per_trans (int) – If the transformations are one-to-many and, nstructures_per_trans structures from each transformation are added to the full list. Defaults to 1, i.e., only best structure.
apply_transformation(structure, return_ranked_list=False)[source]
inverse
is_one_to_many