doped.corrections module
Code to compute finite-size charge corrections for charged defects in periodic systems.
The charge-correction methods implemented are:
Extended FNV (eFNV) / Kumagai correction for isotropic and anistropic systems. This is the recommended default as it works regardless of system (an)isotropy, has a more efficient implementation and tends to be more robust in cases of small supercells. Includes:
anisotropic PC energy
potential alignment by atomic site averaging outside Wigner Seitz radius
Freysoldt (FNV) correction for isotropic systems. Only recommended if eFNV correction fails for some reason. Includes:
point-charge (PC) energy
potential alignment by planar averaging
If you use the corrections implemented in this module, please cite:
Kumagai and Oba, Phys. Rev. B. 89, 195205 (2014) for the eFNV correction or Freysoldt, Neugebauer, and Van de Walle, Phys. Status Solidi B. 248, 1067-1076 (2011) for FNV
Note: Ideally, the “defect site” used for all charge corrections should actually be the centre of the localised charge in the defect supercell. Usually this coincides with the defect site, but not always (e.g. vacancy where the charge localises as a polaron on a neighbouring atom etc.). For sufficiently large supercells this is usually fine, as the defect and centre-of-charge sites are close enough that any resulting quantitative error in the correction is negligible. However, in cases where we have large finite-size correction values, this can be significant. If some efficient methods for determining the centroid of charge difference between bulk/defect LOCPOTs (for FNV) or bulk/defect site potentials (for eFNV) were developed, they should be used here.
- doped.corrections.get_freysoldt_correction(defect_entry, dielectric: float | int | ndarray | list | None = None, defect_locpot: str | Locpot | dict | None = None, bulk_locpot: str | Locpot | dict | None = None, plot: bool = False, filename: str | None = None, axis: int | None = None, verbose: bool = True, style_file: str | None = None, **kwargs) CorrectionResult [source]
Function to compute the isotropic Freysoldt (FNV) correction for the input defect_entry.
This function does not add the correction to
defect_entry.corrections
(but the defect_entry.get_freysoldt_correction method does). If this correction is used, please cite Freysoldt’s original paper; 10.1103/PhysRevLett.102.016402.- Parameters:
defect_entry – DefectEntry object with the following for which to compute the FNV finite-size charge correction.
dielectric (float or int or 3x1 matrix or 3x3 matrix) – Total dielectric constant of the host compound (including both ionic and (high-frequency) electronic contributions). If None, then the dielectric constant is taken from the
defect_entry
calculation_metadata
if available.defect_locpot – Path to the output VASP LOCPOT file from the defect supercell calculation, or the corresponding pymatgen Locpot object, or a dictionary of the planar-averaged potential in the form: {i: Locpot.get_average_along_axis(i) for i in [0,1,2]}. If None, will try to use
defect_locpot_dict
from thedefect_entry
calculation_metadata
if available.bulk_locpot – Path to the output VASP LOCPOT file from the bulk supercell calculation, or the corresponding pymatgen Locpot object, or a dictionary of the planar-averaged potential in the form: {i: Locpot.get_average_along_axis(i) for i in [0,1,2]}. If None, will try to use
bulk_locpot_dict
from thedefect_entry
calculation_metadata
if available.plot (bool) – Whether to plot the FNV electrostatic potential plots (for manually checking the behaviour of the charge correction here).
filename (str) – Filename to save the FNV electrostatic potential plots to. If None, plots are not saved.
axis (int or None) – If int, then the FNV electrostatic potential plot along the specified axis (0, 1, 2 for a, b, c) will be plotted. Note that the output charge correction is still that for all axes. If None, then all three axes are plotted.
verbose (bool) – Whether to print the correction energy (default = True).
style_file (str) – Path to a
.mplstyle
file to use for the plot. IfNone
(default), uses the default doped style (fromdoped/utils/doped.mplstyle
).**kwargs – Additional kwargs to pass to pymatgen.analysis.defects.corrections.freysoldt.get_freysoldt_correction (e.g. energy_cutoff, mad_tol, q_model, step).
- Returns:
CorrectionResults (summary of the corrections applied and metadata), and the matplotlib figure object (or axis object if axis specified) if
plot
is True.
- doped.corrections.get_kumagai_correction(defect_entry, dielectric: float | int | ndarray | list | None = None, defect_region_radius: float | None = None, excluded_indices: list[int] | None = None, defect_outcar: str | Outcar | None = None, bulk_outcar: str | Outcar | None = None, plot: bool = False, filename: str | None = None, verbose: bool = True, style_file: str | None = None, **kwargs) CorrectionResult [source]
Function to compute the Kumagai (eFNV) finite-size charge correction for the input defect_entry. Compatible with both isotropic/cubic and anisotropic systems.
This function does not add the correction to
defect_entry.corrections
(but the defect_entry.get_kumagai_correction method does). If this correction is used, please cite the Kumagai & Oba paper: 10.1103/PhysRevB.89.195205Typically for reasonably well-converged supercell sizes, the default
defect_region_radius
works perfectly well. However, for certain materials at small/intermediate supercell sizes, you may want to adjust this (and/orexcluded_indices
) to ensure the best sampling of the plateau region away from the defect position -doped
should throw a warning in these cases (about the correction error being above the default tolerance (50 meV)). For example, with layered materials, the defect charge is often localised to one layer, so we may want to adjustdefect_region_radius
and/orexcluded_indices
to ensure that only sites in other layers are used for the sampling region (plateau) - see example on doped docs Tips page.- Parameters:
defect_entry (DefectEntry) – DefectEntry object with the following for which to compute the Kumagai finite-size charge correction.
dielectric (float or int or 3x1 matrix or 3x3 matrix) – Total dielectric constant of the host compound (including both ionic and (high-frequency) electronic contributions). If None, then the dielectric constant is taken from the
defect_entry
calculation_metadata
if available.defect_region_radius (float) – Radius of the defect region (in Å). Sites outside the defect region are used for sampling the electrostatic potential far from the defect (to obtain the potential alignment). If None (default), uses the Wigner-Seitz radius of the supercell.
excluded_indices (list) – List of site indices (in the defect supercell) to exclude from the site potential sampling in the correction calculation/plot. If None (default), no sites are excluded.
defect_outcar (str or Outcar) – Path to the output VASP OUTCAR file from the defect supercell calculation, or the corresponding pymatgen Outcar object. If None, will try to use the
defect_site_potentials
from thedefect_entry
calculation_metadata
if available.bulk_outcar (str or Outcar) – Path to the output VASP OUTCAR file from the bulk supercell calculation, or the corresponding pymatgen Outcar object. If None, will try to use the
bulk_site_potentials
from thedefect_entry
calculation_metadata
if available.plot (bool) – Whether to plot the Kumagai site potential plots (for manually checking the behaviour of the charge correction here).
filename (str) – Filename to save the Kumagai site potential plots to. If None, plots are not saved.
verbose (bool) – Whether to print the correction energy (default = True).
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).
**kwargs – Additional kwargs to pass to pydefect.corrections.efnv_correction.ExtendedFnvCorrection (e.g. charge, defect_region_radius, defect_coords).
- Returns:
CorrectionResults (summary of the corrections applied and metadata), and the matplotlib figure object if
plot
is True.
- doped.corrections.plot_FNV(plot_data, title=None, ax=None, style_file=None)[source]
Plots the planar-averaged electrostatic potential against the long range and short range models from the FNV correction method.
Code templated from the original PyCDT and new pymatgen.analysis.defects implementations.
- Parameters:
plot_data (dict) – Dictionary of Freysoldt correction metadata to plot (i.e. defect_entry.corrections_metadata[“plot_data”][axis] where axis is one of [0, 1, 2] specifying which axis to plot along (a, b, c)).
title (str) – Title for the plot. Default is no title.
ax (matplotlib.axes.Axes) – Axes object to plot on. If None, makes new figure.
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).