import pivotpy as pp
ps = pp.download_structure('GaAs',max_sites=8,min_sites=2)
poscar = pp.POSCAR(data = ps[-1].export_poscar())
poscar.write() # writes POSCAR file to stdout if outfile is not specified
If you are not happy with direction of basis, rotate it. Infact you can do many operations, see a list by dir(poscar)
poscar.rotate(45,[0,0,1]).write() # rotate POSCAR by 45 degrees around z-axis
poscar.splot_lat(eqv_sites=True,tol=0.35)
KPOINTS
Provided that you export a POSCAR after relaxation, you can create KPOINTS based on that POSCAR as well (specifically for fermi surface and DOS), KPOINTS for bandstructure are just a relative thing. Good thing about this is KPOINTS labels and indices will be kept in header for post processing.
Tip: Use pivotpy.KPathApp to trace path visually by selecting points on figure.
poscar.get_kmesh(2,2,1,cartesian=True) # KPOINTS are made on a rectangular grid for later easy processing
pp.str2kpath(
'''
0 0 0 ! G
1 0 0 ! X
1/2 1/2 0 ! M 3
3/8 1/8 0 ! K
''',
n = 10
)
# Note 3 after M, Only two kpoints are generated for that patch, first one got 10*distance(G,X) from n = 10, so total 13
INCAR, POTCAR, job.sh and Calculation
These files should be created by user in terminal.
During structure relaxtion, you change CONTCAR -> POSCAR until required accuracy acheived. At end, do not forget to export POSCAR using pp.POSCAR(path) and write KPOINTS file as well based on that before calculating bandstructure/DOS/Fermi Surface etc.
Post Processing
After your calculations, you are back to python for plotting etc. Use your own data, I will show examples from my data only.
Tip: Use pivotpy.VasprunApp to do output analysis visually as a pro.
Bandstructure
Tip : In case you are dealing with plenty of files, your system's memory may not withstand with a lot of data. For such cases, Vasprun and LOCPOT can be used as context managers to release memeory once you exit context manager. e.g.
with pp.Vasprun() as vr:
vr.splot_bands()
# vr is no more
with pp.LOCPOT() as locpot:
locpot.splot()
# locpot is no more
vr = pp.Vasprun(r'E:\Research\graphene_example\ISPIN_1\bands\vasprun.xml')
ax1, ax2 = pp.get_axes((6,3),ncols=2, sharey=True)
vr.splot_bands(ax=ax1,elim=[-15,10])
diff = vr.get_en_diff(3,8)
vr.splot_en_diff(diff.coords,ax1,color='r',lw=0.7)
ax1.add_text(*diff.coords.mean(axis=0),f'ΔE = {diff.de:>6.3f} eV',colors='r',transform=False,rotation=30)
vr.splot_rgb_lines(ax=ax2, elim=[-15,10])
vr1 = pp.Vasprun(r'E:\Research\graphene_example\ISPIN_2\dos\vasprun.xml', dos_only=True)
query_data = {'s':[[0,1],[0]], 'p':[[0,1],[1,2,3]], 'd':[[0,1],range(4,9)]}
ax = vr1.splot_dos_lines(query_data=query_data,elim=[-5,5])
ax.set(ylim=[-0.5,0.5],ylabel='Density of States',xlabel='Energy (eV)')
fig = vr1.iplot_dos_lines(query_data=query_data,elim=[-5,5]) # This appears later
fig.layout.xaxis.title = 'Energy (eV)'
fig.layout.yaxis.title = 'Density of States'
fig.layout.yaxis.range = [-0.5,0.5]
pp.iplot2html(fig,modebar=False) # This appears first
poscar = vr.poscar
ax1, ax2 = pp.get_axes((6,3),ncols=2)
poscar.splot_lat(eqv_sites=True,tol=0.4,plane='xy',ax = ax1)
poscar.splot_bz(ax=ax2, plane='xy',vectors=False).set_axis_off()
poscar.splot_kpath(vertex=5,knn_inds=[0,5,3,0],labels=['','M\n','K\n','G\n'],color='r')
Spin Texture
Although the DOS calculations are not intended to get spin texture/fermi surface, just to have an idea. Default values of parameters for SpinDataFrame are band = 0, elements =[[0],], orbs = [[0],] which you can tweak.
df = pp.SpinDataFrame(path = r'E:\Research\graphene_example\ISPIN_2\dos\vasprun.xml')
df.describe() # This is data for band 0, ion 0 and orbital 0
Since data is in positive y-plane only and we can translate back to equivalent kpoints, we can make it complete.
df_m = df.copy()
df_m[['kx','ky','kz']] = -df[['kx','ky','kz']]
df1 = df.append(df_m)
df.send_metadata(df1) # Enforce metadata to be the same
df1[['kx','ky']] = df1[['kx','ky']].to_numpy() - df1[['kx','ky']].to_numpy().mean(axis=0) # Shift k-points to center, mesh was bad
df2 = df1.masked('eu',-15.5,0.04,n=100,band=1,method='linear')
df2
ax1, ax2 = pp.get_axes((6,3),ncols=2,wspace= 0.5)
df1.poscar.splot_bz(plane='xy',vectors=False,ax=ax1)
df1.poscar.splot_bz(plane='xy',vectors=False,ax=ax2)
df1.splot('kx','ky','eu',ax=ax1,cmap = 'jet',s=3,marker='H')
df1.colorbar(nticks=4) # Only after an splot
ax1.set_title('$E_\\uparrow$')
df1.splot('kx','ky','eu',arrows=['','su_0','eu'], every=15,quiver_kws = dict(cmap='jet',scale=5,linewidth=8),ax=ax2)
#df1.colorbar()
ax2.set_title('$S_\\uparrow (C_1-s)$, color shows $E_\\uparrow$')
df['kz'] = 0
ax = df1.splot3d('kx','ky','eu',s=3)
df1.splot3d('kx','ky','eu',arrows=['','','su_0'],norm=5,every=11,quiver_kws= dict(arrowstyle='wedge',cmap='inferno',zorder=-5),ax=ax)
df1.colorbar()
ax.set_title('$E_\\uparrow $ with $S_\\uparrow (C_1-s)$ in 3D')
ax = df2.poscar.splot_bz(plane='xy',vectors=False)
df2.splot('kx','ky','eu',ax=ax,color='b',s=0.5)
df2.splot('kx','ky','eu',arrows=['','su_0','su_0'], every=7,quiver_kws = dict(cmap='inferno'),ax=ax)
df2.colorbar(digits=3)
ax.set_title('$S_\\uparrow (C_1-s), -15.46 < E_\\uparrow < 15.54$')
There is plenty of options for plotting in matplotlib, plotly, myavi etc, so having a SpinDataFrame ready, you can plot 2D/3D surfaces/spin textures the way you want. Making a way fixed is not a good idea. For example, try matplotlib.pyplot.tricontourf instead of above scattered plot. For 3D surfaces, try plotly Mesh, Isosurface etc and whatever suits you best.
import matplotlib.tri as tri
import numpy as np
ax = pp.get_axes()
data = df1.get_data().band_1
KPTS = np.array([data.kx,data.ky,data.kz]).T
KPTS[:,2] = 0
df1.poscar.splot_bz(plane='xy',ax=ax)
tri1 = tri.Triangulation(*KPTS.T[:2])
ax.tricontourf(tri1, data.eu, levels=np.linspace(-20, -10, 10), cmap='plasma')
ax.add_colorbar('plasma',ticks=np.linspace(-20,-10,4))
sdf = pp.SpinDataFrame(path = r'E:\Research\cart\vasprun.xml')
ax = sdf.splot('kx','ky','e',s=200,marker='s')
sdf.poscar.splot_bz(plane='xy',ax =ax)
import matplotlib.pyplot as plt
image = sdf['e'].to_numpy().reshape((8,8))
ax = sdf.poscar.splot_bz(plane='xy')
XYZ = sdf.poscar.bring_in_bz(sdf[['kx','ky','kz']].to_numpy(),sys_info = sdf.sys_info)
x1, y1, _ = XYZ.min(axis=0) - 0.03 # imshow do not scale properly with coordinates
x2, y2, _ = XYZ.max(axis=0) + 0.03
print(x1,x2,y1,y2)
# Infact limits for hexagonal lattice should be [-0.25,0.25] and [-0.3333,0.3333] in units of 2pi
# But grid is not exactly terminated at the edges of the hexagonal lattice
ax.imshow(image,cmap='hot',interpolation='nearest', extent=[x1,x2,y1,y2])