TABLE OF CONTENTS
ABINIT/getcgqphase [ Functions ]
NAME
getcgqphase
FUNCTION
extract phases from wave functions, to cancel contributions to gkk matrix elements
COPYRIGHT
Copyright (C) 2011-2018 ABINIT group (MJV) This file is distributed under the terms of the GNU General Public License, see ~abinit/COPYING or http://www.gnu.org/copyleft/gpl.txt . For the initials of contributors, see ~abinit/doc/developers/contributors.txt.
INPUTS
dtset <type(dataset_type)>=all input variables for this dataset timrev = flag for use of time reversal symmetry cg = input wavefunctions mcg = dimension of cg = nspinor*mband*mpw*mkmem cgq = input wavefunctions at k+q mcgq = dimension of cgq = nspinor*mband*mpw*mkmem mpi_enreg = datastructure for mpi communication nkpt_rbz = number of k-points in reduced zone for present q point npwarr = array of numbers of plane waves for each k-point npwar1 = array of numbers of plane waves for each k+q point
OUTPUT
phasecg = phase of different wavefunction products <k,n | k+q,n'>
PARENTS
dfpt_looppert
CHILDREN
smatrix,wrtout,xmpi_barrier,xmpi_sum_master
SOURCE
39 #if defined HAVE_CONFIG_H 40 #include "config.h" 41 #endif 42 43 #include "abi_common.h" 44 45 46 47 subroutine getcgqphase(dtset, timrev, cg, mcg, cgq, mcgq, mpi_enreg, & 48 & nkpt_rbz, npwarr, npwar1, phasecg) 49 50 use m_profiling_abi 51 use defs_basis 52 use defs_abitypes 53 use m_xmpi 54 55 !This section has been created automatically by the script Abilint (TD). 56 !Do not modify the following lines by hand. 57 #undef ABI_FUNC 58 #define ABI_FUNC 'getcgqphase' 59 use interfaces_14_hidewrite 60 use interfaces_32_util 61 !End of the abilint section 62 63 implicit none 64 65 !Arguments ------------------------------- 66 ! scalars 67 integer, intent(in) :: mcg, mcgq, timrev 68 integer, intent(in) :: nkpt_rbz 69 type(dataset_type), intent(in) :: dtset 70 71 ! arrays 72 integer, intent(in) :: npwarr(nkpt_rbz) 73 integer, intent(in) :: npwar1(nkpt_rbz) 74 real(dp), intent(in) :: cg(2,mcg) 75 real(dp), intent(in) :: cgq(2,mcgq) 76 type(MPI_type), intent(in) :: mpi_enreg 77 real(dp),intent(out) :: phasecg(2, dtset%mband*dtset%mband*nkpt_rbz& 78 & *dtset%nsppol) 79 80 !Local variables ------------------------- 81 ! local vars 82 integer :: icg, icgq, isppol, ikpt, ipw 83 integer :: istate, iband1, iband2 84 integer :: npw_k, npw_q 85 integer :: me, ierr, master, spaceComm, nprocs 86 integer :: usepaw 87 integer :: ddkflag, itrs, job, maxbd, mcg1_k, minbd, shiftbd 88 real(dp) :: normsmat 89 90 integer, allocatable :: sflag_k(:) 91 integer, allocatable :: pwind_k(:) 92 93 real(dp) :: cg1_dummy(1,1) 94 real(dp) :: smat_inv_dummy(1,1,1) 95 real(dp) :: smat_k_paw_dummy(1,1,1) 96 real(dp) :: dtm_k_dummy(2) 97 real(dp), allocatable :: smat_k(:,:,:) 98 real(dp), allocatable :: pwnsfac_k(:,:) 99 logical, allocatable :: my_kpt(:,:) 100 character(len=500) :: message 101 102 ! ********************************************************************* 103 104 ABI_ALLOCATE(smat_k,(2,dtset%mband,dtset%mband)) 105 ABI_ALLOCATE(sflag_k,(dtset%mband)) 106 107 !dummy use of timrev so abirules stops complaining. 108 icg = timrev 109 110 !!MPI data for future use 111 spaceComm=mpi_enreg%comm_cell 112 nprocs=xmpi_comm_size(spaceComm) 113 master=0 114 me=mpi_enreg%me_kpt 115 116 ABI_ALLOCATE(my_kpt, (nkpt_rbz, dtset%nsppol)) 117 my_kpt = .true. 118 if (mpi_enreg%nproc_kpt > 1) then 119 do isppol = 1, dtset%nsppol 120 do ikpt = 1, nkpt_rbz 121 my_kpt(ikpt, isppol) = .not.(proc_distrb_cycle(mpi_enreg%proc_distrb,ikpt,1,& 122 & dtset%nband(ikpt),isppol,me)) 123 end do 124 end do 125 end if 126 127 128 !make trivial association of G vectors: we just want <psi_k| psi_k+q> 129 !TODO: check this is correct wrt arrangement of kg vectors for k+q 130 !looks ok : usually made in initberry, from the translations associated 131 !to the symops, scalar product with the G vectors. The symop is the one 132 !used to go from the irreducible k to the full zone k. In present context 133 !we should be using only the reduced zone, and anyhow have the same k-grid 134 !for the gkk matrix elements and for the cg here... 135 ABI_ALLOCATE(pwind_k,(dtset%mpw)) 136 ABI_ALLOCATE(pwnsfac_k,(4,dtset%mpw)) 137 do ipw = 1, dtset%mpw 138 pwind_k(ipw) = ipw 139 pwnsfac_k(1,ipw) = one 140 pwnsfac_k(2,ipw) = zero 141 pwnsfac_k(3,ipw) = one 142 pwnsfac_k(4,ipw) = zero 143 end do 144 145 !flags for call to smatrix 146 usepaw = 0 ! for now 147 ddkflag = 0 148 itrs = 0 149 job = 0 150 maxbd = 1 151 mcg1_k = 1 152 minbd = 1 153 shiftbd = 1 154 155 !from overlap matrix for each wavefunction, extract phase 156 icg = 0 157 icgq = 0 158 istate = 0 159 160 phasecg = zero 161 do isppol = 1, dtset%nsppol 162 do ikpt = 1, nkpt_rbz 163 ! each proc only has certain k 164 if (.not. my_kpt(ikpt, isppol)) then 165 istate = istate + dtset%nband(ikpt)*dtset%nband(ikpt) 166 cycle 167 end if 168 169 npw_k = npwarr(ikpt) 170 npw_q= npwar1(ikpt) 171 172 ! TODO: question: are the k-points in the ibz correctly ordered in cg and cgq? if not the icg below have to be adapted. 173 sflag_k = 0 ! make sure all elements are calculated 174 smat_k = zero 175 176 call smatrix(cg, cgq, cg1_dummy, ddkflag, dtm_k_dummy, icg, icgq,& 177 & itrs, job, maxbd, mcg, mcgq, mcg1_k, minbd,dtset%mpw, dtset%mband, dtset%mband,& 178 & npw_k, npw_q, dtset%nspinor, pwind_k, pwnsfac_k, sflag_k, shiftbd,& 179 & smat_inv_dummy, smat_k, smat_k_paw_dummy, usepaw) 180 181 icg = icg + npw_k*dtset%nspinor*dtset%nband(ikpt) 182 icgq = icgq + npw_q*dtset%nspinor*dtset%nband(ikpt) 183 184 do iband1 = 1, dtset%nband(ikpt) 185 do iband2 = 1, dtset%nband(ikpt) 186 istate = istate + 1 187 ! normalise the overlap matrix element to get just the phase difference phi_k - phi_k+q 188 normsmat = sqrt(smat_k(1,iband2, iband1)**2 & 189 & + smat_k(2,iband2, iband1)**2) 190 if (normsmat > tol12) then 191 phasecg(:,istate) = smat_k(:,iband2, iband1) / normsmat 192 ! NOTE: 21/9/2011 these appear to be always 1, i, or -i, to within 1.e-5 at worst! 193 end if 194 end do 195 end do 196 end do 197 end do 198 199 !eventually do an mpi allreduce over the k-points for phasecg 200 if (nprocs>1) then 201 call xmpi_barrier(spaceComm) 202 call xmpi_sum_master(phasecg,master,spaceComm,ierr) 203 call xmpi_barrier(spaceComm) 204 if (1==1) then 205 write(message,'(a)') ' In getcgqphase - contributions to phasecg collected' 206 call wrtout(std_out,message,'PERS') 207 end if 208 end if 209 210 ABI_DEALLOCATE(sflag_k) 211 ABI_DEALLOCATE(smat_k) 212 ABI_DEALLOCATE(pwind_k) 213 ABI_DEALLOCATE(pwnsfac_k) 214 ABI_DEALLOCATE(my_kpt) 215 216 end subroutine getcgqphase