应《网络安全法》要求，自2017年10月1日起，未进行实名认证将不得使用互联网跟帖服务。为保障您的帐号能够正常使用，请尽快对帐号进行手机号验证，感谢您的理解与支持！

24小时热门版块排行榜    

小木虫论坛-学术科研互动平台 » 计算模拟区 » 第一性原理 » Vasp&MedeA » 用dos_procar处理dos数据不准
51/1返回列表
查看: 1604  |  回复: 3
只看楼主 @他人 存档 新回复提醒 (忽略) 收藏    在APP中查看
当前只显示满足指定条件的回帖，点击这里查看本话题的所有回帖

dx0620

木虫 (正式写手)

[求助] 用dos_procar处理dos数据不准 已有1人参与

我的vasp版本是5.3，用的dos_procar处理的PROCAR，想分别得到s，p，d的态密度，但是最后得到的态密度曲线不对，总态密度跟s，p，d的对不上，如下图

用dos_procar处理dos数据不准

我的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,*) gaussian
      write(6,'("ISPIN(1/2) = ?  "')                ;  read(5,*) ispin
      write(6,'("How many atoms are integrated? "') ;  read(5,*) Nna
      allocate(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 do
          do 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
回复此楼

» 猜你喜欢

» 本主题相关价值贴推荐，对您同样有帮助:
1楼 2014-05-02 10:38:18
已阅   关注TA 给TA发消息 送TA红花 TA的回帖

dx0620

木虫 (正式写手)
引用回帖:
2楼: Originally posted by mywai520 at 2014-05-02 12:24:18
换一个处理脚本吧，最好用p4vasp,直观图形化。若是不想也可以用vaspkit，所有的处理打包，是个不错的选择。

您好，用p4vasp的话是怎么处理啊？我用vaspkit处理的话，ipdos是什么啊？为什么感觉tdos跟pdos里面s、p、d加起来差好多，这几个都要除以原子数吗？！谢谢~
赞一下 回复此楼
3楼2014-05-02 15:03:31
已阅   关注TA 给TA发消息 送TA红花 TA的回帖
查看全部 4 个回答

mywai520

铁杆木虫 (著名写手)

【答案】应助回帖

★ ★ ★ ★ ★
感谢参与，应助指数 +1
dx0620: 金币+5, ★有帮助 2014-05-04 09:29:41
换一个处理脚本吧，最好用p4vasp,直观图形化。若是不想也可以用vaspkit，所有的处理打包，是个不错的选择。
赞一下 回复此楼
2楼2014-05-02 12:24:18
已阅   关注TA 给TA发消息 送TA红花 TA的回帖

mywai520

铁杆木虫 (著名写手)
引用回帖:
3楼: Originally posted by dx0620 at 2014-05-02 15:03:31
您好，用p4vasp的话是怎么处理啊？我用vaspkit处理的话，ipdos是什么啊？为什么感觉tdos跟pdos里面s、p、d加起来差好多，这几个都要除以原子数吗？！谢谢~...

你给出的是每一个原子的spd态，要全部原子的加起来才是总的
赞一下 回复此楼
4楼2014-05-02 16:50:24
已阅   关注TA 给TA发消息 送TA红花 TA的回帖
查看全部 4 个回答
最具人气热帖推荐 [查看全部] 作者 回/看 最后发表
[考研] 材料与化工304求B区调剂 +6 邱gl 2026-03-11 7/350 2026-03-17 01:27 by 学员FnSWZj
[考研] 326求调剂 +4 诺贝尔化学奖觊� 2026-03-15 7/350 2026-03-16 17:11 by 诺贝尔化学奖觊�
[考研] 304求调剂 +5 素年祭语 2026-03-15 5/250 2026-03-16 17:00 by 我的船我的海
[考研] 085600材料与化工 求调剂 +13 enenenhui 2026-03-13 14/700 2026-03-16 15:19 by 了了了了。。
[考研] 材料与化工求调剂 +3 为学666 2026-03-16 3/150 2026-03-16 15:09 by 加号+
[基金申请] NSFC申报书里申请人简历中代表性论著还需要在申报书最后的附件里面再上传一遍吗 20+5 NSFC2026我来了 2026-03-10 14/700 2026-03-15 23:53 by 不负韶华的虎
[考研] 268求调剂 +5 一定有学上- 2026-03-14 6/300 2026-03-14 22:20 by 运气yunqi
[考研] 材料与化工 323 英一+数二+物化，一志愿：哈工大 本人本科双一流 +4 自由的_飞翔 2026-03-13 5/250 2026-03-14 19:39 by hmn_wj
[考研] 26考研调剂 +3 ying123. 2026-03-10 3/150 2026-03-14 00:18 by JourneyLucky
[考研] 279求调剂 +3 抓着星星的女孩 2026-03-10 3/150 2026-03-13 23:47 by userper
[考研] 材料与化工（0856）304求B区调剂 +6 邱gl 2026-03-12 7/350 2026-03-13 23:24 by 邱gl
[考研] 一志愿中科院，化学方向，295求调剂 +4 一氧二氮 2026-03-11 4/200 2026-03-13 22:35 by JourneyLucky
[考研] 求材料调剂 085600英一数二总分302 前三科235 精通机器学习 一志愿哈工大 +4 林yaxin 2026-03-12 4/200 2026-03-13 22:04 by 星空星月
[考研] 26调剂/材料科学与工程/总分295/求收留 +9 2026调剂侠 2026-03-12 9/450 2026-03-13 20:46 by 18595523086
[考研] 070303一志愿西北大学学硕310找调剂 +3 d如愿上岸 2026-03-13 3/150 2026-03-13 10:43 by houyaoxu
[考研] 296求调剂 +3 大口吃饭 身体健 2026-03-13 3/150 2026-03-13 10:31 by 学员8dgXkO
[考研] 283求调剂，材料、化工皆可 +8 苏打水7777 2026-03-11 10/500 2026-03-13 09:06 by Linda Hu
[考博] 2026年博士申请 +3 QwQwQW10 2026-03-11 3/150 2026-03-12 17:58 by gxch43
[考研] 341求调剂 +4 捣蛋猪猪 2026-03-11 4/200 2026-03-12 14:47 by ruiyingmiao
[考研] 研究生招生 +3 徐海涛11 2026-03-10 7/350 2026-03-12 14:26 by 徐海涛11
信息提示
关闭
请填处理意见
关闭