TABLE OF CONTENTS
ABINIT/WffReadEigK [ Functions ]
NAME
WffReadEigK
FUNCTION
(Wavefunction file, read eigenvalues at one k point) This subroutine reads the block of records related to one k point, and one spin-polarization, that contains the wavefunctions as well as the eigenvalues and occupations . It outputs only the eigenvalues. The wavefunction file should have been prepared outside of this routine, in order to read or write the correct records.
COPYRIGHT
Copyright (C) 2004-2018 ABINIT group (XG) This file is distributed under the terms of the GNU General Public License, see ~abinit/COPYING or http://www.gnu.org/copyleft/gpl.txt .
INPUTS
formeig=format of the eigenvalues 0 => vector of eigenvalues (Ground-State case) 1 => hermitian matrix of eigenvalues (Response Function case) headform=format of the header of the wf file, also governing the k block format in case headform=0, use the default (current) format and headform ikpt=index of current k point (only needed for error message) isppol=spin polarization currently treated (only needed for error message) mband=maximum number of bands (governs the dimension of eigen) mpi_enreg=information about MPI parallelization nband=number of bands actually in eigen can be equal, larger or smaller than nband_disk, but eigen will not be completely filled if nband>nband_disk) must be less or equal to mband tim_rwwf=timing code of the calling routine (set to 0 if not attributed) wff=structured info for the wavefunction file
OUTPUT
eigen((2*mband)**formeig *mband)=array for holding eigenvalues (hartree)
PARENTS
CHILDREN
rwwf
SOURCE
48 #if defined HAVE_CONFIG_H 49 #include "config.h" 50 #endif 51 52 #include "abi_common.h" 53 54 subroutine WffReadEigK(eigen,formeig,headform,ikpt,isppol,mband,mpi_enreg,nband,tim_rwwf,wff) 55 56 use defs_basis 57 use defs_abitypes 58 use m_wffile 59 use m_profiling_abi 60 61 !This section has been created automatically by the script Abilint (TD). 62 !Do not modify the following lines by hand. 63 #undef ABI_FUNC 64 #define ABI_FUNC 'WffReadEigK' 65 use interfaces_56_io_mpi 66 !End of the abilint section 67 68 implicit none 69 70 !Arguments ------------------------------------ 71 !scalars 72 integer,intent(in) :: formeig,headform,ikpt,isppol,mband,nband,tim_rwwf 73 type(MPI_type),intent(in) :: mpi_enreg 74 type(wffile_type),intent(inout) :: wff 75 !arrays 76 real(dp),intent(out) :: eigen((2*mband)**formeig*mband) 77 78 !Local variables------------------------------- 79 !scalars 80 integer,parameter :: nspinor1=1 81 integer :: icg,mcg,nband_disk,npw,option,optkg 82 !arrays 83 integer,allocatable :: kg_dum(:,:) 84 real(dp) :: cg_dum(2,1) 85 real(dp),allocatable :: occ_dum(:) 86 87 ! ************************************************************************* 88 89 !DEBUG 90 !write(std_out,*)' WffReadEigK : enter ' 91 !ENDDEBUG 92 93 option=3 94 mcg=1 ; icg=0 ; optkg=0 95 nband_disk=mband; npw=0 96 ABI_ALLOCATE(kg_dum,(3,optkg)) 97 ABI_ALLOCATE(occ_dum,(mband)) 98 99 !DEBUG 100 !write(std_out,*)' WffReadEigK : formeig,option=',formeig,option 101 !ENDDEBUG 102 103 call rwwf(cg_dum,eigen,formeig,headform,icg,ikpt,isppol,kg_dum,mband,mcg,mpi_enreg,& 104 & nband,nband_disk,& 105 & npw,nspinor1,occ_dum,option,optkg,tim_rwwf,wff) 106 107 ABI_DEALLOCATE(kg_dum) 108 ABI_DEALLOCATE(occ_dum) 109 110 end subroutine WffReadEigK