TABLE OF CONTENTS
ABINIT/wrtloctens [ Functions ]
NAME
wrtloctens
FUNCTION
Output of the localisation tensor
COPYRIGHT
Copyright (C) 1999-2018 ABINIT group (mkv) 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
blkflg = flags for each element of the 2DTE (=1 if computed) d2bbb = band by band decomposition of second order derivatives d2nl = non-local contributions to the 2DTEs mband = maximum number of bands mpert = maximum number of ipert natom = number of atoms in unit cell prtbbb = if = 1, write the band by band decomposition of the localization tensor rprimd = dimensional primitive translations for real space (bohr) usepaw = flag for PAW
OUTPUT
(only writing)
TODO
The localization tensor cannot be defined in the metallic case. It should not be computed.
PARENTS
respfn
CHILDREN
wrtout
SOURCE
42 #if defined HAVE_CONFIG_H 43 #include "config.h" 44 #endif 45 46 #include "abi_common.h" 47 48 49 subroutine wrtloctens(blkflg,d2bbb,d2nl,mband,mpert,natom,prtbbb,rprimd,usepaw) 50 51 use defs_basis 52 use m_profiling_abi 53 54 !This section has been created automatically by the script Abilint (TD). 55 !Do not modify the following lines by hand. 56 #undef ABI_FUNC 57 #define ABI_FUNC 'wrtloctens' 58 use interfaces_14_hidewrite 59 !End of the abilint section 60 61 implicit none 62 63 !Arguments ------------------------------------ 64 !scalars 65 integer,intent(in) :: mband,mpert,natom,prtbbb,usepaw 66 !arrays 67 integer,intent(in) :: blkflg(3,mpert,3,mpert) 68 real(dp),intent(in) :: rprimd(3,3) 69 real(dp),intent(inout) :: d2bbb(2,3,3,mpert,mband,mband*prtbbb) 70 real(dp),intent(inout) :: d2nl(2,3,mpert,3,mpert) 71 72 !Local variables ------------------------------ 73 !scalars 74 integer :: flag,iband,idir,idir2,jband 75 character(len=500) :: message 76 !arrays 77 real(dp) :: loctenscart(2,3,3) 78 real(dp),allocatable :: loctenscart_bbb(:,:,:,:,:) 79 80 ! ********************************************************************* 81 82 !This feature is disabled in the PAW case 83 if (usepaw==1) return 84 85 if(prtbbb==1)then 86 ABI_ALLOCATE(loctenscart_bbb,(2,3,3,mband,mband*prtbbb)) 87 loctenscart_bbb(:,:,:,:,:)=zero 88 end if 89 90 !complete missing elements 91 flag = 0 92 do idir2 = 1,3 93 do idir = 1,3 94 95 if (blkflg(idir2,natom+1,idir,natom+1)==0) then 96 if (blkflg(idir,natom+1,idir2,natom+1)==0) then 97 flag = 1 98 else 99 d2nl(1,idir2,natom+1,idir,natom+1) = d2nl(1,idir,natom+1,idir2,natom+1) 100 d2nl(2,idir2,natom+1,idir,natom+1) =-d2nl(2,idir,natom+1,idir2,natom+1) 101 if(prtbbb==1)then 102 d2bbb(1,idir2,idir,natom+1,:,:) = d2bbb(1,idir,idir2,natom+1,:,:) 103 d2bbb(2,idir2,idir,natom+1,:,:) =-d2bbb(2,idir,idir2,natom+1,:,:) 104 end if 105 106 end if 107 end if 108 109 end do ! idir=1,3 110 end do ! idir2=1,3 111 112 !Transform the tensor to cartesian coordinates 113 114 loctenscart(1,:,:) = matmul(rprimd,d2nl(1,:,natom+1,:,natom+1)) 115 loctenscart(2,:,:) = matmul(rprimd,d2nl(2,:,natom+1,:,natom+1)) 116 117 loctenscart(1,:,:) = matmul(loctenscart(1,:,:),transpose(rprimd)) 118 loctenscart(2,:,:) = matmul(loctenscart(2,:,:),transpose(rprimd)) 119 120 loctenscart(:,:,:) = loctenscart(:,:,:)/(two_pi**2) 121 122 if (prtbbb == 1) then 123 124 write(message,'(a,a)')ch10, & 125 & ' Band by band decomposition of the localisation tensor (bohr^2)' 126 call wrtout(std_out,message,'COLL') 127 call wrtout(ab_out,message,'COLL') 128 129 do iband = 1,mband 130 do jband = 1,mband 131 132 loctenscart_bbb(1,:,:,iband,jband) = matmul(rprimd,d2bbb(1,:,:,natom+1,iband,jband)) 133 loctenscart_bbb(2,:,:,iband,jband) = matmul(rprimd,d2bbb(2,:,:,natom+1,iband,jband)) 134 135 loctenscart_bbb(1,:,:,iband,jband) = matmul(loctenscart_bbb(1,:,:,iband,jband),transpose(rprimd)) 136 loctenscart_bbb(2,:,:,iband,jband) = matmul(loctenscart_bbb(2,:,:,iband,jband),transpose(rprimd)) 137 138 loctenscart_bbb(:,:,:,iband,jband) = loctenscart_bbb(:,:,:,iband,jband)/(two_pi**2) 139 140 141 write(message,'(a,a,i5,a,i5,a)')ch10, & 142 & ' Localisation tensor (bohr^2) for band ',iband,',',jband, & 143 & ' in cartesian coordinates' 144 call wrtout(std_out,message,'COLL') 145 call wrtout(ab_out,message,'COLL') 146 147 write(ab_out,*)' direction matrix element' 148 write(ab_out,*)' alpha beta real part imaginary part' 149 do idir2 = 1,3 150 do idir = 1,3 151 write(ab_out,'(5x,i1,8x,i1,3x,2f16.10)')idir2,idir,& 152 & loctenscart_bbb(1,idir2,idir,iband,jband),& 153 & loctenscart_bbb(2,idir2,idir,iband,jband) 154 end do 155 end do 156 157 end do !jband 158 end do !iband 159 160 end if !prtbbb 161 162 if (usepaw==0) then 163 write(message,'(a,a,a,a)')ch10, & 164 & ' Total localisation tensor (bohr^2) in cartesian coordinates',ch10,& 165 & ' WARNING : still subject to testing - especially symmetries.' 166 else 167 write(message,'(a,a,a,a,a,a)')ch10, & 168 & ' Total localisation tensor (bohr^2) in cartesian coordinates',ch10,& 169 & ' WARNING : probably wrong for PAW (printing for testing purpose)',ch10,& 170 & ' WARNING : still subject to testing - especially symmetries.' 171 end if 172 call wrtout(std_out,message,'COLL') 173 call wrtout(ab_out,message,'COLL') 174 175 write(ab_out,*)' direction matrix element' 176 write(ab_out,*)' alpha beta real part imaginary part' 177 do idir2 = 1,3 178 do idir = 1,3 179 write(ab_out,'(5x,i1,8x,i1,3x,2f16.10)')idir2,idir,& 180 & loctenscart(1,idir2,idir),& 181 & loctenscart(2,idir2,idir) 182 end do 183 end do 184 185 if (flag == 1) then 186 write(message,'(6a)')ch10,& 187 & ' WARNING : Localization tensor calculation (this does not apply to other properties).',ch10,& 188 & ' Not all d/dk perturbations were computed. So the localization tensor in reciprocal space is incomplete,',ch10,& 189 & ' and transformation to cartesian coordinates may be wrong.' 190 call wrtout(std_out,message,'COLL') 191 call wrtout(ab_out,message,'COLL') 192 end if 193 194 if (prtbbb == 1) then 195 ABI_DEALLOCATE(loctenscart_bbb) 196 end if 197 198 end subroutine wrtloctens