Plot band structures of Pd intercalated Bi2Te3¶

%load_ext autoreload
%autoreload 2

from aiida import load_profile
from aiida import orm, engine
import numpy as np
import matplotlib.pyplot as plt
from masci_tools.vis.plot_methods import plot_bands

profile = load_profile()

def get_bs_data(bs):
    from masci_tools.io.parsers.hdf5 import HDF5Reader
    from masci_tools.io.parsers.hdf5.recipes import FleurBands

    with bs.outputs.banddos_calc.retrieved.open('banddos.hdf', 'rb') as _f:
        #Read in data
        with HDF5Reader(_f) as h5reader:
            data, attributes = h5reader.read(recipe=FleurBands)
            
    return data, attributes


def get_data(bs, relax, return_ticks=False):
    data, attributes = get_bs_data(bs)
    sites = relax.outputs.optimized_structure.sites
    Te_sites = [i for i,s in enumerate(sites) if 'Te' in s.kind_name]
    Bi_sites = [i for i,s in enumerate(sites) if 'Bi' in s.kind_name]
    Pd_sites = [i for i,s in enumerate(sites) if 'Pd' in s.kind_name]
    for iat in Pd_sites:
        if 'Pd' not in data:
            data['Pd:d'] = data[f'MT:{iat+1}d_up']
            data['Pd'] = data[f'MT:{iat+1}_up']
        else:
            data['Pd:d'] += data[f'MT:{iat+1}d_up']
            data['Pd'] += data[f'MT:{iat+1}_up']
    for iat in Te_sites:
        if 'Te' not in data:
            data['Te:p'] = data[f'MT:{iat+1}p_up']
            data['Te'] = data[f'MT:{iat+1}_up']
        else:
            data['Te:p'] += data[f'MT:{iat+1}p_up']
            data['Te'] += data[f'MT:{iat+1}_up']
    for iat in Bi_sites:
        if 'Bi' not in data:
            data['Bi:p'] = data[f'MT:{iat+1}p_up']
            data['Bi'] = data[f'MT:{iat+1}_up']
        else:
            data['Bi:p'] += data[f'MT:{iat+1}p_up']
            data['Bi'] += data[f'MT:{iat+1}_up']
            
    if return_ticks:
        nbands = attributes['nbands']
        ticks = []
        for k_index, label in zip(attributes['special_kpoint_indices'], attributes['special_kpoint_labels']):
            ticks.append((data['kpath'][(k_index * nbands) + 1], label.replace('GAMMA', '$\Gamma$')))
        
        return data, ticks
    else:
        return data
    
def plot_orb_nosoc_soc(data, data_soc, ticks, suptitle=None, fig=True, show=True, subplots=(2,2,0), Te_p_noSOC=True, Te_p_SOC=True, Pd_d_noSOC=True, Pd_d_SOC=True,
                       clims_TepSOC=[0, 0.5], clims_TepnoSOC=[0, 0.6], clims_PddSOC=[0, 0.1], clims_PddnoSOC=[0, 0.25]):
    if fig:
        plt.figure(figsize=(16,9))

    kwargs = dict(
        special_kpoints=[(i[1], i[0]) for i in ticks],
        cmap='seismic_r',
        colorbar=True,
        markersize_min=2,
        markersize_scaling=2,
        show=False,
    )

    i0 = subplots[2]

    if Te_p_noSOC:
        i0 += 1
        ax = plt.subplot(subplots[0],subplots[1],i0)
        plot_bands(
            kpath=data['kpath'],
            bands=data['eigenvalues_up'],
            size_data=data['Te:p'],
            axis=ax,
            title='Te:p, no SOC',
            limits = {'y': [-1.5, 1.5], 'color': clims_TepnoSOC},
            **kwargs
        )

    if Pd_d_noSOC:
        i0 += 1
        ax = plt.subplot(subplots[0],subplots[1],i0)
        plot_bands(
            kpath=data['kpath'],
            bands=data['eigenvalues_up'],
            size_data=data['Pd:d'],
            axis=ax,
            title='Pd:d, no SOC',
            limits = {'y': [-1.5, 1.5], 'color': clims_PddnoSOC},
            **kwargs
        )

    if Te_p_SOC:
        i0 += 1
        ax = plt.subplot(subplots[0],subplots[1],i0)
        plot_bands(
            kpath=data_soc['kpath'],
            bands=data_soc['eigenvalues_up'],
            size_data=data_soc['Te:p'],
            axis=ax,
            title='Te:p, with SOC',
            limits = {'y': [-1.5, 1.5], 'color': clims_TepSOC},
            **kwargs
        )

    if Pd_d_SOC:
        i0 += 1
        ax = plt.subplot(subplots[0],subplots[1],i0)
        plot_bands(
            kpath=data_soc['kpath'],
            bands=data_soc['eigenvalues_up'],
            size_data=data_soc['Pd:d'],
            axis=ax,
            title='Pd:d, with SOC',
            limits = {'y': [-1.5, 1.5], 'color': clims_PddSOC},
            **kwargs
        )
                         
    if suptitle is not None:
        plt.suptitle(suptitle, fontsize='xx-large')
    plt.tight_layout()
    if show:
        plt.show()

# for data import uncomment:
# !verdi archive import fleur_relax_Pd_Bi2Te3.aiida

# 1 Pd in 3x3 cell
Relax_workchain = orm.load_node(uuid = '024e7d89-473e-46dd-bee4-ba364439fc81')
bs_1Pd_3x3_soc = orm.load_node('2492ce29-62fb-4b37-8d02-a79b5819498a') #run with scf, primitive structure, with SOC
bs_1Pd_3x3 = orm.load_node('1d78cbc9-ace0-48d1-a952-a5b22332a4b0') #run with scf, primitive structure

# 2Pd atoms in 3x3 unit cell
Relax_workchain2 = orm.load_node(uuid = 'dcc785ae-fee5-40a2-a568-da301b4f3dd4')
bs_2Pd_3x3_soc = orm.load_node('ef559526-5014-47d7-8a04-dbddfa8b40eb') #run with scf, primitive structure, with SOC
bs_2Pd_3x3 = orm.load_node('1c5b9e72-4cb3-4381-858a-c9cb862050e6') #run with scf, primitive structure

# 3Pd atoms in 3x3 unit cell
Relax_workchain3 = orm.load_node(uuid = '6551418a-39fc-4006-bf2d-9f41b1860d76')
bs_3Pd_3x3_soc = orm.load_node('d3eeafe7-d744-4b9c-b98a-184b108f75fc') #run with scf, primitive structure, with SOC
bs_3Pd_3x3 = orm.load_node('1866ab35-fc11-402d-bc63-8f30143859b3') #run with scf, primitive structure

# 2Pd atoms in 3x3 unit cell, vertically above each other
Relax_workchain4 = orm.load_node(uuid = '214e2068-4aea-4d87-84d3-8596ff27bfc8')
bs_3Pd_3x3_vertical_soc = orm.load_node('ff18b9b1-3fb2-4f4a-87ac-b995754ab7b6') #run with scf, primitive structure, with SOC
bs_3Pd_3x3_vertical = orm.load_node('9c687a96-a581-4f13-9bb2-6faa2a983a5b') #run with scf, primitive structure

data2 = get_data(bs_2Pd_3x3, Relax_workchain2)
data2_soc = get_data(bs_2Pd_3x3_soc, Relax_workchain2)

data3 = get_data(bs_3Pd_3x3, Relax_workchain3)
data3_soc = get_data(bs_3Pd_3x3_soc, Relax_workchain3)

data4 = get_data(bs_3Pd_3x3_vertical, Relax_workchain4)
data4_soc = get_data(bs_3Pd_3x3_vertical_soc, Relax_workchain4)