diff --git a/doped/interface/fermi_solver.py b/doped/interface/fermi_solver.py index a35b46b67..070c940a6 100644 --- a/doped/interface/fermi_solver.py +++ b/doped/interface/fermi_solver.py @@ -1775,6 +1775,392 @@ def _parse_and_check_grid_like_chempots(self, chempots: Optional[dict] = None) - return chempots, el_refs def min_max_X( + self, + target: str, + min_or_max: str, + chempots: dict, + annealing_temperature: Optional[float] = None, + quenched_temperature: float = 300, + temperature: float = 300, + tolerance: float = 0.01, + n_points: int = 10, + effective_dopant_concentration: Optional[float] = None, + fix_charge_states: bool = False, + fixed_defects: Optional[dict[str, float]] = None, + free_defects: Optional[list[str]] = None, + ) -> pd.DataFrame: + r""" + Search for the chemical potentials that minimize or maximize a target + variable, such as electron concentration, within a specified tolerance. + + This function iterates over a grid of chemical potentials and "zooms in" on + the chemical potential that either minimizes or maximizes the target variable. + The process continues until the change in the target variable is less than + the specified tolerance. + + If ``annealing_temperature`` (and ``quenched_temperature``; 300 K by + default) are specified, then the frozen defect approximation is + employed, whereby total defect concentrations are calculated at the + elevated annealing temperature, then fixed at these values (unless + ``free_defects`` or ``fix_charge_states`` are specified) and the Fermi + level and relative charge state populations are recalculated at the + quenched temperature. Otherwise, if only ``temperature`` is specified, + then the Fermi level and defect/carrier concentrations are calculated + assuming thermodynamic equilibrium at that temperature. + + Args: + target (str): + The target variable to minimize or maximize, e.g., "Electrons (cm^-3)". + # TODO: Allow this to match a substring of the column name. + min_or_max (str): + Specify whether to "minimize" ("min") or "maximize" ("max"; default) + the target variable. + chempots (Optional[dict]): + Dictionary of chemical potentials to scan over, in the ``doped`` + format (i.e. ``{"limits": [{'limit': [chempot_dict]}], ...}``) + -- the format generated by ``doped``\'s chemical potential parsing + functions (see tutorials). + + If ``None`` (default), will use ``self.defect_thermodynamics.chempots``. + + Note that you can also set + ``FermiSolver.defect_thermodynamics.chempots = ...`` or + ``DefectThermodynamics.chempots = ...`` to set the default chemical + potentials for all calculations, and you can set + ``FermiSolver.defect_thermodynamics.el_refs = ...`` or + ``DefectThermodynamics.el_refs = ...`` if you want to update the + elemental reference energies for any reason. + annealing_temperature (Optional[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 ``None`` (uses ``temperature`` + under thermodynamic equilibrium). + quenched_temperature (float): + Temperature in Kelvin at which to calculate the self-consistent + (constrained equilibrium) Fermi level and carrier concentrations, + given the fixed total concentrations, which should correspond to + operating temperature of the material (typically room temperature). + Default is 300 K. + temperature (float): + The temperature at which to solve for defect concentrations + and Fermi level, under thermodynamic equilibrium (if + ``annealing_temperature`` is not specified). + Defaults to 300 K. + tolerance (float): + The convergence criterion for the target variable. The search + stops when the target value change is less than this value. + Defaults to ``0.01``. + n_points (int): + The number of points to generate along each axis of the grid for + the initial search. Defaults to ``10``. + effective_dopant_concentration (Optional[float]): + The fixed concentration (in cm^-3) of an arbitrary dopant or + impurity in the material. This value is included in the charge + neutrality condition to analyze the Fermi level and doping + response under hypothetical doping conditions. + A positive value corresponds to donor doping, while a negative + value corresponds to acceptor doping. + Defaults to ``None``, corresponding to no additional extrinsic + dopant. + fix_charge_states (bool): + Whether to fix the concentrations of individual defect charge states + (``True``) or allow charge states to vary while keeping total defect + concentrations fixed (``False``) upon quenching. Not expected to be + physically sensible in most cases. Defaults to ``False``. + fixed_defects (Optional[dict[str, float]]): + A dictionary of defect concentrations to fix at the quenched temperature, + in the format: ``{defect_name: concentration}``. Concentrations should be + given in cm^-3. The this can be used to fix the concentrations of specific + defects regardless of the chemical potentials, or anneal-quench procedure + (e.g. to simulate the effect of a fixed impurity concentration). + If a fixed-concentration of a specific charge state is desired, + the defect name should be formatted as ``"defect_name_charge"``. + I.e. ``"v_O_+2"`` for a doubly positively charged oxygen vacancy. + Defaults to ``None``. + free_defects (Optional[list[str]]): + A list of defects to be excluded from high-temperature concentration + fixing, useful for highly mobile defects that are not expected + to be "frozen-in" upon quenching. Defaults to ``None``. + + Returns: + pd.DataFrame: + A ``DataFrame`` containing the results of the minimization or + maximization process, including the optimal chemical potentials and + the corresponding values of the target variable. + + Raises: + ValueError: + If neither ``chempots`` nor ``self.chempots`` is provided, or if + ``min_or_max`` is not ``"minimize"/"min"`` or ``"maximize"/"max"``. + """ + # Determine the dimensionality of the chemical potential space + chempots, el_refs = self._parse_and_check_grid_like_chempots(chempots) + n_chempots = len(chempots["elemental_refs"]) + + # Call the appropriate method based on dimensionality + if n_chempots == 2: + return self._min_max_X_line( + target, + min_or_max, + chempots, + annealing_temperature, + quenched_temperature, + temperature, + tolerance, + n_points, + effective_dopant_concentration, + fix_charge_states, + fixed_defects, + free_defects, + ) + return self._min_max_X_grid( + target, + min_or_max, + chempots, + annealing_temperature, + quenched_temperature, + temperature, + tolerance, + n_points, + effective_dopant_concentration, + fix_charge_states, + fixed_defects, + free_defects, + ) + + def _min_max_X_line( + self, + target: str, + min_or_max: str, + chempots: dict, + annealing_temperature: Optional[float] = None, + quenched_temperature: float = 300, + temperature: float = 300, + tolerance: float = 0.01, + n_points: int = 100, + effective_dopant_concentration: Optional[float] = None, + fix_charge_states: bool = False, + fixed_defects: Optional[dict[str, float]] = None, + free_defects: Optional[list[str]] = None, + ) -> pd.DataFrame: + r""" + Search for the chemical potentials that minimize or maximize a target + variable, such as electron concentration, within a specified tolerance. + + This function iterates over a grid of chemical potentials and "zooms in" on + the chemical potential that either minimizes or maximizes the target variable. + The process continues until the change in the target variable is less than + the specified tolerance. + + If ``annealing_temperature`` (and ``quenched_temperature``; 300 K by + default) are specified, then the frozen defect approximation is + employed, whereby total defect concentrations are calculated at the + elevated annealing temperature, then fixed at these values (unless + ``free_defects`` or ``fix_charge_states`` are specified) and the Fermi + level and relative charge state populations are recalculated at the + quenched temperature. Otherwise, if only ``temperature`` is specified, + then the Fermi level and defect/carrier concentrations are calculated + assuming thermodynamic equilibrium at that temperature. + + Args: + target (str): + The target variable to minimize or maximize, e.g., "Electrons (cm^-3)". + # TODO: Allow this to match a substring of the column name. + min_or_max (str): + Specify whether to "minimize" ("min") or "maximize" ("max"; default) + the target variable. + chempots (dict): + Dictionary of chemical potentials to scan over, in the ``doped`` + format (i.e. ``{"limits": [{'limit': [chempot_dict]}], ...}``) + -- the format generated by ``doped``\'s chemical potential parsing + functions (see tutorials). + + If ``None`` (default), will use ``self.defect_thermodynamics.chempots``. + + Note that you can also set + ``FermiSolver.defect_thermodynamics.chempots = ...`` or + ``DefectThermodynamics.chempots = ...`` to set the default chemical + potentials for all calculations, and you can set + ``FermiSolver.defect_thermodynamics.el_refs = ...`` or + ``DefectThermodynamics.el_refs = ...`` if you want to update the + elemental reference energies for any reason. + annealing_temperature (Optional[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 ``None`` (uses ``temperature`` + under thermodynamic equilibrium). + quenched_temperature (float): + Temperature in Kelvin at which to calculate the self-consistent + (constrained equilibrium) Fermi level and carrier concentrations, + given the fixed total concentrations, which should correspond to + operating temperature of the material (typically room temperature). + Default is 300 K. + temperature (float): + The temperature at which to solve for defect concentrations + and Fermi level, under thermodynamic equilibrium (if + ``annealing_temperature`` is not specified). + Defaults to 300 K. + tolerance (float): + The convergence criterion for the target variable. The search + stops when the target value change is less than this value. + Defaults to ``0.01``. + n_points (int): + The number of points to generate along each axis of the grid for + the initial search. Defaults to ``10``. + effective_dopant_concentration (Optional[float]): + The fixed concentration (in cm^-3) of an arbitrary dopant or + impurity in the material. This value is included in the charge + neutrality condition to analyze the Fermi level and doping + response under hypothetical doping conditions. + A positive value corresponds to donor doping, while a negative + value corresponds to acceptor doping. + Defaults to ``None``, corresponding to no additional extrinsic + dopant. + fix_charge_states (bool): + Whether to fix the concentrations of individual defect charge states + (``True``) or allow charge states to vary while keeping total defect + concentrations fixed (``False``) upon quenching. Not expected to be + physically sensible in most cases. Defaults to ``False``. + fixed_defects (Optional[dict[str, float]]): + A dictionary of defect concentrations to fix at the quenched temperature, + in the format: ``{defect_name: concentration}``. Concentrations should be + given in cm^-3. The this can be used to fix the concentrations of specific + defects regardless of the chemical potentials, or anneal-quench procedure + (e.g. to simulate the effect of a fixed impurity concentration). + If a fixed-concentration of a specific charge state is desired, + the defect name should be formatted as ``"defect_name_charge"``. + I.e. ``"v_O_+2"`` for a doubly positively charged oxygen vacancy. + Defaults to ``None``. + free_defects (Optional[list[str]]): + A list of defects to be excluded from high-temperature concentration + fixing, useful for highly mobile defects that are not expected + to be "frozen-in" upon quenching. Defaults to ``None``. + + Returns: + pd.DataFrame: + A ``DataFrame`` containing the results of the minimization or + maximization process, including the optimal chemical potentials and + the corresponding values of the target variable. + + Raises: + ValueError: + If neither ``chempots`` nor ``self.chempots`` is provided, or if + ``min_or_max`` is not ``"minimize"/"min"`` or ``"maximize"/"max"``. + """ + el_refs = chempots["elemental_refs"] + chempots_label = list(el_refs.keys()) # Assuming 1D space, focus on one label. + chempots_labels = [f"μ_{label}" for label in chempots_label] + + # Initial line scan: Rich to Poor limit + rich = self._get_single_chempot_dict(f"{chempots_label[0]}-rich") + poor = self._get_single_chempot_dict(f"{chempots_label[0]}-poor") + starting_line = self._get_interpolated_chempots(rich[0], poor[0], n_points) + + previous_value = None + + while True: + # Calculate results based on the given temperature conditions + if annealing_temperature is not None: + results_df = pd.concat( + [ + self.pseudo_equilibrium_solve( + single_chempot_dict={ + k.replace("μ_", ""): v for k, v in chempot_series.items() + }, + el_refs=el_refs, + annealing_temperature=annealing_temperature, + quenched_temperature=quenched_temperature, + effective_dopant_concentration=effective_dopant_concentration, + fix_charge_states=fix_charge_states, + free_defects=free_defects, + fixed_defects=fixed_defects, + ) + for chempot_series in tqdm(starting_line) + ] + ) + else: + results_df = pd.concat( + [ + self.equilibrium_solve( + single_chempot_dict={ + k.replace("μ_", ""): v for k, v in chempot_series.items() + }, + el_refs=el_refs, + temperature=temperature, + effective_dopant_concentration=effective_dopant_concentration, + ) + for chempot_series in tqdm(starting_line) + ] + ) + + if target in results_df.columns: + if "min" in min_or_max: + target_chempot = results_df[results_df[target] == results_df[target].min()][ + chempots_labels + ] + target_dataframe = results_df[results_df[target] == results_df[target].min()] + elif "max" in min_or_max: + target_chempot = results_df[results_df[target] == results_df[target].max()][ + chempots_labels + ] + target_dataframe = results_df[results_df[target] == results_df[target].max()] + current_value = ( + results_df[target].min() if "min" in min_or_max else results_df[target].max() + ) + + else: + # Filter the DataFrame for the specific defect + filtered_df = results_df[results_df.index == target] + # Find the row where "Concentration (cm^-3)" is at its minimum or maximum + if "min" in min_or_max: + min_value = filtered_df["Concentration (cm^-3)"].min() + target_chempot = results_df.loc[results_df["Concentration (cm^-3)"] == min_value][ + chempots_labels + ] + target_dataframe = results_df[ + results_df[chempots_labels].eq(target_chempot.iloc[0]).all(axis=1) + ] + + elif "max" in min_or_max: + max_value = filtered_df["Concentration (cm^-3)"].max() + target_chempot = results_df.loc[results_df["Concentration (cm^-3)"] == max_value][ + chempots_labels + ] + target_dataframe = results_df[ + results_df[chempots_labels].eq(target_chempot.iloc[0]).all(axis=1) + ] + current_value = ( + filtered_df["Concentration (cm^-3)"].min() + if "min" in min_or_max + else filtered_df["Concentration (cm^-3)"].max() + ) + + # Check if the change in the target value is less than the tolerance + if ( + previous_value is not None + and abs((current_value - previous_value) / previous_value) < tolerance + ): + break + previous_value = current_value + target_chempot = target_chempot.drop_duplicates(ignore_index=True) + + starting_line = self._get_interpolated_chempots( + chempot_start={ + k: v - 0.1 for k, v in target_chempot.iloc[0].items() if k in chempots_labels + }, + chempot_end={ + k: v + 0.1 for k, v in target_chempot.iloc[0].items() if k in chempots_labels + }, + n_points=100, + ) + + return target_dataframe + + def _min_max_X_grid( self, target: str, min_or_max: str = "max",