Commit 0e024d79 authored by dalcorso's avatar dalcorso

Bands.x can now write the bands in files with a format readable by gnuplot

that can make a contour plot of the bands.


git-svn-id: http://qeforge.qe-forge.org/svn/q-e/trunk/espresso@9100 c92efa57-630b-4861-b058-cf58834340f0
parent eec16a1f
...@@ -81,6 +81,20 @@ input_description -distribution {Quantum Espresso} -package PWscf -program bands ...@@ -81,6 +81,20 @@ input_description -distribution {Quantum Espresso} -package PWscf -program bands
} }
} }
var plot_2d -type LOGICAL {
default { .false. }
info {
If .true. writes the eigenvalues in the output file
in a 2D format readable by gnuplot. Band ordering is not
changed. Each band is written in a different file called
filband.# with the format:
xk, yk, energy
xk, yk, energy
.. .. ..
enegies are written in eV and xk in units 2\pi/a.
}
vargroup -type INTEGER { vargroup -type INTEGER {
var firstk var firstk
......
...@@ -39,6 +39,7 @@ sym_band.o \ ...@@ -39,6 +39,7 @@ sym_band.o \
work_function.o \ work_function.o \
write_p_avg.o \ write_p_avg.o \
xsf.o \ xsf.o \
xk_et_collect.o \
wannier_hamiltonians.o wannier_hamiltonians.o
PWOBJS = ../../PW/src/libpw.a PWOBJS = ../../PW/src/libpw.a
......
! !
! Copyright (C) 2001-2009 Quantum ESPRESSO group ! Copyright (C) 2001-2012 Quantum ESPRESSO group
! This file is distributed under the terms of the ! This file is distributed under the terms of the
! GNU General Public License. See the file `License' ! GNU General Public License. See the file `License'
! in the root directory of the present distribution, ! in the root directory of the present distribution,
...@@ -29,12 +29,12 @@ PROGRAM do_bands ...@@ -29,12 +29,12 @@ PROGRAM do_bands
CHARACTER(LEN=256), EXTERNAL :: trimcheck CHARACTER(LEN=256), EXTERNAL :: trimcheck
! !
CHARACTER (len=256) :: filband, filp, outdir CHARACTER (len=256) :: filband, filp, outdir
LOGICAL :: lsigma(4), lsym, lp, no_overlap LOGICAL :: lsigma(4), lsym, lp, no_overlap, plot_2d
INTEGER :: spin_component, firstk, lastk INTEGER :: spin_component, firstk, lastk
INTEGER :: ios INTEGER :: ios
! !
NAMELIST / bands / outdir, prefix, filband, filp, spin_component, lsigma,& NAMELIST / bands / outdir, prefix, filband, filp, spin_component, lsigma,&
lsym, lp, filp, firstk, lastk, no_overlap lsym, lp, filp, firstk, lastk, no_overlap, plot_2d
! !
! initialise environment ! initialise environment
! !
...@@ -56,6 +56,7 @@ PROGRAM do_bands ...@@ -56,6 +56,7 @@ PROGRAM do_bands
firstk=0 firstk=0
lastk=10000000 lastk=10000000
spin_component = 1 spin_component = 1
plot_2d=.false.
no_overlap=.false. no_overlap=.false.
! !
ios = 0 ios = 0
...@@ -90,6 +91,13 @@ PROGRAM do_bands ...@@ -90,6 +91,13 @@ PROGRAM do_bands
CALL mp_bcast( lsym, ionode_id ) CALL mp_bcast( lsym, ionode_id )
CALL mp_bcast( lsigma, ionode_id ) CALL mp_bcast( lsigma, ionode_id )
CALL mp_bcast( no_overlap, ionode_id ) CALL mp_bcast( no_overlap, ionode_id )
CALL mp_bcast( plot_2d, ionode_id )
IF (plot_2d) THEN
lsym=.false.
lp=.false.
no_overlap=.true.
ENDIF
IF ( npool > 1 .and..not.(lsym.or.no_overlap)) CALL errore('bands', & IF ( npool > 1 .and..not.(lsym.or.no_overlap)) CALL errore('bands', &
'pools not implemented',npool) 'pools not implemented',npool)
...@@ -109,9 +117,13 @@ PROGRAM do_bands ...@@ -109,9 +117,13 @@ PROGRAM do_bands
CALL openfil_pp() CALL openfil_pp()
! !
IF (lsym) no_overlap=.true. IF (lsym) no_overlap=.true.
CALL punch_band(filband,spin_component,lsigma,no_overlap) IF (plot_2d) THEN
IF (lsym) CALL sym_band(filband,spin_component,firstk,lastk) CALL punch_band_2d(filband,spin_component)
IF (lp) CALL write_p_avg(filp,spin_component,firstk,lastk) ELSE
CALL punch_band(filband,spin_component,lsigma,no_overlap)
IF (lsym) CALL sym_band(filband,spin_component,firstk,lastk)
IF (lp) CALL write_p_avg(filp,spin_component,firstk,lastk)
END IF
! !
CALL stop_pp CALL stop_pp
STOP STOP
...@@ -489,3 +501,108 @@ SUBROUTINE punch_band (filband, spin_component, lsigma, no_overlap) ...@@ -489,3 +501,108 @@ SUBROUTINE punch_band (filband, spin_component, lsigma, no_overlap)
RETURN RETURN
! !
END SUBROUTINE punch_band END SUBROUTINE punch_band
SUBROUTINE punch_band_2d(filband,spin_component)
!
! This routine opens a file for each band and writes on output
! kx, ky, energy,
! kx, ky, energy
! .., .., ..
! where kx and ky are proportional to the length
! of the vectors k_1 and k_2 specified in the input of the 2d plot.
!
! The k points are supposed to be in the form
! xk(i,j) = xk_0 + dkx *(i-1) + dky * (j-1) 1<i<n1, 1<j<n2
!
! kx(i,j) = (i-1) |dkx|
! ky(i,j) = (j-1) |dky|
!
USE kinds, ONLY : DP
USE constants, ONLY : eps8, rytoev
USE lsda_mod, ONLY : nspin
USE klist, ONLY : xk, nkstot, nks
USE wvfct, ONLY : et, nbnd
USE io_files, ONLY : iuntmp
USE io_global, ONLY : ionode, ionode_id
USE mp, ONLY : mp_bcast
IMPLICIT NONE
CHARACTER(LEN=256),INTENT(IN) :: filband
INTEGER, INTENT(IN) :: spin_component
REAL(DP) :: xk0(3), xk1(3), xk2(3), dkx(3), dky(3), xkdum(3), mdkx, mdky
INTEGER :: n1, n2
INTEGER :: ik, i, i1, i2, ibnd, ijk, start_k, last_k, nks_eff, j, ios
CHARACTER(LEN=256) :: filename
CHARACTER(LEN=6), EXTERNAL :: int_to_char
REAL(DP), ALLOCATABLE :: xk_collect(:,:), et_collect(:,:)
ALLOCATE(xk_collect(3,nkstot))
ALLOCATE(et_collect(nbnd,nkstot))
CALL xk_et_collect( xk_collect, et_collect, xk, et, nkstot, nks, nbnd )
start_k=1
last_k=nkstot
nks_eff=nkstot
IF (nspin==2) THEN
nks_eff=nkstot/2
IF (spin_component==1) THEN
start_k=1
last_k=nks_eff
ELSE
start_k=nks_eff+1
last_k=nkstot
ENDIF
ENDIF
!
! Determine xk0
!
xk0(:)=xk_collect(:,start_k)
!
! Determine dkx
!
dky(:)=xk_collect(:,start_k+1)-xk0(:)
!
! Determine n2 and dky
!
loop_k: DO j=start_k+2, nkstot
xkdum(:)=xk0(:)+(j-1)*dky(:)
IF (ABS(xk_collect(1,j)-xkdum(1))>eps8.OR. &
ABS(xk_collect(2,j)-xkdum(2))>eps8.OR. &
ABS(xk_collect(3,j)-xkdum(3))>eps8) THEN
n2=j-1
dkx(:)=xk_collect(:,j)-xk0(:)
EXIT loop_k
ENDIF
ENDDO loop_k
n1=nks_eff/n2
IF (n1*n2 /= nks_eff) CALL errore('punch_band_2d',&
'Problems with k points',1)
mdkx = sqrt( dkx(1)**2 + dkx(2)**2 + dkx(3)**2 )
mdky = sqrt( dky(1)**2 + dky(2)**2 + dky(3)**2 )
!
! write the output, a band per file
!
DO ibnd=1,nbnd
filename=TRIM(filband) // '.' // TRIM(int_to_char(ibnd))
IF (ionode) &
open(unit=iuntmp,file=filename,status='unknown', err=100, iostat=ios)
CALL mp_bcast(ios,ionode_id)
100 CALL errore('punch_band_2d','Problem opening outputfile',ios)
ijk=0
DO i1=1,n1
DO i2=1,n2
ijk=ijk+1
IF (ionode) &
WRITE(iuntmp,'(3f16.6)') mdkx*(i1-1), mdky*(i2-1), &
et_collect(ibnd,ijk)*rytoev
ENDDO
ENDDO
IF (ionode) CLOSE(unit=iuntmp,status='KEEP')
ENDDO
DEALLOCATE(xk_collect)
DEALLOCATE(et_collect)
RETURN
END
!
! Copyright (C) 2012 Quantum ESPRESSO group
! This file is distributed under the terms of the
! GNU General Public License. See the file `License'
! in the root directory of the present distribution,
! or http://www.gnu.org/copyleft/gpl.txt .
!
!
!----------------------------------------------------------------------------
SUBROUTINE xk_et_collect( xk_collect, et_collect, xk, et, nkstot, nks, nbnd )
!----------------------------------------------------------------------------
!
! ... This routine collects the k points (with granularity kunit) among
! ... nodes and sets the variable xk_collect and wk_collect with the total
! ... number of k-points
!
USE kinds, ONLY : DP
USE mp_global, ONLY : my_pool_id, npool, kunit
USE mp_global, ONLY : inter_pool_comm
USE mp, ONLY : mp_sum
!
IMPLICIT NONE
!
INTEGER :: nkstot, nks, nbnd
! total number of k-points
! number of k-points per pool
REAL (DP) :: xk(3,nks), et(nbnd,nks)
REAL (DP) :: xk_collect(3,nkstot), et_collect(nbnd,nkstot)
! k-points
! k-point weights
!
#if defined (__MPI)
!
INTEGER :: nbase, rest, nks1
!
xk_collect=0.d0
!
et_collect=0.d0
!
nks1 = kunit * ( nkstot / kunit / npool )
!
rest = ( nkstot - nks1 * npool ) / kunit
!
IF ( ( my_pool_id + 1 ) <= rest ) nks1 = nks1 + kunit
!
IF (nks1.ne.nks) &
call errore('xk_wk_collect','problems with nks1',1)
!
! ... calculates nbase = the position in the list of the first point that
! ... belong to this npool - 1
!
nbase = nks * my_pool_id
!
IF ( ( my_pool_id + 1 ) > rest ) nbase = nbase + rest * kunit
!
! copy the original points in the correct position of the list
!
xk_collect(:,nbase+1:nbase+nks) = xk(:,1:nks)
!
et_collect(:,nbase+1:nbase+nks)=et(:,1:nks)
!
CALL mp_sum( xk_collect, inter_pool_comm )
!
CALL mp_sum( et_collect, inter_pool_comm )
!
#else
xk_collect=xk
et_collect=et
#endif
!
RETURN
!
END SUBROUTINE xk_et_collect
!
Markdown is supported
0% or
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment