各位老铁们好,相信很多人对布里渊区都不是特别的了解,因此呢,今天就来为大家分享下关于布里渊区以及布里渊区及能带路径的可视化的问题知识,还望可以帮助大家,解决大家的一些困惑,下面一起来看看吧!
先前本公众号已经介绍"五种方法生成能带结构计算的高对称点",使用里面的方法即可获取相应晶体对应的布里渊区图。在VASPKIT1.20新版本中即将实现布里渊区及能带路径的可视化。相应的脚本文件(见文末代码)已经开放内测,感兴趣的可以到VASPKITFAQsQQ群331895604群文件下载。相比其它软件,VASPKIT可以更加方便实现自动标记对应能带计算的k点路径。
晶体结构的第一布里渊区可通过调用Voronoi图算法实现,特别感谢中科大郑奇靖老师提供了关于实现算法的建议。接下来我们演示如何调用VASPKIT结合python脚本实现布里渊区及能带路径的可视化:
第一步:我们以Cu和graphene为例,准备好POSCAR文件,注意一定是原胞结构(primtivecell)。然后分别运行vaspkit-303(Cu体相)和vaspkit-302(二维体系graphene)生成包含推荐能带路径的KPATH.in文件和高对称点HIGH_SYMMETRY_POINTS文件;
=====================StructuralOptions========================01)VASPInputFilesGenerator02)Elastic-Properties03)K-PathGenerator04)StructureEditor05)Catalysis-ElectroChemKit06)SymmetrySearch=====================ElectronicOptions========================11)Density-of-States21)DFTBand-Structure23)3DBand-Structure25)Hybrid-DFTBand-Structure26)Fermi-Surface28)Band-StructureUnfolding===========Charge&Potential&WavefunctionOptions===========31)Charge&SpinDensity42)Potential-Related51)Wave-FunctionAnalysis======================MiscUtilities===========================71)Optical-Properties72)Molecular-DynamicsKit73)VASP2otherInterface91)SemiconductorCalculator92)2D-MaterialsKit0)Quit------------>>302+--------------------------WarmTips--------------------------+SeeAnExampleinvaspkit/examples/seek_kpath/graphene_2D.MoredetailsaregiveninarXiv:1806.04285(2018)andrefs.Thisfeatureisstillexperimental&checkPRIMCELL.vaspfile.+---------------------------------------------------------------+-->>(01)ReadingStructuralParametersfromPOSCARFile...+--------------------------Summary----------------------------+ThevacuumslabissupposedtobealongzaxisPrototype:AB2TotalAtomsinInputCell:3LatticeConstantsinInputCell:3.1903.19012.000LatticeAnglesinInputCell:90.00090.000120.000TotalAtomsinPrimitiveCell:3LatticeConstantsinPrimitiveCell:3.1903.19012.000LatticeAnglesinPrimitiveCell:90.00090.000120.0002DBravaisLattice:HexagonalPointGroup:26[D3h]International:P-6m2SymmetryOperations:12SuggestedK-Path:(showninthenextline)[GAMMA-M-K-GAMMA]+---------------------------------------------------------------+-->>(02)WrittenPRIMCELL.vaspfile.-->>(03)WrittenHIGH_SYMMETRY_POINTSFileforReference.-->>(04)WrittenKPATH.inFileforBand-StructureCalculation.+----------------------------WARNING----------------------------+|DoNOTforgettocopyPRIMCELL.vasptoPOSCARunlessyouknow||whatyouaredoing.Otherwiseyoumightgetwrongresults!|+---------------------------------------------------------------+`
第二步,终端运行python3visualize_BZ.py(代码见下面)得到下图:
importnumpyasnp,numpy.linalgasnplfromscipy.spatialimportVoronoifromitertoolsimportproductfrommatplotlib.patchesimportFancyArrowPatchfrommpl_toolkits.mplot3dimportproj3dimportmatplotlibasmplmpl.rcParams['font.size']=12.
defread_poscar(poscar,species=None):poscar=open(poscar,'r')title=poscar.readline.stripscale=float(poscar.readline.strip)s=float(scale)lattice_vectors=[[float(v)forvinposcar.readline().split()],[float(v)forvinposcar.readline().split()],[float(v)forvinposcar.readline().split()]]lattice_vectors=np.array(lattice_vectors)reciprocal_lattice_vectors=np.linalg.inv(lattice_vectors).Treciprocal_lattice_vectors=reciprocal_lattice_vectors*np.pi*2returnreciprocal_lattice_vectors
defread_high_symmetry_points:f=open("HIGH_SYMMETRY_POINTS",'r')fa=f.readlinen=0hpts=whileTrue:fa=f.readlinehpts.append(fa)n=n+1iffa.find('use')>0:breakf.closekpts=klabels=foriinrange(0,n-3):kpts.append(hpts[i].split[0:3])klabels.append(hpts[i].split[3])kpts=np.array(kpts,dtype=np.float64)returnkpts,klabels
defis_greek_alphabets(klabels):Greek_alphabets=['Alpha','Beta','Gamma','Delta','Epsilon','Zeta','Eta','Theta','Iota','Kappa','Lambda','Mu','Nu','Xi','Omicron','Pi','Rho','Sigma','Tau','Upsilon','Phi','Chi','Psi','Pega']group_labels=foriinrange(len(klabels)):klabel=klabels[i]forjinrange(len(Greek_alphabets)):if(klabel.find(Greek_alphabets[j].upper)>=0):latex_exp=r''+'$\\'+str(Greek_alphabets[j])+'$'klabel=klabel.replace(str(Greek_alphabets[j].upper),str(latex_exp))if(klabel.find('_')>0):n=klabel.find('_')klabel=klabel[:n]+'$'+klabel[n:n+2]+'$'+klabel[n+2:]group_labels.append(klabel)klabels=group_labelsreturnklabels
defread_kpath:kpath=np.loadtxt("KPATH.in",dtype=np.string_,skiprows=4)#print(kpath)kpath_labels=kpath[:,3].tolistkpath_labels=[i.decode('utf-8','ignore')foriinkpath_labels]foriinrange(len(kpath_labels)):ifkpath_labels[i]=="Gamma":kpath_labels[i]=u"Γ"kpaths=np.zeros((len(kpath_labels),3),dtype=float)foriinrange(len(kpath_labels)):kpaths[i,:]=\[float(x)forxinkpath[i][0:3]]returnkpath_labels,kpaths
defget_Wigner_Seitz_BZ(lattice_vectors):#Inspiredbyhttp://www.thp.uni-koeln.de/trebst/Lectures/SolidState-2016/wigner_seitz_3d.py#Inspiredbyhttps://github.com/QijingZheng/VASP_FermiSurface/blob/master/fs.pylatt=prefactors=[0.,-1.,1.]forpinprefactors:foruinlattice_vectors:latt.append(p*u)lattice=forvsinproduct(latt,latt,latt):a=vs[0]+vs[1]+vs[2]ifnotany((a==x).allforxinlattice):lattice.append(a)voronoi=Voronoi(lattice)bz_facets=bz_ridges=bz_vertices=forpid,ridinzip(voronoi.ridge_points,voronoi.ridge_vertices):if(pid[0]==0orpid[1]==0):bz_ridges.append(voronoi.vertices[np.r_[rid,[rid[0]]]])bz_facets.append(voronoi.vertices[rid])bz_vertices+=ridbz_vertices=list(set(bz_vertices))returnvoronoi.vertices[bz_vertices],bz_ridges,bz_facets
classArrow3D(FancyArrowPatch):def__init__(self,xs,ys,zs,*args,**kwargs):FancyArrowPatch.__init__(self,(0,0),(0,0),*args,**kwargs)self._verts3d=xs,ys,zsdefdraw(self,renderer):xs3d,ys3d,zs3d=self._verts3dxs,ys,zs=proj3d.proj_transform(xs3d,ys3d,zs3d,renderer.M)self.set_positions((xs[0],ys[0]),(xs[1],ys[1]))FancyArrowPatch.draw(self,renderer)
defvisualize_BZ_matplotlib(points,ridges,facets,reciprocal_lattice_vectors,kpts,klabels,kpaths):importmatplotlibasmplimportmatplotlib.pyplotaspltfrommpl_toolkits.mplot3dimportAxes3Dfrommpl_toolkits.mplot3d.art3dimportPoly3DCollectionfig=plt.figure(figsize=(8,8))ax=fig.add_subplot(111,projection='3d')basis_vector_clrs=['r','g','b']basis_vector_labs=['x','y','z']foriiinrange(3):arrow=Arrow3D([0,reciprocal_lattice_vectors[ii,0]],[0,reciprocal_lattice_vectors[ii,1]],[0,reciprocal_lattice_vectors[ii,2]],color=basis_vector_clrs[ii],mutation_scale=20,lw=1,arrowstyle="->")ax.add_artist(arrow)ax.text(reciprocal_lattice_vectors[ii,0],reciprocal_lattice_vectors[ii,1],reciprocal_lattice_vectors[ii,2],basis_vector_labs[ii])foririnridges:ax.plot(ir[:,0],ir[:,1],ir[:,2],color='k',lw=1.5,alpha=0.5)foriinrange(len(klabels)):kpt=np.dot(kpts[i,:],reciprocal_lattice_vectors)ax.scatter(kpt[0],kpt[1],kpt[2],c='b',marker='o',s=20,alpha=0.8)ax.text(kpt[0],kpt[1],kpt[2],klabels[i],c='b')foriinrange(kpaths.shape[0]):kpaths[i,:]=np.dot(kpaths[i,:],reciprocal_lattice_vectors)foriinrange(0,kpaths.shape[0],2):arrow=Arrow3D([kpaths[i,0],kpaths[i+1,0]],[kpaths[i,1],kpaths[i+1,1]],[kpaths[i,2],kpaths[i+1,2]],mutation_scale=20,lw=1.5,arrowstyle="->",color="magenta")ax.add_artist(arrow)ax.set_axis_offax.view_init(elev=12,azim=23)plt.savefig('Brillouin_Zone.pdf',dpi=100)plt.show
defwelcome:print('')print('+---------------------------------------------------------------+')print('|AVASPKITPlugintoVisualizeBrillouinZoneUsingMatplotlib|')print('|WrittenbyVeiWANG(wangvei@icloud.com)|')print('+---------------------------------------------------------------+')print('')
if__name__=="__main__":welcomereciprocal_lattice_vectors=read_poscar('POSCAR')lattice_vectors=[np.array(reciprocal_lattice_vectors[0,:]),np.array(reciprocal_lattice_vectors[1,:]),np.array(reciprocal_lattice_vectors[2,:])]kpts,klabels=read_high_symmetry_pointsklabels=is_greek_alphabets(klabels)kpath_labels,kpaths=read_kpathpoints,ridges,facets=get_Wigner_Seitz_BZ(lattice_vectors)visualize_BZ_matplotlib(points,ridges,facets,reciprocal_lattice_vectors,kpts,klabels,kpaths)
布里渊区和布里渊区及能带路径的可视化的问题分享结束啦,以上的文章解决了您的问题吗?欢迎您下次再来哦!
还没有评论,来说两句吧...