Bug Computation of the decomposition pathway

Materials Project
, ,
1

Dear Materials Project team,

Today, I tried to check the stability of a structure using the API. According to my calculations, my material would not be stable and would decompose as follows:

decomp, e_above_hull = pd.get_decomp_and_e_above_hull(corrected_entry)  # "pd" is the Cs-I-S PhaseDiagram, and "corrected_entry" is the CsS3I custom entry for which I am checking the formation energy.

print(f"Energy above hull: {e_above_hull:.4f} eV/atom\n")
print(f"Decomposition pathway of {corrected_entry.composition.reduced_formula}:")
for phase, amount in decomp.items():
    print(f"  {amount:.4f} -> {phase.composition.reduced_formula}")

Energy above hull: 0.4885 eV/atom
Decomposition pathway of CsS3I:
  0.6000 → CsS3
  0.2500 → CsI4
  0.1500 → S

However, as you can see, there is fundamentally something wrong with the suggested decomposition pathway: the computed pathway does not respect the stoichiometry of the compound.

This issue is also observable in the decomposition pathways of other materials, such as GaAs3 and GeBiO3.

Could you please confirm if this is indeed a bug?

2

I think the issue I encountered might actually stem from the number of atoms per unit cell in the structures involved.

For example, in the decomposition reaction of CsS₃I:

  • The CsS₃I structure has 8 formula units per unit cell: CsS₃I ≡ Cs₈S₂₄I₈
  • The CsS₃ structure has 4 formula units per unit cell: CsS₃ ≡ Cs₄S₁₂
  • The CsI₄ structure has 4 formula units per unit cell: CsI₄ ≡ Cs₄I₁₆
  • The S structure has 32 formula units per unit cell: S ≡ S₃₂

If I write the chemical decomposition equation using these values: 0.6 Cs₄S₁₂ + 0.25 Cs₈S₂₄I₈ + 0.15 S₃₂ → 2.4 CsS₃ + CsI₄ + 4.8 S ≈ 4 Cs₀.₈₅SI₃ ≈ 4 CsSI₃

I suspect this is what the code is doing. Am I correct?

If so, it would be helpful to have access to:

  1. The stoichiometric formula coefficients (based on the actual unit cell composition).
  2. The “reduced” stoichiometric coefficients (using the reduced formula units).

Additionally, it would be useful to control the threshold for determining when a decomposition pathway is considered stoichiometrically valid.

3

Hi @gjoalland13, in the pymatgen objects, the formula tied to a Composition or Structure always respects the number of atoms in that object, regardless if they represent multiple formula units. The reduced_formula will then reduce a given composition by the greatest common divisor of the species’ multiplicity, with the exception of a few special formulas defined here :

from pymatgen.core import Composition

for f in ('Li2Fe2','H2SO4','O3'):
  print(f'Original={f}, full formula={Composition(f).formula}, reduced formula={Composition(f).reduced_formula}')
>>> Original=Li2Fe2, full formula=Li2 Fe2, reduced formula=LiFe
>>> Original=H2SO4, full formula=H2 S1 O4, reduced formula=H2SO4
>>> Original=O3, full formula=O3, reduced formula=O2

Which can be confusing. When tracing the decomposition products, pymatgen indeed looks only at the reduced formulas

If so, it would be helpful to have access to:

  1. The stoichiometric formula coefficients (based on the actual unit cell composition).
  2. The “reduced” stoichiometric coefficients (using the reduced formula units).

For these points, you may want to post a feature-request issue on pymatgen-core