How to apply corrections to SCAN energy calculations?

Materials Project Materials Project Data/API
1

Hello all, I am currently calculating the decomposition energy using the code below. At the moment, I am simply using the suggested decomposition fraction and I input self-calculated, uncorrected SCAN energies to obtain the decomposition energy. I understand that this is not necessarily correct, and I’d like to know to what extent if someone with more experience could share their insight.

I have two questions:

  1. Where can I look for the parameters used for the R2SCAN calculations/the INCAR file format? I have checked the R2SCAN calculation parameters site, but the github page is deprecated.
  2. How do I apply the R2SCAN/SCAN correcting scheme to the Python code below so that I’d only need to input the target structure’s energies?

I’d appreciate some suggestions. Thank you for your time!

!pip install pymatgen
from pymatgen.ext.matproj import MPRester
from pymatgen.analysis.phase_diagram import PhaseDiagram, PDPlotter, PDEntry
from pymatgen.entries.computed_entries import ComputedEntry # Import ComputedEntry
from pymatgen.core import Composition # Import Composition

# Set your API key
api_key = "G6pdzeKR62X3CIJvfqRbNyAzZUStiYOI"
m = MPRester(api_key)

entries = m.get_entries_in_chemsys(["La", "Ti", "O", "N"])

my_entry = ComputedEntry(Composition("LaTiO2N"), -18.83607524)
entries.append(my_entry)

decomp, e_above_hull = pd.get_decomp_and_e_above_hull(my_entry)

if e_above_hull == 0:
    print(f"{my_entry.composition.reduced_formula} is stable on the convex hull.")
else:
    print(f"{my_entry.composition.reduced_formula} is metastable with an energy above the hull of {e_above_hull:.4f} eV/atom.")
    print("It decomposes into:")
    for comp, frac in decomp.items():
        print(f"  - {comp.composition.reduced_formula}: {frac:.4f}")

# Build the phase diagram
pd = PhaseDiagram(entries)
plotter = PDPlotter(pd)
plotter.show()
2

Hey @rputra17, first thing: please treat your API key like a password and don’t share it. You can regenerate an API key on the website. If that’s not working, let me know

  1. Thanks for finding the broken link in the docs, that’s been fixed! For compatibility with current MP r2SCAN data, use the MPScanRelaxSet for relaxations, and MPScanStaticSet for single points / statics

  2. We don’t apply corrections to r2SCAN or SCAN data (these are separate functionals, keep that in mind too - only r2SCAN data is used in the MP hull). However we do adjust r2SCAN energies to merge the r2SCAN and PBE GGA hulls using the MaterialsProject2020Compatibility compatibility scheme

2 Likes
3

Dear @Aaron_Kaplan,

Right, thank you for the warning! I’ll regenerate it as soon as I can.

  1. Thank you for fixing the links! I found the MPScanRelaxSet.yaml, but I didn’t find the MPScanStaticSet, only MatPESStaticSet.yaml. Is that the one to use, but I replace/adjust the inputs for an R2SCAN run? Also, the yaml file shows the usage of PBE_.54. Should I use PBE_54 as well, or is it okay to use PBE_64 as long as I use the same pseudopotential type?
  2. I see. Is the adjustment here to directly input the r2SCAN-calculated energy into the Materials Project convex hull scheme? For example, if I would like to know how much energy-above-hull my added input has? If possible, is there a YouTube tutorial on how to apply adjustments on the r2SCAN energies?

Thank you for your time!

4
  1. MPScanStaticSet is a class from pymatgen, you can import them like this:
from pymatgen.io.vasp.sets import MPScanRelaxSet, MPScanStaticSet

You can explore how these sets are configured by giving them a structure:

vasp_input_set = MPScanRelaxSet(structure=< some pymatgen structure>)
print(vasp_input_set.incar)
print(vasp_input_set.kpoints)
print(vasp_input_set.potcar)

The MatPES static set is for generating ML training datasets

  1. No youtube tutorial for this, but you would apply it to a ComputedEntry or ComputedStructureEntry (search pymatgen for these):
from mp_api.client import MPRester

from pymatgen.entries.mixing_scheme import MaterialsProjectDFTMixingScheme

chemsys = "La-Ti-O-N"
with MPRester() as mpr:
    entries = mpr.get_entries_in_chemsys(chemsys,additional_criteria={"thermo_types": ["GGA_GGA+U"]})

corrected_entries = MaterialsProjectDFTMixingScheme(run_type_2="r2SCAN").process_entries(entries)

After that runs, you can create a PhaseDiagram with corrected_entries. To supply your own calculations, you would just add their representation as ComputedEntry or ComputedStructureEntry obejcts to entries

5

Dear @Aaron_Kaplan,

Thank you for your explanation! I managed to find the parameters used. I have two questions.

  1. In the INCAR, some of the tags used such as LMAXMIX = 6 are said to be used only for DFT+U based on the VASP wiki (no mention of SCAN or R2SCAN), and another such as ENAUG is already considered deprecated. Should it still be used for R2SCAN?
  2. How are the MAGMOM in a given entry is decided? For example, in an entry for LaTiO3 (mp-22103), for the entries:
    La4 Ti4 O12
    The MAGMOM is: 1 * 0.0 2 * 0.001 1 * 0.0 1 * 0.008 3 * 0.009000000000000001 12 * 0.0
    Is it decided based on trial and error and/or matching with experimental results?

Thank you for your assistance.

6

Dear @Aaron_Kaplan,

I have calculated using the same parameters as the code. Specifically, here is what I used to get the INCAR, POSCAR and POTCAR parameters:

# Retrieve the structure using the material ID
structure = m.get_structure_by_material_id("mp-277")

# Pass the structure object to MPScanRelaxSet and specify the POTCAR functional
vasp_input_set = MPScanRelaxSet(structure=structure, user_potcar_functional='PBE_64')

print(vasp_input_set.incar)
print(vasp_input_set.kpoints)
print(vasp_input_set.poscar)
print("K points used: ", vasp_input_set.kpoints)

# Iterate through the individual POTCARs and print the TITEL line
for potcar_single in vasp_input_set.potcar:
    print(f"TITEL: {potcar_single.keywords.get('TITEL')}")

The calculated formation energies and decomposition energies I get are different despite using the same parameters (some I commented out because I am only running ionic relaxation). I am using VASP 6.3.1 and PBE.64 POTCAR files. Why am I getting different formation energy values? Does these differences matter if I want to calculate decomposition energies and/or want to append my calculated energies into the Materials Project database to create the convex hull?

Example:
LaTiO2N (POTCARS: La, Ti_sv, O, N (I am still calculating using Ti_pv, will update once I have the numbers)
Calculated formation energy: -3.043 eV/atom; reported formation energy: -3.045 eV/atom
Calculated decomposition energy: -0.025 eV/atom; reported energy above hull: 0.035 eV/atom
LaZrO2N
Calculated formation energy: -3.101 eV/atom; reported formation energy: -3.063 eV/atom

Could I be missing the energy corrections?

Here are the recommended VASP parameters I tried:

SYSTEM = matproj-alltags

ALGO = All
EDIFF = 1e-05
EDIFFG = -0.02
ENAUG = 1360
ENCUT = 680.0
IBRION = 2
# ICHARG = 1
ISIF = 3
ISMEAR = 0
ISPIN = 2
KSPACING = 0.22
# LAECHG = True
LASPH = True
LCHARG = True
LELF = False
LMAXMIX = 6
LMIXTAU = True
# LORBIT = 11
LREAL = Auto
# LVTOT = True
LWAVE = False
MAGMOM = 20*-0.0
METAGGA = R2scan
NELM = 2000
NSW = 999
PREC = Accurate
SIGMA = 0.2

Thank you!

7

  1. The VASP documentation isn’t always up to date. LMAXMIX is important for any DFT calculation from our own benchmarks. ENAUG is fine to include, it’s ignored anyways because of the PREC setting

  2. Read the docstring for the input sets to see how the MAGMOM tag is determined

The MPScanRelaxSet uses the PBE.54 POTCARs but you should get similar results. To calculate formation energy, are you using your own calculations or the energies from MP?

1 Like
8

I used my own calculations to calculate the formation energy using the parameters above.
By energies from MP, which energies are you referring to? Because calculation of formation energy uses total energy which is not in MP, no?

Thank you for your time!

9

The total energies for all unaries (structures containing one element but possibly many sites) are in MP or we wouldn’t be able to calculate the formation energy. See our docs

You can retrieve these unary structures from our API:

from mp_api.client import MPRester

with MPRester() as mpr:
    docs = mpr.materials.summary.search(num_elements=1,energy_above_hull=(0,1e-6))

and then filter the structures by composition / energy above the hull. You can either then use the r2SCAN total energies without correction from MP, or recalculate them.