doped.thermodynamics module

Code for analysing the thermodynamics of defect formation in solids, including calculation of formation energies as functions of Fermi level and chemical potentials, charge transition levels, defect/carrier concentrations etc.

class doped.thermodynamics.DefectThermodynamics(defect_entries: list[DefectEntry] | dict[str, DefectEntry], chempots: dict | None = None, el_refs: dict | None = None, vbm: float | None = None, band_gap: float | None = None, dist_tol: float = 1.5, check_compatibility: bool = True)[source]

Bases: MSONable

Class for analysing the calculated thermodynamics of defects in solids. Similar to a pymatgen PhaseDiagram object, allowing the analysis of formation energies as functions of Fermi level and chemical potentials (i.e. defect formation energy / transition level diagrams), charge transition levels, defect/carrier concentrations etc.

This class is able to get:

  1. stability of charge states for a given defect,

  2. list of all formation energies,

  3. transition levels in the gap,

  4. used as input to doped plotting/analysis functions

Create a DefectThermodynamics object, which can be used to analyse the calculated thermodynamics of defects in solids (formation energies, transition levels, concentrations etc).

Usually initialised using DefectsParser.get_defect_thermodynamics(), but can also be initialised with a list or dict of DefectEntry objects (e.g. from DefectsParser.defect_dict).

Note that the DefectEntry.name attributes are used to label the defects in plots.

Parameters:

  • {str (defect_entries ([DefectEntry] or) – DefectEntry}): A list or dict of DefectEntry objects. Note that DefectEntry.name attributes are used for grouping and plotting purposes! These should be in the format “{defect_name}_{optional_site_info}_{charge_state}”. If the DefectEntry.name attribute is not defined or does not end with the charge state, then the entry will be renamed with the doped default name.

  • chempots (dict) –

    Dictionary of chemical potentials to use for calculating the defect formation energies. This can have the form of {"limits": [{'limit': [chempot_dict]}]} (the format generated by doped's chemical potential parsing functions (see tutorials)) which allows easy analysis over a range of chemical potentials - where limit(s) (chemical potential limit(s)) to analyse/plot can later be chosen using the limits argument.

    Alternatively this can be a dictionary of chemical potentials for a single limit (limit), in the format: {element symbol: chemical potential}. If manually specifying chemical potentials this way, you can set the el_refs option with the DFT reference energies of the elemental phases in order to show the formal (relative) chemical potentials above the formation energy plot, in which case it is the formal chemical potentials (i.e. relative to the elemental references) that should be given here, otherwise the absolute (DFT) chemical potentials should be given.

    If None (default), sets all chemical potentials to zero. Chemical potentials can also be supplied later in each analysis function. (Default: None)

  • el_refs (dict) – Dictionary of elemental reference energies for the chemical potentials in the format: {element symbol: reference energy} (to determine the formal chemical potentials, when chempots has been manually specified as {element symbol: chemical potential}). Unnecessary if chempots is provided in format generated by doped (see tutorials). (Default: None)

  • vbm (float) – VBM eigenvalue in the bulk supercell, to use as Fermi level reference point for analysis. If None (default), will use “vbm” from the calculation_metadata dict attributes of the DefectEntry objects in defect_entries.

  • band_gap (float) – Band gap of the host, to use for analysis. If None (default), will use “gap” from the calculation_metadata dict attributes of the DefectEntry objects in defect_entries.

  • dist_tol (float) – Threshold for the closest distance (in Å) between equivalent defect sites, for different species of the same defect type, to be grouped together (for plotting and transition level analysis). If the minimum distance between equivalent defect sites is less than dist_tol, then they will be grouped together, otherwise treated as separate defects. (Default: 1.5)

  • check_compatibility (bool) – Whether to check the compatibility of the bulk entry for each defect entry (i.e. that all reference bulk energies are the same). (Default: True)

Key Attributes:
defect_entries (list):

List of DefectEntry objects included in the DefectThermodynamics set.

chempots (dict):

Dictionary of chemical potentials to use for calculating the defect formation energies (and hence concentrations etc), in the doped format.

el_refs (dict):

Dictionary of elemental reference energies for the chemical potentials.

vbm (float):

VBM eigenvalue in the bulk supercell, to use as Fermi level reference.

band_gap (float):

Band gap of the host, to use for analysis.

dist_tol (float):

Threshold for the closest distance (in Å) between equivalent defect sites, for different species of the same defect type, to be grouped together (for plotting and transition level analysis).

transition_levels (dict):

Dictionary of charge transition levels for each defect entry. (e.g. {defect_name: {charge: transition_level}}).

check_compatibility (bool):

Whether to check the compatibility of the bulk entry for each defect entry (i.e. that all reference bulk energies are the same).

bulk_formula (str):

The reduced formula of the bulk structure (e.g. “CdTe”).

add_entries(defect_entries: list[DefectEntry] | dict[str, DefectEntry], check_compatibility: bool = True)[source]

Add additional defect entries to the DefectThermodynamics object.

Parameters:

  • {str (defect_entries ([DefectEntry] or) – DefectEntry}): A list or dict of DefectEntry objects, to add to the DefectThermodynamics.defect_entries list. Note that DefectEntry.name attributes are used for grouping and plotting purposes! These should be in the format “{defect_name}_{optional_site_info}_{charge_state}”. If the DefectEntry.name attribute is not defined or does not end with the charge state, then the entry will be renamed with the doped default name.

  • check_compatibility (bool) – Whether to check the compatibility of the bulk entry for each defect entry (i.e. that all reference bulk energies are the same). (Default: True)

property all_stable_entries

List all stable entries (defect + charge) in the DefectThermodynamics set.

property all_unstable_entries

List all unstable entries (defect + charge) in the DefectThermodynamics set.

as_dict()[source]
Returns:

JSON-serializable dict representation of DefectThermodynamics.

property chempots

Get the chemical potentials dictionary used for calculating the defect formation energies.

chempots is a dictionary of chemical potentials to use for calculating the defect formation energies, in the form of: {"limits": [{'limit': [chempot_dict]}]} (the format generated by doped's chemical potential parsing functions (see tutorials)) which allows easy analysis over a range of chemical potentials - where limit(s) (chemical potential limit(s)) to analyse/plot can later be chosen using the limits argument.

property defect_entries

Get the list of parsed DefectEntry objects in the DefectThermodynamics analysis object.

property defect_names

List of names of defects in the DefectThermodynamics set.

property dist_tol

Get the distance tolerance (in Å) used for grouping (equivalent) defects together (for plotting and transition level analysis).

DefectEntrys of the same type and with a closest distance between equivalent defect sites less than dist_tol (1.5 Å by default) are grouped together. If a DefectEntry’s site has a closest distance less than dist_tol to multiple sets of equivalent sites, then it is matched to the one with the lowest minimum distance.

property el_refs

Get the elemental reference energies for the chemical potentials.

This is in the form of a dictionary of elemental reference energies for the chemical potentials, in the format: {element symbol: reference energy}.

classmethod from_dict(d)[source]

Reconstitute a DefectThermodynamics object from a dict representation created using as_dict().

Parameters:

d (dict) – dict representation of DefectThermodynamics.

Returns:

DefectThermodynamics object

classmethod from_json(filename: str)[source]

Load a DefectThermodynamics object from a json file.

Parameters:

  • filename (str) – Filename of json file to load DefectThermodynamics

  • from. (object)

Returns:

DefectThermodynamics object

get_dopability_limits(chempots: dict | None = None, limit: str | None = None, el_refs: dict | None = None) → DataFrame[source]

Find the dopability limits of the defect system, searching over all limits (chemical potential limits) in chempots and returning the most p/n-type conditions, or for a given chemical potential limit (if limit is set or chempots corresponds to a single chemical potential limit; i.e. {element symbol: chemical potential}).

The dopability limites are defined by the (first) Fermi level positions at which defect formation energies become negative as the Fermi level moves towards/beyond the band edges, thus determining the maximum possible Fermi level range upon doping for this chemical potential limit.

This is computed by obtaining the formation energy for every stable defect with non-zero charge, and then finding the highest Fermi level position at which a donor defect (positive charge) has zero formation energy (crosses the x-axis) - giving the lower dopability limit, and the lowest Fermi level position at which an acceptor defect (negative charge) has zero formation energy - giving the upper dopability limit.

Parameters:

  • chempots (dict) –

    Dictionary of chemical potentials to use for calculating the defect formation energies (and thus dopability limits). If None (default), will use self.chempots. This can have the form of {"limits": [{'limit': [chempot_dict]}]} (the format generated by doped's chemical potential parsing functions (see tutorials)) and specific limits (chemical potential limits) can then be chosen using limit.

    Alternatively this can be a dictionary of chemical potentials for a single limit (limit), in the format: {element symbol: chemical potential}. If manually specifying chemical potentials this way, you can set the el_refs option with the DFT reference energies of the elemental phases, in which case it is the formal chemical potentials (i.e. relative to the elemental references) that should be given here, otherwise the absolute (DFT) chemical potentials should be given.

  • limit (str) –

    The chemical potential limit for which to calculate formation energies (and thus dopability limits). Can be either:

    • None, in which case we search over all limits (chemical potential limits) in chempots and return the most n/p-type conditions, unless chempots corresponds to a single chemical potential limit.

    • "X-rich"/"X-poor" where X is an element in the system, in which case the most X-rich/poor limit will be used (e.g. “Li-rich”).

    • A key in the (self.)chempots["limits"] dictionary.

    The latter two options can only be used if chempots is in the doped format (see chemical potentials tutorial). (Default: None)

  • el_refs (dict) – Dictionary of elemental reference energies for the chemical potentials in the format: {element symbol: reference energy} (to determine the formal chemical potentials, when chempots has been manually specified as {element symbol: chemical potential}). Unnecessary if chempots is provided/present in format generated by doped (see tutorials). (Default: None)

Returns:

“limit”, “Compensating Defect”, “Dopability Limit” for both p/n-type where ‘Dopability limit’ are the corresponding Fermi level positions in eV, relative to the VBM.

Return type:

pandas DataFrame of dopability limits, with columns

get_doping_windows(chempots: dict | None = None, limit: str | None = None, el_refs: dict | None = None) → DataFrame[source]

Find the doping windows of the defect system, searching over all limits (chemical potential limits) in chempots and returning the most p/n-type conditions, or for a given chemical potential limit (if limit is set or chempots corresponds to a single chemical potential limit; i.e. {element symbol: chemical potential}).

Doping window is defined by the formation energy of the lowest energy compensating defect species at the corresponding band edge (i.e. VBM for hole doping and CBM for electron doping), as these set the upper limit to the formation energy of dopants which could push the Fermi level close to the band edge without being negated by defect charge compensation.

If a dopant has a higher formation energy than the doping window at the band edge, then its charge will be compensated by formation of the corresponding limiting defect species (rather than free carrier populations).

Parameters:

  • chempots (dict) –

    Dictionary of chemical potentials to use for calculating the defect formation energies (and thus doping windows). If None (default), will use self.chempots. This can have the form of {"limits": [{'limit': [chempot_dict]}]} (the format generated by doped's chemical potential parsing functions (see tutorials)) and specific limits (chemical potential limits) can then be chosen using limit.

    Alternatively this can be a dictionary of chemical potentials for a single limit (limit), in the format: {element symbol: chemical potential}. If manually specifying chemical potentials this way, you can set the el_refs option with the DFT reference energies of the elemental phases, in which case it is the formal chemical potentials (i.e. relative to the elemental references) that should be given here, otherwise the absolute (DFT) chemical potentials should be given.

  • limit (str) –

    The chemical potential limit for which to calculate formation energies (and thus doping windows). Can be either:

    • None, in which case we search over all limits (chemical potential limits) in chempots and return the most n/p-type conditions, unless chempots corresponds to a single chemical potential limit.

    • "X-rich"/"X-poor" where X is an element in the system, in which case the most X-rich/poor limit will be used (e.g. “Li-rich”).

    • A key in the (self.)chempots["limits"] dictionary.

    The latter two options can only be used if chempots is in the doped format (see chemical potentials tutorial). (Default: None)

  • el_refs (dict) – Dictionary of elemental reference energies for the chemical potentials in the format: {element symbol: reference energy} (to determine the formal chemical potentials, when chempots has been manually specified as {element symbol: chemical potential}). Unnecessary if chempots is provided/present in format generated by doped (see tutorials). (Default: None)

Returns:

“limit”, “Compensating Defect”, “Doping Window” for both p/n-type where ‘Doping Window’ are the corresponding doping windows in eV.

Return type:

pandas DataFrame of doping windows, with columns

get_equilibrium_concentrations(chempots: dict | None = None, limit: str | None = None, el_refs: dict | None = None, fermi_level: float | None = None, temperature: float = 300, per_charge: bool = True, per_site: bool = False, skip_formatting: bool = False, lean: bool = False) → DataFrame[source]

Compute the equilibrium concentrations (in cm^-3) for all DefectEntrys in the DefectThermodynamics object, at a given chemical potential limit, fermi_level and temperature, assuming the dilute limit approximation.

Note that these are the equilibrium defect concentrations! DefectThermodynamics.get_quenched_fermi_level_and_concentrations() can instead be used to calculate the Fermi level and defect concentrations for a material grown/annealed at higher temperatures and then cooled (quenched) to room/operating temperature (where defect concentrations are assumed to remain fixed) - this is known as the frozen defect approach and is typically the most valid approximation (see its docstring for more information).

The degeneracy/multiplicity factor “g” is an important parameter in the defect concentration equation (see discussion in https://doi.org/10.1039/D2FD00043A and https://doi.org/10.1039/D3CS00432E), affecting the final concentration by up to ~2 orders of magnitude. This factor is taken from the product of the defect_entry.defect.multiplicity and defect_entry.degeneracy_factors attributes.

Parameters:

  • chempots (dict) –

    Dictionary of chemical potentials to use for calculating the defect formation energies (and thus concentrations). If None (default), will use self.chempots (= 0 for all chemical potentials by default). This can have the form of {"limits": [{'limit': [chempot_dict]}]} (the format generated by doped's chemical potential parsing functions (see tutorials)) and specific limits (chemical potential limits) can then be chosen using limit.

    Alternatively this can be a dictionary of chemical potentials for a single limit (limit), in the format: {element symbol: chemical potential}. If manually specifying chemical potentials this way, you can set the el_refs option with the DFT reference energies of the elemental phases, in which case it is the formal chemical potentials (i.e. relative to the elemental references) that should be given here, otherwise the absolute (DFT) chemical potentials should be given.

  • limit (str) –

    The chemical potential limit for which to obtain the equilibrium concentrations. Can be either:

    • None, if chempots corresponds to a single chemical potential limit - otherwise will use the first chemical potential limit in the chempots dict.

    • "X-rich"/"X-poor" where X is an element in the system, in which case the most X-rich/poor limit will be used (e.g. “Li-rich”).

    • A key in the (self.)chempots["limits"] dictionary.

    The latter two options can only be used if chempots is in the doped format (see chemical potentials tutorial). (Default: None)

  • el_refs (dict) – Dictionary of elemental reference energies for the chemical potentials in the format: {element symbol: reference energy} (to determine the formal chemical potentials, when chempots has been manually specified as {element symbol: chemical potential}). Unnecessary if chempots is provided/present in format generated by doped (see tutorials). (Default: None)

  • fermi_level (float) – Value corresponding to the electron chemical potential, referenced to the VBM. If None (default), set to the mid-gap Fermi level (E_g/2).

  • temperature (float) – Temperature in Kelvin at which to calculate the equilibrium concentrations. Default is 300 K.

  • per_charge (bool) – Whether to break down the defect concentrations into individual defect charge states (e.g. v_Cd_0, v_Cd_-1, v_Cd_-2 instead of v_Cd). (default: True)

  • per_site (bool) – Whether to return the concentrations as percent concentrations per site, rather than the default of per cm^3. (default: False)

  • skip_formatting (bool) – Whether to skip formatting the defect charge states and concentrations as strings (and keep as ints and floats instead). (default: False)

  • lean (bool) – Whether to return a leaner DataFrame with only the defect name, charge state, and concentration in cm^-3 (assumes skip_formatting=True). Only really intended for internal doped usage, to reduce compute times when calculating defect concentrations repeatedly. (default: False)

Returns:

pandas DataFrame of defect concentrations (and formation energies) for each DefectEntry in the DefectThermodynamics object.

get_equilibrium_fermi_level(bulk_dos_vr: str | Vasprun | FermiDos, chempots: dict | None = None, limit: str | None = None, el_refs: dict | None = None, temperature: float = 300, return_concs: bool = False) → float | tuple[float, float, float][source]

Calculate the self-consistent Fermi level, at a given chemical potential limit and temperature, assuming equilibrium defect concentrations (i.e. under annealing) and the dilute limit approximation, by self-consistently solving for the Fermi level which yields charge neutrality.

Note that this assumes equilibrium defect concentrations! DefectThermodynamics.get_quenched_fermi_level_and_concentrations() can instead be used to calculate the Fermi level and defect concentrations for a material grown/annealed at higher temperatures and then cooled (quenched) to room/operating temperature (where defect concentrations are assumed to remain fixed) - this is known as the frozen defect approach and is typically the most valid approximation (see its docstring for more information).

Note that the bulk DOS calculation should be well-converged with respect to k-points for accurate Fermi level predictions! Note: For dense k-point DOS calculations, loading the vasprun.xml(.gz) file can be the most time-consuming part of this function, so if you are running this function multiple times, it is faster to initialise the FermiDos object yourself (and set nelecs appropriately - important!), and set bulk_dos_vr to this object instead.

The degeneracy/multiplicity factor “g” is an important parameter in the defect concentration equation and thus Fermi level calculation (see discussion in https://doi.org/10.1039/D2FD00043A and https://doi.org/10.1039/D3CS00432E), affecting the final concentration by up to ~2 orders of magnitude. This factor is taken from the product of the defect_entry.defect.multiplicity and defect_entry.degeneracy_factors attributes.

Parameters:

  • bulk_dos_vr (str or Vasprun or FermiDos) – Path to the vasprun.xml(.gz) output of a bulk electronic density of states (DOS) calculation, or the corresponding pymatgen Vasprun object. Usually this is a static calculation with the primitive cell of the bulk, with relatively dense k-point sampling (especially for materials with disperse band edges) to ensure an accurately-converged DOS and thus Fermi level. ISMEAR = -5 (tetrahedron smearing) is usually recommended for best convergence. Consistent functional settings should be used for the bulk DOS and defect supercell calculations. Alternatively, a pymatgen FermiDos object can be supplied directly (e.g. in case you are using a DFT code other than VASP).

  • chempots (dict) –

    Dictionary of chemical potentials to use for calculating the defect formation energies (and thus concentrations and Fermi level). If None (default), will use self.chempots (= 0 for all chemical potentials by default). This can have the form of {"limits": [{'limit': [chempot_dict]}]} (the format generated by doped's chemical potential parsing functions (see tutorials)) and specific limits (chemical potential limits) can then be chosen using limit.

    Alternatively this can be a dictionary of chemical potentials for a single limit (limit), in the format: {element symbol: chemical potential}. If manually specifying chemical potentials this way, you can set the el_refs option with the DFT reference energies of the elemental phases, in which case it is the formal chemical potentials (i.e. relative to the elemental references) that should be given here, otherwise the absolute (DFT) chemical potentials should be given.

  • limit (str) –

    The chemical potential limit for which to determine the equilibrium Fermi level. Can be either:

    • None, if chempots corresponds to a single chemical potential limit - otherwise will use the first chemical potential limit in the chempots dict.

    • "X-rich"/"X-poor" where X is an element in the system, in which case the most X-rich/poor limit will be used (e.g. “Li-rich”).

    • A key in the (self.)chempots["limits"] dictionary.

    The latter two options can only be used if chempots is in the doped format (see chemical potentials tutorial). (Default: None)

  • el_refs (dict) – Dictionary of elemental reference energies for the chemical potentials in the format: {element symbol: reference energy} (to determine the formal chemical potentials, when chempots has been manually specified as {element symbol: chemical potential}). Unnecessary if chempots is provided/present in format generated by doped (see tutorials). (Default: None)

  • temperature (float) – Temperature in Kelvin at which to calculate the equilibrium Fermi level. Default is 300 K.

  • return_concs (bool) – Whether to return the corresponding electron and hole concentrations (in cm^-3) as well as the Fermi level. (default: False)

Returns:

Self consistent Fermi level (in eV from the VBM), and the corresponding electron and hole concentrations (in cm^-3) if return_concs=True.

get_formation_energies(chempots: dict | None = None, limit: str | None = None, el_refs: dict | None = None, fermi_level: float | None = None, skip_formatting: bool = False) → DataFrame | list[DataFrame][source]

Generates defect formation energy tables (DataFrames) for either a single chemical potential limit (i.e. phase diagram limit) or each limit in the phase diagram (chempots dict), depending on input limit and chempots.

Table Key: (all energies in eV):

  • ‘Defect’:

    Defect name (without charge)

  • ‘q’:

    Defect charge state.

  • ‘ΔEʳᵃʷ’:

    Raw DFT energy difference between defect and host supercell (E_defect - E_host).

  • ‘qE_VBM’:

    Defect charge times the VBM eigenvalue (to reference the Fermi level to the VBM)

  • ‘qE_F’:

    Defect charge times the Fermi level (referenced to the VBM if qE_VBM is not 0 (if “vbm” in DefectEntry.calculation_metadata)

  • ‘Σμ_ref’:

    Sum of reference energies of the elemental phases in the chemical potentials sum.

  • ‘Σμ_formal’:

    Sum of formal atomic chemical potential terms (Σμ_DFT = Σμ_ref + Σμ_formal).

  • ‘E_corr’:

    Finite-size supercell charge correction.

  • ‘ΔEᶠᵒʳᵐ’:

    Defect formation energy, with the specified chemical potentials and Fermi level. Equals the sum of all other terms.

Parameters:

  • chempots (dict) –

    Dictionary of chemical potentials to use for calculating the defect formation energies. If None (default), will use self.chempots. This can have the form of {"limits": [{'limit': [chempot_dict]}]} (the format generated by doped's chemical potential parsing functions (see tutorials)) and specific limits (chemical potential limits) can then be chosen using limit.

    Alternatively this can be a dictionary of chemical potentials for a single limit (limit), in the format: {element symbol: chemical potential}. If manually specifying chemical potentials this way, you can set the el_refs option with the DFT reference energies of the elemental phases, in which case it is the formal chemical potentials (i.e. relative to the elemental references) that should be given here, otherwise the absolute (DFT) chemical potentials should be given.

    If None (default), sets all chemical potentials to zero. (Default: None)

  • limit (str) –

    The chemical potential limit for which to tabulate formation energies. Can be either:

    • None, in which case tables are generated for all limits in chempots.

    • ”X-rich”/”X-poor” where X is an element in the system, in which case the most X-rich/poor limit will be used (e.g. “Li-rich”).

    • A key in the (self.)chempots["limits"] dictionary.

    The latter two options can only be used if chempots is in the doped format (see chemical potentials tutorial). (Default: None)

  • el_refs (dict) – Dictionary of elemental reference energies for the chemical potentials in the format: {element symbol: reference energy} (to determine the formal chemical potentials, when chempots has been manually specified as {element symbol: chemical potential}). Unnecessary if chempots is provided/present in format generated by doped (see tutorials). (Default: None)

  • fermi_level (float) – Value corresponding to the electron chemical potential, referenced to the VBM. The VBM value is taken from the calculation_metadata dict attributes of DefectEntrys in self.defect_entries if present, otherwise self.vbm. If None (default), set to the mid-gap Fermi level (E_g/2).

  • skip_formatting (bool) – Whether to skip formatting the defect charge states as strings (and keep as ints and floats instead). (default: False)

Returns:

pandas DataFrame or list of DataFrames

get_formation_energy(defect_entry: str | DefectEntry, chempots: dict | None = None, limit: str | None = None, el_refs: dict | None = None, fermi_level: float | None = None) → float[source]

Compute the formation energy for a DefectEntry at a given chemical potential limit and fermi_level. defect_entry can be a string of the defect name, of the DefectEntry object itself.

Parameters:

  • defect_entry (str or DefectEntry) – Either a string of the defect entry name (in DefectThermodynamics.defect_entries), or a DefectEntry object. If the defect name is given without the charge state, then the formation energy of the lowest energy (stable) charge state at the chosen Fermi level will be given.

  • chempots (dict) –

    Dictionary of chemical potentials to use for calculating the defect formation energy. If None (default), will use self.chempots. This can have the form of {"limits": [{'limit': [chempot_dict]}]} (the format generated by doped's chemical potential parsing functions (see tutorials)) and specific limits (chemical potential limits) can then be chosen using limit.

    Alternatively this can be a dictionary of chemical potentials for a single limit (limit), in the format: {element symbol: chemical potential}. If manually specifying chemical potentials this way, you can set the el_refs option with the DFT reference energies of the elemental phases, in which case it is the formal chemical potentials (i.e. relative to the elemental references) that should be given here, otherwise the absolute (DFT) chemical potentials should be given.

    If None (default), sets all chemical potentials to zero. (Default: None)

  • limit (str) –

    The chemical potential limit for which to calculate the formation energy. Can be either:

    • None (default), if chempots corresponds to a single chemical potential limit - otherwise will use the first chemical potential limit in the chempots dict.

    • ”X-rich”/”X-poor” where X is an element in the system, in which case the most X-rich/poor limit will be used (e.g. “Li-rich”).

    • A key in the (self.)chempots["limits"] dictionary.

    The latter two options can only be used if chempots is in the doped format (see chemical potentials tutorial). (Default: None)

  • el_refs (dict) – Dictionary of elemental reference energies for the chemical potentials in the format: {element symbol: reference energy} (to determine the formal chemical potentials, when chempots has been manually specified as {element symbol: chemical potential}). Unnecessary if chempots is provided/present in format generated by doped (see tutorials). (Default: None)

  • fermi_level (float) – Value corresponding to the electron chemical potential, referenced to the VBM. The VBM value is taken from DefectEntry.calculation_metadata if present, otherwise from self.vbm. If None (default), set to the mid-gap Fermi level (E_g/2).

Returns:

Formation energy value (float)

get_quenched_fermi_level_and_concentrations(bulk_dos_vr: str | Vasprun | FermiDos, chempots: dict | None = None, limit: str | None = None, el_refs: dict | None = None, annealing_temperature: float = 1000, quenched_temperature: float = 300, delta_gap: float = 0, per_charge: bool = True, per_site: bool = False, skip_formatting: bool = False, **kwargs) → tuple[float, float, float, DataFrame][source]

Calculate the self-consistent Fermi level and corresponding carrier/defect calculations, for a given chemical potential limit, annealing temperature and quenched/operating temperature, using the frozen defect and dilute limit approximations under the constraint of charge neutrality.

According to the ‘frozen defect’ approximation, we typically expect defect concentrations to reach equilibrium during annealing/crystal growth (at elevated temperatures), but not upon quenching (i.e. at room/operating temperature) where we expect kinetic inhibition of defect annhiliation and hence non-equilibrium defect concentrations / Fermi level. Typically this is approximated by computing the equilibrium Fermi level and defect concentrations at the annealing temperature, and then assuming the total concentration of each defect is fixed to this value, but that the relative populations of defect charge states (and the Fermi level) can re-equilibrate at the lower (room) temperature. See discussion in https://doi.org/10.1039/D3CS00432E (brief), https://doi.org/10.1016/j.cpc.2019.06.017 (detailed) and doped/py-sc-fermi tutorials for more information. In certain cases (such as Li-ion battery materials or extremely slow charge capture/emission), these approximations may have to be adjusted such that some defects/charge states are considered fixed and some are allowed to re-equilibrate (e.g. highly mobile Li vacancies/interstitials). Modelling these specific cases is demonstrated in: py-sc-fermi.readthedocs.io/en/latest/tutorial.html#3.-Applying-concentration-constraints

This function works by calculating the self-consistent Fermi level and total concentration of each defect at the annealing temperature, then fixing the total concentrations to these values and re-calculating the self-consistent (constrained equilibrium) Fermi level and relative charge state concentrations under this constraint at the quenched/operating temperature.

Note that the bulk DOS calculation should be well-converged with respect to k-points for accurate Fermi level predictions! Note: For dense k-point DOS calculations, loading the vasprun.xml(.gz) file can be the most time-consuming part of this function, so if you are running this function multiple times, it is faster to initialise the FermiDos object yourself (and set nelecs appropriately - important!), and set bulk_dos_vr to this object instead.

The degeneracy/multiplicity factor “g” is an important parameter in the defect concentration equation and thus Fermi level calculation (see discussion in https://doi.org/10.1039/D2FD00043A and https://doi.org/10.1039/D3CS00432E), affecting the final concentration by up to 2 orders of magnitude. This factor is taken from the product of the defect_entry.defect.multiplicity and defect_entry.degeneracy_factors attributes.

If you use this code in your work, please also cite: Squires et al., (2023). Journal of Open Source Software, 8(82), 4962 https://doi.org/10.21105/joss.04962

Parameters:

  • bulk_dos_vr (str or Vasprun or FermiDos) – Path to the vasprun.xml(.gz) output of a bulk electronic density of states (DOS) calculation, or the corresponding pymatgen Vasprun object. Usually this is a static calculation with the primitive cell of the bulk, with relatively dense k-point sampling (especially for materials with disperse band edges) to ensure an accurately-converged DOS and thus Fermi level. ISMEAR = -5 (tetrahedron smearing) is usually recommended for best convergence. Consistent functional settings should be used for the bulk DOS and defect supercell calculations. Alternatively, a pymatgen FermiDos object can be supplied directly (e.g. in case you are using a DFT code other than VASP).

  • chempots (dict) –

    Dictionary of chemical potentials to use for calculating the defect formation energies (and thus concentrations and Fermi level). If None (default), will use self.chempots (= 0 for all chemical potentials by default). This can have the form of {"limits": [{'limit': [chempot_dict]}]} (the format generated by doped's chemical potential parsing functions (see tutorials)) and specific limits (chemical potential limits) can then be chosen using limit.

    Alternatively this can be a dictionary of chemical potentials for a single limit (limit), in the format: {element symbol: chemical potential}. If manually specifying chemical potentials this way, you can set the el_refs option with the DFT reference energies of the elemental phases, in which case it is the formal chemical potentials (i.e. relative to the elemental references) that should be given here, otherwise the absolute (DFT) chemical potentials should be given.

  • limit (str) –

    The chemical potential limit for which to determine the Fermi level and concentrations. Can be either:

    • None, if chempots corresponds to a single chemical potential limit - otherwise will use the first chemical potential limit in the chempots dict.

    • "X-rich"/"X-poor" where X is an element in the system, in which case the most X-rich/poor limit will be used (e.g. “Li-rich”).

    • A key in the (self.)chempots["limits"] dictionary.

    The latter two options can only be used if chempots is in the doped format (see chemical potentials tutorial). (Default: None)

  • el_refs (dict) – Dictionary of elemental reference energies for the chemical potentials in the format: {element symbol: reference energy} (to determine the formal chemical potentials, when chempots has been manually specified as {element symbol: chemical potential}). Unnecessary if chempots is provided/present in format generated by doped (see tutorials). (Default: None)

  • annealing_temperature (float) – Temperature in Kelvin at which to calculate the high temperature (fixed) total defect concentrations, which should correspond to the highest temperature during annealing/synthesis of the material (at which we assume equilibrium defect concentrations) within the frozen defect approach. Default is 1000 K.

  • quenched_temperature (float) – Temperature in Kelvin at which to calculate the self-consistent (constrained equilibrium) Fermi level, given the fixed total concentrations, which should correspond to operating temperature of the material (typically room temperature). Default is 300 K.

  • delta_gap (float) – Change in band gap (in eV) of the host material at the annealing temperature (e.g. due to thermal renormalisation), relative to the original band gap of the FermiDos object. If set, applies a scissor correction to fermi_dos which renormalises the band gap symmetrically about the VBM and CBM (i.e. assuming equal up/downshifts of the band-edges around their original eigenvalues) while the defect levels remain fixed. (Default: 0)

  • per_charge (bool) – Whether to break down the defect concentrations into individual defect charge states (e.g. v_Cd_0, v_Cd_-1, v_Cd_-2 instead of v_Cd). (default: True)

  • per_site (bool) – Whether to return the concentrations as percent concentrations per site, rather than the default of per cm^3. (default: False)

  • skip_formatting (bool) – Whether to skip formatting the defect charge states and concentrations as strings (and keep as ints and floats instead). (default: False)

  • **kwargs – Additional keyword arguments to pass to scissor_dos (if delta_gap is not 0).

Returns:

Predicted quenched Fermi level (in eV from the VBM), the corresponding electron and hole concentrations (in cm^-3) and a dataframe of the quenched defect concentrations (in cm^-3). (fermi_level, e_conc, h_conc, conc_df)

get_symmetries_and_degeneracies(skip_formatting: bool = False, symprec: float | None = None) → DataFrame[source]

Generates a table of the bulk-site & relaxed defect point group symmetries, spin/orientational/total degeneracies and (bulk-)site multiplicities for each defect in the DefectThermodynamics object.

Table Key:

  • ‘Defect’: Defect name (without charge)

  • ‘q’: Defect charge state.

  • ‘Site_Symm’: Point group symmetry of the defect site in the bulk cell.

  • ‘Defect_Symm’: Point group symmetry of the relaxed defect.

  • ‘g_Orient’: Orientational degeneracy of the defect.

  • ‘g_Spin’: Spin degeneracy of the defect.

  • ‘g_Total’: Total degeneracy of the defect.

  • ‘Mult’: Multiplicity of the defect site in the bulk cell, per primitive unit cell.

For interstitials, the bulk site symmetry corresponds to the point symmetry of the interstitial site with no relaxation of the host structure, while for vacancies/substitutions it is simply the symmetry of their corresponding bulk site. This corresponds to the point symmetry of DefectEntry.defect, or calculation_metadata["bulk_site"]/["unrelaxed_defect_structure"].

Point group symmetries are taken from the calculation_metadata (“relaxed point symmetry” and “bulk site symmetry”) if present (should be, if parsed with doped and defect supercell doesn’t break host periodicity), otherwise are attempted to be recalculated.

Note: doped tries to use the defect_entry.defect_supercell to determine the relaxed site symmetry. However, it should be noted that this is not guaranteed to work in all cases; namely for non-diagonal supercell expansions, or sometimes for non-scalar supercell expansion matrices (e.g. a 2x1x2 expansion)(particularly with high-symmetry materials) which can mess up the periodicity of the cell. doped tries to automatically check if this is the case, and will warn you if so.

This can also be checked by using this function on your doped generated defects:

from doped.generation import get_defect_name_from_entry
for defect_name, defect_entry in defect_gen.items():
    print(defect_name, get_defect_name_from_entry(defect_entry, relaxed=False),
          get_defect_name_from_entry(defect_entry), "\n")

And if the point symmetries match in each case, then doped should be able to correctly determine the final relaxed defect symmetry (and orientational degeneracy) - otherwise periodicity-breaking prevents this.

If periodicity-breaking prevents auto-symmetry determination, you can manually determine the relaxed defect and bulk-site point symmetries, and/or orientational degeneracy, from visualising the structures (e.g. using VESTA)(can use get_orientational_degeneracy to obtain the corresponding orientational degeneracy factor for given defect/bulk-site point symmetries) and setting the corresponding values in the calculation_metadata['relaxed point symmetry']/['bulk site symmetry'] and/or degeneracy_factors['orientational degeneracy'] attributes. Note that the bulk-site point symmetry corresponds to that of DefectEntry.defect, or equivalently calculation_metadata["bulk_site"]/["unrelaxed_defect_structure"], which for vacancies/substitutions is the symmetry of the corresponding bulk site, while for interstitials it is the point symmetry of the final relaxed interstitial site when placed in the (unrelaxed) bulk structure. The degeneracy factor is used in the calculation of defect/carrier concentrations and Fermi level behaviour (see e.g. https://doi.org/10.1039/D2FD00043A & https://doi.org/10.1039/D3CS00432E).

Parameters:

  • skip_formatting (bool) – Whether to skip formatting the defect charge states as strings (and keep as ints and floats instead). (default: False)

  • symprec (float) – Symmetry tolerance for spglib to use when determining relaxed defect point symmetries and thus orientational degeneracies. Default is 0.1 which matches that used by the Materials Project and is larger than the pymatgen default of 0.01 (which is used by doped for unrelaxed/bulk structures) to account for residual structural noise in relaxed defect supercells. You may want to adjust for your system (e.g. if there are very slight octahedral distortions etc.). If symprec is set, then the point symmetries and corresponding orientational degeneracy will be re-parsed/computed even if already present in the DefectEntry object calculation_metadata.

Returns:

pandas DataFrame

get_transition_levels(all: bool = False, format_charges: bool = True) → DataFrame[source]

Return a DataFrame of the charge transition levels for the defects in the DefectThermodynamics object (stored in the transition_level_map attribute).

By default, only returns the thermodynamic ground-state transition levels (i.e. those visible on the defect formation energy diagram), not including metastable defect states (which can be important for recombination, migration, degeneracy/concentrations etc, see e.g. https://doi.org/10.1039/D2FD00043A & https://doi.org/10.1039/D3CS00432E). e.g. negative-U defects will show the 2-electron transition level (N+1/N-1) rather than (N+1/N) and (N/N-1). If instead all single-electron transition levels are desired, set all = True.

Returns a DataFrame with columns:

  • “Defect”: Defect name

  • “Charges”: Defect charge states which make up the transition level

    (as a string if format_charges=True, otherwise as a list of integers)

  • “eV from VBM”: Transition level position in eV from the VBM

  • “In Band Gap?”: Whether the transition level is within the host band gap

  • “N(Metastable)”: Number of metastable states involved in the transition level

    (0, 1 or 2). Only included if all = True.

Parameters:

  • all (bool) – Whether to print all single-electron transition levels (i.e. including metastable defect states), or just the thermodynamic ground-state transition levels (default).

  • format_charges (bool) – Whether to format the transition level charge states as strings (e.g. “ε(+1/+2)”) or keep in list format (e.g. [1,2]). (Default: True)

plot(chempots: dict | None = None, limit: str | None = None, el_refs: dict | None = None, chempot_table: bool = True, all_entries: bool | str = False, style_file: str | None = None, xlim: tuple | None = None, ylim: tuple | None = None, fermi_level: float | None = None, colormap: str | Colormap | None = None, auto_labels: bool = False, filename: str | None = None) → Figure | list[Figure][source]

Produce a defect formation energy vs Fermi level plot (a.k.a. a defect formation energy / transition level diagram). Returns the Matplotlib Figure object to allow further plot customisation.

Parameters:

  • chempots (dict) –

    Dictionary of chemical potentials to use for calculating the defect formation energies. If None (default), will use self.chempots. This can have the form of {"limits": [{'limit': [chempot_dict]}]} (the format generated by doped's chemical potential parsing functions (see tutorials)) and specific limits (chemical potential limits) can then be chosen using limit.

    Alternatively this can be a dictionary of chemical potentials for a single limit (limit), in the format: {element symbol: chemical potential}. If manually specifying chemical potentials this way, you can set the el_refs option with the DFT reference energies of the elemental phases in order to show the formal (relative) chemical potentials above the formation energy plot, in which case it is the formal chemical potentials (i.e. relative to the elemental references) that should be given here, otherwise the absolute (DFT) chemical potentials should be given.

    If None (default), sets all chemical potentials to zero. (Default: None)

  • limit (str) –

    The chemical potential limit for which to plot formation energies. Can be either:

    • None, in which case plots are generated for all limits in chempots.

    • ”X-rich”/”X-poor” where X is an element in the system, in which case the most X-rich/poor limit will be used (e.g. “Li-rich”).

    • A key in the (self.)chempots["limits"] dictionary.

    The latter two options can only be used if chempots is in the doped format (see chemical potentials tutorial). (Default: None)

  • el_refs (dict) – Dictionary of elemental reference energies for the chemical potentials in the format: {element symbol: reference energy} (to determine the formal chemical potentials, when chempots has been manually specified as {element symbol: chemical potential}). Unnecessary if chempots is provided/present in format generated by doped (see tutorials). (Default: None)

  • chempot_table (bool) – Whether to print the chemical potential table above the plot. (Default: True)

  • all_entries (bool, str) – Whether to plot the formation energy lines of all defect entries, rather than the default of showing only the equilibrium states at each Fermi level position (traditional). If instead set to “faded”, will plot the equilibrium states in bold, and all unstable states in faded grey (Default: False)

  • style_file (str) – Path to a mplstyle file to use for the plot. If None (default), uses the default doped style (from doped/utils/doped.mplstyle).

  • xlim – Tuple (min,max) giving the range of the x-axis (Fermi level). May want to set manually when including transition level labels, to avoid crossing the axes. Default is to plot from -0.3 to +0.3 eV above the band gap.

  • ylim – Tuple (min,max) giving the range for the y-axis (formation energy). May want to set manually when including transition level labels, to avoid crossing the axes. Default is from 0 to just above the maximum formation energy value in the band gap.

  • fermi_level (float) – If set, plots a dashed vertical line at this Fermi level value, typically used to indicate the equilibrium Fermi level position (e.g. calculated with py-sc-fermi). (Default: None)

  • colormap (str, matplotlib.colors.Colormap) – Colormap to use for the formation energy lines, either as a string (i.e. name from https://matplotlib.org/stable/users/explain/colors/colormaps.html) or a Colormap / ListedColormap object. If None (default), uses Dark2 (if 8 or less lines) or tab20 (if more than 8 lines being plotted).

  • auto_labels (bool) – Whether to automatically label the transition levels with their charge states. If there are many transition levels, this can be quite ugly. (Default: False)

  • filename (str) – Filename to save the plot to. (Default: None (not saved))

Returns:

Matplotlib Figure object, or list of Figure objects if multiple limits chosen.

print_transition_levels(all: bool = False)[source]

Iteratively prints the charge transition levels for the defects in the DefectThermodynamics object (stored in the transition_level_map attribute).

By default, only returns the thermodynamic ground-state transition levels (i.e. those visible on the defect formation energy diagram), not including metastable defect states (which can be important for recombination, migration, degeneracy/concentrations etc, see e.g. https://doi.org/10.1039/D2FD00043A & https://doi.org/10.1039/D3CS00432E). e.g. negative-U defects will show the 2-electron transition level (N+1/N-1) rather than (N+1/N) and (N/N-1). If instead all single-electron transition levels are desired, set all = True.

Parameters:

all (bool) – Whether to print all single-electron transition levels (i.e. including metastable defect states), or just the thermodynamic ground-state transition levels (default).

to_json(filename: str | None = None)[source]

Save the DefectThermodynamics object as a json file, which can be reloaded with the DefectThermodynamics.from_json() class method.

Parameters:

filename (str) – Filename to save json file as. If None, the filename will be set as “{Chemical Formula}_defect_thermodynamics.json” where {Chemical Formula} is the chemical formula of the host material.

doped.thermodynamics.bold_print(string: str) → None[source]

Does what it says on the tin.

Prints the input string in bold.

doped.thermodynamics.get_e_h_concs(fermi_dos: FermiDos, fermi_level: float, temperature: float) → tuple[float, float][source]

Get the corresponding electron and hole concentrations (in cm^-3) for a given Fermi level (in eV) and temperature (in K), for a FermiDos object.

Note that the Fermi level here is NOT referenced to the VBM! So the Fermi level should be the corresponding eigenvalue within the calculation (or in other words, the Fermi level relative to the VBM plus the VBM eigenvalue).

doped.thermodynamics.get_rich_poor_limit_dict(chempots: dict) → dict[source]

Get a dictionary of {“X-rich”: limit, “X-poor”: limit…} for each element X in the chempots phase diagram.

Parameters:

chempots (dict) – The chempots dict, in the doped format.

doped.thermodynamics.group_defects_by_distance(entry_list: list[DefectEntry], dist_tol: float = 1.5) → dict[str, dict[tuple, list[DefectEntry]]][source]

Given an input list of DefectEntry objects, returns a dictionary of {simple defect name: {(equivalent defect sites): [DefectEntry]}, where ‘simple defect name’ is the nominal defect type (e.g. Te_i for Te_i_Td_Cd2.83_+2), (equivalent defect sites) is a tuple of the equivalent defect sites (in the bulk supercell), and [DefectEntry] is a list of DefectEntry objects with the same simple defect name and a closest distance between symmetry-equivalent sites less than dist_tol (1.5 Å by default).

If a DefectEntry’s site has a closest distance less than dist_tol to multiple sets of equivalent sites, then it is matched to the one with the lowest minimum distance.

Parameters:

  • entry_list ([DefectEntry]) – A list of DefectEntry objects to group together.

  • dist_tol (float) – Threshold for the closest distance (in Å) between equivalent defect sites, for different species of the same defect type, to be grouped together (for plotting and transition level analysis). If the minimum distance between equivalent defect sites is less than dist_tol, then they will be grouped together, otherwise treated as separate defects. (Default: 1.5)

Returns:

{simple defect name: {(equivalent defect sites): [DefectEntry]}

Return type:

dict

doped.thermodynamics.group_defects_by_name(entry_list: list[DefectEntry]) → dict[str, list[DefectEntry]][source]

Given an input list of DefectEntry objects, returns a dictionary of {defect name without charge: [DefectEntry]}, where the values are lists of DefectEntry objects with the same defect name (excluding charge state).

The DefectEntry.name attributes are used to get the defect names. These should be in the format: “{defect_name}_{optional_site_info}_{charge_state}”. If the DefectEntry.name attribute is not defined or does not end with the charge state, then the entry will be renamed with the doped default name for the unrelaxed defect (i.e. using the point symmetry of the defect site in the bulk cell).

For example, v_Cd_C3v_+1, v_Cd_Td_+1 and v_Cd_C3v_+2 will be grouped as {v_Cd_C3v: [v_Cd_C3v_+1, v_Cd_C3v_+2], v_Cd_Td: [v_Cd_Td_+1]}.

Parameters:

entry_list ([DefectEntry]) – A list of DefectEntry objects to group together by defect name (without charge).

Returns:

Dictionary of {defect name without charge: [DefectEntry]}.

Return type:

dict

doped.thermodynamics.scissor_dos(delta_gap: float, dos: Dos | FermiDos, tol: float = 1e-08, verbose: bool = True)[source]

Given an input Dos/FermiDos object, rigidly shifts the valence and conduction bands of the DOS object to give a band gap that is now increased/decreased by delta_gap eV, where this rigid scissor shift is applied symmetrically around the original gap (i.e. the VBM is downshifted by delta_gap/2 and the CBM is upshifted by delta_gap/2).

Note this assumes a non-spin-polarised (i.e. non-magnetic) density of states!

Parameters:

  • delta_gap (float) – The amount by which to increase/decrease the band gap (in eV).

  • dos (Dos/FermiDos) – The input DOS object to scissor.

  • tol (float) – The tolerance to use for determining the VBM and CBM (used in Dos.get_gep(tol=tol)). Default: 1e-8.

  • verbose (bool) – Whether to print information about the original and new band gaps.

Returns:

The scissored DOS object.

Return type:

FermiDos