| 查看: 1607 | 回复: 3 | |||
dx0620木虫 (正式写手)
|
[求助]
用dos_procar处理dos数据不准 已有1人参与
|
我的vasp版本是5.3，用的dos_procar处理的PROCAR，想分别得到s，p，d的态密度，但是最后得到的态密度曲线不对，总态密度跟s，p，d的对不上，如下图 我的dos_procar.f90的代码如下 !=========================================================================== ! INPUT FILES : OUTCAR , PROCAR !=========================================================================== implicit real*8(a-h,o-z) character*64 tmp dimension ddosP(4),ddosT(4) allocatable wt(  ,oc(:,:,:,,eigen(:,,ee(,na(allocatable gdosP(:, ,gdosT(:,pi=4.0d0*atan(1.0d0) ; netot = 1001 allocate(ee(netot),gdosP(4,netot),gdosT(4,netot)) open(1,file='OUTCAR') !================================================================= nline = 0 ; open(99,status='scratch') 111 read(1,'(a64)',end=222) tmp if(tmp(2:8) .eq. 'E-fermi') nline=nline+1 if(tmp(2:8) .eq. 'E-fermi') write(99,'(a64)') tmp go to 111 222 rewind(99) do i=1,nline ; read(99,*) tmp,tmp,fermi ; end do close(99) !================================================================= write(6,'("fermi energy = ",f10.5)') fermi write(6,'("enter the gaussian "  ') ; read(5,*) gaussianwrite(6,'("ISPIN(1/2) = ? " ') ; read(5,*) ispinwrite(6,'("How many atoms are integrated? " ') ; read(5,*) Nnaallocate(na(Nna)) write(6,*) ' atom number :' ; read(5,*) na open( 7,file='PROCAR') if(ispin .eq. 1) then open(11,file='gdosP-s.dat') ; open(51,file='gdosT-s.dat') open(12,file='gdosP-p.dat') ; open(52,file='gdosT-p.dat') open(13,file='gdosP-d.dat') ; open(53,file='gdosT-d.dat') open(14,file='gdosP-A.dat') ; open(54,file='gdosT-A.dat') else if(ispin .eq. 2) then open(11,file='gdosP-up-s.dat') ; open(51,file='gdosT-up-s.dat') open(12,file='gdosP-up-p.dat') ; open(52,file='gdosT-up-p.dat') open(13,file='gdosP-up-d.dat') ; open(53,file='gdosT-up-d.dat') open(14,file='gdosP-up-A.dat') ; open(54,file='gdosT-up-A.dat') open(21,file='gdosP-dn-s.dat') ; open(61,file='gdosT-dn-s.dat') open(22,file='gdosP-dn-p.dat') ; open(62,file='gdosT-dn-p.dat') open(23,file='gdosP-dn-d.dat') ; open(63,file='gdosT-dn-d.dat') open(24,file='gdosP-dn-A.dat') ; open(64,file='gdosT-dn-A.dat') else write(6,*) 'Input ERROR, you must input 1 or 2 ' ; stop end if !=========================================================================== do is=1,ispin if(is .eq. 1) read(7,'(a64)') tmp read(7,'(16x,i3,20x,i4,19x,i4)') nk,nband,nion if(is .eq. 1) then write(6,'("number of kpoints: ",i4)') nk write(6,'("number of bands: ",i4)') nband write(6,'("number of ions: ",i4)') nion end if emin = 9999.0 ; emax= -9999.0 allocate(wt(nk),eigen(nk,nband),oc(nk,nband,nion+1,4)) do k=1,nk read(7,'(a64)') tmp read(7,'(10x,i3,5x,3f11.8,13x,f11.8)') kp,pt1,pt2,pt3,wt(k) read(7,'(a64)') tmp do nb=1,nband read(7,'(8x,9x,f14.8,7x,f12.8)') eigen(k,nb),occ eigen(k,nb) = eigen(k,nb) - fermi if(eigen(k,nb) .gt. emax) emax=eigen(k,nb) if(eigen(k,nb) .lt. emin) emin=eigen(k,nb) read(7,'(a64)') tmp ; read(7,'(a64)') tmp niont = nion +1 if(nion .eq. 1) niont = 1 do ion = 1,niont read(7,'(3x,4f7.3)') (oc(k,nb,ion,j),j=1,4) end do read(7,'(a64)') tmp end do end do ! ############################################################## weight = 0.0 ; gdosP=0.0d0 ; gdosT=0.0d0 do k=1,nk ; weight = weight + wt(k) ; end do do k=1,nk ; wt(k) = wt(k) / weight ; end do estart_ev = int(emin -5.0) eend_ev = int(emax +5.0) de_ev = (eend_ev - estart_ev)/(netot-1) do ne = 1,netot ; ee(ne) = estart_ev + (ne-1) * de_ev ; end do ascal = 1.0/(gaussian*sqrt(pi)) ! ############################################################### do k=1,nk ; do nb=1,nband do i=1,4 ; ddosT(i) = oc(k,nb,niont,i)*wt(k) ; end do do i=1,4 ; ddosP(i) = sum(oc(k,nb,na( ,i)*wt(k)) ; end dodo ne=1,netot dij = ( eigen(k,nb)-ee(ne) )**2 / (gaussian**2) do i=1,4 gdosT(i,ne) = gdosT(i,ne) + ascal*ddosT(i)*exp(-dij) gdosP(i,ne) = gdosP(i,ne) + ascal*ddosP(i)*exp(-dij) end do end do end do ; end do ! ############################################################### if(is .eq. 1) ns=+1 ; if(is .eq. 2) ns=-1 do i=1,4 ; write(is*10+i ,1) (ee(j),ns*gdosP(i,j),j=1,netot) ; end do do i=1,4 ; write(is*10+i+40,1) (ee(j),ns*gdosT(i,j),j=1,netot) ; end do deallocate(wt,eigen,oc) end do !=========================================================================== 1 format(f10.5,f12.4) end ############################################################################ INCAR 如下： ENCUT = 400 ISTART = 1 ICHARG = 11 ISMEAR = 1 LORBIT = 11 #EMIN = -25 #EMAX = 22.5 NBANDS = 64 #SIGMA = 0.1 #EDIFF = 1E-6 #EDIFFG = -1E-3 #NSW = 200 #IBRION = 2 #POTIM = 0.2 #ISIF = 3 #ALGO = Normal #LREAL = Auto PREC = Accurate PSTRESS = 500.00 |
回复此楼» 猜你喜欢
293求调剂
已经有12人回复
-
070300化学319求调剂
已经有4人回复
-
268求调剂
已经有9人回复
-
301求调剂
已经有4人回复
-
0703化学336分求调剂
已经有5人回复
-
301求调剂
已经有4人回复
-
材料专硕306英一数二
已经有9人回复
-
304求调剂
已经有8人回复
-
材料工程专硕调剂
已经有5人回复
-
344求调剂
已经有4人回复
» 本主题相关价值贴推荐，对您同样有帮助:
- dos-procar.f90
已经有94人回复
- 求解split_dos分割文件该怎么理解！？
已经有11人回复
- vasp计算，DOS图的分析，分波态密度的数据分析~·
已经有4人回复
- 请教split_dos处理DOS的相关问题
已经有4人回复
- 【分享】 最新VASP脚本大全 支持VASP 5.3 (有split_dos; VP; Sumdos; sumup_dos 等)
已经有615人回复
- split_dos分割以后，如何找到对应的原子？
已经有11人回复
- 【求助成功】split_dos分割DOSCAR的问题
已经有14人回复
- 【求助】那位强人能在c++程序中调用DOS命令行啊！指点一下呗（悬赏给二十金币）
已经有5人回复
- 【求助成功】怎么使用DOSCAR画出DOS图形
已经有23人回复
- 【求助】murn.f gk.f和pbnd.f split_dos和vp程序
已经有17人回复
mywai520
铁杆木虫 (著名写手)
- 应助: 339 (大学生)
- 金币: 10425.6
- 红花: 45
- 帖子: 1664
- 在线: 296.9小时
- 虫号: 265142
- 注册: 2006-07-09
- 专业: 凝聚态物性 II ：电子结构
2楼2014-05-02 12:24:18
dx0620
木虫 (正式写手)
- 应助: 14 (小学生)
- 金币: 4053.6
- 散金: 5
- 红花: 6
- 帖子: 393
- 在线: 910小时
- 虫号: 1303311
- 注册: 2011-05-22
- 性别: GG
- 专业: 凝聚态物性 II ：电子结构
3楼2014-05-02 15:03:31
mywai520
铁杆木虫 (著名写手)
- 应助: 339 (大学生)
- 金币: 10425.6
- 红花: 45
- 帖子: 1664
- 在线: 296.9小时
- 虫号: 265142
- 注册: 2006-07-09
- 专业: 凝聚态物性 II ：电子结构
4楼2014-05-02 16:50:24
20