VaspUtilities
Yuan Fang
Tue, 11/03/2020 - 05:05 am
#pragma TextEncoding = "UTF-8"
#pragma rtGlobals=3 // Use modern global access method and strict wave access.
function vaspbands(path1,fname1,path2,fname2,fermilevel)
string path1,fname1,path2,fname2
variable fermilevel //fermilevel
variable a,i,t,j,k
string str,fword,nkpts,nbnds
Open/R a path1 + ":" + fname1
for(i=0;i<5;i++)
FReadLine a, str
endfor
FReadLine a, str
i=0
t=0
do
fword=stringfromlist(i,str," ")
if (stringmatch(fword,"")!=1)
t=t+1
if(t==2)
nkpts=fword
endif
if(t==3)
nbnds=fword
endif
endif
i=i+1
while(t<3)
make/d/n=(str2num(nkpts),str2num(nbnds)) w
for(i=0;i<str2num(nkpts);i++)
FReadLine a, str
FReadLine a, str
for(j=0;j<str2num(nbnds);j++)
FReadLine a, str
k=0
t=0
do
fword=stringfromlist(k,str," ")
if (stringmatch(fword,"")!=1)
t=t+1
if(t==2)
w[i][j]=str2num(fword)
endif
endif
k=k+1
while(t<2)
endfor
endfor
for(i=0;i<str2num(nbnds);i++)
make/d/n=(str2num(nkpts)) tmp
tmp=w[p][i]-fermilevel
duplicate tmp $"wave"+num2str(i)
killwaves tmp
endfor
killwaves w
close a
Open/R a path2 + ":" + fname2
make/d/n=3 kold,knew
make/d/n=(str2num(nkpts)) kpt
do
FReadLine a, str
if(strlen(str)==0)
break
endif
if(stringmatch(str,"*k-points in units of 2pi/SCALE*")==1)
FReadLine a, str
i=0
t=0
do
fword=stringfromlist(i,str," ")
if (stringmatch(fword,"")!=1)
t=t+1
if(t==1)
kold[0]=str2num(fword)
elseif(t==2)
kold[1]=str2num(fword)
elseif(t==3)
kold[2]=str2num(fword)
endif
endif
i=i+1
while(t<3)
kpt[0]=0
for(i=1;i<str2num(nkpts);i++)
FReadLine a, str
k=0
t=0
do
fword=stringfromlist(k,str," ")
if (stringmatch(fword,"")!=1)
t=t+1
if(t==1)
knew[0]=str2num(fword)
elseif(t==2)
knew[1]=str2num(fword)
elseif(t==3)
knew[2]=str2num(fword)
endif
endif
k=k+1
while(t<3)
if(mod(i,100)==0)
kpt[i]=kpt[i-1]
else
kpt[i]=kpt[i-1]+sqrt((knew[0]-kold[0])^2+(knew[1]-kold[1])^2+(knew[2]-kold[2])^2)
endif
kold=knew
endfor
endif
while(1)
killwaves kold,knew
close a
for(i=0;i<str2num(nbnds);i++)
if(i==0)
display $"wave"+num2str(i) vs kpt
else
appendtograph $"wave"+num2str(i) vs kpt
endif
endfor
end
function vaspbands_mbj(path1,fname1,path2,fname2,fermilevel)
string path1,fname1,path2,fname2
variable fermilevel //fermilevel
variable a,i,t,j,k
variable cc=120
string str,fword,nkpts,nbnds
Open/R a path1 + ":" + fname1
for(i=0;i<5;i++)
FReadLine a, str
endfor
FReadLine a, str
i=0
t=0
do
fword=stringfromlist(i,str," ")
if (stringmatch(fword,"")!=1)
t=t+1
if(t==2)
nkpts=fword
endif
if(t==3)
nbnds=fword
endif
endif
i=i+1
while(t<3)
make/d/n=(str2num(nkpts),str2num(nbnds)) w
for(i=0;i<str2num(nkpts);i++)
FReadLine a, str
FReadLine a, str
for(j=0;j<str2num(nbnds);j++)
FReadLine a, str
k=0
t=0
do
fword=stringfromlist(k,str," ")
if (stringmatch(fword,"")!=1)
t=t+1
if(t==2)
w[i][j]=str2num(fword)
endif
endif
k=k+1
while(t<2)
endfor
endfor
for(i=0;i<str2num(nbnds);i++)
make/d/n=(str2num(nkpts)) tmp
tmp=w[p][i]-fermilevel
duplicate tmp $"wave"+num2str(i)
killwaves tmp
endfor
killwaves w
close a
Open/R a path2 + ":" + fname2
make/d/n=3 kold,knew
make/d/n=(str2num(nkpts)) kpt
do
FReadLine a, str
if(strlen(str)==0)
break
endif
if(stringmatch(str,"*k-points in units of 2pi/SCALE*")==1)
FReadLine a, str
i=0
t=0
do
fword=stringfromlist(i,str," ")
if (stringmatch(fword,"")!=1)
t=t+1
if(t==1)
kold[0]=str2num(fword)
elseif(t==2)
kold[1]=str2num(fword)
elseif(t==3)
kold[2]=str2num(fword)
endif
endif
i=i+1
while(t<3)
kpt[0]=0
for(i=1;i<str2num(nkpts);i++)
FReadLine a, str
k=0
t=0
do
fword=stringfromlist(k,str," ")
if (stringmatch(fword,"")!=1)
t=t+1
if(t==1)
knew[0]=str2num(fword)
elseif(t==2)
knew[1]=str2num(fword)
elseif(t==3)
knew[2]=str2num(fword)
endif
endif
k=k+1
while(t<3)
if(mod(i,100)==0)
kpt[i]=kpt[i-1]
else
kpt[i]=kpt[i-1]+sqrt((knew[0]-kold[0])^2+(knew[1]-kold[1])^2+(knew[2]-kold[2])^2)
endif
kold=knew
endfor
endif
while(1)
killwaves kold,knew
close a
for(i=0;i<str2num(nbnds);i++)
DeletePoints 0,cc, $"wave"+num2str(i)
endfor
DeletePoints 0,cc,kpt
for(i=0;i<str2num(nbnds);i++)
if(i==0)
display $"wave"+num2str(i) vs kpt
else
appendtograph $"wave"+num2str(i) vs kpt
endif
endfor
end
function splitdos(path1,fname1,path2,fname2) //fname1=DOSCAR,fname2=POSCAR
string path1,fname1,path2,fname2
variable a,b,i,t,j,k
string str,fword,natms,nspin
Open/R a path1 + ":" + fname1
FReadLine a, str
i=0
t=0
do
fword=stringfromlist(i,str," ")
if (stringmatch(fword,"")!=1)
t=t+1
if(t==1)
natms=fword
endif
if(t==3)
nspin=fword
endif
endif
i=i+1
while(t<3)
for(i=0;i<4;i++)
FReadLine a, str
endfor
printf "Number of atoms: %d\r Number of spin: %d\r" str2num(natms),str2num(nspin)
Open/R b path2 + ":" + fname2
for(i=0;i<5;i++)
FReadLine b, str
endfor
FReadLine b, str
i=0
variable/g nspec=0
do
fword=stringfromlist(i,str," ")
if (stringmatch(fword,"")!=1)
nspec=nspec+1
endif
i=i+1
while(i<strlen(str))
make/t/n=(nspec)/o atmnam
i=0
t=0
do
fword=stringfromlist(i,str," ")
if (stringmatch(fword,"")!=1)
t=t+1
if(t==nspec)
variable memory=strlen(fword)
string lastnm
sscanf fword[0,memory-1],"%s",lastnm
atmnam[t-1]=lastnm
else
atmnam[t-1]=cleanupname(fword,0)
endif
endif
i=i+1
while(t<nspec)
make/d/n=(nspec)/o atmnum
FReadLine b, str
i=0
t=0
do
fword=stringfromlist(i,str," ")
if (stringmatch(fword,"")!=1)
t=t+1
atmnum[t-1]=str2num(fword)
endif
i=i+1
while(t<nspec)
close b
printf "Number of species: %d\r" nspec
variable temp=0
for(i=0;i<nspec;i++)
temp+=atmnum[i]
endfor
if(str2num(natms)!=temp)
printf "DOSCAR info & POSCAR does not match!\r"
close a
return -1
endif
FReadLine a, str
i=0
t=0
do
fword=stringfromlist(i,str," ",32)
if (stringmatch(fword,"")!=1)
t=t+1
if(t==1)
variable/g nedos=str2num(fword)
endif
if(t==2)
variable/g efermi=str2num(fword)
endif
endif
i=i+1
while(t<2)
for(i=0;i<nedos;i++)
FReadLine a, str
if(i==0)
j=0
variable tl=0
do
fword=stringfromlist(j,str," ")
if (stringmatch(fword,"")!=1)
tl=tl+1
endif
j=j+1
while(j<strlen(str))
make/o/n=(nedos,tl)/d pdos_tot
endif
j=0
t=0
do
fword=stringfromlist(j,str," ")
if (stringmatch(fword,"")!=1)
t=t+1
pdos_tot[i][t-1]=str2num(fword)
endif
j=j+1
while(t<tl)
endfor
for(i=0;i<nspec;i++)
for(j=0;j<atmnum[i];j++)
FReadLine a, str
for(k=0;k<nedos;k++)
FReadLine a, str
if(k==0)
variable in=0
tl=0
do
fword=stringfromlist(in,str," ")
if (stringmatch(fword,"")!=1)
tl=tl+1
endif
in=in+1
while(in<strlen(str))
make/o/n=(nedos,tl)/d tempwave
endif
in=0
t=0
do
fword=stringfromlist(in,str," ")
if (stringmatch(fword,"")!=1)
t=t+1
tempwave[k][t-1]=str2num(fword)
endif
in=in+1
while(t<tl)
endfor
duplicate/o tempwave, $"pdos_"+atmnam[i]+"_"+num2str(j+1)
endfor
endfor
killwaves tempwave
end
function splitdosif(path1,fname1,path2,fname2)
string path1,fname1,path2,fname2
variable a,b,i,t,j,k
string str,fword,natms,nspin
Open/R a path1 + ":" + fname1
FReadLine a, str
i=0
t=0
do
fword=stringfromlist(i,str," ")
if (stringmatch(fword,"")!=1)
t=t+1
if(t==1)
natms=fword
endif
if(t==3)
nspin=fword
endif
endif
i=i+1
while(t<3)
for(i=0;i<4;i++)
FReadLine a, str
endfor
printf "Number of atoms: %d\r Number of spin: %d\r" str2num(natms),str2num(nspin)
Open/R b path2 + ":" + fname2
for(i=0;i<5;i++)
FReadLine b, str
endfor
FReadLine b, str
i=0
variable/g nspec=0
do
fword=stringfromlist(i,str," ")
if (stringmatch(fword,"")!=1)
nspec=nspec+1
endif
i=i+1
while(i<strlen(str))
make/t/n=(nspec)/o atmnam
i=0
t=0
do
fword=stringfromlist(i,str," ")
if (stringmatch(fword,"")!=1)
t=t+1
if(t==nspec)
variable memory=strlen(fword)
string lastnm
sscanf fword[0,memory-1],"%s",lastnm
atmnam[t-1]=lastnm
else
atmnam[t-1]=cleanupname(fword,0)
endif
endif
i=i+1
while(t<nspec)
make/d/n=(nspec)/o atmnum
FReadLine b, str
i=0
t=0
do
fword=stringfromlist(i,str," ")
if (stringmatch(fword,"")!=1)
t=t+1
atmnum[t-1]=str2num(fword)
endif
i=i+1
while(t<nspec)
close b
printf "Number of species: %d\r" nspec
variable temp=0
for(i=0;i<nspec;i++)
temp+=atmnum[i]
endfor
if(str2num(natms)!=temp)
printf "DOSCAR info & POSCAR does not match!\r"
close a
return -1
endif
FReadLine a, str
i=0
t=0
do
fword=stringfromlist(i,str," ",32)
if (stringmatch(fword,"")!=1)
t=t+1
if(t==1)
variable/g nedos=str2num(fword)
endif
if(t==2)
variable/g efermi=str2num(fword)
endif
endif
i=i+1
while(t<2)
for(i=0;i<nedos;i++)
FReadLine a, str
if(i==0)
j=0
variable tl=0
do
fword=stringfromlist(j,str," ")
if (stringmatch(fword,"")!=1)
tl=tl+1
endif
j=j+1
while(j<strlen(str))
make/o/n=(nedos,tl)/d pdos_tot
endif
j=0
t=0
do
fword=stringfromlist(j,str," ")
if (stringmatch(fword,"")!=1)
t=t+1
pdos_tot[i][t-1]=str2num(fword)
endif
j=j+1
while(t<tl)
endfor
for(i=0;i<nspec;i++)
for(j=0;j<atmnum[i];j++)
FReadLine a, str
for(k=0;k<(2*nedos);k++)
FReadLine a, str
if(k==0)
variable in=0
tl=0
do
fword=stringfromlist(in,str," ")
if (stringmatch(fword,"")!=1)
tl=tl+1
endif
in=in+1
while(in<strlen(str))
make/o/n=(nedos,tl)/d tempwave
endif
if(k==1)
in=0
variable tll=0
do
fword=stringfromlist(in,str," ")
if (stringmatch(fword,"")!=1)
tll=tll+1
endif
in=in+1
while(in<strlen(str))
redimension/n=(nedos,(tl+tll))/d tempwave
endif
in=0
t=0
if(mod(k,2)==0)
do
fword=stringfromlist(in,str," ")
if (stringmatch(fword,"")!=1)
t=t+1
tempwave[(k/2)][t-1]=str2num(fword)
endif
in=in+1
while(t<tl)
else
do
fword=stringfromlist(in,str," ")
if (stringmatch(fword,"")!=1)
t=t+1
tempwave[floor(k/2)][tl+t-1]=str2num(fword)
endif
in=in+1
while(t<tll)
endif
endfor
duplicate/o tempwave, $"pdos_"+atmnam[i]+"_"+num2str(j+1)
endfor
endfor
killwaves tempwave
end
function sumdos()
variable i,j
variable/g nspec,nedos
wave atmnum
wave/t atmnam
for(i=0;i<nspec;i++)
for(j=0;j<atmnum[i];j++)
if(j==0)
duplicate/o $"pdos_"+atmnam[i]+"_"+num2str(j+1),tempwave
duplicate/o/r=[][0] $"pdos_"+atmnam[i]+"_"+num2str(j+1),energywave
else
duplicate/o $"pdos_"+atmnam[i]+"_"+num2str(j+1),tempwave2
tempwave+=tempwave2
killwaves tempwave2
endif
endfor
variable k
for(k=0;k<nedos;k++)
tempwave[k][0]=energywave[k]
endfor
duplicate/o tempwave $"pdos_"+atmnam[i]+"_tot"
endfor
killwaves tempwave,energywave
end
function plotdosif()
variable i
variable/g nspec
wave/t atmnam
for(i=0;i<nspec;i++)
variable j
for(j=0;j<65;j++)
if(j==1)
duplicate/o/r=[][j] $"pdos_"+atmnam[i]+"_tot" $"pdos_"+atmnam[i]+"_tot_s"
endif
if(j==5)
duplicate/o/r=[][j] $"pdos_"+atmnam[i]+"_tot" temp_p
endif
if(j==9)
duplicate/o/r=[][j] $"pdos_"+atmnam[i]+"_tot" temp_p_t
temp_p+=temp_p_t
endif
if(j==13)
duplicate/o/r=[][j] $"pdos_"+atmnam[i]+"_tot" temp_p_t
temp_p+=temp_p_t
duplicate/o temp_p,$"pdos_"+atmnam[i]+"_tot_p"
killwaves temp_p_t,temp_p
endif
if(j==17)
duplicate/o/r=[][j] $"pdos_"+atmnam[i]+"_tot" temp_d
endif
if(j==21)
duplicate/o/r=[][j] $"pdos_"+atmnam[i]+"_tot" temp_d_t
temp_d+=temp_d_t
endif
if(j==25)
duplicate/o/r=[][j] $"pdos_"+atmnam[i]+"_tot" temp_d_t
temp_d+=temp_d_t
endif
if(j==29)
duplicate/o/r=[][j] $"pdos_"+atmnam[i]+"_tot" temp_d_t
temp_d+=temp_d_t
endif
if(j==33)
duplicate/o/r=[][j] $"pdos_"+atmnam[i]+"_tot" temp_d_t
temp_d+=temp_d_t
duplicate/o temp_d,$"pdos_"+atmnam[i]+"_tot_d"
killwaves temp_d_t,temp_d
endif
if(j==37)
duplicate/o/r=[][j] $"pdos_"+atmnam[i]+"_tot" temp_f
endif
if(j==41)
duplicate/o/r=[][j] $"pdos_"+atmnam[i]+"_tot" temp_f_t
temp_f+=temp_f_t
endif
if(j==45)
duplicate/o/r=[][j] $"pdos_"+atmnam[i]+"_tot" temp_f_t
temp_f+=temp_f_t
endif
if(j==49)
duplicate/o/r=[][j] $"pdos_"+atmnam[i]+"_tot" temp_f_t
temp_f+=temp_f_t
endif
if(j==53)
duplicate/o/r=[][j] $"pdos_"+atmnam[i]+"_tot" temp_f_t
temp_f+=temp_f_t
endif
if(j==57)
duplicate/o/r=[][j] $"pdos_"+atmnam[i]+"_tot" temp_f_t
temp_f+=temp_f_t
endif
if(j==61)
duplicate/o/r=[][j] $"pdos_"+atmnam[i]+"_tot" temp_f_t
temp_f+=temp_f_t
duplicate/o temp_f,$"pdos_"+atmnam[i]+"_tot_f"
killwaves temp_f_t,temp_f
endif
endfor
endfor
for(i=0;i<nspec;i++)
if(i==0)
display $"pdos_"+atmnam[i]+"_tot_s" vs $"pdos_"+atmnam[i]+"_tot"[][0]
else
appendtograph $"pdos_"+atmnam[i]+"_tot_s" vs $"pdos_"+atmnam[i]+"_tot"[][0]
endif
appendtograph $"pdos_"+atmnam[i]+"_tot_p" vs $"pdos_"+atmnam[i]+"_tot"[][0]
appendtograph $"pdos_"+atmnam[i]+"_tot_d" vs $"pdos_"+atmnam[i]+"_tot"[][0]
appendtograph $"pdos_"+atmnam[i]+"_tot_f" vs $"pdos_"+atmnam[i]+"_tot"[][0]
endfor
end
function plotdoslf()
variable i
variable/g nspec
wave/t atmnam
for(i=0;i<nspec;i++)
variable j
for(j=0;j<65;j++)
if(j==1)
duplicate/o/r=[][j] $"pdos_"+atmnam[i]+"_tot" $"pdos_"+atmnam[i]+"_tot_s"
endif
if(j==5)
duplicate/o/r=[][j] $"pdos_"+atmnam[i]+"_tot" temp_p
endif
if(j==9)
duplicate/o/r=[][j] $"pdos_"+atmnam[i]+"_tot" temp_p_t
temp_p+=temp_p_t
endif
if(j==13)
duplicate/o/r=[][j] $"pdos_"+atmnam[i]+"_tot" temp_p_t
temp_p+=temp_p_t
duplicate/o temp_p,$"pdos_"+atmnam[i]+"_tot_p"
killwaves temp_p_t,temp_p
endif
if(j==17)
duplicate/o/r=[][j] $"pdos_"+atmnam[i]+"_tot" temp_d
endif
if(j==21)
duplicate/o/r=[][j] $"pdos_"+atmnam[i]+"_tot" temp_d_t
temp_d+=temp_d_t
endif
if(j==25)
duplicate/o/r=[][j] $"pdos_"+atmnam[i]+"_tot" temp_d_t
temp_d+=temp_d_t
endif
if(j==29)
duplicate/o/r=[][j] $"pdos_"+atmnam[i]+"_tot" temp_d_t
temp_d+=temp_d_t
endif
if(j==33)
duplicate/o/r=[][j] $"pdos_"+atmnam[i]+"_tot" temp_d_t
temp_d+=temp_d_t
duplicate/o temp_d,$"pdos_"+atmnam[i]+"_tot_d"
killwaves temp_d_t,temp_d
endif
endfor
endfor
for(i=0;i<nspec;i++)
if(i==0)
display $"pdos_"+atmnam[i]+"_tot_s" vs $"pdos_"+atmnam[i]+"_tot"[][0]
else
appendtograph $"pdos_"+atmnam[i]+"_tot_s" vs $"pdos_"+atmnam[i]+"_tot"[][0]
endif
appendtograph $"pdos_"+atmnam[i]+"_tot_p" vs $"pdos_"+atmnam[i]+"_tot"[][0]
appendtograph $"pdos_"+atmnam[i]+"_tot_d" vs $"pdos_"+atmnam[i]+"_tot"[][0]
endfor
end
function Findbands(numbnds,efermi)
variable numbnds,efermi
variable i
make/o/n=(numbnds) mini,maxm
for(i=0;i<numbnds;i++)
loadwave/a/o/d/g/L={0,21+1030302*i,1030301,0,0} "C:\Users\Lenovo\Desktop\wannier90.bxsf"
endfor
for(i=0;i<numbnds;i++)
wavestats $"wave"+num2str(i)
mini[i]=V_min
maxm[i]=V_max
killwaves $"wave"+num2str(i)
endfor
for(i=0;i<numbnds;i++)
if((efermi>mini[i])&&(efermi<maxm[i]))
print i+1
endif
endfor
end
function extract_fs(path1,fname1,path2,fname2,nstart,nend) //fname1 read;fname2 write
string path1,fname1,path2,fname2
variable nstart,nend
variable a,b,i,t,j,k
string str,temp,fword
variable nx,ny,nz
Open/R a path1 + ":" + fname1
Open b path2 + ":" + fname2
for(i=0;i<14;i++)
FReadLine a, str
temp=str[0,strlen(str)-2]
fprintf b,temp
fprintf b,"\n"
endfor
FReadLine a, str
i=0
t=0
do
fword=stringfromlist(i,str," ")
if (stringmatch(fword,"")!=1)
t=t+1
if(t==1)
temp=fword[0,strlen(fword)-2]
endif
endif
i=i+1
while(t<1)
if((nend<nstart)||(str2num(temp)<nend))
fprintf b,"Wrong indices,exit...\n"
close a
close b
return -1
endif
fprintf b,"%15d\n",nend-nstart+1
FReadLine a, str
temp=str[0,strlen(str)-2]
fprintf b,temp
fprintf b,"\n"
i=0
t=0
do
fword=stringfromlist(i,str," ")
if (stringmatch(fword,"")!=1)
t=t+1
if(t==1)
temp=fword[0,strlen(fword)-1]
nx=str2num(temp)
endif
if(t==2)
temp=fword[0,strlen(fword)-1]
ny=str2num(temp)
endif
if(t==3)
temp=fword[0,strlen(fword)-2]
nz=str2num(temp)
endif
endif
i=i+1
while(t<3)
for(i=0;i<4;i++)
FReadLine a, str
temp=str[0,strlen(str)-2]
fprintf b,temp
fprintf b,"\n"
endfor
for(i=1;i<nstart;i++)
for(j=0;j<(nx*ny*nz+1);j++)
FReadLine a, str
endfor
endfor
for(i=0;i<(nend-nstart+1);i++)
FReadLine a, str
fprintf b,"BAND:%12d\n",i+1
for(j=0;j<(nx*ny*nz);j++)
FReadLine a, str
temp=str[0,strlen(str)-2]
fprintf b,temp
fprintf b,"\n"
endfor
endfor
fprintf b," END_BANDGRID_3D\n"
fprintf b," END_BLOCK_BANDGRID_3D\n"
close a
close b
end
function vaspbands_proj(path1,fname1,atoms,orbitals) //read procar non-soc
wave atoms
wave orbitals
string path1,fname1
variable a,i,t,j,k,u
variable s
variable inner,np,ccc
string str,fword,nkpts,nbnds,nions
Open/R a path1 + ":" + fname1
FReadLine a, str
FReadLine a, str
i=0
t=0
do
fword=stringfromlist(i,str," ")
if (stringmatch(fword,"")!=1)
t=t+1
if(t==4)
nkpts=fword
endif
if(t==8)
nbnds=fword
endif
if(t==12)
nions=fword
endif
endif
i=i+1
while(t<12)
make/o/d/n=(str2num(nkpts),str2num(nbnds)) tem
for(u=0;u<str2num(nkpts);u++)
for(i=0;i<3;i++)
FReadLine a, str
endfor
for(k=0;k<str2num(nbnds);k++)
for(i=0;i<3;i++)
FReadLine a, str
endfor
s=0
for(np=0;np<numpnts(atoms);np++)
if(np==0)
if(atoms[0]==1)
FReadLine a, str
j=0
t=0
do
fword=stringfromlist(j,str," ")
if (stringmatch(fword,"")!=1)
t=t+1
for(ccc=0;ccc<numpnts(orbitals);ccc++)
if(t==orbitals[ccc])
s+=str2num(fword)
endif
endfor
endif
j=j+1
while(t<orbitals[numpnts(orbitals)-1])
endif
if(atoms[0]>1)
for(inner=0;inner<(atoms[0]-1);inner++)
FReadLine a, str
endfor
FReadLine a, str
j=0
t=0
do
fword=stringfromlist(j,str," ")
if (stringmatch(fword,"")!=1)
t=t+1
for(ccc=0;ccc<numpnts(orbitals);ccc++)
if(t==orbitals[ccc])
s+=str2num(fword)
endif
endfor
endif
j=j+1
while(t<orbitals[numpnts(orbitals)-1])
endif
else
if((atoms[np]-atoms[np-1])==1)
FReadLine a, str
j=0
t=0
do
fword=stringfromlist(j,str," ")
if (stringmatch(fword,"")!=1)
t=t+1
for(ccc=0;ccc<numpnts(orbitals);ccc++)
if(t==orbitals[ccc])
s+=str2num(fword)
endif
endfor
endif
j=j+1
while(t<orbitals[numpnts(orbitals)-1])
endif
if((atoms[np]-atoms[np-1])>1)
for(inner=0;inner<(atoms[np]-atoms[np-1]-1);inner++)
FReadLine a, str
endfor
FReadLine a, str
j=0
t=0
do
fword=stringfromlist(j,str," ")
if (stringmatch(fword,"")!=1)
t=t+1
for(ccc=0;ccc<numpnts(orbitals);ccc++)
if(t==orbitals[ccc])
s+=str2num(fword)
endif
endfor
endif
j=j+1
while(t<orbitals[numpnts(orbitals)-1])
endif
endif
endfor
if(str2num(nions)==atoms[numpnts(atoms)-1])
FReadLine a, str
FReadLine a, str
else
for(i=0;i<(str2num(nions)-atoms[numpnts(atoms)-1]);i++)
FReadLine a, str
endfor
FReadLine a, str
FReadLine a, str
endif
tem[u][k]=s
endfor
endfor
for(i=0;i<str2num(nbnds);i++)
duplicate/o/r=[0,str2num(nkpts)-1][i,i] tem $"weight"+num2str(i)
endfor
close a
end
function vaspbands_projsoc(path1,fname1,atoms,orbitals) //read procar soc
wave atoms
wave orbitals
string path1,fname1
variable a,i,t,j,k,u
variable s
variable inner,np,ccc
string str,fword,nkpts,nbnds,nions
Open/R a path1 + ":" + fname1
FReadLine a, str
FReadLine a, str
i=0
t=0
do
fword=stringfromlist(i,str," ")
if (stringmatch(fword,"")!=1)
t=t+1
if(t==4)
nkpts=fword
endif
if(t==8)
nbnds=fword
endif
if(t==12)
nions=fword
endif
endif
i=i+1
while(t<12)
make/o/d/n=(str2num(nkpts),str2num(nbnds)) tem
for(u=0;u<str2num(nkpts);u++)
for(i=0;i<3;i++)
FReadLine a, str
endfor
for(k=0;k<str2num(nbnds);k++)
for(i=0;i<3;i++)
FReadLine a, str
endfor
s=0
for(np=0;np<numpnts(atoms);np++)
if(np==0)
if(atoms[0]==1)
FReadLine a, str
j=0
t=0
do
fword=stringfromlist(j,str," ")
if (stringmatch(fword,"")!=1)
t=t+1
for(ccc=0;ccc<numpnts(orbitals);ccc++)
if(t==orbitals[ccc])
s+=str2num(fword)
endif
endfor
endif
j=j+1
while(t<orbitals[numpnts(orbitals)-1])
endif
if(atoms[0]>1)
for(inner=0;inner<(atoms[0]-1);inner++)
FReadLine a, str
endfor
FReadLine a, str
j=0
t=0
do
fword=stringfromlist(j,str," ")
if (stringmatch(fword,"")!=1)
t=t+1
for(ccc=0;ccc<numpnts(orbitals);ccc++)
if(t==orbitals[ccc])
s+=str2num(fword)
endif
endfor
endif
j=j+1
while(t<orbitals[numpnts(orbitals)-1])
endif
else
if((atoms[np]-atoms[np-1])==1)
FReadLine a, str
j=0
t=0
do
fword=stringfromlist(j,str," ")
if (stringmatch(fword,"")!=1)
t=t+1
for(ccc=0;ccc<numpnts(orbitals);ccc++)
if(t==orbitals[ccc])
s+=str2num(fword)
endif
endfor
endif
j=j+1
while(t<orbitals[numpnts(orbitals)-1])
endif
if((atoms[np]-atoms[np-1])>1)
for(inner=0;inner<(atoms[np]-atoms[np-1]-1);inner++)
FReadLine a, str
endfor
FReadLine a, str
j=0
t=0
do
fword=stringfromlist(j,str," ")
if (stringmatch(fword,"")!=1)
t=t+1
for(ccc=0;ccc<numpnts(orbitals);ccc++)
if(t==orbitals[ccc])
s+=str2num(fword)
endif
endfor
endif
j=j+1
while(t<orbitals[numpnts(orbitals)-1])
endif
endif
endfor
variable jjj
for(jjj=0;jjj<(3*(str2num(nions)+1+1));jjj++)
FReadLine a, str
endfor
if(str2num(nions)==atoms[numpnts(atoms)-1])
FReadLine a, str
FReadLine a, str
FReadLine a, str
else
for(i=0;i<(str2num(nions)-atoms[numpnts(atoms)-1]);i++)
FReadLine a, str
endfor
FReadLine a, str
FReadLine a, str
FReadLine a, str
endif
tem[u][k]=s
endfor
endfor
for(i=0;i<str2num(nbnds);i++)
duplicate/o/r=[0,str2num(nkpts)-1][i,i] tem $"weight"+num2str(i)
endfor
close a
end
function projbulkbnds(m,n,p) //projected bulk bands plot.m:the number of points between two kpoints;
//n:number of different kzs;p:numberof bands.
variable m,n,p
variable i,j,k
wave kpt
make/n=(2*m) refer
for(i=0;i<(2*m);i++)
refer[i]=kpt[i]
endfor
make/n=(2*m) tem1
for(i=0;i<p;i++)
duplicate $"wave"+num2str(i),tem2
for(j=i;j<i+p*(n-1)+1;j=j+p)
for(k=0;k<(2*m);k++)
tem1[k]=tem2[k]
endfor
duplicate tem1,$"w"+num2str(j)
if(j==0)
display $"w"+num2str(j) vs refer
else
appendtograph $"w"+num2str(j) vs refer
endif
DeletePoints 0,3*m, tem2
endfor
killwaves tem2
endfor
killwaves tem1
end
function difkzplot(nbnds,nkz,ymin,ymax,j)
variable nbnds
variable nkz
variable ymin,ymax //display range
variable j //j can take value 0-20,corresponding to kz=j/(nkz-1)*(gamma-z)
variable kz=j/(nkz-1)
wave xwave
variable i
for(i=0;i<nbnds;i++)
if(i==0)
display $"w"+num2str(j*nbnds) vs xwave
TextBox/C/N=text0/F=0/A=MC num2str(kz)+"*G-Z"
SetAxis left ymin,ymax
else
appendtograph $"w"+num2str(j*nbnds+i) vs xwave
endif
endfor
end
function fcc111() //for writing fcc 111 KPOINTS
variable a,b,c
a= 5.429
make/o/d/n=3 a1={0,a/2,a/2}
make/o/d/n=3 a2={a/2,0,a/2}
make/o/d/n=3 a3={a/2,a/2,0}
cross a2,a3
variable vprim
wave W_Cross
vprim=matrixDot(a1,W_Cross)
make/o/d/n=3 b1,b2,b3
b1=2*pi*W_Cross/vprim
cross a3,a1
b2=2*pi*W_Cross/vprim
cross a1,a2
b3=2*pi*W_Cross/vprim
Concatenate/O {b1,b2,b3},primbasis
make/o/d/n=3 aa1={-a/2,a/2,0}
make/o/d/n=3 aa2={0,-a/2,a/2}
make/o/d/n=3 aa3={a,a,a}
cross aa2,aa3
variable vconv
vconv=matrixDot(aa1,W_Cross)
make/o/d/n=3 bb1,bb2,bb3
bb1=2*pi*W_Cross/vconv
cross aa3,aa1
bb2=2*pi*W_Cross/vconv
bb3=1/2*(b1+b2+b3)
Concatenate/O {bb1,bb2,bb3},convbasis
matrixinverse primbasis
wave M_inverse
matrixop/o matC=M_inverse x convbasis
print vprim
print vprim/vconv
cross b2,b3
variable kvprim
kvprim=matrixDot(b1,W_Cross)
cross bb2,bb3
variable kvconv
kvconv=matrixDot(bb1,W_Cross)
print kvprim/kvconv
end
function simplecubic111() //for writing simple cubic 111 KPOINTS
variable a,b,c
a= 2.3537117117709037*2
make/o/d/n=3 a1={a,0,0}
make/o/d/n=3 a2={0,a,0}
make/o/d/n=3 a3={0,0,a}
cross a2,a3
variable vprim
wave W_Cross
vprim=matrixDot(a1,W_Cross)
make/o/d/n=3 b1,b2,b3
b1=2*pi*W_Cross/vprim
cross a3,a1
b2=2*pi*W_Cross/vprim
cross a1,a2
b3=2*pi*W_Cross/vprim
Concatenate/O {b1,b2,b3},primbasis
make/o/d/n=3 aa1={a,-a,0}
make/o/d/n=3 aa2={0,a,-a}
make/o/d/n=3 aa3={a,a,a}
cross aa2,aa3
variable vconv
vconv=matrixDot(aa1,W_Cross)
make/o/d/n=3 bb1,bb2,bb3
bb1=2*pi*W_Cross/vconv
cross aa3,aa1
bb2=2*pi*W_Cross/vconv
bb3=1/2*(b1+b2+b3)
Concatenate/O {bb1,bb2,bb3},convbasis
matrixinverse primbasis
wave M_inverse
matrixop/o matC=M_inverse x convbasis
print vprim/vconv
cross b2,b3
variable kvprim
kvprim=matrixDot(b1,W_Cross)
cross bb2,bb3
variable kvconv
kvconv=matrixDot(bb1,W_Cross)
print kvprim/kvconv
end
function hr111() //for writing hr 111 KPOINTS
variable a,c
a= 3.6281/sqrt(3)
c=26.233/3
make/o/d/n=3 a1={-a/2,sqrt(3)/2*a,c}
make/o/d/n=3 a2={-a/2,-sqrt(3)/2*a,c}
make/o/d/n=3 a3={a,0,c}
cross a2,a3
variable vprim
wave W_Cross
vprim=matrixDot(a1,W_Cross)
make/o/d/n=3 b1,b2,b3
b1=2*pi*W_Cross/vprim
cross a3,a1
b2=2*pi*W_Cross/vprim
cross a1,a2
b3=2*pi*W_Cross/vprim
Concatenate/O {b1,b2,b3},primbasis
make/o/d/n=3 aa1=a2-a3
make/o/d/n=3 aa2=a3-a1
make/o/d/n=3 aa3=a1+a2+a3
cross aa2,aa3
variable vconv
wave W_cross
vconv=matrixDot(aa1,W_Cross)
make/o/d/n=3 bb1,bb2,bb3
bb1=2*pi*W_Cross/vconv
cross aa3,aa1
bb2=2*pi*W_Cross/vconv
bb3=1/2*(b1+b2+b3)
Concatenate/O {bb1,bb2,bb3},convbasis
matrixinverse primbasis
wave M_inverse
matrixop/o matC=M_inverse x convbasis
print vprim
print vprim/vconv
cross b2,b3
variable kvprim
kvprim=matrixDot(b1,W_Cross)
cross bb2,bb3
variable kvconv
kvconv=matrixDot(bb1,W_Cross)
print kvprim/kvconv
end
function changecoord() //for KPOINTS coordinate transformation
variable x,y,z
variable m,n,l
wave matC
variable i
variable k=0
for(i=0;i<21;i++)
x=1/2;y=0;z=i/20 // M bar point
m=matC[0][0]*x+matC[0][1]*y+matC[0][2]*z
n=matC[1][0]*x+matC[1][1]*y+matC[1][2]*z
l=matC[2][0]*x+matC[2][1]*y+matC[2][2]*z
if(i==0)
printf "%.12f %.12f %.12f\r" m, n, l
else
printf "%.12f %.12f %.12f\r" m, n, l
printf "%.12f %.12f %.12f\r" m, n, l
endif
x=0;y=0;z=i/20 //gamma bar point
m=matC[0][0]*x+matC[0][1]*y+matC[0][2]*z
n=matC[1][0]*x+matC[1][1]*y+matC[1][2]*z
l=matC[2][0]*x+matC[2][1]*y+matC[2][2]*z
printf "%.12f %.12f %.12f\r" m, n, l
printf "%.12f %.12f %.12f\r" m, n, l
x=1/3;y=1/3;z=i/20 //K bar point
m=matC[0][0]*x+matC[0][1]*y+matC[0][2]*z
n=matC[1][0]*x+matC[1][1]*y+matC[1][2]*z
l=matC[2][0]*x+matC[2][1]*y+matC[2][2]*z
if(i==20)
printf "%.12f %.12f %.12f\r" m, n, l
else
printf "%.12f %.12f %.12f\r" m, n, l
printf "%.12f %.12f %.12f\r" m, n, l
endif
endfor
end
function fcc111slab(nlayer,a) //for writing fcc 111 slab POSCAR
variable nlayer,a
make/o/d matA={{-1/2,1/2,0},{0,-1/2,1/2},{1,1,1}}
matrixinverse matA
wave M_inverse
make/o/d matB={{0,1/2,1/2},{1/2,0,1/2},{1/2,1/2,0}}
matrixop/o matC=M_inverse x matB
make/o/d/n=(3,3) Cefilm
variable x,y,z,m,n,l,i
for(i=0;i<3;i++)
x=0;y=0;z=i
m=matC[0][0]*x+matC[0][1]*y+matC[0][2]*z
n=matC[1][0]*x+matC[1][1]*y+matC[1][2]*z
l=matC[2][0]*x+matC[2][1]*y+matC[2][2]*z
do
if(m<0)
m=m+1
elseif((m>1)||(m==1))
m=m-1
endif
while((m<0)||(m>1)||(m==1))
do
if(n<0)
n=n+1
elseif((n>1)||(n==1))
n=n-1
endif
while((n<0)||(n>1)||(n==1))
do
if(l<0)
l=l+1
elseif((l>1)||(l==1))
l=l-1
endif
while((l<0)||(l>1)||(l==1))
Cefilm[i][0]=m ;Cefilm[i][1]=n ;Cefilm[i][2]=l
endfor
bubble(Cefilm,3)
redimension/n=(nlayer,3) Cefilm
for(i=3;i<nlayer;i++)
Cefilm[i][0]=Cefilm[i-3][0]
Cefilm[i][1]=Cefilm[i-3][1]
Cefilm[i][2]=Cefilm[i-3][2]+1
endfor
for(i=0;i<nlayer;i++)
printf "%.12f %.12f %.12f\r",Cefilm[i][0],Cefilm[i][1], (Cefilm[i][2]*a*sqrt(3)/((Cefilm[nlayer-1][2]-Cefilm[0][2])*a*sqrt(3)+18))
endfor
print "\r"
printf "%.12f\r" (Cefilm[nlayer-1][2]-Cefilm[0][2])*a*sqrt(3)+18
end
function bubble(w,n)
wave w
variable n
variable i,j,k,t1,t2,t3
for(i=1;i<n;i++)
k=0
for(j=0;j<(n-i);j++)
if(w[j][2]>w[j+1][2])
t1=w[j][0];w[j][0]=w[j+1][0];w[j+1][0]=t1
t2=w[j][1];w[j][1]=w[j+1][1];w[j+1][1]=t2
t3=w[j][2];w[j][2]=w[j+1][2];w[j+1][2]=t3
else
k++
endif
endfor
if(k==(n-i))
break
endif
endfor
end
function xcoordi()
wave b1,b2
make/o/d/n=200 xwave
variable i
for(i=0;i<100;i++)
xwave[i]=1/2*sqrt(b1[0]^2+b1[1]^2+b1[2]^2)/99*(i-99)
endfor
for(i=100;i<200;i++)
xwave[i]=1/3*sqrt((b1[0]+b2[0])^2+(b1[1]+b2[1])^2+(b1[2]+b2[2])^2)/99*(i-100)
endfor
end
#pragma rtGlobals=3 // Use modern global access method and strict wave access.
function vaspbands(path1,fname1,path2,fname2,fermilevel)
string path1,fname1,path2,fname2
variable fermilevel //fermilevel
variable a,i,t,j,k
string str,fword,nkpts,nbnds
Open/R a path1 + ":" + fname1
for(i=0;i<5;i++)
FReadLine a, str
endfor
FReadLine a, str
i=0
t=0
do
fword=stringfromlist(i,str," ")
if (stringmatch(fword,"")!=1)
t=t+1
if(t==2)
nkpts=fword
endif
if(t==3)
nbnds=fword
endif
endif
i=i+1
while(t<3)
make/d/n=(str2num(nkpts),str2num(nbnds)) w
for(i=0;i<str2num(nkpts);i++)
FReadLine a, str
FReadLine a, str
for(j=0;j<str2num(nbnds);j++)
FReadLine a, str
k=0
t=0
do
fword=stringfromlist(k,str," ")
if (stringmatch(fword,"")!=1)
t=t+1
if(t==2)
w[i][j]=str2num(fword)
endif
endif
k=k+1
while(t<2)
endfor
endfor
for(i=0;i<str2num(nbnds);i++)
make/d/n=(str2num(nkpts)) tmp
tmp=w[p][i]-fermilevel
duplicate tmp $"wave"+num2str(i)
killwaves tmp
endfor
killwaves w
close a
Open/R a path2 + ":" + fname2
make/d/n=3 kold,knew
make/d/n=(str2num(nkpts)) kpt
do
FReadLine a, str
if(strlen(str)==0)
break
endif
if(stringmatch(str,"*k-points in units of 2pi/SCALE*")==1)
FReadLine a, str
i=0
t=0
do
fword=stringfromlist(i,str," ")
if (stringmatch(fword,"")!=1)
t=t+1
if(t==1)
kold[0]=str2num(fword)
elseif(t==2)
kold[1]=str2num(fword)
elseif(t==3)
kold[2]=str2num(fword)
endif
endif
i=i+1
while(t<3)
kpt[0]=0
for(i=1;i<str2num(nkpts);i++)
FReadLine a, str
k=0
t=0
do
fword=stringfromlist(k,str," ")
if (stringmatch(fword,"")!=1)
t=t+1
if(t==1)
knew[0]=str2num(fword)
elseif(t==2)
knew[1]=str2num(fword)
elseif(t==3)
knew[2]=str2num(fword)
endif
endif
k=k+1
while(t<3)
if(mod(i,100)==0)
kpt[i]=kpt[i-1]
else
kpt[i]=kpt[i-1]+sqrt((knew[0]-kold[0])^2+(knew[1]-kold[1])^2+(knew[2]-kold[2])^2)
endif
kold=knew
endfor
endif
while(1)
killwaves kold,knew
close a
for(i=0;i<str2num(nbnds);i++)
if(i==0)
display $"wave"+num2str(i) vs kpt
else
appendtograph $"wave"+num2str(i) vs kpt
endif
endfor
end
function vaspbands_mbj(path1,fname1,path2,fname2,fermilevel)
string path1,fname1,path2,fname2
variable fermilevel //fermilevel
variable a,i,t,j,k
variable cc=120
string str,fword,nkpts,nbnds
Open/R a path1 + ":" + fname1
for(i=0;i<5;i++)
FReadLine a, str
endfor
FReadLine a, str
i=0
t=0
do
fword=stringfromlist(i,str," ")
if (stringmatch(fword,"")!=1)
t=t+1
if(t==2)
nkpts=fword
endif
if(t==3)
nbnds=fword
endif
endif
i=i+1
while(t<3)
make/d/n=(str2num(nkpts),str2num(nbnds)) w
for(i=0;i<str2num(nkpts);i++)
FReadLine a, str
FReadLine a, str
for(j=0;j<str2num(nbnds);j++)
FReadLine a, str
k=0
t=0
do
fword=stringfromlist(k,str," ")
if (stringmatch(fword,"")!=1)
t=t+1
if(t==2)
w[i][j]=str2num(fword)
endif
endif
k=k+1
while(t<2)
endfor
endfor
for(i=0;i<str2num(nbnds);i++)
make/d/n=(str2num(nkpts)) tmp
tmp=w[p][i]-fermilevel
duplicate tmp $"wave"+num2str(i)
killwaves tmp
endfor
killwaves w
close a
Open/R a path2 + ":" + fname2
make/d/n=3 kold,knew
make/d/n=(str2num(nkpts)) kpt
do
FReadLine a, str
if(strlen(str)==0)
break
endif
if(stringmatch(str,"*k-points in units of 2pi/SCALE*")==1)
FReadLine a, str
i=0
t=0
do
fword=stringfromlist(i,str," ")
if (stringmatch(fword,"")!=1)
t=t+1
if(t==1)
kold[0]=str2num(fword)
elseif(t==2)
kold[1]=str2num(fword)
elseif(t==3)
kold[2]=str2num(fword)
endif
endif
i=i+1
while(t<3)
kpt[0]=0
for(i=1;i<str2num(nkpts);i++)
FReadLine a, str
k=0
t=0
do
fword=stringfromlist(k,str," ")
if (stringmatch(fword,"")!=1)
t=t+1
if(t==1)
knew[0]=str2num(fword)
elseif(t==2)
knew[1]=str2num(fword)
elseif(t==3)
knew[2]=str2num(fword)
endif
endif
k=k+1
while(t<3)
if(mod(i,100)==0)
kpt[i]=kpt[i-1]
else
kpt[i]=kpt[i-1]+sqrt((knew[0]-kold[0])^2+(knew[1]-kold[1])^2+(knew[2]-kold[2])^2)
endif
kold=knew
endfor
endif
while(1)
killwaves kold,knew
close a
for(i=0;i<str2num(nbnds);i++)
DeletePoints 0,cc, $"wave"+num2str(i)
endfor
DeletePoints 0,cc,kpt
for(i=0;i<str2num(nbnds);i++)
if(i==0)
display $"wave"+num2str(i) vs kpt
else
appendtograph $"wave"+num2str(i) vs kpt
endif
endfor
end
function splitdos(path1,fname1,path2,fname2) //fname1=DOSCAR,fname2=POSCAR
string path1,fname1,path2,fname2
variable a,b,i,t,j,k
string str,fword,natms,nspin
Open/R a path1 + ":" + fname1
FReadLine a, str
i=0
t=0
do
fword=stringfromlist(i,str," ")
if (stringmatch(fword,"")!=1)
t=t+1
if(t==1)
natms=fword
endif
if(t==3)
nspin=fword
endif
endif
i=i+1
while(t<3)
for(i=0;i<4;i++)
FReadLine a, str
endfor
printf "Number of atoms: %d\r Number of spin: %d\r" str2num(natms),str2num(nspin)
Open/R b path2 + ":" + fname2
for(i=0;i<5;i++)
FReadLine b, str
endfor
FReadLine b, str
i=0
variable/g nspec=0
do
fword=stringfromlist(i,str," ")
if (stringmatch(fword,"")!=1)
nspec=nspec+1
endif
i=i+1
while(i<strlen(str))
make/t/n=(nspec)/o atmnam
i=0
t=0
do
fword=stringfromlist(i,str," ")
if (stringmatch(fword,"")!=1)
t=t+1
if(t==nspec)
variable memory=strlen(fword)
string lastnm
sscanf fword[0,memory-1],"%s",lastnm
atmnam[t-1]=lastnm
else
atmnam[t-1]=cleanupname(fword,0)
endif
endif
i=i+1
while(t<nspec)
make/d/n=(nspec)/o atmnum
FReadLine b, str
i=0
t=0
do
fword=stringfromlist(i,str," ")
if (stringmatch(fword,"")!=1)
t=t+1
atmnum[t-1]=str2num(fword)
endif
i=i+1
while(t<nspec)
close b
printf "Number of species: %d\r" nspec
variable temp=0
for(i=0;i<nspec;i++)
temp+=atmnum[i]
endfor
if(str2num(natms)!=temp)
printf "DOSCAR info & POSCAR does not match!\r"
close a
return -1
endif
FReadLine a, str
i=0
t=0
do
fword=stringfromlist(i,str," ",32)
if (stringmatch(fword,"")!=1)
t=t+1
if(t==1)
variable/g nedos=str2num(fword)
endif
if(t==2)
variable/g efermi=str2num(fword)
endif
endif
i=i+1
while(t<2)
for(i=0;i<nedos;i++)
FReadLine a, str
if(i==0)
j=0
variable tl=0
do
fword=stringfromlist(j,str," ")
if (stringmatch(fword,"")!=1)
tl=tl+1
endif
j=j+1
while(j<strlen(str))
make/o/n=(nedos,tl)/d pdos_tot
endif
j=0
t=0
do
fword=stringfromlist(j,str," ")
if (stringmatch(fword,"")!=1)
t=t+1
pdos_tot[i][t-1]=str2num(fword)
endif
j=j+1
while(t<tl)
endfor
for(i=0;i<nspec;i++)
for(j=0;j<atmnum[i];j++)
FReadLine a, str
for(k=0;k<nedos;k++)
FReadLine a, str
if(k==0)
variable in=0
tl=0
do
fword=stringfromlist(in,str," ")
if (stringmatch(fword,"")!=1)
tl=tl+1
endif
in=in+1
while(in<strlen(str))
make/o/n=(nedos,tl)/d tempwave
endif
in=0
t=0
do
fword=stringfromlist(in,str," ")
if (stringmatch(fword,"")!=1)
t=t+1
tempwave[k][t-1]=str2num(fword)
endif
in=in+1
while(t<tl)
endfor
duplicate/o tempwave, $"pdos_"+atmnam[i]+"_"+num2str(j+1)
endfor
endfor
killwaves tempwave
end
function splitdosif(path1,fname1,path2,fname2)
string path1,fname1,path2,fname2
variable a,b,i,t,j,k
string str,fword,natms,nspin
Open/R a path1 + ":" + fname1
FReadLine a, str
i=0
t=0
do
fword=stringfromlist(i,str," ")
if (stringmatch(fword,"")!=1)
t=t+1
if(t==1)
natms=fword
endif
if(t==3)
nspin=fword
endif
endif
i=i+1
while(t<3)
for(i=0;i<4;i++)
FReadLine a, str
endfor
printf "Number of atoms: %d\r Number of spin: %d\r" str2num(natms),str2num(nspin)
Open/R b path2 + ":" + fname2
for(i=0;i<5;i++)
FReadLine b, str
endfor
FReadLine b, str
i=0
variable/g nspec=0
do
fword=stringfromlist(i,str," ")
if (stringmatch(fword,"")!=1)
nspec=nspec+1
endif
i=i+1
while(i<strlen(str))
make/t/n=(nspec)/o atmnam
i=0
t=0
do
fword=stringfromlist(i,str," ")
if (stringmatch(fword,"")!=1)
t=t+1
if(t==nspec)
variable memory=strlen(fword)
string lastnm
sscanf fword[0,memory-1],"%s",lastnm
atmnam[t-1]=lastnm
else
atmnam[t-1]=cleanupname(fword,0)
endif
endif
i=i+1
while(t<nspec)
make/d/n=(nspec)/o atmnum
FReadLine b, str
i=0
t=0
do
fword=stringfromlist(i,str," ")
if (stringmatch(fword,"")!=1)
t=t+1
atmnum[t-1]=str2num(fword)
endif
i=i+1
while(t<nspec)
close b
printf "Number of species: %d\r" nspec
variable temp=0
for(i=0;i<nspec;i++)
temp+=atmnum[i]
endfor
if(str2num(natms)!=temp)
printf "DOSCAR info & POSCAR does not match!\r"
close a
return -1
endif
FReadLine a, str
i=0
t=0
do
fword=stringfromlist(i,str," ",32)
if (stringmatch(fword,"")!=1)
t=t+1
if(t==1)
variable/g nedos=str2num(fword)
endif
if(t==2)
variable/g efermi=str2num(fword)
endif
endif
i=i+1
while(t<2)
for(i=0;i<nedos;i++)
FReadLine a, str
if(i==0)
j=0
variable tl=0
do
fword=stringfromlist(j,str," ")
if (stringmatch(fword,"")!=1)
tl=tl+1
endif
j=j+1
while(j<strlen(str))
make/o/n=(nedos,tl)/d pdos_tot
endif
j=0
t=0
do
fword=stringfromlist(j,str," ")
if (stringmatch(fword,"")!=1)
t=t+1
pdos_tot[i][t-1]=str2num(fword)
endif
j=j+1
while(t<tl)
endfor
for(i=0;i<nspec;i++)
for(j=0;j<atmnum[i];j++)
FReadLine a, str
for(k=0;k<(2*nedos);k++)
FReadLine a, str
if(k==0)
variable in=0
tl=0
do
fword=stringfromlist(in,str," ")
if (stringmatch(fword,"")!=1)
tl=tl+1
endif
in=in+1
while(in<strlen(str))
make/o/n=(nedos,tl)/d tempwave
endif
if(k==1)
in=0
variable tll=0
do
fword=stringfromlist(in,str," ")
if (stringmatch(fword,"")!=1)
tll=tll+1
endif
in=in+1
while(in<strlen(str))
redimension/n=(nedos,(tl+tll))/d tempwave
endif
in=0
t=0
if(mod(k,2)==0)
do
fword=stringfromlist(in,str," ")
if (stringmatch(fword,"")!=1)
t=t+1
tempwave[(k/2)][t-1]=str2num(fword)
endif
in=in+1
while(t<tl)
else
do
fword=stringfromlist(in,str," ")
if (stringmatch(fword,"")!=1)
t=t+1
tempwave[floor(k/2)][tl+t-1]=str2num(fword)
endif
in=in+1
while(t<tll)
endif
endfor
duplicate/o tempwave, $"pdos_"+atmnam[i]+"_"+num2str(j+1)
endfor
endfor
killwaves tempwave
end
function sumdos()
variable i,j
variable/g nspec,nedos
wave atmnum
wave/t atmnam
for(i=0;i<nspec;i++)
for(j=0;j<atmnum[i];j++)
if(j==0)
duplicate/o $"pdos_"+atmnam[i]+"_"+num2str(j+1),tempwave
duplicate/o/r=[][0] $"pdos_"+atmnam[i]+"_"+num2str(j+1),energywave
else
duplicate/o $"pdos_"+atmnam[i]+"_"+num2str(j+1),tempwave2
tempwave+=tempwave2
killwaves tempwave2
endif
endfor
variable k
for(k=0;k<nedos;k++)
tempwave[k][0]=energywave[k]
endfor
duplicate/o tempwave $"pdos_"+atmnam[i]+"_tot"
endfor
killwaves tempwave,energywave
end
function plotdosif()
variable i
variable/g nspec
wave/t atmnam
for(i=0;i<nspec;i++)
variable j
for(j=0;j<65;j++)
if(j==1)
duplicate/o/r=[][j] $"pdos_"+atmnam[i]+"_tot" $"pdos_"+atmnam[i]+"_tot_s"
endif
if(j==5)
duplicate/o/r=[][j] $"pdos_"+atmnam[i]+"_tot" temp_p
endif
if(j==9)
duplicate/o/r=[][j] $"pdos_"+atmnam[i]+"_tot" temp_p_t
temp_p+=temp_p_t
endif
if(j==13)
duplicate/o/r=[][j] $"pdos_"+atmnam[i]+"_tot" temp_p_t
temp_p+=temp_p_t
duplicate/o temp_p,$"pdos_"+atmnam[i]+"_tot_p"
killwaves temp_p_t,temp_p
endif
if(j==17)
duplicate/o/r=[][j] $"pdos_"+atmnam[i]+"_tot" temp_d
endif
if(j==21)
duplicate/o/r=[][j] $"pdos_"+atmnam[i]+"_tot" temp_d_t
temp_d+=temp_d_t
endif
if(j==25)
duplicate/o/r=[][j] $"pdos_"+atmnam[i]+"_tot" temp_d_t
temp_d+=temp_d_t
endif
if(j==29)
duplicate/o/r=[][j] $"pdos_"+atmnam[i]+"_tot" temp_d_t
temp_d+=temp_d_t
endif
if(j==33)
duplicate/o/r=[][j] $"pdos_"+atmnam[i]+"_tot" temp_d_t
temp_d+=temp_d_t
duplicate/o temp_d,$"pdos_"+atmnam[i]+"_tot_d"
killwaves temp_d_t,temp_d
endif
if(j==37)
duplicate/o/r=[][j] $"pdos_"+atmnam[i]+"_tot" temp_f
endif
if(j==41)
duplicate/o/r=[][j] $"pdos_"+atmnam[i]+"_tot" temp_f_t
temp_f+=temp_f_t
endif
if(j==45)
duplicate/o/r=[][j] $"pdos_"+atmnam[i]+"_tot" temp_f_t
temp_f+=temp_f_t
endif
if(j==49)
duplicate/o/r=[][j] $"pdos_"+atmnam[i]+"_tot" temp_f_t
temp_f+=temp_f_t
endif
if(j==53)
duplicate/o/r=[][j] $"pdos_"+atmnam[i]+"_tot" temp_f_t
temp_f+=temp_f_t
endif
if(j==57)
duplicate/o/r=[][j] $"pdos_"+atmnam[i]+"_tot" temp_f_t
temp_f+=temp_f_t
endif
if(j==61)
duplicate/o/r=[][j] $"pdos_"+atmnam[i]+"_tot" temp_f_t
temp_f+=temp_f_t
duplicate/o temp_f,$"pdos_"+atmnam[i]+"_tot_f"
killwaves temp_f_t,temp_f
endif
endfor
endfor
for(i=0;i<nspec;i++)
if(i==0)
display $"pdos_"+atmnam[i]+"_tot_s" vs $"pdos_"+atmnam[i]+"_tot"[][0]
else
appendtograph $"pdos_"+atmnam[i]+"_tot_s" vs $"pdos_"+atmnam[i]+"_tot"[][0]
endif
appendtograph $"pdos_"+atmnam[i]+"_tot_p" vs $"pdos_"+atmnam[i]+"_tot"[][0]
appendtograph $"pdos_"+atmnam[i]+"_tot_d" vs $"pdos_"+atmnam[i]+"_tot"[][0]
appendtograph $"pdos_"+atmnam[i]+"_tot_f" vs $"pdos_"+atmnam[i]+"_tot"[][0]
endfor
end
function plotdoslf()
variable i
variable/g nspec
wave/t atmnam
for(i=0;i<nspec;i++)
variable j
for(j=0;j<65;j++)
if(j==1)
duplicate/o/r=[][j] $"pdos_"+atmnam[i]+"_tot" $"pdos_"+atmnam[i]+"_tot_s"
endif
if(j==5)
duplicate/o/r=[][j] $"pdos_"+atmnam[i]+"_tot" temp_p
endif
if(j==9)
duplicate/o/r=[][j] $"pdos_"+atmnam[i]+"_tot" temp_p_t
temp_p+=temp_p_t
endif
if(j==13)
duplicate/o/r=[][j] $"pdos_"+atmnam[i]+"_tot" temp_p_t
temp_p+=temp_p_t
duplicate/o temp_p,$"pdos_"+atmnam[i]+"_tot_p"
killwaves temp_p_t,temp_p
endif
if(j==17)
duplicate/o/r=[][j] $"pdos_"+atmnam[i]+"_tot" temp_d
endif
if(j==21)
duplicate/o/r=[][j] $"pdos_"+atmnam[i]+"_tot" temp_d_t
temp_d+=temp_d_t
endif
if(j==25)
duplicate/o/r=[][j] $"pdos_"+atmnam[i]+"_tot" temp_d_t
temp_d+=temp_d_t
endif
if(j==29)
duplicate/o/r=[][j] $"pdos_"+atmnam[i]+"_tot" temp_d_t
temp_d+=temp_d_t
endif
if(j==33)
duplicate/o/r=[][j] $"pdos_"+atmnam[i]+"_tot" temp_d_t
temp_d+=temp_d_t
duplicate/o temp_d,$"pdos_"+atmnam[i]+"_tot_d"
killwaves temp_d_t,temp_d
endif
endfor
endfor
for(i=0;i<nspec;i++)
if(i==0)
display $"pdos_"+atmnam[i]+"_tot_s" vs $"pdos_"+atmnam[i]+"_tot"[][0]
else
appendtograph $"pdos_"+atmnam[i]+"_tot_s" vs $"pdos_"+atmnam[i]+"_tot"[][0]
endif
appendtograph $"pdos_"+atmnam[i]+"_tot_p" vs $"pdos_"+atmnam[i]+"_tot"[][0]
appendtograph $"pdos_"+atmnam[i]+"_tot_d" vs $"pdos_"+atmnam[i]+"_tot"[][0]
endfor
end
function Findbands(numbnds,efermi)
variable numbnds,efermi
variable i
make/o/n=(numbnds) mini,maxm
for(i=0;i<numbnds;i++)
loadwave/a/o/d/g/L={0,21+1030302*i,1030301,0,0} "C:\Users\Lenovo\Desktop\wannier90.bxsf"
endfor
for(i=0;i<numbnds;i++)
wavestats $"wave"+num2str(i)
mini[i]=V_min
maxm[i]=V_max
killwaves $"wave"+num2str(i)
endfor
for(i=0;i<numbnds;i++)
if((efermi>mini[i])&&(efermi<maxm[i]))
print i+1
endif
endfor
end
function extract_fs(path1,fname1,path2,fname2,nstart,nend) //fname1 read;fname2 write
string path1,fname1,path2,fname2
variable nstart,nend
variable a,b,i,t,j,k
string str,temp,fword
variable nx,ny,nz
Open/R a path1 + ":" + fname1
Open b path2 + ":" + fname2
for(i=0;i<14;i++)
FReadLine a, str
temp=str[0,strlen(str)-2]
fprintf b,temp
fprintf b,"\n"
endfor
FReadLine a, str
i=0
t=0
do
fword=stringfromlist(i,str," ")
if (stringmatch(fword,"")!=1)
t=t+1
if(t==1)
temp=fword[0,strlen(fword)-2]
endif
endif
i=i+1
while(t<1)
if((nend<nstart)||(str2num(temp)<nend))
fprintf b,"Wrong indices,exit...\n"
close a
close b
return -1
endif
fprintf b,"%15d\n",nend-nstart+1
FReadLine a, str
temp=str[0,strlen(str)-2]
fprintf b,temp
fprintf b,"\n"
i=0
t=0
do
fword=stringfromlist(i,str," ")
if (stringmatch(fword,"")!=1)
t=t+1
if(t==1)
temp=fword[0,strlen(fword)-1]
nx=str2num(temp)
endif
if(t==2)
temp=fword[0,strlen(fword)-1]
ny=str2num(temp)
endif
if(t==3)
temp=fword[0,strlen(fword)-2]
nz=str2num(temp)
endif
endif
i=i+1
while(t<3)
for(i=0;i<4;i++)
FReadLine a, str
temp=str[0,strlen(str)-2]
fprintf b,temp
fprintf b,"\n"
endfor
for(i=1;i<nstart;i++)
for(j=0;j<(nx*ny*nz+1);j++)
FReadLine a, str
endfor
endfor
for(i=0;i<(nend-nstart+1);i++)
FReadLine a, str
fprintf b,"BAND:%12d\n",i+1
for(j=0;j<(nx*ny*nz);j++)
FReadLine a, str
temp=str[0,strlen(str)-2]
fprintf b,temp
fprintf b,"\n"
endfor
endfor
fprintf b," END_BANDGRID_3D\n"
fprintf b," END_BLOCK_BANDGRID_3D\n"
close a
close b
end
function vaspbands_proj(path1,fname1,atoms,orbitals) //read procar non-soc
wave atoms
wave orbitals
string path1,fname1
variable a,i,t,j,k,u
variable s
variable inner,np,ccc
string str,fword,nkpts,nbnds,nions
Open/R a path1 + ":" + fname1
FReadLine a, str
FReadLine a, str
i=0
t=0
do
fword=stringfromlist(i,str," ")
if (stringmatch(fword,"")!=1)
t=t+1
if(t==4)
nkpts=fword
endif
if(t==8)
nbnds=fword
endif
if(t==12)
nions=fword
endif
endif
i=i+1
while(t<12)
make/o/d/n=(str2num(nkpts),str2num(nbnds)) tem
for(u=0;u<str2num(nkpts);u++)
for(i=0;i<3;i++)
FReadLine a, str
endfor
for(k=0;k<str2num(nbnds);k++)
for(i=0;i<3;i++)
FReadLine a, str
endfor
s=0
for(np=0;np<numpnts(atoms);np++)
if(np==0)
if(atoms[0]==1)
FReadLine a, str
j=0
t=0
do
fword=stringfromlist(j,str," ")
if (stringmatch(fword,"")!=1)
t=t+1
for(ccc=0;ccc<numpnts(orbitals);ccc++)
if(t==orbitals[ccc])
s+=str2num(fword)
endif
endfor
endif
j=j+1
while(t<orbitals[numpnts(orbitals)-1])
endif
if(atoms[0]>1)
for(inner=0;inner<(atoms[0]-1);inner++)
FReadLine a, str
endfor
FReadLine a, str
j=0
t=0
do
fword=stringfromlist(j,str," ")
if (stringmatch(fword,"")!=1)
t=t+1
for(ccc=0;ccc<numpnts(orbitals);ccc++)
if(t==orbitals[ccc])
s+=str2num(fword)
endif
endfor
endif
j=j+1
while(t<orbitals[numpnts(orbitals)-1])
endif
else
if((atoms[np]-atoms[np-1])==1)
FReadLine a, str
j=0
t=0
do
fword=stringfromlist(j,str," ")
if (stringmatch(fword,"")!=1)
t=t+1
for(ccc=0;ccc<numpnts(orbitals);ccc++)
if(t==orbitals[ccc])
s+=str2num(fword)
endif
endfor
endif
j=j+1
while(t<orbitals[numpnts(orbitals)-1])
endif
if((atoms[np]-atoms[np-1])>1)
for(inner=0;inner<(atoms[np]-atoms[np-1]-1);inner++)
FReadLine a, str
endfor
FReadLine a, str
j=0
t=0
do
fword=stringfromlist(j,str," ")
if (stringmatch(fword,"")!=1)
t=t+1
for(ccc=0;ccc<numpnts(orbitals);ccc++)
if(t==orbitals[ccc])
s+=str2num(fword)
endif
endfor
endif
j=j+1
while(t<orbitals[numpnts(orbitals)-1])
endif
endif
endfor
if(str2num(nions)==atoms[numpnts(atoms)-1])
FReadLine a, str
FReadLine a, str
else
for(i=0;i<(str2num(nions)-atoms[numpnts(atoms)-1]);i++)
FReadLine a, str
endfor
FReadLine a, str
FReadLine a, str
endif
tem[u][k]=s
endfor
endfor
for(i=0;i<str2num(nbnds);i++)
duplicate/o/r=[0,str2num(nkpts)-1][i,i] tem $"weight"+num2str(i)
endfor
close a
end
function vaspbands_projsoc(path1,fname1,atoms,orbitals) //read procar soc
wave atoms
wave orbitals
string path1,fname1
variable a,i,t,j,k,u
variable s
variable inner,np,ccc
string str,fword,nkpts,nbnds,nions
Open/R a path1 + ":" + fname1
FReadLine a, str
FReadLine a, str
i=0
t=0
do
fword=stringfromlist(i,str," ")
if (stringmatch(fword,"")!=1)
t=t+1
if(t==4)
nkpts=fword
endif
if(t==8)
nbnds=fword
endif
if(t==12)
nions=fword
endif
endif
i=i+1
while(t<12)
make/o/d/n=(str2num(nkpts),str2num(nbnds)) tem
for(u=0;u<str2num(nkpts);u++)
for(i=0;i<3;i++)
FReadLine a, str
endfor
for(k=0;k<str2num(nbnds);k++)
for(i=0;i<3;i++)
FReadLine a, str
endfor
s=0
for(np=0;np<numpnts(atoms);np++)
if(np==0)
if(atoms[0]==1)
FReadLine a, str
j=0
t=0
do
fword=stringfromlist(j,str," ")
if (stringmatch(fword,"")!=1)
t=t+1
for(ccc=0;ccc<numpnts(orbitals);ccc++)
if(t==orbitals[ccc])
s+=str2num(fword)
endif
endfor
endif
j=j+1
while(t<orbitals[numpnts(orbitals)-1])
endif
if(atoms[0]>1)
for(inner=0;inner<(atoms[0]-1);inner++)
FReadLine a, str
endfor
FReadLine a, str
j=0
t=0
do
fword=stringfromlist(j,str," ")
if (stringmatch(fword,"")!=1)
t=t+1
for(ccc=0;ccc<numpnts(orbitals);ccc++)
if(t==orbitals[ccc])
s+=str2num(fword)
endif
endfor
endif
j=j+1
while(t<orbitals[numpnts(orbitals)-1])
endif
else
if((atoms[np]-atoms[np-1])==1)
FReadLine a, str
j=0
t=0
do
fword=stringfromlist(j,str," ")
if (stringmatch(fword,"")!=1)
t=t+1
for(ccc=0;ccc<numpnts(orbitals);ccc++)
if(t==orbitals[ccc])
s+=str2num(fword)
endif
endfor
endif
j=j+1
while(t<orbitals[numpnts(orbitals)-1])
endif
if((atoms[np]-atoms[np-1])>1)
for(inner=0;inner<(atoms[np]-atoms[np-1]-1);inner++)
FReadLine a, str
endfor
FReadLine a, str
j=0
t=0
do
fword=stringfromlist(j,str," ")
if (stringmatch(fword,"")!=1)
t=t+1
for(ccc=0;ccc<numpnts(orbitals);ccc++)
if(t==orbitals[ccc])
s+=str2num(fword)
endif
endfor
endif
j=j+1
while(t<orbitals[numpnts(orbitals)-1])
endif
endif
endfor
variable jjj
for(jjj=0;jjj<(3*(str2num(nions)+1+1));jjj++)
FReadLine a, str
endfor
if(str2num(nions)==atoms[numpnts(atoms)-1])
FReadLine a, str
FReadLine a, str
FReadLine a, str
else
for(i=0;i<(str2num(nions)-atoms[numpnts(atoms)-1]);i++)
FReadLine a, str
endfor
FReadLine a, str
FReadLine a, str
FReadLine a, str
endif
tem[u][k]=s
endfor
endfor
for(i=0;i<str2num(nbnds);i++)
duplicate/o/r=[0,str2num(nkpts)-1][i,i] tem $"weight"+num2str(i)
endfor
close a
end
function projbulkbnds(m,n,p) //projected bulk bands plot.m:the number of points between two kpoints;
//n:number of different kzs;p:numberof bands.
variable m,n,p
variable i,j,k
wave kpt
make/n=(2*m) refer
for(i=0;i<(2*m);i++)
refer[i]=kpt[i]
endfor
make/n=(2*m) tem1
for(i=0;i<p;i++)
duplicate $"wave"+num2str(i),tem2
for(j=i;j<i+p*(n-1)+1;j=j+p)
for(k=0;k<(2*m);k++)
tem1[k]=tem2[k]
endfor
duplicate tem1,$"w"+num2str(j)
if(j==0)
display $"w"+num2str(j) vs refer
else
appendtograph $"w"+num2str(j) vs refer
endif
DeletePoints 0,3*m, tem2
endfor
killwaves tem2
endfor
killwaves tem1
end
function difkzplot(nbnds,nkz,ymin,ymax,j)
variable nbnds
variable nkz
variable ymin,ymax //display range
variable j //j can take value 0-20,corresponding to kz=j/(nkz-1)*(gamma-z)
variable kz=j/(nkz-1)
wave xwave
variable i
for(i=0;i<nbnds;i++)
if(i==0)
display $"w"+num2str(j*nbnds) vs xwave
TextBox/C/N=text0/F=0/A=MC num2str(kz)+"*G-Z"
SetAxis left ymin,ymax
else
appendtograph $"w"+num2str(j*nbnds+i) vs xwave
endif
endfor
end
function fcc111() //for writing fcc 111 KPOINTS
variable a,b,c
a= 5.429
make/o/d/n=3 a1={0,a/2,a/2}
make/o/d/n=3 a2={a/2,0,a/2}
make/o/d/n=3 a3={a/2,a/2,0}
cross a2,a3
variable vprim
wave W_Cross
vprim=matrixDot(a1,W_Cross)
make/o/d/n=3 b1,b2,b3
b1=2*pi*W_Cross/vprim
cross a3,a1
b2=2*pi*W_Cross/vprim
cross a1,a2
b3=2*pi*W_Cross/vprim
Concatenate/O {b1,b2,b3},primbasis
make/o/d/n=3 aa1={-a/2,a/2,0}
make/o/d/n=3 aa2={0,-a/2,a/2}
make/o/d/n=3 aa3={a,a,a}
cross aa2,aa3
variable vconv
vconv=matrixDot(aa1,W_Cross)
make/o/d/n=3 bb1,bb2,bb3
bb1=2*pi*W_Cross/vconv
cross aa3,aa1
bb2=2*pi*W_Cross/vconv
bb3=1/2*(b1+b2+b3)
Concatenate/O {bb1,bb2,bb3},convbasis
matrixinverse primbasis
wave M_inverse
matrixop/o matC=M_inverse x convbasis
print vprim
print vprim/vconv
cross b2,b3
variable kvprim
kvprim=matrixDot(b1,W_Cross)
cross bb2,bb3
variable kvconv
kvconv=matrixDot(bb1,W_Cross)
print kvprim/kvconv
end
function simplecubic111() //for writing simple cubic 111 KPOINTS
variable a,b,c
a= 2.3537117117709037*2
make/o/d/n=3 a1={a,0,0}
make/o/d/n=3 a2={0,a,0}
make/o/d/n=3 a3={0,0,a}
cross a2,a3
variable vprim
wave W_Cross
vprim=matrixDot(a1,W_Cross)
make/o/d/n=3 b1,b2,b3
b1=2*pi*W_Cross/vprim
cross a3,a1
b2=2*pi*W_Cross/vprim
cross a1,a2
b3=2*pi*W_Cross/vprim
Concatenate/O {b1,b2,b3},primbasis
make/o/d/n=3 aa1={a,-a,0}
make/o/d/n=3 aa2={0,a,-a}
make/o/d/n=3 aa3={a,a,a}
cross aa2,aa3
variable vconv
vconv=matrixDot(aa1,W_Cross)
make/o/d/n=3 bb1,bb2,bb3
bb1=2*pi*W_Cross/vconv
cross aa3,aa1
bb2=2*pi*W_Cross/vconv
bb3=1/2*(b1+b2+b3)
Concatenate/O {bb1,bb2,bb3},convbasis
matrixinverse primbasis
wave M_inverse
matrixop/o matC=M_inverse x convbasis
print vprim/vconv
cross b2,b3
variable kvprim
kvprim=matrixDot(b1,W_Cross)
cross bb2,bb3
variable kvconv
kvconv=matrixDot(bb1,W_Cross)
print kvprim/kvconv
end
function hr111() //for writing hr 111 KPOINTS
variable a,c
a= 3.6281/sqrt(3)
c=26.233/3
make/o/d/n=3 a1={-a/2,sqrt(3)/2*a,c}
make/o/d/n=3 a2={-a/2,-sqrt(3)/2*a,c}
make/o/d/n=3 a3={a,0,c}
cross a2,a3
variable vprim
wave W_Cross
vprim=matrixDot(a1,W_Cross)
make/o/d/n=3 b1,b2,b3
b1=2*pi*W_Cross/vprim
cross a3,a1
b2=2*pi*W_Cross/vprim
cross a1,a2
b3=2*pi*W_Cross/vprim
Concatenate/O {b1,b2,b3},primbasis
make/o/d/n=3 aa1=a2-a3
make/o/d/n=3 aa2=a3-a1
make/o/d/n=3 aa3=a1+a2+a3
cross aa2,aa3
variable vconv
wave W_cross
vconv=matrixDot(aa1,W_Cross)
make/o/d/n=3 bb1,bb2,bb3
bb1=2*pi*W_Cross/vconv
cross aa3,aa1
bb2=2*pi*W_Cross/vconv
bb3=1/2*(b1+b2+b3)
Concatenate/O {bb1,bb2,bb3},convbasis
matrixinverse primbasis
wave M_inverse
matrixop/o matC=M_inverse x convbasis
print vprim
print vprim/vconv
cross b2,b3
variable kvprim
kvprim=matrixDot(b1,W_Cross)
cross bb2,bb3
variable kvconv
kvconv=matrixDot(bb1,W_Cross)
print kvprim/kvconv
end
function changecoord() //for KPOINTS coordinate transformation
variable x,y,z
variable m,n,l
wave matC
variable i
variable k=0
for(i=0;i<21;i++)
x=1/2;y=0;z=i/20 // M bar point
m=matC[0][0]*x+matC[0][1]*y+matC[0][2]*z
n=matC[1][0]*x+matC[1][1]*y+matC[1][2]*z
l=matC[2][0]*x+matC[2][1]*y+matC[2][2]*z
if(i==0)
printf "%.12f %.12f %.12f\r" m, n, l
else
printf "%.12f %.12f %.12f\r" m, n, l
printf "%.12f %.12f %.12f\r" m, n, l
endif
x=0;y=0;z=i/20 //gamma bar point
m=matC[0][0]*x+matC[0][1]*y+matC[0][2]*z
n=matC[1][0]*x+matC[1][1]*y+matC[1][2]*z
l=matC[2][0]*x+matC[2][1]*y+matC[2][2]*z
printf "%.12f %.12f %.12f\r" m, n, l
printf "%.12f %.12f %.12f\r" m, n, l
x=1/3;y=1/3;z=i/20 //K bar point
m=matC[0][0]*x+matC[0][1]*y+matC[0][2]*z
n=matC[1][0]*x+matC[1][1]*y+matC[1][2]*z
l=matC[2][0]*x+matC[2][1]*y+matC[2][2]*z
if(i==20)
printf "%.12f %.12f %.12f\r" m, n, l
else
printf "%.12f %.12f %.12f\r" m, n, l
printf "%.12f %.12f %.12f\r" m, n, l
endif
endfor
end
function fcc111slab(nlayer,a) //for writing fcc 111 slab POSCAR
variable nlayer,a
make/o/d matA={{-1/2,1/2,0},{0,-1/2,1/2},{1,1,1}}
matrixinverse matA
wave M_inverse
make/o/d matB={{0,1/2,1/2},{1/2,0,1/2},{1/2,1/2,0}}
matrixop/o matC=M_inverse x matB
make/o/d/n=(3,3) Cefilm
variable x,y,z,m,n,l,i
for(i=0;i<3;i++)
x=0;y=0;z=i
m=matC[0][0]*x+matC[0][1]*y+matC[0][2]*z
n=matC[1][0]*x+matC[1][1]*y+matC[1][2]*z
l=matC[2][0]*x+matC[2][1]*y+matC[2][2]*z
do
if(m<0)
m=m+1
elseif((m>1)||(m==1))
m=m-1
endif
while((m<0)||(m>1)||(m==1))
do
if(n<0)
n=n+1
elseif((n>1)||(n==1))
n=n-1
endif
while((n<0)||(n>1)||(n==1))
do
if(l<0)
l=l+1
elseif((l>1)||(l==1))
l=l-1
endif
while((l<0)||(l>1)||(l==1))
Cefilm[i][0]=m ;Cefilm[i][1]=n ;Cefilm[i][2]=l
endfor
bubble(Cefilm,3)
redimension/n=(nlayer,3) Cefilm
for(i=3;i<nlayer;i++)
Cefilm[i][0]=Cefilm[i-3][0]
Cefilm[i][1]=Cefilm[i-3][1]
Cefilm[i][2]=Cefilm[i-3][2]+1
endfor
for(i=0;i<nlayer;i++)
printf "%.12f %.12f %.12f\r",Cefilm[i][0],Cefilm[i][1], (Cefilm[i][2]*a*sqrt(3)/((Cefilm[nlayer-1][2]-Cefilm[0][2])*a*sqrt(3)+18))
endfor
print "\r"
printf "%.12f\r" (Cefilm[nlayer-1][2]-Cefilm[0][2])*a*sqrt(3)+18
end
function bubble(w,n)
wave w
variable n
variable i,j,k,t1,t2,t3
for(i=1;i<n;i++)
k=0
for(j=0;j<(n-i);j++)
if(w[j][2]>w[j+1][2])
t1=w[j][0];w[j][0]=w[j+1][0];w[j+1][0]=t1
t2=w[j][1];w[j][1]=w[j+1][1];w[j+1][1]=t2
t3=w[j][2];w[j][2]=w[j+1][2];w[j+1][2]=t3
else
k++
endif
endfor
if(k==(n-i))
break
endif
endfor
end
function xcoordi()
wave b1,b2
make/o/d/n=200 xwave
variable i
for(i=0;i<100;i++)
xwave[i]=1/2*sqrt(b1[0]^2+b1[1]^2+b1[2]^2)/99*(i-99)
endfor
for(i=100;i<200;i++)
xwave[i]=1/3*sqrt((b1[0]+b2[0])^2+(b1[1]+b2[1])^2+(b1[2]+b2[2])^2)/99*(i-100)
endfor
end
Vasp utilities written in Igor. It can be used to plot bands, projected dos, weighted bands and find bands cross the fermi level. It can also be used to write POSCAR and KPOINTS for several crystal structures.
November 3, 2020 at 05:47 am - Permalink
Would you be so kind to do two things: Package this as an Igor Pro pxp file and provide an example of the plots that it can produce.
November 3, 2020 at 09:25 am - Permalink
In reply to Would you be so kind to do… by jjweimer
I don't quite understand what is meant by "Package this as an Igor Pro pxp file". Could you explain it to me? If you want to do some plot, you can just run one of the functions, e.g. vaspbands("D","EIGENVAL","D","OUTCAR",1.2). As for the plots it can give, I can show some graphs I have done later.
November 3, 2020 at 05:45 pm - Permalink
In reply to Would you be so kind to do… by jjweimer
slab bands plot
projected bulk bands plot
dos
weighted bands read from EIGENVAL and PROCAR
Fermi surface extracted using Igor and plotted using Fermisurfer
November 3, 2020 at 06:20 pm - Permalink
In reply to I don't quite understand… by Yuan Fang
Well, it seems that some files are required which many of the forum users will not have readily available.
November 4, 2020 at 12:04 am - Permalink
In reply to Yuan Fang wrote: If you… by ChrLie
Yes. It requires output files from Vienna Ab initio Simulation Package (VASP). And for function extract_fs(), it requires output file from Wannier90.
November 4, 2020 at 12:12 am - Permalink
> I don't quite understand what is meant by "Package this as an Igor Pro pxp file".
* Open Igor Pro.
* Copy the code into the procedure window.
* Save the Igor Pro experiment (e.g. as VASP Utitlities Code.pxp).
* Post the pxp file.
The better approach for long code that does something specific such as your example is to create an Igor Pro project rather than post the raw code as forum discussion. See the brief information and longer list here. https://www.wavemetrics.com/projects
Finally, it would be helpful when the pxp project file contained sample raw data so that the functions could be run to create example plots directly. Then, it could be called VASP Utilities Demo.pxp by example.
November 4, 2020 at 06:07 am - Permalink
In reply to > I don't quite understand… by jjweimer
Thank you! I think if I create a project for this code maybe I will also need to upload some demo input files generated by VASP onto this site. Is it possible? Maybe I can do this when I have time later.
November 10, 2020 at 09:58 pm - Permalink
Here is a synopsis of the steps to make what you have posted into a project on Igor Exchange.
Create an Igor Pro experiment (pxp file) with the code. Include example data in the experiment. Make plots using the example data. Save that file as VASP Demo. Compress that file as a ZIP archive. Create a NEW PROJECT on this site. Upload the ZIP archive.
November 11, 2020 at 08:20 am - Permalink
In reply to Here is a synopsis of the… by jjweimer
I have created a new project based on your suggestion. Now the input and output of the files containing band information can be realized through a panel.
January 8, 2022 at 08:40 pm - Permalink