TABLE OF CONTENTS
ABINIT/mksubham [ Functions ]
NAME
mksubham
FUNCTION
Build the Hamiltonian matrix in the eigenfunctions subspace, for one given band (or for one given block of bands)
COPYRIGHT
Copyright (C) 1998-2017 ABINIT group (DCA, XG, GMR, MT) 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
cg(2,mcg)=wavefunctions gsc(2,mgsc)=<g|S|c> matrix elements (S=overlap) iblock=index of block of bands icg=shift to be applied on the location of data in the array cg igsc=shift to be applied on the location of data in the array cg istwf_k=input parameter that describes the storage of wfs mcg=second dimension of the cg array mgsc=second dimension of the gsc array nband_k=number of bands at this k point for that spin polarization nbdblock=number of bands in a block npw_k=number of plane waves at this k point nspinor=number of spinorial components of the wavefunctions use_subovl=1 if the overlap matrix is not identity in WFs subspace use_vnl= 1 if <C band,k|H|C band_prime,k> has to be computed me_g0=1 if this processors has G=0, 0 otherwise
OUTPUT
SIDE EFFECTS
ghc(2,npw_k*nspinor)=<G|H|C band,k> for the current state This is an input in non-blocked algorithm an output in blocked algorithm gvnlc(2,npw_k*nspinor)=<G|Vnl|C band,k> for the current state This is an input in non-blocked algorithm an output in blocked algorithm isubh=index of current state in array subham isubo=index of current state in array subovl subham(nband_k*(nband_k+1))=Hamiltonian expressed in the WFs subspace subovl(nband_k*(nband_k+1)*use_subovl)=overlap matrix expressed in the WFs subspace subvnl(nband_k*(nband_k+1)*use_vnl)=non-local Hamiltonian expressed in the WFs subspace
PARENTS
cgwf
CHILDREN
SOURCE
56 #if defined HAVE_CONFIG_H 57 #include "config.h" 58 #endif 59 60 #include "abi_common.h" 61 62 63 subroutine mksubham(cg,ghc,gsc,gvnlc,iblock,icg,igsc,istwf_k,& 64 & isubh,isubo,mcg,mgsc,nband_k,nbdblock,npw_k,& 65 & nspinor,subham,subovl,subvnl,use_subovl,use_vnl,me_g0) 66 67 use defs_basis 68 use m_profiling_abi 69 use m_cgtools 70 71 !This section has been created automatically by the script Abilint (TD). 72 !Do not modify the following lines by hand. 73 #undef ABI_FUNC 74 #define ABI_FUNC 'mksubham' 75 !End of the abilint section 76 77 implicit none 78 79 !Arguments ------------------------------------ 80 !scalars 81 integer,intent(in) :: iblock,icg,igsc,istwf_k,mcg,mgsc,nband_k 82 integer,intent(in) :: nbdblock,npw_k,nspinor,use_subovl,use_vnl,me_g0 83 integer,intent(inout) :: isubh,isubo 84 !arrays 85 real(dp),intent(in) :: cg(2,mcg) 86 real(dp),intent(in) :: gsc(2,mgsc) 87 real(dp),intent(inout) :: ghc(2,npw_k*nspinor),gvnlc(2,npw_k*nspinor) 88 real(dp),intent(inout) :: subham(nband_k*(nband_k+1)) 89 real(dp),intent(inout) :: subovl(nband_k*(nband_k+1)*use_subovl) 90 real(dp),intent(inout) :: subvnl(nband_k*(nband_k+1)*use_vnl) 91 92 !Local variables------------------------------- 93 !scalars 94 integer :: iband,ibdblock,ii,ipw,ipw1,isp,iwavef,jwavef 95 real(dp) :: cgimipw,cgreipw,chcim,chcre,cscim,cscre,cvcim,cvcre 96 !real(dp) :: chc(2),cvc(2),csc(2) 97 98 ! ********************************************************************* 99 100 !Loop over bands in a block This loop can be parallelized 101 do iband=1+(iblock-1)*nbdblock,min(iblock*nbdblock,nband_k) 102 ibdblock=iband-(iblock-1)*nbdblock 103 104 ! Compute elements of subspace Hamiltonian <C(i)|H|C(n)> and <C(i)|Vnl|C(n)> 105 if(istwf_k==1)then 106 107 do ii=1,iband 108 iwavef=(ii-1)*npw_k*nspinor+icg 109 chcre=zero ; chcim=zero 110 if (use_vnl==0) then 111 do ipw=1,npw_k*nspinor 112 cgreipw=cg(1,ipw+iwavef) 113 cgimipw=cg(2,ipw+iwavef) 114 chcre=chcre+cgreipw*ghc(1,ipw)+cgimipw*ghc(2,ipw) 115 chcim=chcim+cgreipw*ghc(2,ipw)-cgimipw*ghc(1,ipw) 116 end do 117 ! chc = cg_zdotc(npw_k*nspinor,cg(1,1+iwavef),ghc) 118 else 119 #if 1 120 do ipw=1,npw_k*nspinor 121 cgreipw=cg(1,ipw+iwavef) 122 cgimipw=cg(2,ipw+iwavef) 123 chcre=chcre+cgreipw*ghc(1,ipw)+cgimipw*ghc(2,ipw) 124 chcim=chcim+cgreipw*ghc(2,ipw)-cgimipw*ghc(1,ipw) 125 end do 126 cvcre=zero ; cvcim=zero 127 do ipw=1,npw_k*nspinor 128 cgreipw=cg(1,ipw+iwavef) 129 cgimipw=cg(2,ipw+iwavef) 130 cvcre=cvcre+cgreipw*gvnlc(1,ipw)+cgimipw*gvnlc(2,ipw) 131 cvcim=cvcim+cgreipw*gvnlc(2,ipw)-cgimipw*gvnlc(1,ipw) 132 end do 133 subvnl(isubh )=cvcre 134 subvnl(isubh+1)=cvcim 135 #else 136 ! New version with BLAS1, will require some update of the refs. 137 cvc = cg_zdotc(npw_k*nspinor,cg(1,1+iwavef),gvnlc) 138 subvnl(isubh )=cvc(1) 139 subvnl(isubh+1)=cvc(2) 140 chc = cg_zdotc(npw_k*nspinor,cg(1,1+iwavef),ghc) 141 chcre = chc(1) 142 chcim = chc(2) 143 #endif 144 ! Store real and imag parts in Hermitian storage mode: 145 end if 146 subham(isubh )=chcre 147 subham(isubh+1)=chcim 148 ! subham(isubh )=chc(1) 149 ! subham(isubh+1)=chc(2) 150 isubh=isubh+2 151 end do 152 153 else if(istwf_k>=2)then 154 do ii=1,iband 155 iwavef=(ii-1)*npw_k+icg 156 ! Use the time-reversal symmetry, but should not double-count G=0 157 if(istwf_k==2 .and. me_g0==1) then 158 chcre = half*cg(1,1+iwavef)*ghc(1,1) 159 if (use_vnl==1) cvcre=half*cg(1,1+iwavef)*gvnlc(1,1) 160 ipw1=2 161 else 162 chcre=zero; ipw1=1 163 if (use_vnl==1) cvcre=zero 164 end if 165 if (use_vnl==0) then 166 do isp=1,nspinor 167 do ipw=ipw1+(isp-1)*npw_k,npw_k*isp 168 cgreipw=cg(1,ipw+iwavef) 169 cgimipw=cg(2,ipw+iwavef) 170 chcre=chcre+cgreipw*ghc(1,ipw)+cgimipw*ghc(2,ipw) 171 end do 172 end do 173 chcre=two*chcre 174 else 175 do isp=1,nspinor 176 do ipw=ipw1+(isp-1)*npw_k,npw_k*isp 177 cgreipw=cg(1,ipw+iwavef) 178 cgimipw=cg(2,ipw+iwavef) 179 chcre=chcre+cgreipw*ghc(1,ipw)+cgimipw*ghc(2,ipw) 180 cvcre=cvcre+cgreipw*gvnlc(1,ipw)+cgimipw*gvnlc(2,ipw) 181 end do 182 end do 183 chcre=two*chcre 184 cvcre=two*cvcre 185 ! Store real and imag parts in Hermitian storage mode: 186 subvnl(isubh )=cvcre 187 subvnl(isubh+1)=zero 188 end if 189 subham(isubh )=chcre 190 subham(isubh+1)=zero 191 isubh=isubh+2 192 end do 193 end if 194 195 ! Compute elements of subspace <C(i)|S|C(n)> (S=overlap matrix) 196 ! <C(i)|S|C(n)> should be closed to Identity. 197 if (use_subovl==1) then 198 jwavef=(iband-1)*npw_k*nspinor+igsc 199 if(istwf_k==1)then 200 do ii=1,iband 201 iwavef=(ii-1)*npw_k*nspinor+icg 202 cscre=zero ; cscim=zero 203 do ipw=1,npw_k*nspinor 204 cgreipw=cg(1,ipw+iwavef) 205 cgimipw=cg(2,ipw+iwavef) 206 cscre=cscre+cgreipw*gsc(1,ipw+jwavef)+cgimipw*gsc(2,ipw+jwavef) 207 cscim=cscim+cgreipw*gsc(2,ipw+jwavef)-cgimipw*gsc(1,ipw+jwavef) 208 end do 209 ! csc = cg_zdotc(npw_k*nspinor,cg(1,1+iwavef),gsc) 210 ! subovl(isubo )=csc(1) 211 ! subovl(isubo+1)=csc(2) 212 ! Store real and imag parts in Hermitian storage mode: 213 subovl(isubo )=cscre 214 subovl(isubo+1)=cscim 215 isubo=isubo+2 216 end do 217 else if(istwf_k>=2)then 218 do ii=1,iband 219 iwavef=(ii-1)*npw_k*nspinor+icg 220 if(istwf_k==2 .and. me_g0==1)then 221 cscre=half*cg(1,1+iwavef)*gsc(1,1+jwavef) 222 ipw1=2 223 else 224 cscre=zero; ipw1=1 225 end if 226 do isp=1,nspinor 227 do ipw=ipw1+(isp-1)*npw_k,npw_k*isp 228 cgreipw=cg(1,ipw+iwavef) 229 cgimipw=cg(2,ipw+iwavef) 230 cscre=cscre+cg(1,ipw+iwavef)*gsc(1,ipw+jwavef)+cg(2,ipw+iwavef)*gsc(2,ipw+jwavef) 231 end do 232 end do 233 cscre=two*cscre 234 ! Store real and imag parts in Hermitian storage mode: 235 subovl(isubo )=cscre 236 subovl(isubo+1)=zero 237 isubo=isubo+2 238 end do 239 end if 240 end if 241 242 end do ! iband in a block 243 244 end subroutine mksubham