def get_structure_from_mp(formula):
    """
    Convenience method to get a crystal from the Materials Project database via
    the API. Requires PMG_MAPI_KEY to be set.

    Args:
        formula (str): A formula

    Returns:
        (Structure) The lowest energy structure in Materials Project with that
            formula.
    """
    if not SETTINGS.get("PMG_MAPI_KEY"):
        raise RuntimeError("PMG_MAPI_KEY must be set in .pmgrc.yaml to use this "
                           "function.")

    from pymatgen.matproj.rest import MPRester
    m = MPRester()
    entries = m.get_entries(formula, inc_structure=True)
    if len(entries) == 0:
        raise ValueError("No structure with formula %s in Materials Project!" %
                         formula)
    elif len(entries) > 1:
        warnings.warn("%d structures with formula %s found in Materials Project."
                      "The lowest energy structure will be returned." %
                      (len(entries), formula))
    return min(entries, key=lambda e: e.energy_per_atom).structure

def main(api="", queryid=""):
    """Get VASP inputs for Materials Project structure
        Args:
            api <str>: Materials Project API key
            queryid <str>: Materials Project ID of the structure
        Returns:
            creates a folder named with that mpid, and including
            some VASP input files.
    """
    if api == "":
        print "Must have an API key from materialsproject.org"
        return None
    if queryid == "":
        print "No MP structure ID given. Exiting."
        return None
    rest_adapter = MPRester(api)
    entries=list()
    proplist=list()
    proplist.append('pretty_formula')
    proplist.append('structure')
    proplist.append('potcar')
    proplist.append('material_id')
    
    myentry = rest_adapter.mpquery(criteria={'material_id':queryid}, properties=proplist)
    if len(myentry) == 0:
        print "Could not find entry for %s as material_id. Trying entry_id." % queryid
        myentry = rest_adapter.mpquery(criteria={'entry_id':queryid}, properties=proplist)
    if len(myentry) == 0:
        print "Could not find entry for %s" % queryid
        return None
    entries.extend(myentry)

    workdir = os.getcwd()
    from pymatgen.io.vaspio_set import MITVaspInputSet, MPVaspInputSet
    for entry in entries: 
        mpvis = MPVaspInputSet()
        myname = str(entry['pretty_formula'])
        #print entry['structure'].composition
        #print entry['structure'].formula
        #myname = entry['pretty_formula']
        myname = myname.replace("(","_").replace(")","_")
        myname = myname + "_" + entry['material_id']
        os.mkdir(myname)
        os.chdir(myname)
        mystructure = entry['structure']
        if mystructure.num_sites <= 10:
            mystructure.make_supercell([2,2,2])
        #mystructure.perturb(0.01)
        incar = mpvis.get_incar(mystructure)
        incar.write_file("INCAR")
        potcar = mpvis.get_potcar(mystructure)
        potcar.write_file("POTCAR")
        #potcar_symbols = mpvis.get_potcar_symbols(mystructure) 
        myposcar=Poscar(mystructure)
        mykpoints=mpvis.get_kpoints(mystructure)
        mykpoints.write_file("KPOINTS")


        myposcar.write_file("POSCAR")
        os.chdir(workdir)

 def _test_from_pm(self):
     """
     MPのデータをプロット
     """
     ent_list = []
     mpr = MPRester("WTxsDhRV7g2Mcbqw")
     composition = ['Fe', 'Ni', 'Si']
     entries = mpr.get_entries_in_chemsys(composition)
     for entry in entries:
         formula = entry.as_dict()['data']['unit_cell_formula']
         formation_e = entry.as_dict()['data']['formation_energy_per_atom']
         if formation_e <= 0:
             single_data = []
             sum_atoms = 0
             for elements, index in formula.items():
                 sum_atoms += index
             for element in composition:
                 try:
                     single_data.append(formula[element] / sum_atoms)
                 except KeyError:
                     single_data.append(0)
             single_data.append(formation_e * 96.485344520851)
             ent_list.append(single_data)
     initial_base = [[x[0], x[1], x[3]] for x in ent_list[0:3]]
     not_base = [[x[0], x[1], x[3]] for x in ent_list[3:]]
     fig = pylab.figure()
     ax = Axes3D(fig)
     convex_hull.draw_convex_hull(ax, initial_base, not_base,
                                  ['Fe', 'Ni', 'Si'], [-60, 5],
                                  color='magenta')
     #pint(not_base)
     pylab.show()

def main(comp="La0.5Sr0.5MnO3", energy=-43.3610, ostart="", oend="", ostep=""):
    """Get energy above hull for a composition
        Args:
            comp <str>: Composition in string form
            energy <float>: Energy PER FORMULA UNIT of composition given
            (Leave the following arguments blank for a non-grand potential
                phase diagram.)
            ostart <float>: Starting oxygen chemical potential. 
            oend <float>: Ending oxygen chemical potential. 
            ostep <float>: Step for oxygen chemical potential
        Returns:
            Prints to screen
    """        
    #a = MPRester("<YOUR_MPREST_API_KEY_HERE>")
    a = MPRester("wfmUu5VSsDCvIrhz")
    
    mycomp=Composition(comp)
    print "Composition: ", mycomp
    myenergy=energy
    print "Energy: ", myenergy
    myPDEntry = PDEntry(mycomp, myenergy)

    elements = mycomp.elements
    ellist = map(str, elements)
    
    chemsys_entries = a.get_entries_in_chemsys(ellist)
    #For reference: other ways of getting entries
    #entries = a.mpquery(criteria={'elements':{'$in':['La','Mn'],'$all':['O']},'nelements':3})
    #entries = a.mpquery(criteria={'elements':{'$in':['La','Mn','O'],'$all':['O']}},properties=['pretty_formula'])
    #entries = a.get_entries_in_chemsys(['La', 'Mn', 'O', 'Sr'])
   
    if ostart=="": #Regular phase diagram
        entries = list(chemsys_entries)
        entries.append(myPDEntry)
        pd = PhaseDiagram(entries)
        #plotter = PDPlotter(gppd)
        #plotter.show()
        ppda = PDAnalyzer(pd)
        eabove=ppda.get_decomp_and_e_above_hull(myPDEntry)
        print "Energy above hull: ", eabove[1]
        print "Decomposition: ", eabove[0]
        return eabove
    else: #Grand potential phase diagram
        orange = np.arange(ostart, oend+ostep, ostep) #add ostep because otherwise the range ends before oend
        for o_chem_pot in orange:
            entries = list(chemsys_entries)
            myGrandPDEntry = GrandPotPDEntry(myPDEntry,{Element('O'): float(o_chem_pot)}) #need grand pot pd entry for GPPD
            entries.append(myGrandPDEntry)
            gppd = GrandPotentialPhaseDiagram(entries,{Element('O'): float(o_chem_pot)})
            gppda = PDAnalyzer(gppd)
            geabove=gppda.get_decomp_and_e_above_hull(myGrandPDEntry, True)
            print "******** Decomposition for mu_O = %s eV ********" % o_chem_pot
            print "%30s%1.4f" % ("mu_O: ",o_chem_pot)
            print "%30s%1.4f" % ("Energy above hull (eV): ",geabove[1])
            decomp=geabove[0]
            #print "Decomp: ", decomp
            print "%30s" % "Decomposition: "
            for dkey in decomp.keys():
                print "%30s:%1.4f" % (dkey.composition,decomp[dkey])
    return

def mp_id(mpid):
    """
    materials project の ID から POSCAR を作成
    """
    mpr = MPRester("WTxsDhRV7g2Mcbqw")
    strctr = mpr.get_structure_by_material_id(mpid)
    poscar = Poscar(strctr)
    poscar.write_file('POSCAR_{0}'.format(mpid))

 def test_find_structure(self):
     # nosetests pymatgen/matproj/tests/test_rest.py:MPResterTest.test_find_structure
     m = MPRester()
     ciffile = os.path.join(test_dir, 'Fe3O4.cif')
     data = m.find_structure(ciffile)
     self.assertTrue(len(data) > 1)
     s = CifParser(ciffile).get_structures()[0]
     data = m.find_structure(s)
     self.assertTrue(len(data) > 1)

def get_INCARs(comp_name, mp_api_key, SOC):
    from pymatgen.matproj.rest import MPRester
    matproj = MPRester(mp_api_key)

    structure = matproj.get_structure_by_material_id(comp_name)
    from pymatgen.io.vasp.sets import MPVaspInputSet
    from pymatgen.io.vasp.sets import MPGGAVaspInputSet   # +U is turned off (if any)
    from pymatgen.io.vasp.sets import MPHSEVaspInputSet
    from pymatgen.io.vasp.sets import MPStaticVaspInputSet
    from pymatgen.io.vasp.sets import MPNonSCFVaspInputSet
    from pymatgen.io.vasp.sets import MPStaticDielectricDFPTVaspInputSet
    from pymatgen.io.vasp.sets import MPBSHSEVaspInputSet
    from pymatgen.io.vasp.sets import MPOpticsNonSCFVaspInputSet

    v = MPVaspInputSet()
    indic = v.get_incar(structure)
    k = v.get_all_vasp_input(structure,generate_potcar = False)
    incar_file = open('INCAR','w')
    if SOC in ['TRUE', 'True', 'true']:
	indic['LMAXMIX'] = 4
	indic['NBANDS'] = 128
	indic['LSORBIT'] = '.TRUE.'
	indic.pop('MAGMOM', None)
    incar_file.write(str(indic))

    self = MPStaticVaspInputSet()
    indic = self.get_incar(structure)
    indic['NEDOS'] = 9001
    incar_file = open('INCAR_SELF','w')
    if SOC in ['TRUE', 'True', 'true']:
	indic['LMAXMIX'] = 4
	indic['NBANDS'] = 128
	indic['LSORBIT'] = '.TRUE.'
	indic.pop('MAGMOM', None)
    incar_file.write(str(indic))

    user_incar_settings = {"NBANDS" : 100}
    nself = MPNonSCFVaspInputSet(user_incar_settings)
    indic = nself.get_incar(structure)
    incar_file = open('INCAR_NSELF','w')
    if SOC in ['TRUE', 'True', 'true']:
	indic['LMAXMIX'] = 4
	indic['NBANDS'] = 128
	indic['LSORBIT'] = '.TRUE.'
	indic.pop('MAGMOM', None)
    incar_file.write(str(indic))

    HSE = MPHSEVaspInputSet()
    indic = HSE.get_incar(structure)
    incar_file = open('INCAR_HSE','w')
    if SOC in ['TRUE', 'True', 'true']:
	indic['LMAXMIX'] = 4
	indic['NBANDS'] = 128
	indic['LSORBIT'] = '.TRUE.'
	indic.pop('MAGMOM', None)
    incar_file.write(str(indic))

 def test_find_structure(self):
     # nosetests pymatgen/matproj/tests/test_rest.py:MPResterTest.test_find_structure
     # self.rester points to rest/v2 by default which doesn't have the find_structure endpoint
     m = MPRester(endpoint="https://www.materialsproject.org/rest")
     ciffile = os.path.join(test_dir, 'Fe3O4.cif')
     data = m.find_structure(ciffile)
     self.assertTrue(len(data) > 1)
     s = CifParser(ciffile).get_structures()[0]
     data = m.find_structure(s)
     self.assertTrue(len(data) > 1)

def benchmark_input_scf(request):
    pseudos = abidata.pseudos("14si.pspnc", "6c.pspnc", "3li.pspnc", "9f.pspnc",
                              "12mg.pspnc", "8o.pspnc", "31ga.pspnc", "7n.pspnc")
    rest = MPRester()
    structure = rest.get_structure_by_material_id(request.param)
    try:
        return ebands_input(structure, pseudos, kppa=100, ecut=6).split_datasets()[0]
    except:
        #to deal with missing pseudos
        pytest.skip('Cannot create input for material {}.'.format(request.param))

 def calculate_stability(self, d):
     m = MPRester(self.mapi_key)
     functional = d["pseudo_potential"]["functional"]
     syms = ["{} {}".format(functional, l) for l in d["pseudo_potential"]["labels"]]
     entry = ComputedEntry(
         Composition(d["unit_cell_formula"]),
         d["output"]["final_energy"],
         parameters={"hubbards": d["hubbards"], "potcar_symbols": syms},
     )
     data = m.get_stability([entry])[0]
     for k in ("e_above_hull", "decomposes_to"):
         d["analysis"][k] = data[k]

    def _test_from_pm_volume(self):
        """
        MPのデータ volume をプリント
        MD 計算の初期値用
        """
        ent_list = []
        mpr = MPRester("WTxsDhRV7g2Mcbqw")
        composition = ['Fe', 'P']
        entries = mpr.get_entries_in_chemsys(composition)

        for entry in entries:
            if entry.as_dict()['data']['e_above_hull'] == 0:
                print(entry.as_dict()['data']['pretty_formula'])
                print(entry.as_dict()['data']['volume']/entry.as_dict()['data']['nsites'])

    def _test_from_pm_volume(self):
        """
        MPのデータ volume をプリント
        MD 計算の初期値用
        """
        ent_list = []
        mpr = MPRester("WTxsDhRV7g2Mcbqw")
        composition = ["Fe", "P"]
        entries = mpr.get_entries_in_chemsys(composition)

        for entry in entries:
            if entry.as_dict()["data"]["e_above_hull"] == 0:
                print(entry.as_dict()["data"]["pretty_formula"])
                print(entry.as_dict()["data"]["volume"] / entry.as_dict()["data"]["nsites"])

    def test_parse_criteria(self):
        crit = MPRester.parse_criteria("mp-1234 Li-*")
        self.assertIn("Li-O", crit["$or"][1]["chemsys"]["$in"])
        self.assertIn({"task_id": "mp-1234"}, crit["$or"])

        crit = MPRester.parse_criteria("Li2*")
        self.assertIn("Li2O", crit["pretty_formula"]["$in"])
        self.assertIn("Li2I", crit["pretty_formula"]["$in"])
        self.assertIn("CsLi2", crit["pretty_formula"]["$in"])

        crit = MPRester.parse_criteria("Li-*-*")
        self.assertIn("Li-Re-Ru", crit["chemsys"]["$in"])
        self.assertNotIn("Li-Li", crit["chemsys"]["$in"])

        comps = MPRester.parse_criteria("**O3")["pretty_formula"]["$in"]
        for c in comps:
            self.assertEqual(len(Composition(c)), 3, "Failed in %s" % c)

        chemsys = MPRester.parse_criteria("{Fe,Mn}-O")["chemsys"]["$in"]
        self.assertEqual(len(chemsys), 2)
        comps = MPRester.parse_criteria("{Fe,Mn,Co}O")["pretty_formula"]["$in"]
        self.assertEqual(len(comps), 3, comps)

        #Let's test some invalid symbols

        self.assertRaises(ValueError, MPRester.parse_criteria, "li-fe")
        self.assertRaises(ValueError, MPRester.parse_criteria, "LO2")

        crit = MPRester.parse_criteria("POPO2")
        self.assertIn("P2O3", crit["pretty_formula"]["$in"])

    def __init__(self, api_key=None):
        """
        Args:
            api_key: (str) Your Citrine API key, or None if you've set the MAPI_KEY environment variable

        Returns: None
        """
        api_key = api_key if api_key else os.environ["MAPI_KEY"]
        self.mprester = MPRester(api_key=api_key)

def get_decomp(o_chem_pot, mycomp, verbose=1):
    """Get decomposition from open phase diagram
        Args:
            o_chem_pot <float>: Oxygen chemical potential
            mycomp <pymatgen Composition>: Composition
            verbose <int>: 1 - verbose (default)
                           0 - silent
        Returns:
            decomposition string
    """        
    a = MPRester("<YOUR_MPREST_API_KEY_HERE>")
    elements = mycomp.elements
    ellist = map(str, elements)
    entries = a.get_entries_in_chemsys(ellist)
    #entries = a.get_entries_in_chemsys(['La', 'Mn', 'O', 'Fe'])
    pd = PhaseDiagram(entries)
    gppd = GrandPotentialPhaseDiagram(entries,{Element('O'): float(o_chem_pot)})
    print gppd
    #plotter = PDPlotter(gppd)
    #plotter.show()

    gppda = PDAnalyzer(gppd)
    #mychempots = gppda.get_composition_chempots(mycomp)
    #print "My chem pots:"
    #print mychempots
    mydecompgppd = gppda.get_decomposition(mycomp)
    #pdentry = PDEntry(mycomp, 0)
    #print "Decomp and energy:"
    #decompandenergy = gppda.get_decomp_and_e_above_hull(pdentry)
    #print decompandenergy
    #mydecomppd = pda.get_decomposition(mycomp)
    #print "Mn profile:"
    #mnprof= gppda.get_element_profile(Element('Mn'),mycomp)
    #print mnprof

    if verbose:
        for (entry,amount) in mydecompgppd.iteritems():
            print "%s: %3.3f" % (entry.name, amount)
            #mymurangegppd = gppda.getmu_range_stability_phase(Composition(entry.name),Element('O'))
            #print mymurangegppd
        #for (entry,amount) in mydecomppd.iteritems():
        #    print "%s: %3.3f" % (entry.name, amount)
        print ""
    return mydecompgppd

def test():
    """for test"""
    #This initializes the REST adaptor. Put your own API key in.
    mpr = MPRester("WTxsDhRV7g2Mcbqw")
    #Entries are the basic unit for thermodynamic
    #and other analyses in pymatgen.
    #This gets all entries belonging to the Ca-C-O system.
    composition = ['Fe', 'Ni', 'Si']
    entries = mpr.get_entries_in_chemsys(composition)

    # entryの内容を確認
    entry = entries[0]
    print(entry.as_dict())
    print(entry.as_dict()['data']['volume']/6)

    #With entries, you can do many sophisticated analyses,
    #like creating phase diagrams.
    pd = PhaseDiagram(entries) #pylint: disable=C0103

    #print(dir(pd))
    #print(pd.get_form_energy_per_atom(list(pd.stable_entries)[0]))

    plotter = PDPlotter(pd)

def calculate_phase_stability(args):
    #This initializes the REST adaptor.
    a = MPRester(args.api_key)

    drone = VaspToComputedEntryDrone()
    entry = drone.assimilate(args.directory)

    compat = MaterialsProjectCompatibility()
    entry = compat.process_entry(entry)

    if not entry:
        print "Calculation parameters are not consistent with Materials " + \
            "Project parameters."
        sys.exit()

    syms = [el.symbol for el in entry.composition.elements]
    #This gets all entries belonging to the relevant system.
    entries = a.get_entries_in_chemsys(syms)
    entries.append(entry)

    #Process entries with Materials Project compatibility.
    entries = compat.process_entries(entries)

    print [e.composition.reduced_formula for e in entries]

    comps_cids = None
    try:
        from bson import ObjectId
        cids = open(cids_file).read().splitlines()
        comps_cids = [(None, ObjectId(cid_str)) for cid_str in cids]
        if args.dev:
            raise ValueError('need full MPFile for on-the-fly analysis!')
    except (InvalidId, TypeError):
        from mpcontribs.io.mpfile import MPFile
        mpfile = MPFile.from_file(cids_file) # TODO cid-only read mode?
        comps_cids = [(x[0], ObjectId(x[1])) for x in mpfile.get_identifiers()]

    # set up rester and PDF
    SITE = 'http://localhost:8000'
    ENDPOINT, API_KEY = "{}/rest".format(SITE), os.environ.get('MAPI_KEY_LOC')
    mpr = MPRester(API_KEY, endpoint=ENDPOINT)
    merger = PdfFileMerger()

    # get all XAS/XMCD spectra
    if args.dev:
        port = 5000 + random.randint(0, 999)
        url = "http://127.0.0.1:{0}/pdf?mpfile={1}".format(port, cids_file)
        #threading.Timer(1.25, lambda: webbrowser.open(url)).start()
        #app.run(port=port, debug=False)
        # TODO find way to automatically append to PDF
        #pdf = pdfkit.from_url(url, False)
        #merger.append(pdf)
    else:
        doc = {}
        for comp,cid in comps_cids:
            print cid

class MPResterTest(unittest.TestCase):

    def setUp(self):
        self.rester = MPRester()

    def test_get_data(self):
        props = ["energy", "energy_per_atom", "formation_energy_per_atom",
                 "nsites", "unit_cell_formula", "pretty_formula", "is_hubbard",
                 "elements", "nelements", "e_above_hull", "hubbards",
                 "is_compatible", "task_ids",
                 "density", "icsd_ids", "total_magnetization"]
        # unicode literals have been reintroduced in py>3.2
        expected_vals = [-191.33812137, -6.833504334642858, -2.551358929370749,
                         28, {k: v for k, v in {'P': 4, 'Fe': 4, 'O': 16, 'Li': 4}.items()},
                         "LiFePO4", True, ['Li', 'O', 'P', 'Fe'], 4, 0.0,
                         {k: v for k, v in {'Fe': 5.3, 'Li': 0.0, 'O': 0.0, 'P': 0.0}.items()}, True,
                         [u'mp-601412', u'mp-19017', u'mp-796535', u'mp-797820',
                          u'mp-540081', u'mp-797269'],
                         3.4662026991351147,
                         [159107, 154117, 160776, 99860, 181272, 166815,
                          260571, 92198, 165000, 155580, 38209, 161479, 153699,
                          260569, 260570, 200155, 260572, 181341, 181342,
                          72545, 56291, 97764, 162282, 155635],
                         16.0002716]

        for (i, prop) in enumerate(props):
            if prop not in ['hubbards', 'unit_cell_formula', 'elements',
                            'icsd_ids', 'task_ids']:
                val = self.rester.get_data("mp-19017", prop=prop)[0][prop]
                self.assertAlmostEqual(expected_vals[i], val)
            elif prop in ["elements", "icsd_ids", "task_ids"]:
                self.assertEqual(set(expected_vals[i]),
                                 set(self.rester.get_data("mp-19017",
                                                          prop=prop)[0][prop]))
            else:
                self.assertEqual(expected_vals[i],
                                 self.rester.get_data("mp-19017",
                                                      prop=prop)[0][prop])

        props = ['structure', 'initial_structure', 'final_structure', 'entry']
        for prop in props:
            print(prop)
            obj = self.rester.get_data("mp-19017", prop=prop)[0][prop]
            if prop.endswith("structure"):
                self.assertIsInstance(obj, Structure)
            elif prop == "entry":
                obj = self.rester.get_data("mp-19017", prop=prop)[0][prop]
                self.assertIsInstance(obj, ComputedEntry)

        #Test chemsys search
        data = self.rester.get_data('Fe-Li-O', prop='unit_cell_formula')
        self.assertTrue(len(data) > 1)
        elements = {Element("Li"), Element("Fe"), Element("O")}
        for d in data:
            self.assertTrue(
                set(Composition(d['unit_cell_formula']).elements).issubset(
                    elements))

        self.assertRaises(MPRestError, self.rester.get_data, "Fe2O3",
                          "badmethod")

    def test_get_materials_id_from_task_id(self):
        self.assertEqual(self.rester.get_materials_id_from_task_id(
            "mp-540081"), "mp-19017")

    def test_get_materials_id_references(self):
        # nosetests pymatgen/matproj/tests/test_rest.py:MPResterTest.test_get_materials_id_references
        m = MPRester()
        data = m.get_materials_id_references('mp-123')
        self.assertTrue(len(data) > 1000)

    def test_find_structure(self):
        # nosetests pymatgen/matproj/tests/test_rest.py:MPResterTest.test_find_structure
        m = MPRester()
        ciffile = os.path.join(test_dir, 'Fe3O4.cif')
        data = m.find_structure(ciffile)
        self.assertTrue(len(data) > 1)
        s = CifParser(ciffile).get_structures()[0]
        data = m.find_structure(s)
        self.assertTrue(len(data) > 1)

    def test_get_entries_in_chemsys(self):
        syms = ["Li", "Fe", "O"]
        entries = self.rester.get_entries_in_chemsys(syms)
        elements = set([Element(sym) for sym in syms])
        for e in entries:
            self.assertIsInstance(e, ComputedEntry)
            self.assertTrue(set(e.composition.elements).issubset(elements))

    def test_get_structure_by_material_id(self):
        s1 = self.rester.get_structure_by_material_id("mp-1")
        self.assertEqual(s1.formula, "Cs1")

    def test_get_entry_by_material_id(self):
        e = self.rester.get_entry_by_material_id("mp-19017")
        self.assertIsInstance(e, ComputedEntry)
        self.assertTrue(e.composition.reduced_formula, "LiFePO4")

    def test_query(self):
        criteria = {'elements': {'$in': ['Li', 'Na', 'K'], '$all': ['O']}}
#.........这里部分代码省略.........

from pymatgen import Composition, Element
from pymatgen.matproj.rest import MPRester

__author__ = 'Anubhav Jain'
__copyright__ = 'Copyright 2014, The Materials Project'
__version__ = '0.1'
__maintainer__ = 'Anubhav Jain'
__email__ = '[email protected]'
__date__ = 'Oct 01, 2014'


def get_band_center(form):
    c = Composition(str(form))
    prod = 1.0
    for el, amt in c.get_el_amt_dict().iteritems():
        prod = prod * (Element(el).X ** amt)

    return -prod ** (1 / sum(c.get_el_amt_dict().values()))


if __name__ == '__main__':
    mpr = MPRester()
    for i in mpr.query({"elements": {"$in": ["O", "S", "Se", "Te"]}, "nelements": 2, "e_above_hull": {"$lte": 0.05}},
                       properties=['band_gap', 'formula', 'material_id', 'icsd_id']):
        c = Composition(i['formula'])
        center = get_band_center(c)
        print ','.join([str(x) for x in
                        c.formula, c.reduced_formula, center, center + i['band_gap'] / 2, center - i['band_gap'] / 2, i[
                            'material_id'], i['icsd_id']])