# Import necessary tools from pymatgen
from pymatgen import MPRester
from pymatgen.analysis.pourbaix.maker import PourbaixDiagram
from pymatgen.analysis.pourbaix.analyzer import PourbaixAnalyzer
from pymatgen.analysis.pourbaix.plotter import PourbaixPlotter
%matplotlib inline

# Initialize the MP Rester
mpr = MPRester()

"""
Explanation: Plotting a pourbaix diagram using the Materials Project API
If you use this infrastructure, please consider citing the following work:
K. A. Persson, B. Waldwick, P. Lazic and G. Ceder, Phys. Rev. B, 2012, 85, 235438.
A. K. Singh, L. Zhou, A. Shinde, S. K. Suram, J. H. Montoya, D. Winston, J. M. Gregoire and K. A. Persson, Chem. Mater., 2017, DOI: acs.chemmater.7b03980.
Notebook Author: Joseph Montoya

The Materials Project REST interface includes functionality to construct pourbaix diagrams from computed entries.  Note that the Pourbaix infrastructure is still undergoing revision, but now includes a simplified interface that enables MP and pymatgen users to fetch entries that have been processed according to the Materials Project Aqueous Compatibility scheme.  Thus, users can reproduce web Pourbaix diagrams in two or three steps in pymatgen.
End of explanation
"""


# Get all pourbaix entries corresponding to the Cu-O-H chemical system.
entries = mpr.get_pourbaix_entries(["Cu"])

"""
Explanation: To retrieve entries necessary to construct a Pourbaix Diagram use MPRester.get_pourbaix_entries(LIST_OF_ENTRIES) with a list of entries comprising your chemical system.  It is not necessary to include 'O' and 'H' in your list, as they are added automatically.  This function also contains all of necessary preprocessing to ensure the PourbaixEntries are compatible with the pymatgen PourbaixDiagram constructor.
End of explanation
"""


# Construct the PourbaixDiagram object
pd = PourbaixDiagram(entries)

"""
Explanation: Pourbaix diagrams can be constructed using PourbaixDiagram(RETRIEVED_ENTRIES) as below.  Note that a comp_dict keyword argument may also be supplied to the PourbaixDiagram constructor if a fixed composition for a multi-element pourbaix diagram is desired.
End of explanation
"""


# Construct the pourbaix analyzer
pa = PourbaixAnalyzer(pd)

# Get an entry stability as a function of pH and V
entry = [e for e in entries if e.entry_id == 'mp-1692'][0]
print("CuO's potential energy per atom relative to the most",
      "stable decomposition product is {} eV/atom".format(
          pa.get_entry_stability(entry, pH=7, V=-0.2)))

"""
Explanation: The PourbaixAnalyzer includes a number of useful functions for determining stable species and stability of entries relative to a given pourbaix facet (i.e. as a function of pH and V).
End of explanation
"""


plotter = PourbaixPlotter(pd)

plotter.get_pourbaix_plot().show()

"""
Explanation: This suggests that CuO, for example, has a large driving force for decomposition at neutral pH and mildly reducing conditions.
To see this in more detail, we can plot the pourbaix diagram.  The PourbaixPlotter object is also initialized using an instance of the PourbaixDiagram object.
End of explanation
"""


plt = plotter.plot_entry_stability(entry)
plt.show()

"""
Explanation: The PourbaixPlotter object can also plot the relative stability of an entry across the pourbaix diagram.  To do this, use the PourbaixPlotter.plot_entry_stability method.
End of explanation
"""


# Get all pourbaix entries corresponding to the Fe-O-H chemical system.
entries = mpr.get_pourbaix_entries(["Bi", "V"])
# Construct the PourbaixDiagram object
pd = PourbaixDiagram(entries, comp_dict={"Bi": 0.5, "V": 0.5},
                     conc_dict={"Bi": 1e-8, "V": 1e-8})

# Construct the pourbaix analyzer
pa = PourbaixAnalyzer(pd)
plotter = PourbaixPlotter(pd)
plt = plotter.get_pourbaix_plot()
plt.show()

"""
Explanation: Plotting k-nary systems
Pymatgen also supports binary/ternary pourbaix diagrams with fixed compositions of the non-H or O elements.  This is achieved by finding all possible combinations of entries that fulfill the composition constraint and treating them as individual entries in pourbaix space.  Note that you can supply a composition dictionary and to further tune the pourbaix diagram.
End of explanation
"""


bivo4_entry = [entry for entry in entries if entry.entry_id=="mp-613172"][0]
plt = plotter.plot_entry_stability(bivo4_entry)

"""
Explanation: Getting the heatmaps for a solid entry in these cases is a bit more involved, because many of the regions of the pourbaix diagram include multiphase entries.  Here's an example for this case.
End of explanation
"""


bio2_entry = [entry for entry in entries if entry.entry_id=="mp-557993"][0]
plt = plotter.plot_entry_stability(bio2_entry)

"""
Explanation: If a compound is provided that doesn't meet the constraints of the composition (in this case, equal parts Bi and V), the plotter will attempt to find the most stable entry containing that solid and any spectating ions.  In this case, it estimates the stability of BiO$_2$(s) + VO$_4$-.
End of explanation
"""


