diff --git a/src/mp_api/client.py b/src/mp_api/client.py index 731bb709..2ee6bddb 100644 --- a/src/mp_api/client.py +++ b/src/mp_api/client.py @@ -15,7 +15,7 @@ from pymatgen.analysis.magnetism import Ordering from pymatgen.analysis.phase_diagram import PhaseDiagram from pymatgen.analysis.pourbaix_diagram import IonEntry -from pymatgen.core import Element, Structure +from pymatgen.core import Element, Structure, Composition from pymatgen.core.ion import Ion from pymatgen.entries.computed_entries import ComputedEntry from pymatgen.io.vasp import Chgcar @@ -513,14 +513,14 @@ def get_pourbaix_entries( ion_ref_comps = [ Ion.from_formula(d["data"]["RefSolid"]).composition for d in ion_data ] - ion_ref_elts = list( + ion_ref_elts = set( itertools.chain.from_iterable(i.elements for i in ion_ref_comps) ) # TODO - would be great if the commented line below would work # However for some reason you cannot process GibbsComputedStructureEntry with # MaterialsProjectAqueousCompatibility ion_ref_entries = self.get_entries_in_chemsys( - list(set([str(e) for e in ion_ref_elts] + ["O", "H"])), + list([str(e) for e in ion_ref_elts] + ["O", "H"]), # use_gibbs=use_gibbs ) @@ -553,10 +553,6 @@ def get_pourbaix_entries( pbx_entries = [PourbaixEntry(e, f"ion-{n}") for n, e in enumerate(ion_entries)] # Construct the solid pourbaix entries from filtered ion_ref entries - ion_ref_comps = [e.composition for e in ion_entries] - ion_ref_elts = list( - itertools.chain.from_iterable(i.elements for i in ion_ref_comps) - ) extra_elts = ( set(ion_ref_elts) - {Element(s) for s in chemsys} diff --git a/tests/test_mprester.py b/tests/test_mprester.py index 5f250938..aed28f14 100644 --- a/tests/test_mprester.py +++ b/tests/test_mprester.py @@ -1,6 +1,8 @@ import os import random import typing +import numpy as np +import itertools import pytest from emmet.core.summary import HasProps @@ -14,12 +16,14 @@ from pymatgen.analysis.pourbaix_diagram import IonEntry, PourbaixDiagram, PourbaixEntry from pymatgen.analysis.wulff import WulffShape from pymatgen.core.periodic_table import Element +from pymatgen.core.ion import Ion from pymatgen.electronic_structure.bandstructure import ( BandStructure, BandStructureSymmLine, ) from pymatgen.electronic_structure.dos import CompleteDos from pymatgen.entries.computed_entries import ComputedEntry, GibbsComputedStructureEntry +from pymatgen.entries.compatibility import MaterialsProjectAqueousCompatibility from pymatgen.io.cif import CifParser from pymatgen.io.vasp import Chgcar from pymatgen.phonon.bandstructure import PhononBandStructureSymmLine @@ -170,7 +174,15 @@ def test_get_pourbaix_entries(self, mpr): with pytest.raises(ValueError, match="Solid compatibility can only be"): mpr.get_pourbaix_entries("Ti-O", solid_compat=None) - # TODO - old tests copied from pymatgen. Update or delete + # test removal of extra elements from reference solids + # Li-Zn-S has Na in reference solids + pbx_entries = mpr.get_pourbaix_entries("Li-Zn-S") + assert not any([e for e in pbx_entries if 'Na' in e.composition]) + + # Ensure entries are pourbaix compatible + PourbaixDiagram(pbx_entries) + + # TODO - old tests copied from pymatgen with specific energy values. Update or delete # fe_two_plus = [e for e in pbx_entries if e.entry_id == "ion-0"][0] # self.assertAlmostEqual(fe_two_plus.energy, -1.12369, places=3) # @@ -182,8 +194,6 @@ def test_get_pourbaix_entries(self, mpr): # so4_two_minus = pbx_entries[9] # self.assertAlmostEqual(so4_two_minus.energy, 0.301511, places=3) - # Ensure entries are pourbaix compatible - PourbaixDiagram(pbx_entries) def test_get_ion_entries(self, mpr): entries = mpr.get_entries_in_chemsys("Ti-O-H") @@ -191,8 +201,7 @@ def test_get_ion_entries(self, mpr): ion_entry_data = mpr.get_ion_reference_data_for_chemsys("Ti-O-H") ion_entries = mpr.get_ion_entries(pd, ion_entry_data) assert len(ion_entries) == 5 - - assert isinstance(ion_entries[0], IonEntry) + assert all([isinstance(i, IonEntry) for i in ion_entries]) # test an incomplete phase diagram entries = mpr.get_entries_in_chemsys("Ti-O") @@ -200,6 +209,36 @@ def test_get_ion_entries(self, mpr): with pytest.raises(ValueError, match="The phase diagram chemical system"): mpr.get_ion_entries(pd) + # test ion energy calculation + ion_data = mpr.get_ion_reference_data_for_chemsys('S') + ion_ref_comps = [ + Ion.from_formula(d["data"]["RefSolid"]).composition for d in ion_data + ] + ion_ref_elts = set( + itertools.chain.from_iterable(i.elements for i in ion_ref_comps) + ) + ion_ref_entries = mpr.get_entries_in_chemsys( + list([str(e) for e in ion_ref_elts] + ["O", "H"]) + ) + mpc = MaterialsProjectAqueousCompatibility() + ion_ref_entries = mpc.process_entries(ion_ref_entries) + ion_ref_pd = PhaseDiagram(ion_ref_entries) + ion_entries = mpr.get_ion_entries(ion_ref_pd, ion_ref_data=ion_data) + + # In ion ref data, SO4-2 is -744.27 kJ/mol; ref solid is -1,279.0 kJ/mol + # so the ion entry should have an energy (-744.27 +1279) = 534.73 kJ/mol + # or 5.542 eV/f.u. above the energy of Na2SO4 + so4_two_minus = [e for e in ion_entries if e.ion.reduced_formula == "SO4[-2]"][0] + + # the ref solid is Na2SO4, ground state mp-4770 + # the rf factor correction is necessary to make sure the composition + # of the reference solid is normalized to a single formula unit + ref_solid_entry = [e for e in ion_ref_entries if e.entry_id == 'mp-4770'][0] + rf = ref_solid_entry.composition.get_reduced_composition_and_factor()[1] + solid_energy = ion_ref_pd.get_form_energy(ref_solid_entry) / rf + + assert np.allclose(so4_two_minus.energy, solid_energy + 5.542, atol=1e-3) + def test_get_phonon_data_by_material_id(self, mpr): bs = mpr.get_phonon_bandstructure_by_material_id("mp-11659") assert isinstance(bs, PhononBandStructureSymmLine)