TABLE OF CONTENTS
ABINIT/prtefield [ Functions ]
NAME
prtefield
FUNCTION
Print components of electric field, displacement field and polarization in nice format
COPYRIGHT
Copyright (C) 1998-2018 ABINIT group (DCA, XG, GMR, LBoeri, 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
dtset <type(dataset_type)>=all input variables in this dataset | berryopt | efield | dfield | red_efield | red_efieldbar | red_dfield dtefield <type(efield_type)> | efield2 | red_ptot1 iunit = unit number to which the data is printed rprimd
OUTPUT
(only writing)
PARENTS
gstate,update_e_field_vars
CHILDREN
metric,wrtout
SOURCE
42 #if defined HAVE_CONFIG_H 43 #include "config.h" 44 #endif 45 46 #include "abi_common.h" 47 48 49 subroutine prtefield(dtset,dtefield,iunit,rprimd) 50 51 use defs_basis 52 use defs_abitypes 53 use m_efield 54 ! use m_profiling_abi 55 56 use m_geometry, only : metric 57 58 !This section has been created automatically by the script Abilint (TD). 59 !Do not modify the following lines by hand. 60 #undef ABI_FUNC 61 #define ABI_FUNC 'prtefield' 62 use interfaces_14_hidewrite 63 !End of the abilint section 64 65 implicit none 66 67 !Arguments ------------------------------------ 68 integer :: iunit 69 real(dp),intent(in) :: rprimd(3,3) 70 type(efield_type),intent(in) :: dtefield 71 type(dataset_type),intent(inout) :: dtset 72 73 !Local variables------------------------------- 74 ! Do not modify the length of this string 75 !scalars 76 integer :: idir,ii 77 character(len=1500) :: message 78 character(len=7) :: flag_field(3) 79 80 real(dp) :: ucvol 81 ! arrays 82 real(dp) :: ptot_cart(3),gmet(3,3),gprimd(3,3),rmet(3,3),red_pbar(3),red_dbar(3),red_dfieldbar(3) 83 real(dp) :: red_efieldbar_lc(3),red_efield_lc(3) 84 85 86 ! ************************************************************************* 87 88 !DEBUG 89 !write(iout, '(a)') ' prtefield : enter ' 90 !ENDDEBUG 91 92 !write here 93 94 call metric(gmet,gprimd,-1,rmet,rprimd,ucvol) 95 96 97 ptot_cart(:)=zero 98 do idir = 1,3 99 ptot_cart(idir)=rprimd(idir,1)*dtefield%red_ptot1(1) + rprimd(idir,2)*dtefield%red_ptot1(2) + & 100 & rprimd(idir,3)*dtefield%red_ptot1(3) 101 end do 102 ptot_cart(:)=ptot_cart(:)/ucvol 103 104 if (dtset%berryopt == 4) then 105 106 ! to calculate e Eq.(25) 107 do idir=1,3 108 dtset%red_efield(idir) =(ucvol/(4*pi))*dot_product(dtset%efield(:),gprimd(:,idir)) 109 end do 110 111 ! to calculate ebar Eq.(25) 112 do idir=1,3 113 dtset%red_efieldbar(idir) =dot_product(dtset%efield(:),rprimd(:,idir)) 114 end do 115 116 117 ! to calculate pbar 118 do idir=1,3 119 red_pbar(idir) = (4*pi/ucvol)*dot_product(dtefield%red_ptot1(:),rmet(:,idir)) 120 end do 121 122 !MGNAG: This message is to long and causes 123 ! Runtime Error: wrtout_cpp.f90, line 893: Buffer overflow on output 124 ! with NAG in test seq_tsv6_125 where we write to std_out! 125 ! I cannot change the RECLEN of std_out! 126 127 write(message,'(a,a,a,3(es16.9,2x),a,a,3(es16.9,2x))') ' (a. u.)', ch10,& 128 & ' E: ', (dtset%efield(ii), ii=1,3), ch10, & 129 & ' P: ', (ptot_cart(ii), ii=1,3) 130 call wrtout(iunit,message,'COLL') 131 132 write(message,'(a,a,3(es16.9,2x),a,a,3(es16.9,2x))') ch10,& 133 & ' ebar: ', (dtset%red_efieldbar(ii),ii=1,3), ch10, & !!HONG need to change 134 & ' pbar: ', (red_pbar(ii),ii=1,3) 135 call wrtout(iunit,message,'COLL') 136 137 write(message,'(a,a,3(es16.9,2x),a,a,3(es16.9,2x))') ch10,& 138 & ' e: ', (dtset%red_efield(ii),ii=1,3), ch10, & 139 & ' p: ', (dtefield%red_ptot1(ii), ii=1,3) 140 call wrtout(iunit,message,'COLL') 141 142 write(message,'(a,a,a,a,3(es16.9,2x),a,a,3(es16.9,2x),a)') ch10,& 143 & ' (S.I.), that is V/m for E, and C/m^2 for P', ch10, & 144 & '- E: ', (dtset%efield(ii)*(Ha_J/(e_Cb*Bohr_Ang*1d-10)), ii=1,3), ch10, & !(Ha_J/(e_Cb*Bohr_Ang*1d-10))= 5.14220652*1d+11 145 & ' P: ', (ptot_cart(ii)*(e_Cb)/(Bohr_Ang*1d-10)**2, ii=1,3),ch10 146 call wrtout(iunit,message,'COLL') 147 148 end if ! berryopt ==4 149 150 151 if (dtset%berryopt == 6) then 152 do idir=1,3 153 dtset%red_efield(idir) =(ucvol/(4*pi))*dot_product(dtset%efield(:),gprimd(:,idir)) 154 end do 155 156 ! to calculate ebar !! Need to be changed 157 do idir=1,3 158 dtset%red_efieldbar(idir) = dot_product(dtset%efield(:),rprimd(:,idir)) 159 end do 160 161 ! to calculate red_pbar 162 do idir=1,3 163 red_pbar(idir) = (4*pi/ucvol)*dot_product(dtefield%red_ptot1(:),rmet(:,idir)) 164 end do 165 166 ! to calculate red_dbar 167 do idir=1,3 168 red_dbar(idir) = dot_product(dtset%dfield(:),rprimd(:,idir)) 169 end do 170 171 ! to calculate d 172 do idir=1,3 173 dtset%red_dfield(idir) =(ucvol/(4*pi))*dot_product(dtset%dfield(:),gprimd(:,idir)) 174 end do 175 176 write(message,'(a,a,a,3(es16.9,2x),a,a,3(es16.9,2x),a,a,3(es16.9,2x),a,a,3(es16.9,2x))') ' (a. u.)', ch10,& 177 & ' e: ', (dtset%red_efield(ii),ii=1,3), ch10, & 178 & ' p: ', (dtefield%red_ptot1(ii), ii=1,3), ch10, & 179 & ' d: ', (dtset%red_dfield(ii),ii = 1, 3), ch10, & 180 & ' e + p: ', (dtset%red_efield(ii)+dtefield%red_ptot1(ii),ii=1,3) 181 call wrtout(iunit,message,'COLL') 182 183 write(message,'(a,a,3(es16.9,2x),a,a,3(es16.9,2x),a,a,3(es16.9,2x),a,a,3(es16.9,2x))') ch10,& 184 & ' ebar: ', (dtset%red_efieldbar(ii),ii=1,3), ch10, & !!HONG need to change 185 & ' pbar: ', (red_pbar(ii),ii=1,3), ch10, & 186 & ' dbar: ', (red_dbar(ii),ii=1,3), ch10, & 187 & ' eba+pba: ', (dtset%red_efieldbar(ii)+red_pbar(ii),ii=1,3) 188 call wrtout(iunit,message,'COLL') 189 190 write(message,'(a,a,3(es16.9,2x),a,a,3(es16.9,2x),a,a,3(es16.9,2x),a,a,3(es16.9,2x))') ch10, & 191 & ' E: ', (dtset%efield(ii), ii=1,3), ch10, & 192 & ' P: ', (ptot_cart(ii), ii=1,3), ch10, & 193 & ' D: ', (dtset%dfield(ii),ii = 1, 3), ch10, & 194 & 'E+4*pi*P: ', (dtset%efield(ii)+4.0d0*pi*ptot_cart(ii),ii=1,3) 195 call wrtout(iunit,message,'COLL') 196 197 write(message,'(a,a,a,a,3(es16.9,2x),a,a,3(es16.9,2x),a,a,3(es16.9,2x),a,a,3(es16.9,2x),a)') ch10,& 198 & ' (S.I.), that is V/m for E, and C/m^2 for P and D', ch10, & 199 & '- E: ', (dtset%efield(ii)*(Ha_J/(e_Cb*Bohr_Ang*1d-10)), ii=1,3), ch10, & !(Ha_J/(e_Cb*Bohr_Ang*1d-10))= 5.14220652*1d+11 200 & ' P: ', (ptot_cart(ii)*(e_Cb)/(Bohr_Ang*1d-10)**2, ii=1,3), ch10, & 201 & ' D: ', ((1.0d0/(4*pi))*dtset%dfield(ii)*(e_Cb)/(Bohr_Ang*1d-10)**2,ii = 1, 3), ch10, & 202 & 'eps0*E+P: ', (dtset%efield(ii)*eps0*(Ha_J/(e_Cb*Bohr_Ang*1d-10))+ptot_cart(ii)*(e_Cb)/(Bohr_Ang*1d-10)**2,ii=1,3),ch10 ! eps0*(Ha_J/(e_Cb*Bohr_Ang*1d-10))=8.854187817620*5.14220652*1d-1 203 call wrtout(iunit,message,'COLL') 204 205 !MGNAG Runtime Error: wrtout_cpp.f90, line 896: Buffer overflow on output 206 207 end if ! berryopt ==6 208 209 210 if (dtset%berryopt == 14) then 211 212 do idir=1,3 ! ebar local 213 red_efieldbar_lc(idir)=dot_product(dtefield%efield2(:),rprimd(:,idir)) 214 end do 215 216 217 do idir=1,3 218 dtset%red_efield(idir) =(ucvol/(4*pi))*dot_product(dtefield%efield2(:),gprimd(:,idir)) 219 end do 220 221 ! to calculate pbar 222 do idir=1,3 223 red_pbar(idir) = (4*pi/ucvol)*dot_product(dtefield%red_ptot1(:),rmet(:,idir)) 224 end do 225 226 write(message,'(a,a,a,3(es16.9,2x),a,a,3(es16.9,2x),a,a,3(es16.9,2x))') ' (a. u.)', ch10,& 227 & ' ebar0: ', (dtset%red_efieldbar(ii),ii=1,3), ch10, & 228 & ' ebar: ', (red_efieldbar_lc(ii),ii=1,3), ch10, & 229 & ' pbar: ', (red_pbar(ii),ii=1,3) 230 call wrtout(iunit,message,'COLL') 231 232 write(message,'(a,a,3(es16.9,2x),a,a,3(es16.9,2x))') ch10, & 233 & ' e: ', (dtset%red_efield(ii),ii=1,3), ch10, & 234 & ' p: ', (dtefield%red_ptot1(ii), ii=1,3) 235 call wrtout(iunit,message,'COLL') 236 237 write(message,'(a,a,3(es16.9,2x),a,a,3(es16.9,2x))') ch10, & 238 & ' E: ', (dtefield%efield2(ii), ii=1,3), ch10, & 239 & ' P: ', (ptot_cart(ii), ii=1,3) 240 call wrtout(iunit,message,'COLL') 241 242 write(message,'(a,a,a,a,3(es16.9,2x),a,a,3(es16.9,2x),a)') ch10, & 243 & ' (S.I.), that is V/m for E, and C/m^2 for P', ch10, & 244 & '- E: ', (dtefield%efield2(ii)*(Ha_J/(e_Cb*Bohr_Ang*1d-10)), ii=1,3), ch10, & !(Ha_J/(e_Cb*Bohr_Ang*1d-10))= 5.14220652*1d+11 245 & ' P: ', (ptot_cart(ii)*(e_Cb)/(Bohr_Ang*1d-10)**2, ii=1,3),ch10 246 call wrtout(iunit,message,'COLL') 247 248 249 end if ! berryopt ==14 250 251 252 if (dtset%berryopt == 16) then 253 254 ! to calculate e Eq.(25) 255 do idir=1,3 256 dtset%red_efield(idir) =(ucvol/(4*pi))*dot_product(dtset%efield(:),gprimd(:,idir)) 257 end do 258 259 ! to calculate ebar 260 do idir=1,3 261 dtset%red_efieldbar(idir) = dot_product(dtset%efield(:),rprimd(:,idir)) 262 end do 263 264 ! to calculate pbar 265 do idir=1,3 266 red_pbar(idir) = (4*pi/ucvol)*dot_product(dtefield%red_ptot1(:),rmet(:,idir)) 267 end do 268 269 ! to calculate dbar 270 do idir=1,3 271 red_dfieldbar(idir) = (4*pi/ucvol)*dot_product(dtset%red_dfield(:),rmet(:,idir)) 272 end do 273 274 ! to calculate D 275 do idir=1,3 276 dtset%dfield(idir) =(4*pi/ucvol)*dot_product(dtset%red_dfield(:),rprimd(:,idir)) 277 end do 278 279 write(message,'(a,a,a,3(es16.9,2x),a,a,3(es16.9,2x),a,a,3(es16.9,2x),a,a,3(es16.9,2x))') ' (a. u.)', ch10,& 280 & ' e: ', (dtset%red_efield(ii),ii=1,3), ch10, & 281 & ' p: ', (dtefield%red_ptot1(ii), ii=1,3), ch10, & 282 & ' d: ', (dtset%red_dfield(ii),ii = 1, 3), ch10, & 283 & ' e + p: ', (dtset%red_efield(ii)+dtefield%red_ptot1(ii),ii=1,3) 284 call wrtout(iunit,message,'COLL') 285 286 write(message,'(a,a,3(es16.9,2x),a,a,3(es16.9,2x),a,a,3(es16.9,2x),a,a,3(es16.9,2x))') ch10, & 287 & ' ebar: ', (dtset%red_efieldbar(ii),ii=1,3), ch10, & 288 & ' pbar: ', (red_pbar(ii),ii=1,3), ch10, & 289 & ' dbar: ', (red_dfieldbar(ii),ii=1,3), ch10, & 290 & ' eba+pba: ', (dtset%red_efieldbar(ii)+red_pbar(ii),ii=1,3) 291 call wrtout(iunit,message,'COLL') 292 293 write(message,'(a,a,3(es16.9,2x),a,a,3(es16.9,2x),a,a,3(es16.9,2x),a,a,3(es16.9,2x))') ch10, & 294 & ' E: ', (dtset%efield(ii), ii=1,3), ch10, & 295 & ' P: ', (ptot_cart(ii), ii=1,3), ch10, & 296 & ' D: ', (dtset%dfield(ii),ii = 1, 3), ch10, & 297 & 'E+4*pi*P: ', (dtset%efield(ii)+4.0d0*pi*ptot_cart(ii),ii=1,3) 298 call wrtout(iunit,message,'COLL') 299 300 write(message,'(a,a,a,a,3(es16.9,2x),a,a,3(es16.9,2x),a,a,3(es16.9,2x),a,a,3(es16.9,2x),a)') ch10, & 301 & ' (S.I.), that is V/m for E, and C/m^2 for P and D', ch10, & 302 & '- E: ', (dtset%efield(ii)*(Ha_J/(e_Cb*Bohr_Ang*1d-10)), ii=1,3), ch10, & !(Ha_J/(e_Cb*Bohr_Ang*1d-10))= 5.14220652*1d+11 303 & ' P: ', (ptot_cart(ii)*(e_Cb)/(Bohr_Ang*1d-10)**2, ii=1,3), ch10, & 304 & ' D: ', ((1.0d0/(4*pi))*dtset%dfield(ii)*(e_Cb)/(Bohr_Ang*1d-10)**2,ii = 1, 3), ch10, & 305 & 'eps0*E+P: ', (dtset%efield(ii)*eps0*(Ha_J/(e_Cb*Bohr_Ang*1d-10))+ptot_cart(ii)*(e_Cb)/(Bohr_Ang*1d-10)**2,ii=1,3),ch10 306 call wrtout(iunit,message,'COLL') 307 308 end if ! berryopt ==16 309 310 if (dtset%berryopt == 17) then 311 312 do idir=1,3 ! ebar local 313 red_efieldbar_lc(idir)=dot_product(dtefield%efield2(:),rprimd(:,idir)) 314 end do 315 316 do idir=1,3 317 dtset%red_efield(idir) =(ucvol/(4*pi))*dot_product(dtefield%efield2(:),gprimd(:,idir)) 318 end do 319 320 ! to calculate pbar 321 do idir=1,3 322 red_pbar(idir) = (4*pi/ucvol)*dot_product(dtefield%red_ptot1(:),rmet(:,idir)) 323 end do 324 325 326 ! do idir=1,3 327 ! if (dtset%rfdir(idir)==1) then ! fixed ebar 328 ! red_efieldbar_lc(idir)=dot_product(dtefield%efield2(:),rprimd(:,idir)) ! local efieldbar 329 ! dtset%red_efield(idir) =(ucvol/(4*pi))*dot_product(dtefield%efield2(:),gprimd(:,idir)) 330 ! red_pbar(idir) = (4*pi/ucvol)*dot_product(dtefield%red_ptot1(:),rmet(:,idir)) 331 ! dtset%dfield(idir)=dtefield%efield2(idir)+4*pi*ptot_cart(idir) 332 ! dtset%red_dfield(idir)=dtset%red_efield+dtefield%red_ptot1(idir) 333 ! dtset%red_dfieldbar(idir)=red_efieldbar_lc(idir)+red_pbar(idir) 334 ! E_lc(idir)=dtefield%efield2(idir) 335 ! e_lc(idir)=red_efieldbar_lc(idir) 336 ! ebar_lc(idir)=dtset%red_efieldbar(idir) 337 ! else if (dtset%rfdir(idir)==2) then ! fixed d 338 ! dtset%red_efield(idir) =(ucvol/(4*pi))*dot_product(dtefield%efield2(:),gprimd(:,idir)) 339 ! dtset%red_efieldbar(idir) = dot_product(dtefield%efield2(:),rprimd(:,idir)) 340 ! red_pbar(idir) = (4*pi/ucvol)*dot_product(dtefield%red_ptot1(:),rmet(:,idir)) 341 ! red_dfieldbar(idir) = (4*pi/ucvol)*dot_product(dtset%red_dfield(:),rmet(:,idir)) 342 ! dtset%dfield(idir) =(4*pi/ucvol)*dot_product(dtset%red_dfield(:),rprimd(:,idir)) 343 ! E_lc(idir)=dtefield%efield2(idir) 344 ! e_lc(idir)=dtset%red_efield(idir) 345 ! ebar_lc(idir)=dtset%red_efieldbar(idir) 346 ! end if 347 ! enddo 348 349 350 351 do idir=1,3 352 red_efield_lc(idir)= (ucvol/(4*pi))*dot_product(dtefield%efield2(:),gprimd(:,idir)) 353 red_efieldbar_lc(idir)=dot_product(dtefield%efield2(:),rprimd(:,idir)) ! local efieldbar 354 red_pbar(idir) = (4*pi/ucvol)*dot_product(dtefield%red_ptot1(:),rmet(:,idir)) 355 red_dfieldbar(idir) = (4*pi/ucvol)*dot_product(dtset%red_dfield(:),rmet(:,idir)) 356 dtset%dfield(idir) =(4*pi/ucvol)*dot_product(dtset%red_dfield(:),rprimd(:,idir)) 357 358 end do 359 360 do idir=1,3 361 if(dtset%jfielddir(idir)==1) then 362 flag_field(idir)="E-field" 363 else 364 flag_field(idir)="D-field" 365 end if 366 end do 367 368 write(message,'(a,a,a,6x,a,11x,a,11x,a,a,a,3(es16.9,2x),a,a,3(es16.9,2x),a,a,3(es16.9,2x),a,a,3(es16.9,2x))') & 369 & ' (a. u.)', ch10,& 370 & ' ', (flag_field(ii),ii=1,3),ch10, & 371 & ' ebar0: ', (dtset%red_efieldbar(ii),ii=1,3), ch10, & 372 & ' ebar: ', (red_efieldbar_lc(ii),ii=1,3), ch10, & 373 & ' d: ', (dtset%red_dfield(ii),ii = 1, 3), ch10, & 374 & ' e + p: ', (red_efield_lc(ii)+dtefield%red_ptot1(ii),ii=1,3) 375 call wrtout(iunit,message,'COLL') 376 377 write(message,'(a,a,3(es16.9,2x),a,a,3(es16.9,2x),a,a,3(es16.9,2x),a,a,3(es16.9,2x))') ch10, & 378 & ' e: ', (red_efield_lc(ii),ii=1,3), ch10, & 379 & ' p: ', (dtefield%red_ptot1(ii), ii=1,3), ch10, & 380 & ' d: ', (dtset%red_dfield(ii),ii = 1, 3), ch10, & 381 & ' e + p: ', (red_efield_lc(ii)+dtefield%red_ptot1(ii),ii=1,3) 382 call wrtout(iunit,message,'COLL') 383 384 write(message,'(a,a,3(es16.9,2x),a,a,3(es16.9,2x),a,a,3(es16.9,2x),a,a,3(es16.9,2x))') ch10, & 385 & ' ebar: ', (red_efieldbar_lc(ii),ii=1,3), ch10, & 386 & ' pbar: ', (red_pbar(ii),ii=1,3), ch10, & 387 & ' dbar: ', (red_dfieldbar(ii),ii=1,3), ch10, & 388 & ' eba+pba: ', (red_efieldbar_lc(ii)+red_pbar(ii),ii=1,3) 389 call wrtout(iunit,message,'COLL') 390 391 write(message,'(a,a,3(es16.9,2x),a,a,3(es16.9,2x),a,a,3(es16.9,2x),a,a,3(es16.9,2x))') ch10, & 392 & ' E: ', (dtefield%efield2(ii), ii=1,3), ch10, & 393 & ' P: ', (ptot_cart(ii), ii=1,3), ch10, & 394 & ' D: ', (dtset%dfield(ii),ii = 1, 3), ch10, & 395 & 'E+4*pi*P: ', (dtset%efield(ii)+4.0d0*pi*ptot_cart(ii),ii=1,3) 396 call wrtout(iunit,message,'COLL') 397 398 write(message,'(a,a,a,a,3(es16.9,2x),a,a,3(es16.9,2x),a,a,3(es16.9,2x),a,a,3(es16.9,2x),a)') ch10, & 399 & ' (S.I.), that is V/m for E, and C/m^2 for P and D', ch10, & 400 & ' E: ', (dtefield%efield2(ii)*(Ha_J/(e_Cb*Bohr_Ang*1d-10)), ii=1,3), ch10, & !(Ha_J/(e_Cb*Bohr_Ang*1d-10))= 5.14220652*1d+11 401 & ' P: ', (ptot_cart(ii)*(e_Cb)/(Bohr_Ang*1d-10)**2, ii=1,3), ch10, & 402 & ' D: ', ((1.0d0/(4*pi))*dtset%dfield(ii)*(e_Cb)/(Bohr_Ang*1d-10)**2,ii = 1, 3), ch10, & 403 & 'eps0*E+P: ', (dtefield%efield2(ii)*eps0*(Ha_J/(e_Cb*Bohr_Ang*1d-10))+ptot_cart(ii)*(e_Cb)/(Bohr_Ang*1d-10)**2,ii=1,3),ch10 404 call wrtout(iunit,message,'COLL') 405 406 end if ! berryopt ==17 407 408 409 410 end subroutine prtefield