How to use the abipy.data function in abipy

#!/usr/bin/env python
"""
This example shows how to plot the L-projected fatbands of MgB2 
using the FATBANDS.nc files produced by abinit with prtdos 3.
See also PhysRevLett.86.4656
"""
import abipy.abilab as abilab
import abipy.data as abidata

# Open the file (alternatively one can use the shell and `abiopen.py FILE -nb`
# to open the file in a jupyter notebook
# This file has been produced on a k-path so it's not suitable for DOS calculations.
fbnc_kpath = abilab.abiopen(abidata.ref_file("mgb2_kpath_FATBANDS.nc"))

# Print file info (dimensions, variables ...)
# Note that prtdos = 3, so LM decomposition is not available.
print(fbnc_kpath)


# Plot the k-points belonging to the path.
fbnc_kpath.ebands.kpoints.plot()

# NC files have contributions up to L=4 (g channel)
# but here we are intererested in s,p,d terms only so
# we use the optional argument lmax
lmax = 2

# Plot the electronic fatbands grouped by atomic type.
fbnc_kpath.plot_fatbands_typeview(lmax=lmax, tight_layout=True)

def build_flow(options):

    # Set working directory (default is the name of the script with '.py' removed and "run_" replaced by "flow_")
    if not options.workdir:
        options.workdir = os.path.basename(sys.argv[0]).replace(".py", "").replace("run_", "flow_")

    # Get structure from DDB file.
    ddb_path = abidata.ref_file("refs/mgo_v8t57/mgo_zpr_t57o_DS3_DDB")
    with abilab.abiopen(ddb_path) as ddb:
        structure = ddb.structure

    # Build SCF input using structure from DDB file.
    pseudos = abidata.pseudos("Ca.psp8", "O.psp8")
    scf_input = abilab.AbinitInput(structure=structure, pseudos=pseudos)

    # Set other input variables. These quantities are system-depedent.
    # Here we use parameters similar to https://docs.abinit.org/tests/v8/Input/t57.in
    scf_input.set_vars(
        nband=12,
        nbdbuf=2,
        diemac=6,
        ecut=30,                # Underconverged ecut.
        #ecut=15,
        nstep=100,
        tolvrs=1e-16,
        kptrlatt=[-2,  2,  2,   # In cartesian coordinates, this grid is simple cubic
                   2, -2,  2,
                   2,  2, -2],
    )

Spin-polarized e-bands
======================

This example shows how to plot the band structure of nickel
using the eigenvalues stored in the GSR file produced at the end of the GS run.
"""
from abipy import abilab
import abipy.data as abidata

# Open the GSR file and extract the band structure.
# (alternatively one can use the shell and `abiopen.py OUT_GSR.nc -nb`
# to open the file in a jupyter notebook.
with abilab.abiopen(abidata.ref_file("ni_666k_GSR.nc")) as ncfile:
    ni_ebands_kmesh = ncfile.ebands

with abilab.abiopen(abidata.ref_file("ni_kpath_GSR.nc")) as ncfile:
    ni_ebands_kpath = ncfile.ebands

# Energy limits in eV for plots. The pseudo contains semi-core states but
# we are not interested in this energy region. Fermi level set to zero.
elims = [-10, 2]

# Plot band energies on k-path
ni_ebands_kpath.plot(ylims=elims, title="Ni band structure")

# Compute DOS with Gaussian method.
ni_edos = ni_ebands_kmesh.get_edos()

# Plot energies on k-path + DOS
ni_ebands_kpath.plot_with_edos(ni_edos, ylims=elims,
                               title="Ni band structure + DOS")

import abipy.data as abidata

# VKS = Hartree + XC potential + sum of local part of pseudos.
with abiopen(abidata.ref_file("ni_666k_POT.nc")) as ncfile:
    vks = ncfile.vks

#vks.plot_line(point1=[0, 0, 0], point2=[0, 4, 0], num=400, title="$V_{ks}(r)$")

# Hartree potential.
with abiopen(abidata.ref_file("ni_666k_VHA.nc")) as ncfile:
    vh = ncfile.vh

vh.plot_line(point1=[0, 0, 0], point2=[0, 4, 0], num=400, title="$V_{hartree}(r)$")

# XC potential.
with abiopen(abidata.ref_file("ni_666k_VXC.nc")) as ncfile:
    vxc = ncfile.vxc

vxc.plot_line(point1=[0, 0, 0], point2=[0, 4, 0], num=400, title="$V_{xc}(r)$")

# Hartree + XC potential.
with abiopen(abidata.ref_file("ni_666k_VHXC.nc")) as ncfile:
    vhxc = ncfile.vhxc

vhxc.plot_line(point1=[0, 0, 0], point2=[0, 4, 0], num=400, title="$V_{hxc}(r)$")

vloc = vks - vhxc
vloc.plot_line(point1=[0, 0, 0], point2=[0, 4, 0], num=400, title="$V_{loc}(r)$")

foo = vhxc - vh - vxc
#foo.plot_line(point1=[0, 0, 0], point2=[0, 4, 0], num=400, title="$V_{hxc - h - xc}(r)$")

def build_flow(options):
    # Get the unperturbed structure.
    pseudos = abidata.pseudos("14si.pspnc")
    base_structure = abidata.structure_from_ucell("Si")

    ngkpt = [6, 6, 6]

    etas = [-.001, 0, +.001]
    ph_displ = np.reshape(np.zeros(3*len(base_structure)), (-1,3))
    ph_displ[0,:] = [+1, 0, 0]
    ph_displ[1,:] = [-1, 0, 0]

    # Build new structures by displacing atoms according to the phonon displacement
    # ph_displ (in cartesian coordinates). The Displacement is normalized so that
    # the maximum atomic diplacement is 1 Angstrom and then multiplied by eta.
    modifier = abilab.StructureModifier(base_structure)

    displaced_structures = modifier.displace(ph_displ, etas, frac_coords=False)

    # Initialize flow. Each workflow in the flow defines a complete BSE calculation for given eta.

r"""
Lobster COHPCAR
===============

This example shows how to analyze the COHPCAR file
produced by Lobster code 

Use `abiopen.py FILE` with --expose or --print for a command line interface
and --notebook to generate a jupyter notebook.
"""
import os
import abipy.data as abidata

from abipy.abilab import abiopen

dirpath = os.path.join(abidata.dirpath, "refs", "lobster_gaas")
filename = os.path.join(dirpath, "GaAs_COHPCAR.lobster.gz")

# Open the COHPCAR.lobster file (same API for COOPCAR.lobster)
cohp_file = abiopen(filename)
print(cohp_file)

# Plot COHP.
cohp_file.plot(title="GaAs COHP")

# Plot integrated COHP.
cohp_file.plot(what="i", title="GaAs integrated COHP")

# Plot total overlap for all sites listed in `from_site_index`
cohp_file.plot_site_pairs_total(from_site_index=[0], title="COHP total overlap for site index 0")

# Plot partial crystal orbital projections for all sites listed in `from_site_index`

def build_flow(options):
    # Working directory (default is the name of the script with '.py' removed and "run_" replaced by "flow_")
    if not options.workdir:
        options.workdir = os.path.basename(sys.argv[0]).replace(".py", "").replace("run_", "flow_")

    structure = abidata.structure_from_ucell("MgB2")

    # Get pseudos from a table.
    table = abilab.PseudoTable(abidata.pseudos("12mg.pspnc", "5b.pspnc"))
    pseudos = table.get_pseudos_for_structure(structure)

    flow = flowtk.Flow(workdir=options.workdir)

    # Build work of GS task. Each gs_task uses different (ngkpt, tsmear) values
    # and represent the starting point of the phonon works.
    scf_work = flowtk.Work()
    ngkpt_list = [[4, 4, 4], [8, 8, 8]] #, [12, 12, 12]]
    tsmear_list = [0.01, 0.02] # , 0.04]
    for ngkpt in ngkpt_list:
        for tsmear in tsmear_list:
            scf_input = make_scf_input(structure, ngkpt, tsmear, pseudos)
            scf_work.register_scf_task(scf_input)

#!/usr/bin/env python
r"""
Wavefunction file
=================

This example shows how to analyze the wavefunctions
stored in the WFK.nc file.
"""
from abipy.abilab import abiopen
import abipy.data as abidata

# Open the DEN.nc file
ncfile = abiopen(abidata.ref_file("si_nscf_WFK.nc"))
print(ncfile)

# The DEN file has a `Structure` an `ElectronBands` object and wavefunctions.
#print(ncfile.structure)

# To plot the KS eigenvalues.
#ncfile.ebands.plot()

# Extract the wavefunction for the first spin, the first band and k=[0.5, 0, 0]
wave = ncfile.get_wave(spin=0, kpoint=[0.5, 0, 0], band=0)
print(wave)

# This is equivalent to
#wave = ncfile.get_wave(spin=0, kpoint=0, band=0)
# because [0.5, 0, 0] is the first k-point in the WFK file

def make_ion_inputs():
    cif_file = abidata.cif_file("si.cif")
    structure = abilab.Structure.from_file(cif_file)

    # Perturb the structure (random perturbation of 0.1 Angstrom)
    structure.perturb(distance=0.1)

    pseudos = abidata.pseudos("14si.pspnc")

    global_vars = dict(
        ecut=4,
        ngkpt=[4,4,4],
        shiftk=[0,0,0],
        nshiftk=1,
        chksymbreak=0,
        paral_kgb=0,
    )

def update_gs_input(spin_mode, kppra, gs_type, json_structure):
    if not json_structure: return structure_undefined_error()

    structure = abilab.mjson_loads(json_structure)
    pseudos = os.path.join(abidata.pseudo_dir, "14si.pspnc")

    # Build input for GS calculation
    # ecut must be specified because this pseudopotential does not provide hints for ecut.
    try:
        gs_inp = abilab.gs_input(
            structure, pseudos,
            kppa=kppra, ecut=8, spin_mode=spin_mode, smearing=None)
        gs_inp.pop_vars(("charge", "chksymbreak"))

        if gs_type == "relax":
            gs_inp.set_vars(optcell=2, ionmov=2, ecutsm=0.5, dilatmx=1.05)

        #multi = ebands_input(structure, pseudos,
        #                 kppa=kppra, nscf_nband=None, ndivsm=15,
        #                 ecut=8, pawecutdg=None, scf_nband=None, accuracy="normal", spin_mode=spin_mode,
        #                 smearing="fermi_dirac:0.1 eV", charge=None, dos_kppa=None):