Add get_doped_chempots_from_entries convenience function to `doped.…

…chemical_potentials`

doped/chemical_potentials.py

@@ -9,7 +9,7 @@
 import itertools
 import os
 import warnings
-from collections.abc import Iterable
+from collections.abc import Iterable, Sequence
 from copy import deepcopy
 from pathlib import Path
 from typing import Any, Optional, Union
@@ -1681,6 +1681,77 @@ def __init__(
  _name_entries_and_handle_duplicates(self.entries) # set entry names
 
 
+def get_doped_chempots_from_entries(
+ entries: Sequence[Union[ComputedEntry, ComputedStructureEntry, PDEntry]],
+ composition: Union[str, Composition, ComputedEntry],
+ sort_by: Optional[str] = None,
+ single_chempot_limit: bool = False,
+) -> dict:
+ r"""
+ Given a list of ``ComputedEntry``\s / ``ComputedStructureEntry``\s /
+ ``PDEntry``\s and the bulk ``composition``, returns the chemical potential
+ limits dictionary in the ``doped`` format (i.e. ``{"limits": [{'limit':
+ [chempot_dict]}], ...}``) for the host material.
+
+ Args:
+ entries (list[ComputedEntry]):
+ List of ``ComputedEntry``\s / ``ComputedStructureEntry``\s /
+ ``PDEntry``\s for the chemical system, from which to determine
+ the chemical potential limits for the host material (``composition``).
+ composition (str, Composition, ComputedEntry):
+ Composition of the host material either as a string
+ (e.g. 'LiFePO4') a ``pymatgen`` ``Composition`` object (e.g.
+ ``Composition('LiFePO4')``), or a ``ComputedEntry`` object.
+ sort_by (str):
+ If set, will sort the chemical potential limits in the output
+ ``DataFrame`` according to the chemical potential of the specified
+ element (from element-rich to element-poor conditions).
+ single_chempot_limit (bool):
+ If set to ``True``, only returns the first chemical potential limit
+ in the calculated chemical potentials dictionary. Mainly intended for
+ internal ``doped`` usage when the host material is calculated to be
+ unstable with respect to the competing phases.
+
+ Returns:
+ dict:
+ Dictionary of chemical potential limits in the ``doped`` format.
+ """
+ if isinstance(composition, (str, Composition)):
+ composition = Composition(composition)
+ else:
+ composition = composition.composition
+
+ phase_diagram = PhaseDiagram(
+ entries,
+ list(map(Element, composition.elements)), # preserve bulk comp element ordering
+ )
+ chem_lims = phase_diagram.get_all_chempots(composition.reduced_composition)
+ chem_lims_iterator = list(chem_lims.items())[:1] if single_chempot_limit else chem_lims.items()
+
+ # remove Element to make it JSONable:
+ no_element_chem_lims = {k: {str(kk): vv for kk, vv in v.items()} for k, v in chem_lims_iterator}
+
+ if sort_by is not None:
+ no_element_chem_lims = dict(
+ sorted(no_element_chem_lims.items(), key=lambda x: x[1][sort_by], reverse=True)
+ )
+
+ chempots = {
+ "limits": no_element_chem_lims,
+ "elemental_refs": {str(el): ent.energy_per_atom for el, ent in phase_diagram.el_refs.items()},
+ "limits_wrt_el_refs": {},
+ }
+
+ # relate the limits to the elemental energies
+ for limit, chempot_dict in chempots["limits"].items():
+ relative_chempot_dict = copy.deepcopy(chempot_dict)
+ for e in relative_chempot_dict:
+ relative_chempot_dict[e] -= chempots["elemental_refs"][e]
+ chempots["limits_wrt_el_refs"].update({limit: relative_chempot_dict})
+
+ return chempots
+
+
 class CompetingPhasesAnalyzer:
  # TODO: Allow parsing using pymatgen ComputedEntries as well, to aid interoperability with
  # high-througput architectures like AiiDA or atomate2. See:
@@ -2151,7 +2222,7 @@ def calculate_chempots(
 
  self._intrinsic_phase_diagram = PhaseDiagram(
  intrinsic_phase_diagram_entries,
- map(Element, self.bulk_composition.elements),
+ list(map(Element, self.bulk_composition.elements)), # preserve bulk comp element ordering
  )
 
  # check if it's stable and if not, warn user and downshift to get _least_ unstable point on convex
@@ -2179,38 +2250,15 @@ def calculate_chempots(
  )
  self._intrinsic_phase_diagram = PhaseDiagram(
  [*intrinsic_phase_diagram_entries, renormalised_bulk_pde],
- map(Element, self.bulk_composition.elements),
+ list(map(Element, self.bulk_composition.elements)), # preserve bulk comp element ordering
  )
 
- chem_lims = self._intrinsic_phase_diagram.get_all_chempots(self.bulk_composition)
-
- # remove Element to make it JSONable:
- no_element_chem_lims = {k: {str(kk): vv for kk, vv in v.items()} for k, v in chem_lims.items()}
-
- if unstable_host:
- no_element_chem_lims = {
- k: {str(kk): vv for kk, vv in v.items()} for k, v in list(chem_lims.items())[:1]
- }
-
- if sort_by is not None:
- no_element_chem_lims = dict(
- sorted(no_element_chem_lims.items(), key=lambda x: x[1][sort_by], reverse=True)
- )
-
- self._intrinsic_chempots = {
- "limits": no_element_chem_lims,
- "elemental_refs": {
- str(el): ent.energy_per_atom for el, ent in self._intrinsic_phase_diagram.el_refs.items()
- },
- "limits_wrt_el_refs": {},
- }
-
- # relate the limits to the elemental energies
- for limit, chempot_dict in self._intrinsic_chempots["limits"].items():
- relative_chempot_dict = copy.deepcopy(chempot_dict)
- for e in relative_chempot_dict:
- relative_chempot_dict[e] -= self._intrinsic_chempots["elemental_refs"][e]
- self._intrinsic_chempots["limits_wrt_el_refs"].update({limit: relative_chempot_dict})
+ self._intrinsic_chempots = get_doped_chempots_from_entries(
+ self._intrinsic_phase_diagram.entries,
+ self.bulk_composition,
+ sort_by=sort_by,
+ single_chempot_limit=unstable_host,
+ )
 
  # get chemical potentials as pandas dataframe
  chemical_potentials = []
@@ -2340,7 +2388,7 @@ def intrinsic_chempots(self) -> dict:
  return self._intrinsic_chempots
 
  @property
- def intrinsic_phase_diagram(self) -> dict:
+ def intrinsic_phase_diagram(self) -> PhaseDiagram:
  """
  Returns the calculated intrinsic phase diagram.
  """