TABLE OF CONTENTS
ABINIT/electrooptic [ Functions ]
NAME
electrooptic
FUNCTION
Compute the electrooptic tensor and the raman tensors of zone-center phonons
COPYRIGHT
Copyright (C) 1999-2018 ABINIT group (MVeithen) This file is distributed under the terms of the GNU General Public Licence, see ~abinit/COPYING or http://www.gnu.org/copyleft/gpl.txt . For the initials of contributors, see ~abinit/doc/developers/contributors.txt .
INPUTS
dchide(3,3,3) = non-linear optical coefficients dieflag= dielectric tensor flag. 0=> no dielectric tensor, 1=> frequency-dependent dielectric tensor, 2=> only the electronic dielectric tensor. epsinf=electronic dielectric tensor fact_oscstr(2,3,3*natom)=factors of the oscillator strengths for the different eigenmodes, for different direction of the electric field natom=number of atoms in unit cell phfrq(3*natom)=phonon frequencies (square root of the dynamical matrix eigenvalues, except if these are negative, and in this case, give minus the square root of the absolute value of the matrix eigenvalues). Hartree units. prtmbm= if equal to 1 write out the mode by mode decomposition of the EO tensor rsus = Raman susceptibilities ucvol=unit cell volume
OUTPUT
(to be completed ?)
NOTES
1. The phonon frequencies phfrq should correspond to the wavevector at Gamma, without any non-analyticities. 2. Should clean for no imaginary part ... This routine should be used only by one processor. 3. frdiel(3,3,nfreq)= frequency-dependent dielectric tensor mode effective charges for the different eigenmodes, for different direction of the electric field
PARENTS
anaddb
CHILDREN
matr3inv
SOURCE
54 #if defined HAVE_CONFIG_H 55 #include "config.h" 56 #endif 57 58 #include "abi_common.h" 59 60 61 subroutine electrooptic(dchide,dieflag,epsinf,fact_oscstr,natom,phfrq,prtmbm,rsus,ucvol) 62 63 use defs_basis 64 use m_profiling_abi 65 use m_errors 66 67 !This section has been created automatically by the script Abilint (TD). 68 !Do not modify the following lines by hand. 69 #undef ABI_FUNC 70 #define ABI_FUNC 'electrooptic' 71 use interfaces_32_util 72 !End of the abilint section 73 74 implicit none 75 76 !Arguments ------------------------------- 77 !scalars 78 integer,intent(in) :: dieflag,natom,prtmbm 79 real(dp),intent(in) :: ucvol 80 !arrays 81 real(dp),intent(in) :: dchide(3,3,3),epsinf(3,3),fact_oscstr(2,3,3*natom) 82 real(dp),intent(in) :: phfrq(3*natom),rsus(3*natom,3,3) 83 84 !Local variables ------------------------- 85 !scalars 86 integer :: flag,i1,i2,ii,imode,jj,kk 87 real(dp) :: dtm,fac 88 logical :: iwrite 89 character(len=500) :: message 90 !arrays 91 integer :: voigtindex(6,2) 92 real(dp) :: eta(3,3),rvoigt(6,3),work(3,3,3) 93 real(dp),allocatable :: rijk(:,:,:,:),rijk_tot(:,:,:) 94 95 ! ********************************************************************* 96 97 !rijk(1:3*natom,:,:,:) = mode by mode decomposition of the electrooptic tensor 98 !rijk(3*natom+1,:,:,:) = electronic contribution 99 iwrite = ab_out > 0 100 101 voigtindex(1,1) = 1 ; voigtindex(1,2) = 1 102 voigtindex(2,1) = 2 ; voigtindex(2,2) = 2 103 voigtindex(3,1) = 3 ; voigtindex(3,2) = 3 104 voigtindex(4,1) = 2 ; voigtindex(4,2) = 3 105 voigtindex(5,1) = 1 ; voigtindex(5,2) = 3 106 voigtindex(6,1) = 1 ; voigtindex(6,2) = 2 107 108 ABI_ALLOCATE(rijk,(3*natom+1,3,3,3)) 109 ABI_ALLOCATE(rijk_tot,(3,3,3)) 110 rijk(:,:,:,:) = 0._dp 111 rijk_tot(:,:,:) = 0._dp 112 113 114 !In case there is no mode with truly negative frequency 115 !and the electronic dielectric tensor is available 116 !compute the electro-optic tensor 117 118 flag = 1 119 120 if (abs(phfrq(1)) > abs(phfrq(4))) then 121 flag = 0 122 write(message,'(6a)')& 123 & 'The lowest mode appears to be a "true" negative mode,',ch10,& 124 & 'and not an acoustic mode. This precludes the computation',ch10,& 125 & 'of the EO tensor.',ch10 126 MSG_WARNING(message) 127 end if 128 129 dtm = epsinf(1,1)*epsinf(2,2)*epsinf(3,3) + & 130 & epsinf(1,2)*epsinf(2,3)*epsinf(3,1) + & 131 & epsinf(1,3)*epsinf(2,1)*epsinf(3,2) - & 132 & epsinf(3,1)*epsinf(2,2)*epsinf(1,3) - & 133 & epsinf(3,2)*epsinf(2,3)*epsinf(1,1) - & 134 & epsinf(3,3)*epsinf(2,1)*epsinf(1,2) 135 136 if (abs(dtm) < tol6) then 137 flag = 0 138 write(message,'(a,a,a,a,a,a,a,a)')& 139 & 'The determinant of the electronic dielectric tensor is zero.',ch10,& 140 & 'This preludes the computation fo the EO tensor since',ch10,& 141 & 'this quantity requires the inverse of epsilon.',ch10,& 142 & 'Action : check you database and the value of dieflag in the input file.',ch10 143 MSG_WARNING(message) 144 end if 145 146 !dieflag is required to be one since the EO tensor 147 !requires the oscillator strengths 148 149 if ((flag == 1).and.(dieflag==1)) then 150 151 ! Factor to convert atomic units to MKS units 152 153 fac = -16._dp*pi*pi*eps0*(Bohr_Ang**2)*1.0d-8/(e_Cb*sqrt(ucvol)) 154 155 ! Compute inverse of dielectric tensor 156 ! needed to convert the nonlinear optical susceptibility tensor 157 ! to the electrooptic tensor 158 159 call matr3inv(epsinf,eta) 160 161 if (iwrite) then 162 write(ab_out,*)ch10 163 write(ab_out,*)'Output of the EO tensor (pm/V) in Voigt notations' 164 write(ab_out,*)'=================================================' 165 write(ab_out,*) 166 if (prtmbm == 1) then 167 write(ab_out,*)'Mode by mode decomposition' 168 write(ab_out,*) 169 end if 170 end if 171 172 ! Compute the ionic contribution to the EO tensor 173 174 do imode = 4, 3*natom 175 176 if (prtmbm == 1 .and. iwrite) then 177 write(ab_out,*) 178 write(ab_out,'(a4,i3,2x,a2,f7.2,a6)')'Mode',imode,' (',phfrq(imode)*Ha_cmm1,' cm-1)' 179 end if 180 181 do ii = 1, 3 182 do jj = 1, 3 183 do kk = 1, 3 184 rijk(imode,ii,jj,kk) = rsus(imode,ii,jj)*fact_oscstr(1,kk,imode)/(phfrq(imode)**2) 185 end do 186 end do 187 end do 188 189 work(:,:,:) = 0._dp 190 do ii = 1,3 191 do jj = 1, 3 192 do kk = 1, 3 193 194 do i1 = 1, 3 195 do i2 = 1, 3 196 work(ii,jj,kk) = work(ii,jj,kk) + eta(ii,i1)*rijk(imode,i1,i2,kk)*eta(i2,jj) 197 end do ! i2 198 end do ! i1 199 200 rijk(imode,ii,jj,kk) = fac*work(ii,jj,kk) 201 rijk_tot(ii,jj,kk) = rijk_tot(ii,jj,kk) + rijk(imode,ii,jj,kk) 202 end do 203 204 end do 205 end do 206 207 if (prtmbm == 1) then 208 rvoigt(:,:) = 0._dp 209 do i1 = 1, 6 210 ii = voigtindex(i1,1) 211 jj = voigtindex(i1,2) 212 do kk = 1, 3 213 rvoigt(i1,kk) = (rijk(imode,ii,jj,kk) + rijk(imode,jj,ii,kk))/2._dp 214 end do 215 if (iwrite) write(ab_out,'(5x,3(2x,f16.9))')rvoigt(i1,:) 216 end do 217 end if 218 219 end do ! imode 220 221 ! Compute the electronic contribution to the EO tensor 222 223 if (prtmbm == 1 .and. iwrite) then 224 write(ab_out,*) 225 write(ab_out,*)'Electronic contribution to the EO tensor' 226 end if 227 228 fac = 16*(pi**2)*(Bohr_Ang**2)*1.0d-8*eps0/e_Cb 229 230 do ii = 1,3 231 do jj = 1, 3 232 do kk = 1, 3 233 234 do i1 = 1, 3 235 do i2 = 1, 3 236 rijk(3*natom+1,ii,jj,kk) = rijk(3*natom+1,ii,jj,kk) + & 237 & eta(ii,i1)*dchide(i1,i2,kk)*eta(i2,jj) 238 end do ! i2 239 end do ! i1 240 241 rijk(3*natom+1,ii,jj,kk) = -4._dp*rijk(3*natom+1,ii,jj,kk)*fac 242 rijk_tot(ii,jj,kk) = rijk_tot(ii,jj,kk) + rijk(3*natom+1,ii,jj,kk) 243 244 end do 245 246 end do 247 end do 248 249 if (prtmbm == 1) then 250 rvoigt(:,:) = 0._dp 251 do i1 = 1, 6 252 ii = voigtindex(i1,1) 253 jj = voigtindex(i1,2) 254 do kk = 1, 3 255 rvoigt(i1,kk) = (rijk(3*natom+1,ii,jj,kk) + rijk(3*natom+1,jj,ii,kk))/2._dp 256 end do 257 if (iwrite) write(ab_out,'(5x,3(2x,f16.9))')rvoigt(i1,:) 258 end do 259 if (iwrite) write(ab_out,*)ch10 260 end if 261 262 if (iwrite) write(ab_out,*)'Total EO tensor (pm/V) in Voigt notations' 263 rvoigt(:,:) = 0._dp 264 do i1 = 1, 6 265 ii = voigtindex(i1,1) 266 jj = voigtindex(i1,2) 267 do kk = 1, 3 268 rvoigt(i1,kk) = (rijk_tot(ii,jj,kk) + rijk_tot(jj,ii,kk))/2._dp 269 end do 270 if (iwrite) write(ab_out,'(5x,3(2x,f16.9))')rvoigt(i1,:) 271 end do 272 273 end if ! flag 274 275 ABI_DEALLOCATE(rijk) 276 ABI_DEALLOCATE(rijk_tot) 277 278 end subroutine electrooptic