Problem with finding band structures of the crystal

  • so I have a crystal now I want to perform dft on that crystal and find the energy, formation energy, bandgap and band structures. idea is to validate the generated crystals

so the cif I am using is

# generated using pymatgen
data_Nd3Al11
_symmetry_space_group_name_H-M   Immm
_cell_length_a   4.36307790
_cell_length_b   9.97432892
_cell_length_c   12.95781714
_cell_angle_alpha   90.00000000
_cell_angle_beta   90.00000000
_cell_angle_gamma   90.00000000
_symmetry_Int_Tables_number   71
_chemical_formula_structural   Nd3Al11
_chemical_formula_sum   'Nd6 Al22'
_cell_volume   563.90831666
_cell_formula_units_Z   2
loop_
 _symmetry_equiv_pos_site_id
 _symmetry_equiv_pos_as_xyz
  1  'x, y, z'
  2  '-x, -y, -z'
  3  '-x, -y, z'
  4  'x, y, -z'
  5  'x, -y, -z'
  6  '-x, y, z'
  7  '-x, y, -z'
  8  'x, -y, z'
  9  'x+1/2, y+1/2, z+1/2'
  10  '-x+1/2, -y+1/2, -z+1/2'
  11  '-x+1/2, -y+1/2, z+1/2'
  12  'x+1/2, y+1/2, -z+1/2'
  13  'x+1/2, -y+1/2, -z+1/2'
  14  '-x+1/2, y+1/2, z+1/2'
  15  '-x+1/2, y+1/2, -z+1/2'
  16  'x+1/2, -y+1/2, z+1/2'
loop_
 _atom_site_type_symbol
 _atom_site_label
 _atom_site_symmetry_multiplicity
 _atom_site_fract_x
 _atom_site_fract_y
 _atom_site_fract_z
 _atom_site_occupancy
  Nd  Nd0  4  0.00000000  0.00000000  0.31720041  1
  Nd  Nd1  2  0.00000000  0.00000000  0.00000000  1
  Al  Al2  8  0.00000000  0.27452085  0.13593955  1
  Al  Al3  8  0.00000000  0.36810718  0.33399153  1
  Al  Al4  4  0.00000000  0.21533047  0.50000000  1
  Al  Al5  2  0.00000000  0.50000000  0.00000000  1

now code I am using for energy and formation energy is

# 0_run.py

import os
import shutil
from custodian.vasp.jobs import VaspJob
from pymatgen.io.vasp.sets import MPStaticSet
from pymatgen.core.structure import Structure
from custodian.custodian import Custodian
from custodian.vasp.handlers import VaspErrorHandler, UnconvergedErrorHandler, NonConvergingErrorHandler, \
    PotimErrorHandler, DriftErrorHandler, AliasingErrorHandler, PositiveEnergyErrorHandler, MeshSymmetryErrorHandler
from pymatgen.io.vasp.outputs import Vasprun
import json

def run_dft_calculation(structure_file, vasp_cmd, directory="./dft"):
    """
    Set up and run a DFT calculation for a given structure file.

    Args:
        structure_file (str): Path to the input structure file (CIF).
        vasp_cmd (list): Command to run VASP as a list of args.
        directory (str): Directory to run the calculations.

    Returns:
        dict: A dictionary containing the calculation results.
    """
    structure = Structure.from_file(structure_file)
    if os.path.exists(directory):
        shutil.rmtree(directory)
    os.makedirs(directory)

    # Set up DFT calculation
    dft_set = MPStaticSet(structure, 
                          user_incar_settings={
                              'SYSTEM': 'qvasp',
                              'PREC': 'Normal',
                              'ENCUT': 400,
                              'NELMIN': 5,
                              'LREAL': 'Auto',
                              'EDIFF': 1E-6,
                              'ISMEAR': 1,
                              'SIGMA': 0.2,
                              'ISPIN': 1,
                              'NCORE': 4,
                              'ICHARG': 2,
                              'EDIFFG': 0.01,
                              'NSW': 300,
                              'IBRION': 2,
                              'POTIM': 0.5,
                              'ISIF': 2,
                              'PSTRESS': 0,
                              'LCHARG': False,
                              'LWAVE': False,
                              'ISYM': 0,
                              'SYMPREC': 1E-10
                          },
                          user_kpoints_settings={'reciprocal_density': 100},
                          user_potcar_settings={
                              "POTCAR_FUNCTIONAL": "PBE_54",
                              "POTCAR": {"Nd": "Nd_3", "Al": "Al"}
                          })
    dft_set.write_input(directory)

    # Change to the calculation directory
    original_dir = os.getcwd()
    os.chdir(directory)

    # Set up and run the job
    job = VaspJob(vasp_cmd, final=True, suffix=".static")
    handlers = [
        VaspErrorHandler(),
        UnconvergedErrorHandler(),
        NonConvergingErrorHandler(),
        PotimErrorHandler(),
        DriftErrorHandler(),
        AliasingErrorHandler(),
        PositiveEnergyErrorHandler(),
        MeshSymmetryErrorHandler(),
    ]

    c = Custodian(handlers, [job], max_errors=5)
    c.run()

    # Analyze results
    results = {}
    vasprun_file = "vasprun.xml.static"
    if os.path.exists(vasprun_file):
        try:
            vasprun = Vasprun(vasprun_file)
            results['final_energy'] = vasprun.final_energy
            results['bandgap'] = vasprun.get_band_structure().get_band_gap()['energy']
            print(f"Final energy: {results['final_energy']} eV")
            print(f"Bandgap: {results['bandgap']} eV")
        except Exception as e:
            print(f"Error parsing vasprun.xml: {str(e)}")
            results['error'] = f"Error parsing vasprun.xml: {str(e)}"
    else:
        print(f"{vasprun_file} not found. DFT calculation may have failed.")
        results['error'] = f"{vasprun_file} not found"

    # Copy important output files to the original directory
    output_files = ['INCAR', 'KPOINTS', 'POSCAR', 'CONTCAR', 'OUTCAR', vasprun_file]
    for file in output_files:
        if os.path.exists(file):
            shutil.copy(file, original_dir)

    # Change back to the original directory
    os.chdir(original_dir)

    return results

def main():
    # Set up environment variables
    os.environ['VASP_PP_PATH'] = os.path.abspath('./PMG_VASP_PSP_DIR/POT_GGA_PAW_PBE')

    # Define VASP command
    vasp_cmd = ["mpirun", "-np", "1", "vasp_gpu"]

    # Run DFT calculation on Nd3Al11.cif
    cif_file = "Nd3Al11.cif"
    try:
        results = run_dft_calculation(cif_file, vasp_cmd)

        # Save results to a JSON file
        with open('dft_results.json', 'w') as f:
            json.dump(results, f, indent=2)
        print("Results saved to dft_results.json")
    except Exception as e:
        print(f"Error during DFT calculation: {str(e)}")

if __name__ == "__main__":
    main()

and for band structures I am using this code

import os
import shutil
from custodian.vasp.jobs import VaspJob
from pymatgen.io.vasp.sets import MPStaticSet, MPNonSCFSet
from pymatgen.core.structure import Structure
from custodian.custodian import Custodian
from custodian.vasp.handlers import VaspErrorHandler, UnconvergedErrorHandler, NonConvergingErrorHandler, \
    PotimErrorHandler, DriftErrorHandler, AliasingErrorHandler, PositiveEnergyErrorHandler, MeshSymmetryErrorHandler
from pymatgen.io.vasp.outputs import Vasprun
from pymatgen.electronic_structure.plotter import BSPlotter
import json
import matplotlib.pyplot as plt
from pymatgen.electronic_structure.core import Spin

def run_dft_calculation(structure_file, vasp_cmd, directory="./dft"):
    structure = Structure.from_file(structure_file)
    if os.path.exists(directory):
        shutil.rmtree(directory)
    os.makedirs(directory)

    dft_set = MPStaticSet(structure, 
                          user_incar_settings={
                              'SYSTEM': 'qvasp',
                              'PREC': 'Accurate',
                              'ENCUT': 520,
                              'EDIFF': 1E-6,
                              'ISMEAR': 0,
                              'SIGMA': 0.05,
                              'ISPIN': 2,  # Enable spin polarization
                              'LORBIT': 11,
                              'LCHARG': True,
                              'LWAVE': False,
                          },
                          user_kpoints_settings={'reciprocal_density': 100},
                          user_potcar_settings={
                              "POTCAR_FUNCTIONAL": "PBE_54",
                              "POTCAR": {"Nd": "Nd_3", "Al": "Al"}
                          })
    dft_set.write_input(directory)

    original_dir = os.getcwd()
    os.chdir(directory)

    job = VaspJob(vasp_cmd, final=True, suffix=".static")
    handlers = [
        VaspErrorHandler(),
        UnconvergedErrorHandler(),
        NonConvergingErrorHandler(),
        PotimErrorHandler(),
        DriftErrorHandler(),
        AliasingErrorHandler(),
        PositiveEnergyErrorHandler(),
        MeshSymmetryErrorHandler(),
    ]

    c = Custodian(handlers, [job], max_errors=5)
    c.run()

    results = {}
    vasprun_file = "vasprun.xml.static"
    if os.path.exists(vasprun_file):
        try:
            vasprun = Vasprun(vasprun_file)
            results['final_energy'] = vasprun.final_energy
            bs = vasprun.get_band_structure()
            results['bandgap'] = bs.get_band_gap()['energy']
            print(f"Final energy: {results['final_energy']} eV")
            print(f"Bandgap: {results['bandgap']} eV")
        except Exception as e:
            print(f"Error parsing vasprun.xml: {str(e)}")
            results['error'] = f"Error parsing vasprun.xml: {str(e)}"
    else:
        print(f"{vasprun_file} not found. DFT calculation may have failed.")
        results['error'] = f"{vasprun_file} not found"

    os.chdir(original_dir)
    return results

def run_band_structure_calculation(structure_file, vasp_cmd, directory="./bands"):
    structure = Structure.from_file(structure_file)
    if os.path.exists(directory):
        shutil.rmtree(directory)
    os.makedirs(directory)

    bs_set = MPNonSCFSet(structure, 
                         user_incar_settings={
                             'ICHARG': 11,
                             'LORBIT': 11,
                             'LWAVE': False,
                             'LCHARG': False,
                             'ISMEAR': 0,
                             'SIGMA': 0.05,
                             'ISPIN': 2,  # Enable spin polarization
                         },
                         mode="Line")
    bs_set.write_input(directory)

    original_dir = os.getcwd()
    os.chdir(directory)

    job = VaspJob(vasp_cmd, final=True, suffix=".bands")
    handlers = [
        VaspErrorHandler(),
        UnconvergedErrorHandler(),
        NonConvergingErrorHandler(),
        PotimErrorHandler(),
        DriftErrorHandler(),
        AliasingErrorHandler(),
        PositiveEnergyErrorHandler(),
        MeshSymmetryErrorHandler(),
    ]

    c = Custodian(handlers, [job], max_errors=5)
    c.run()

    results = {}
    vasprun_file = "vasprun.xml.bands"
    if os.path.exists(vasprun_file):
        try:
            vasprun = Vasprun(vasprun_file, parse_projected_eigen=True)
            bs = vasprun.get_band_structure(line_mode=True)
            results['bandgap'] = bs.get_band_gap()['energy']
            print(f"Bandgap from band structure: {results['bandgap']} eV")

            # Plot spin-polarized band structure
            plt.figure(figsize=(10, 6))
            plot_spin_polarized_band_structure(bs)
            plt.tight_layout()
            plt.savefig('spin_polarized_band_structure.png', dpi=300)
            print("Spin-polarized band structure plot saved as 'spin_polarized_band_structure.png'")
        except Exception as e:
            print(f"Error analyzing band structure: {str(e)}")
            results['error'] = f"Error analyzing band structure: {str(e)}"
    else:
        print(f"{vasprun_file} not found. Band structure calculation may have failed.")
        results['error'] = f"{vasprun_file} not found"

    os.chdir(original_dir)
    return results

def plot_spin_polarized_band_structure(bs):
    plt.figure(figsize=(10, 6))

    # Plot spin up bands
    for i in range(bs.nb_bands):
        plt.plot(range(len(bs.kpoints)), [e[i] for e in bs.bands[Spin.up]], 'b-', label='Spin up' if i == 0 else "")

    # Plot spin down bands
    if bs.is_spin_polarized:
        for i in range(bs.nb_bands):
            plt.plot(range(len(bs.kpoints)), [e[i] for e in bs.bands[Spin.down]], 'r--', label='Spin down' if i == 0 else "")

    plt.xlabel("k-points")
    plt.ylabel("E - E$_f$ (eV)")
    plt.title(f"Band Structure of {bs.structure.composition.reduced_formula}")

    # Set x-ticks to high symmetry k-points
    special_k_points = [k.label for k in bs.kpoints if k.label]
    special_k_points_indices = [i for i, k in enumerate(bs.kpoints) if k.label]
    plt.xticks(special_k_points_indices, special_k_points)

    # Draw vertical lines at high symmetry points
    for x in special_k_points_indices:
        plt.axvline(x=x, color='k', linestyle=':', linewidth=0.5)

    plt.axhline(y=0, color='k', linestyle='--', linewidth=0.5)  # Fermi level
    plt.legend()
    plt.grid(True, linestyle=':', alpha=0.5)

def main():
    os.environ['VASP_PP_PATH'] = os.path.abspath('./PMG_VASP_PSP_DIR/POT_GGA_PAW_PBE')
    vasp_cmd = ["mpirun", "-np", "1", "vasp_gpu"]
    cif_file = "Nd3Al11.cif"

    try:
        print("Running DFT calculation...")
        dft_results = run_dft_calculation(cif_file, vasp_cmd)

        print("\nRunning band structure calculation...")
        bs_results = run_band_structure_calculation(cif_file, vasp_cmd)

        results = {**dft_results, **bs_results}

        with open('dft_results.json', 'w') as f:
            json.dump(results, f, indent=2)
        print("Results saved to dft_results.json")
    except Exception as e:
        print(f"Error during calculations: {str(e)}")

if __name__ == "__main__":
    main()

and this is the error I am getting

==========================================
SLURM_JOB_ID = 1196832
SLURM_NODELIST = gnode055
SLURM_JOB_GPUS = 2
==========================================
chdir
0_run.py 1_run.py bands chemPotMP.json dft dft.log dft_results.json INCAR Nd3Al11.cif plotter_mod.py PMG_VASP_PSP_DIR run.sh
/home2/harsha.vasamsetti/miniconda3/envs/slices/lib/python3.9/site-packages/pymatgen/io/vasp/sets.py:379: BadInputSetWarning: Overriding POTCARs is generally not recommended as it significantly affect the results of calculations and compatibility with other calculations done with the same input set. In many instances, it is better to write a subclass of a desired input set and override the POTCAR in the subclass to be explicit on the differences.
  warnings.warn(
/home2/harsha.vasamsetti/miniconda3/envs/slices/lib/python3.9/site-packages/pymatgen/io/vasp/sets.py:567: BadInputSetWarning: Relaxation of likely metal with ISMEAR < 1 detected. See VASP recommendations on ISMEAR for metals.
  warnings.warn(
/home2/harsha.vasamsetti/miniconda3/envs/slices/lib/python3.9/site-packages/pymatgen/io/vasp/sets.py:118: BadInputSetWarning: POTCAR data with symbol Nd_3 is not known by pymatgen to correspond with the selected user_potcar_functional='PBE'. This POTCAR is known to correspond with functionals ['PBE_52', 'PBE_54']. Please verify that you are using the right POTCARs!
  warnings.warn(
/home2/harsha.vasamsetti/miniconda3/envs/slices/lib/python3.9/site-packages/pymatgen/io/vasp/sets.py:118: BadInputSetWarning: POTCAR data with symbol Al is not known by pymatgen to correspond with the selected user_potcar_functional='PBE'. This POTCAR is known to correspond with functionals ['PBE_52', 'PBE_54']. Please verify that you are using the right POTCARs!
  warnings.warn(
/home2/harsha.vasamsetti/miniconda3/envs/slices/lib/python3.9/site-packages/pymatgen/symmetry/kpath.py:177: UserWarning: The input structure does not match the expected standard primitive! The path may be incorrect. Use at your own risk.
  warn(
/home2/harsha.vasamsetti/miniconda3/envs/slices/lib/python3.9/site-packages/pymatgen/io/vasp/sets.py:567: BadInputSetWarning: Relaxation of likely metal with ISMEAR < 1 detected. See VASP recommendations on ISMEAR for metals.
  warnings.warn(
/home2/harsha.vasamsetti/miniconda3/envs/slices/lib/python3.9/site-packages/pymatgen/io/vasp/sets.py:118: BadInputSetWarning: POTCAR data with symbol Nd_3 is not known by pymatgen to correspond with the selected user_potcar_functional='PBE'. This POTCAR is known to correspond with functionals ['PBE_52', 'PBE_54']. Please verify that you are using the right POTCARs!
  warnings.warn(
/home2/harsha.vasamsetti/miniconda3/envs/slices/lib/python3.9/site-packages/pymatgen/io/vasp/sets.py:118: BadInputSetWarning: POTCAR data with symbol Al is not known by pymatgen to correspond with the selected user_potcar_functional='PBE'. This POTCAR is known to correspond with functionals ['PBE_52', 'PBE_54']. Please verify that you are using the right POTCARs!
  warnings.warn(
Running DFT calculation...
Final energy: -121.00338703 eV eV
Bandgap: 0.0 eV

Running band structure calculation...
Error analyzing band structure: no element found: line 780, column 0
Results saved to dft_results.json

the INCAR file I am using is

SYSTEM = qvasp

# Basic settings
PREC = Normal
ENCUT = 400
NELMIN = 5
LREAL = AUTO
EDIFF = 1E-6
ISMEAR = -5
SIGMA = 0.1
NCORE = 4

# Electronic structure calculation
ISTART = 1
ICHARG = 11
ISPIN = 2  # Changed to 2 for spin-polarized calculation
LORBIT = 11
NEDOS = 1000
EMIN = -5
EMAX = 5
# NBANDS will be set automatically by VASP

# Ionic Relaxation (for static calculation, these can be commented out)
# ICHARG = 2
# EDIFFG = 0.01
# NSW = 300
# IBRION = 2
# POTIM = 0.5
# ISIF = 2
# PSTRESS = 0

# Output settings
LCHARG = .TRUE.
LWAVE = .FALSE.

# Symmetry settings
ISYM = 0
SYMPREC = 1E-10

# LDA+U settings (uncomment if needed)
# LDAU = .TRUE.
# LDAUTYPE = 2
# LDAUL = 2 -1
# LDAUU = 3 0
# LDAUJ = 0 0

# Magnetism settings (uncomment if needed)
# MAGMOM = 3*1

can some one suggest me modifications to make this correct as I am new to vasp dft calculations I am not understanding where I am doing wrong