TABLE OF CONTENTS
ABINIT/qmatrix [ Functions ]
NAME
qmatrix
FUNCTION
calculation of the inverse of the overlap matrix
COPYRIGHT
Copyright (C) 2004-2018 ABINIT group (XW). 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,mpw*nspinor*mband*mkmem*nsppol) = planewave coefficients of wavefunctions RF wavefunctions at k,q. dtefield = variables related to response Berry-phase calculation ikpt = the index of the current k point mband = maximum number of bands mkmem = maximum number of k-points in core memory mpw = maximum number of plane waves mpw1 = maximum number of plane waves for response wavefunctions nkpt = number of k points npwarr(nkpt) = number of planewaves in basis and boundary at this k point npwar1(nkpt) = number of planewaves in basis and boundary for response wfs nspinor = 1 for scalar wfs, 2 for spinor wfs nsppol = 1 for unpolarized, 2 for spin-polarized pwindall(max(mpw,mpw1)*mkmem,8,3) = array used to compute the overlap matrices
OUTPUT
qmat(2,dtefield%nband_occ,dtefield%nband_occ,nkpt,2,3) = inverse of the overlap matrix
PARENTS
dfpt_scfcv
CHILDREN
dzgedi,dzgefa,overlap_g
SOURCE
44 #if defined HAVE_CONFIG_H 45 #include "config.h" 46 #endif 47 48 #include "abi_common.h" 49 50 subroutine qmatrix(cg,dtefield,qmat,mpw,mpw1,mkmem,mband,npwarr,nkpt,nspinor,nsppol,pwindall) 51 52 use defs_basis 53 use m_profiling_abi 54 use m_efield 55 56 use m_cgtools, only : overlap_g 57 use m_abilasi, only : dzgedi, dzgefa 58 59 !This section has been created automatically by the script Abilint (TD). 60 !Do not modify the following lines by hand. 61 #undef ABI_FUNC 62 #define ABI_FUNC 'qmatrix' 63 !End of the abilint section 64 65 implicit none 66 67 !Arguments ---------------------------------------- 68 !scalars 69 integer,intent(in) :: mband,mkmem,mpw,mpw1,nkpt,nspinor,nsppol 70 type(efield_type),intent(in) :: dtefield 71 !arrays 72 integer,intent(in) :: npwarr(nkpt),pwindall(max(mpw,mpw1)*mkmem,8,3) 73 real(dp),intent(in) :: cg(2,mpw*nspinor*mband*mkmem*nsppol) 74 real(dp),intent(out) :: qmat(2,dtefield%mband_occ,dtefield%mband_occ,nkpt,2,3) 75 76 !Local variables ------------------------- 77 !scalars 78 integer :: iband,icg,icg1,idir,ifor,ikpt,ikpt2,info,jband,job 79 integer :: npw_k1,npw_k2,pwmax,pwmin 80 integer :: isppol 81 real(dp) :: doti,dotr 82 !arrays 83 integer,allocatable :: ipvt(:),pwind_k(:) 84 real(dp) :: det(2,2) 85 real(dp),allocatable :: sinv(:,:,:),smat_k(:,:,:),vect1(:,:),vect2(:,:) 86 real(dp),allocatable :: zgwork(:,:) 87 88 ! ************************************************************************* 89 90 ABI_ALLOCATE(ipvt,(dtefield%mband_occ)) 91 ABI_ALLOCATE(sinv,(2,dtefield%mband_occ,dtefield%mband_occ)) 92 ABI_ALLOCATE(zgwork,(2,dtefield%mband_occ)) 93 ABI_ALLOCATE(vect1,(2,0:mpw)) 94 ABI_ALLOCATE(vect2,(2,0:mpw)) 95 ABI_ALLOCATE(smat_k,(2,dtefield%mband_occ,dtefield%mband_occ)) 96 ABI_ALLOCATE(pwind_k,(max(mpw,mpw1))) 97 vect1(:,0) = zero ; vect2(:,0) = zero 98 99 job = 11 100 101 !************** 102 !loop over k points 103 do isppol = 1, nsppol 104 do ikpt = 1, nkpt 105 npw_k1 = npwarr(ikpt) 106 icg = dtefield%cgindex(ikpt,1) 107 108 do idir = 1, 3 109 110 do ifor = 1, 2 111 112 ikpt2 = dtefield%ikpt_dk(ikpt,ifor,idir) 113 npw_k2 = npwarr(ikpt2) 114 icg1 = dtefield%cgindex(ikpt2,1) 115 pwind_k(1:npw_k1) = pwindall((ikpt-1)*max(mpw,mpw1)+1:(ikpt-1)*max(mpw,mpw1)+npw_k1,ifor,idir) 116 117 do jband = 1, dtefield%nband_occ(isppol) 118 vect2(:,1:npw_k2) = & 119 & cg(:,icg1 + 1 + (jband-1)*npw_k2*nspinor:icg1 + jband*npw_k2*nspinor) 120 if (npw_k2 < mpw) vect2(:,npw_k2+1:mpw) = zero 121 122 do iband = 1, dtefield%nband_occ(isppol) 123 124 pwmin = (iband-1)*npw_k1*nspinor 125 pwmax = pwmin + npw_k1*nspinor 126 vect1(:,1:npw_k1) = & 127 & cg(:,icg + 1 + pwmin:icg + pwmax) 128 if (npw_k1 < mpw) vect1(:,npw_k1+1:mpw) = zero 129 call overlap_g(doti,dotr,mpw,npw_k1,npw_k2,nspinor,pwind_k,& 130 & vect1,vect2) 131 smat_k(1,iband,jband) = dotr 132 smat_k(2,iband,jband) = doti 133 134 end do ! iband 135 136 end do !jband 137 138 sinv(:,:,:) = smat_k(:,:,:) 139 140 call dzgefa(sinv,dtefield%mband_occ,dtefield%nband_occ(isppol),ipvt,info) 141 call dzgedi(sinv,dtefield%mband_occ,dtefield%nband_occ(isppol),ipvt,det,zgwork,job) 142 143 qmat(:,:,:,ikpt,ifor,idir) = sinv(:,:,:) 144 145 146 end do 147 148 end do 149 150 end do !end loop over k 151 end do 152 153 ABI_DEALLOCATE(ipvt) 154 ABI_DEALLOCATE(sinv) 155 ABI_DEALLOCATE(zgwork) 156 ABI_DEALLOCATE(vect1) 157 ABI_DEALLOCATE(vect2) 158 ABI_DEALLOCATE(smat_k) 159 ABI_DEALLOCATE(pwind_k) 160 161 end subroutine qmatrix