Skip to content

Commit

Permalink
change the plot module
Browse files Browse the repository at this point in the history
  • Loading branch information
@xh125
xh125 committed Sep 17, 2021
1 parent ccfa68c commit 07c4a97
Show file tree
Hide file tree
Showing 2 changed files with 116 additions and 33 deletions.
7 changes: 5 additions & 2 deletions src_complex/lvcsh.f90
View file
Original file line number Diff line number Diff line change
Expand Up @@ -93,6 +93,7 @@ program lvcsh
psit_e_file,psit_h_file,save_psit, read_psit,&
plot_csit,plot_wsit,plot_psit,&
band_e_file,band_h_file,&
e_temp_file,h_temp_file,plot_eh_temp,&
plot_band_occupatin_withtime,plot_ph_temp
implicit none

Expand Down Expand Up @@ -458,7 +459,7 @@ program lvcsh
do ifre=1,nefre
csit_e(ifre,isnap) = csit_e(ifre,isnap)+REAL(c_e(ifre)*CONJG(c_e(ifre)))
wsit_e(ifre,isnap) = wsit_e(ifre,isnap)+REAL(w_e(ifre)*CONJG(w_e(ifre)))
psit_e(ifre,isnap) = psit_e(ifre,isnap)+P_e(ifre,iesurface)**2
psit_e(ifre,isnap) = psit_e(ifre,isnap)+ABS(P_e(ifre,iesurface))**2
pes_e( ifre,isnap) = pes_e( ifre,isnap) + E_e(ifre)
if(iaver == 1) pes_one_e(ifre,isnap) = E_e(ifre)
enddo
Expand All @@ -470,7 +471,7 @@ program lvcsh
do ifre=1,nhfre
csit_h(ifre,isnap) = csit_h(ifre,isnap)+REAL(c_h(ifre)*CONJG(c_h(ifre)))
wsit_h(ifre,isnap) = wsit_h(ifre,isnap)+REAL(w_h(ifre)*CONJG(w_h(ifre)))
psit_h(ifre,isnap) = psit_h(ifre,isnap)+P_h(ifre,iesurface)**2
psit_h(ifre,isnap) = psit_h(ifre,isnap)+ABS(P_h(ifre,ihsurface))**2
pes_h( ifre,isnap) = pes_h( ifre,isnap)-E_h(ifre)
if(iaver == 1) pes_one_h(ifre,isnap) = -E_h(ifre)
enddo
Expand Down Expand Up @@ -604,6 +605,7 @@ program lvcsh
call plot_wsit(nefre,nsnap,naver,wsit_e,wsit_e_file)
call plot_psit(nefre,nsnap,naver,psit_e,psit_e_file)
call plot_band_occupatin_withtime(neband,nktotf,Enk_e,xkf,nsnap,psit_e,csit_e,savedsnap,band_e_file)
call plot_eh_temp(neband,nktotf,Enk_e,xkf,nsnap,psit_e,csit_e,savedsnap,e_temp_file)
deallocate(pes_e,pes_one_e,csit_e,wsit_e,psit_e)

endif
Expand All @@ -615,6 +617,7 @@ program lvcsh
call plot_wsit(nhfre,nsnap,naver,wsit_h,wsit_h_file)
call plot_psit(nhfre,nsnap,naver,psit_h,psit_h_file)
call plot_band_occupatin_withtime(nhband,nktotf,Enk_h,xkf,nsnap,psit_h,csit_h,savedsnap,band_h_file)
call plot_eh_temp(nhband,nktotf,Enk_h,xkf,nsnap,psit_h,csit_h,savedsnap,h_temp_file)
deallocate(pes_h,pes_one_h,csit_h,wsit_h,psit_h)
endif

Expand Down
142 changes: 111 additions & 31 deletions src_complex/saveinf.f90
View file
Original file line number Diff line number Diff line change
Expand Up @@ -22,12 +22,14 @@ module saveinf
character(len=maxlen) :: phU_file = "phUsit.dat"
character(len=maxlen) :: phK_file = "phKsit.dat"
character(len=maxlen) :: phT_file = "phLOtemp.dat"
character(len=maxlen) :: band_e_file = "band_e.dat"
character(len=maxlen) :: band_h_file = "band_h.dat"
character(len=maxlen) :: band_e_file = "band_e"
character(len=maxlen) :: band_h_file = "band_h"
character(len=maxlen) :: e_temp_file = "temp_e"
character(len=maxlen) :: h_temp_file = "temp_h"

character(len=12) :: ctmp1,ctmp2
character(len=12),allocatable :: cphmode(:,:)
character(len=12),allocatable :: cefree(:)
character(len=12),allocatable :: cefree(:),cnk(:,:)

contains

Expand Down Expand Up @@ -121,7 +123,6 @@ subroutine plot_phQ(nmodes,nq,nsnap,phQsit)
write(phq_unit,"(A)") "Average of Normal mode coordinate for one core sample. ((phQ(imode,iq),imode=1,nmodes),iq=1,nq)"
write(phq_unit,"(A12,*(2(1X,A12)))") "time",(("Re[Q"//trim(adjustl(cphmode(imode,iq)))//"]",&
"Im[Q"//trim(adjustl(cphmode(imode,iq)))//"]",imode=1,nmodes),iq=1,nq)
!write(phq_unit,"(A12,*(2(1X,A12)))") "time",(("Re[phQ]","Im[phQ]",imode=1,nmodes),iq=1,nq)
write(phq_unit,"(A12,*(2(1X,A12)))") "fs ",((" a.u. "," a.u. ",imode=1,nmodes),iq=1,nq)
write(phq_unit,"(A12,*(2(1X,F12.5)))") "Omega(meV)",(((wf(imode,iq)*ryd2meV,wf(imode,iq)*ryd2meV),imode=1,nmodes),iq=1,nq)
do isnap=0,nsnap
Expand Down Expand Up @@ -157,7 +158,7 @@ subroutine save_phP(nmodes,nq,nsnap,phPsit)
enddo


write(php_unit,"(5X,A)") "Average of Normal mode verlocity for one core sample.((phP(imode,iq),imode=1,nmodes),iq=1,nq)"
write(php_unit,"(A)") "Average of Normal mode verlocity for one core sample.((phP(imode,iq),imode=1,nmodes),iq=1,nq)"
write(php_unit,"(A12,*(2(1X,A12)))") "time",(("Re[P"//trim(adjustl(cphmode(imode,iq)))//"]",&
"Im[P"//trim(adjustl(cphmode(imode,iq)))//"]",imode=1,nmodes),iq=1,nq)
!write(php_unit,"(A12,*(2(1X,A12)))") "time",(("Re[phP]","Im[phP]",imode=1,nmodes),iq=1,nq)
Expand Down Expand Up @@ -223,7 +224,7 @@ subroutine plot_phP(nmodes,nq,nsnap,phPsit)
enddo
enddo

write(php_unit,"(5X,A)") "Average of Normal mode verlocity for one core sample.((phP(imode,iq),imode=1,nmodes),iq=1,nq)"
write(php_unit,"(A)") "Average of Normal mode verlocity for one core sample.((phP(imode,iq),imode=1,nmodes),iq=1,nq)"
write(php_unit,"(A12,*(2(1X,A12)))") "time",(("Re[P"//trim(adjustl(cphmode(imode,iq)))//"]",&
"Im[P"//trim(adjustl(cphmode(imode,iq)))//"]",imode=1,nmodes),iq=1,nq)
!write(php_unit,"(A12,*(2(1X,A12)))") "time",(("Re[phP]","Im[phP]",imode=1,nmodes),iq=1,nq)
Expand Down Expand Up @@ -261,7 +262,7 @@ subroutine save_phK(nmodes,nq,nsnap,phKsit)
enddo
enddo

write(phK_unit,"(5X,A)") "Average of Normal mode kinetic energy for one core sample.SUM_phK,((phK(imode,iq),imode=1,nmodes),iq=1,nq)"
write(phK_unit,"(A)") "Average of Normal mode kinetic energy for one core sample.SUM_phK,((phK(imode,iq),imode=1,nmodes),iq=1,nq)"
write(phK_unit,"(*(1X,A12))") "time","SUM_phK",(("K"//trim(adjustl(cphmode(imode,iq))),imode=1,nmodes),iq=1,nq)
!write(phK_unit,"(*(1X,A12))") "time ","SUM_phK",(("phK(mode,q)",imode=1,nmodes),iq=1,nq)
write(phK_unit,"(*(1X,A12))") "fs "," meV ",((" meV ",imode=1,nmodes),iq=1,nq)
Expand Down Expand Up @@ -359,9 +360,9 @@ subroutine plot_phK(nmodes,nq,nsnap,phKsit)
endif

if(l_ph_quantum) then
write(phK_unit,"(5X,A)") "Average of Normal mode kinetic energy - 0.5*nqv*hbar*wqv for all trajecotry.SUM_phK,((phK(imode,iq),imode=1,nmodes),iq=1,nq)"
write(phK_unit,"(A)") "Average of Normal mode kinetic energy - 0.5*nqv*hbar*wqv for all trajecotry.SUM_phK,((phK(imode,iq),imode=1,nmodes),iq=1,nq)"
else
write(phK_unit,"(5X,A)") "Average of Normal mode kinetic energy - 0.5*KB*T for all trajecotry.SUM_phK,((phK(imode,iq),imode=1,nmodes),iq=1,nq)"
write(phK_unit,"(A)") "Average of Normal mode kinetic energy - 0.5*KB*T for all trajecotry.SUM_phK,((phK(imode,iq),imode=1,nmodes),iq=1,nq)"
endif
write(phK_unit,"(*(1X,A12))") "time ","phK(mode,q)",(("K"//trim(adjustl(cphmode(imode,iq))),imode=1,nmodes),iq=1,nq)
write(phK_unit,"(*(1X,A12))") "fs "," meV ",((" meV ",imode=1,nmodes),iq=1,nq)
Expand Down Expand Up @@ -400,7 +401,7 @@ subroutine save_phU(nmodes,nq,nsnap,phUsit)
enddo
enddo

write(phU_unit,"(5X,A)") "Average of Normal mode potential energy for one core sample.SUM_phU,((phU(imode,iq),imode=1,nmodes),iq=1,nq)"
write(phU_unit,"(A)") "Average of Normal mode potential energy for one core sample.SUM_phU,((phU(imode,iq),imode=1,nmodes),iq=1,nq)"
write(phU_unit,"(*(1X,A12))") "time ","SUM_phU",(("U"//trim(adjustl(cphmode(imode,iq))),imode=1,nmodes),iq=1,nq)
write(phU_unit,"(*(1X,A12))") "fs "," meV ",((" meV ",imode=1,nmodes),iq=1,nq)
write(phU_unit,"(2(1X,A12),*(1X,F12.5))") "Omega(meV)","SUM_phU",((wf(imode,iq)*ryd2meV,imode=1,nmodes),iq=1,nq)
Expand Down Expand Up @@ -481,9 +482,9 @@ subroutine plot_phU(nmodes,nq,nsnap,phUsit)
endif

if(l_ph_quantum) then
write(phU_unit,"(5X,A)") "Average of Normal mode potential energy - 0.5*nqv*hbar*wqv for all trajecotry.SUM_phU,((phU(imode,iq),imode=1,nmodes),iq=1,nq)"
write(phU_unit,"(A)") "Average of Normal mode potential energy - 0.5*nqv*hbar*wqv for all trajecotry.SUM_phU,((phU(imode,iq),imode=1,nmodes),iq=1,nq)"
else
write(phU_unit,"(5X,A)") "Average of Normal mode potential energy - 0.5*KB*T for all trajecotry.SUM_phU,((phU(imode,iq),imode=1,nmodes),iq=1,nq)"
write(phU_unit,"(A)") "Average of Normal mode potential energy - 0.5*KB*T for all trajecotry.SUM_phU,((phU(imode,iq),imode=1,nmodes),iq=1,nq)"
endif
write(phU_unit,"(*(1X,A12))") "time ","phU(mode,q)",(("U"//trim(adjustl(cphmode(imode,iq))),imode=1,nmodes),iq=1,nq)
write(phU_unit,"(*(1X,A12))") "fs "," meV ",((" meV ",imode=1,nmodes),iq=1,nq)
Expand Down Expand Up @@ -675,7 +676,7 @@ subroutine plot_wsit(nfre,nsnap,naver,wsit,wsit_file)

write(wsit_unit,"(A)") "The average electron(hole) population on different adiabatic wave function for all trajecotry. "
write(wsit_unit,"(*(1X,A12))") "time ",("wsit"//trim(adjustl(cefree(ifre))),ifre=1,nfre)
write(wsit_unit,"((1X,A12),*(1X,I12))") "fs ",(ifre,ifre=1,nfre)
write(wsit_unit,"((1X,A12),*(1X,A12))") "fs ",("|w|2",ifre=1,nfre)
do isnap=0,nsnap
write(wsit_unit,"(*(1X,E12.5))") dt*nstep*isnap*ry_to_fs,(wsit(ifre,isnap),ifre=1,nfre)
enddo
Expand Down Expand Up @@ -860,9 +861,15 @@ subroutine save_csit(nfre,nsnap,naver,csit,csit_file)
csit_unit = io_file_unit()
csit_file_ = trim(outdir)//trim(adjustl(csit_file))
call open_file(csit_file_,csit_unit)

if(.not. allocated(cefree)) allocate(cefree(nfre))
do ifre=1,nfre
write(ctmp1,*) ifre
cefree(ifre) = "("//trim(adjustl(ctmp1))//")"
enddo

write(csit_unit,"(A)") "The average electron(hole) population on different diabatic wave function for one core sample. "
write(csit_unit,"(*(1X,A12))") "time",("csit(ifre)",ifre=1,nfre)
write(csit_unit,"(*(1X,A12))") "time",("csit"//trim(adjustl(cefree(ifre))),ifre=1,nfre)
write(csit_unit,"((1X,A12),*(1X,I12))") "fs ",(ifre,ifre=1,nfre)
do isnap=0,nsnap
write(csit_unit,"(*(1X,E12.5))") dt*nstep*isnap*ry_to_fs,(csit(ifre,isnap),ifre=1,nfre)
Expand All @@ -872,6 +879,7 @@ subroutine save_csit(nfre,nsnap,naver,csit,csit_file)
write(stdout,"(A,A)")"Save the average electron(hole) wavefction population on &
&different diabatic wave function to file:",trim(csit_file_)

deallocate(cefree)
end subroutine save_csit

subroutine read_csit(inode,icore,nfre,nsnap,naver,csit,csit_file)
Expand Down Expand Up @@ -917,10 +925,16 @@ subroutine plot_csit(nfre,nsnap,naver,csit,csit_file)
csit_file_ = trim(outdir)//trim(adjustl(csit_file))
csit_unit = io_file_unit()
call open_file(csit_file_,csit_unit)

if(.not. allocated(cefree)) allocate(cefree(nfre))
do ifre=1,nfre
write(ctmp1,*) ifre
cefree(ifre) = "("//trim(adjustl(ctmp1))//")"
enddo

write(csit_unit,"(A)") "The electron(hole) population on different diabatic wave function. "
write(csit_unit,"(*(1X,A12))") "time",("csit(ifre)",ifre=1,nfre)
write(csit_unit,"((1X,A12),*(1X,I12))") "fs ",(ifre,ifre=1,nfre)
write(csit_unit,"(*(1X,A12))") "time",("csit"//trim(adjustl(cefree(ifre))),ifre=1,nfre)
write(csit_unit,"((1X,A12),*(1X,A12))") "fs ",("|c|2",ifre=1,nfre)
do isnap=0,nsnap
write(csit_unit,"(*(1X,E12.5))") dt*nstep*isnap*ry_to_fs,(csit(ifre,isnap),ifre=1,nfre)
enddo
Expand All @@ -930,11 +944,12 @@ subroutine plot_csit(nfre,nsnap,naver,csit,csit_file)
write(stdout,"(A,A)")"Write the average electron(hole) wavefction population on &
&different diabatic wave function to file:",trim(csit_file_)

deallocate(cefree)
end subroutine plot_csit



!The active adiabatic PES project to diabatic states. psit(ifre)=P(ifre,isurface)**2
!The active adiabatic PES project to diabatic states. psit(ifre)=|P(ifre,isurface)|**2
subroutine save_psit(nfre,nsnap,naver,psit,psit_file)
implicit none
integer,intent(in) :: nfre,nsnap,naver
Expand Down Expand Up @@ -1021,7 +1036,7 @@ subroutine plot_psit(nfre,nsnap,naver,psit,psit_file)

write(psit_unit,"(A)") "The average active PES project on different diabatic wave function for all trajecotry. psit(ifre)=P(ifre,isurface)**2 "
write(psit_unit,"(*(1X,A12))") "time",("psit"//trim(adjustl(cefree(ifre))),ifre=1,nfre)
write(psit_unit,"((1X,A12),*(1X,I12))") "fs ",(ifre,ifre=1,nfre)
write(psit_unit,"((1X,A12),*(1X,A12))") "fs ",("|P|2",ifre=1,nfre)
do isnap=0,nsnap
write(psit_unit,"(*(1X,E12.5))") dt*nstep*isnap*ry_to_fs,(psit(ifre,isnap),ifre=1,nfre)
enddo
Expand All @@ -1047,26 +1062,91 @@ subroutine plot_band_occupatin_withtime(nband,nk,Enk,xk,nsnap,psit,csit,dsnap,ba
character(len=maxlen) :: band_file_
integer :: band_unit
integer :: ipol,iband,ik,ifre

band_file_=trim(outdir)//trim(adjustl(band_file))
band_unit = io_file_unit()
call open_file(band_file_,band_unit)
write(band_unit,"(A)") "Carrier occupation on the band structure at different time"
write(band_unit,"(*(1X,A12))") "iband","ik","kx","ky","kz","E_nk",("psit","wsit",isnap=0,nsnap,dsnap)
write(band_unit,"(6(1X,A12),*(1X,F12.2))") "index_band"," index_k","2pi/a0","2pi/a0", "2pi/a0", " eV" ,&
real(kind=dp) :: kx,ky,kz

do iband=1,nband
write(ctmp1,*) iband
band_file_=trim(outdir)//trim(adjustl(band_file))//"_"//trim(adjustl(ctmp1))//".dat"
band_unit = io_file_unit()
call open_file(band_file_,band_unit)
write(band_unit,"(A32,I12,A)") "Carrier occupation on the iband=",iband," structure at different time"
write(band_unit,"(*(1X,A12))") "iband","ik","kx","ky","kz","E_nk",("psit","csit",isnap=0,nsnap,dsnap)
write(band_unit,"(6(1X,A12),*(1X,F12.2))") "index_band"," index_k","b1","b2", "b3", " eV" ,&
(dt*nstep*isnap*ry_to_fs,dt*nstep*isnap*ry_to_fs,isnap=0,nsnap,dsnap)
do iband=1,nband
write(band_unit,"(A,I12)") "#band=",iband
do ik=1,nk
ifre = (ik-1)*nband+iband
write(band_unit,"(2(1X,I12),*(1X,F12.5))") iband,ik,(xk(ipol,ik),ipol=1,3),Enk(iband,ik)*RYTOEV,(psit(ifre,isnap),csit(ifre,isnap),isnap=0,nsnap,dsnap)
kx = xk(1,ik)
ky = xk(2,ik)
kz = xk(3,ik)
if(kx > 0.5) kx=kx-1.0
if(ky > 0.5) ky=ky-1.0
if(kz > 0.5) kz=kz-1.0
write(band_unit,"(2(1X,I12),*(1X,F12.5))") iband,ik,kx,ky,kz,Enk(iband,ik)*RYTOEV,(psit(ifre,isnap),csit(ifre,isnap),isnap=0,nsnap,dsnap)
if(kx == 0.5 .or. ky==0.5 .or. kz==0.5) then
if(kx==0.5) kx=kx-1.0
if(ky==0.5) ky=ky-1.0
if(kz==0.5) kz=kz-1.0
write(band_unit,"(2(1X,I12),*(1X,F12.5))") iband,ik,kx,ky,kz,Enk(iband,ik)*RYTOEV,(psit(ifre,isnap),csit(ifre,isnap),isnap=0,nsnap,dsnap)
endif
enddo
write(band_unit,*)

call close_file(band_file_,band_unit)
enddo
call close_file(band_file_,band_unit)
write(stdout,"(A,A)") "Write Carrier occupation on the band structure at different time to file:",trim(band_file_)
write(stdout,"(A,A)") "Write Carrier occupation on the band structure at different time to file:",trim(band_file_)

end subroutine plot_band_occupatin_withtime


subroutine plot_eh_temp(nband,nk,Enk,xk,nsnap,psit,csit,dsnap,temp_file)
implicit none
integer , intent(in) :: nband,nk,dsnap
real(kind=dp),intent(in) :: Enk(nband,nk)
real(kind=dp),intent(in) :: xk(3,nk)
integer , intent(in) :: nsnap
real(kind=dp),intent(in) :: psit(1:nband*nk,0:nsnap),csit(1:nband*nk,0:nsnap)
character(len=*),intent(in) :: temp_file

character(len=maxlen) :: temp_file_
integer :: temp_unit
integer :: ipol,iband,ik,ifre
real(kind=dp) :: kx,ky,kz

do iband=1,nband
write(ctmp1,*) iband
temp_file_=trim(outdir)//trim(adjustl(temp_file))//"_"//trim(adjustl(ctmp1))//".dat"
temp_unit = io_file_unit()
call open_file(temp_file_,temp_unit)

write(temp_unit,"(A34,I12,A)") "Carrier occupation on the iband = ",iband," structure at different time."
write(temp_unit,"(9(A12,1X))") " time "," iband "," ik "," kx "," ky "," kz "," E_nk ","psit "," csit "
write(temp_unit,"(9(A12,1X))") " fs "," iband "," ik "," b1 "," b2 ", " b3 "," eV "," psit "," csit "

do isnap =0,nsnap
do ik=1,nk
ifre = (ik-1)*nband+iband
kx = xk(1,ik)
ky = xk(2,ik)
kz = xk(3,ik)
if(kx > 0.5) kx=kx-1.0
if(ky > 0.5) ky=ky-1.0
if(kz > 0.5) kz=kz-1.0
write(temp_unit,"(F12.5,2(1X,I12),*(1X,F12.5))") &
dt*nstep*isnap*ry_to_fs,iband,ik,kx,ky,kz,Enk(iband,ik)*RYTOEV,psit(ifre,isnap),csit(ifre,isnap)
if(kx == 0.5 .or. ky==0.5 .or. kz==0.5) then
if(kx==0.5) kx=kx-1.0
if(ky==0.5) ky=ky-1.0
if(kz==0.5) kz=kz-1.0
write(temp_unit,"(F12.5,2(1X,I12),*(1X,F12.5))") &
dt*nstep*isnap*ry_to_fs,iband,ik,kx,ky,kz,Enk(iband,ik)*RYTOEV,psit(ifre,isnap),csit(ifre,isnap)
endif
enddo
write(temp_unit,*)
enddo
call close_file(temp_file_,temp_unit)
enddo
write(stdout,"(A,A)") "Write Carrier occupation on the band structure at different time to file:",trim(temp_file_)

end subroutine plot_eh_temp


end module saveinf

0 comments on commit 07c4a97

Please sign in to comment.