1! --- 2! Copyright (C) 1996-2016 The SIESTA group 3! This file is distributed under the terms of the 4! GNU General Public License: see COPYING in the top directory 5! or http://www.gnu.org/copyleft/gpl.txt . 6! See Docs/Contributors.txt for a list of contributors. 7! --- 8 9! This code segment has been fully created by: 10! Nick Papior, 2014 11 12module m_tbt_options 13 14 use precision, only : dp 15 16 use units, only: Ang 17 18 use m_ts_tdir, only: ts_tidx 19 20 use m_ts_electype 21 use m_ts_chem_pot 22 23 use dictionary 24 25 implicit none 26 27 ! Common flags for parameters 28 public 29 save 30 31 ! The standard name_prefix 32#ifdef TBT_PHONON 33 character(len=*), parameter :: name_prefix = 'PHT' 34#else 35 character(len=*), parameter :: name_prefix = 'TBT' 36#endif 37 38 ! The temperature 39 real(dp) :: kT 40 41 ! Electrodes and different chemical potentials 42 integer :: N_Elec = 0 43 type(Elec), allocatable, target :: Elecs(:) 44 integer :: N_mu = 0 45 type(ts_mu), allocatable, target :: mus(:) 46 47 ! Whether we should stop right after having created 48 ! the Green's function files 49 logical :: stop_after_GS = .false. 50 51 ! Dictionary to contain the data saving methods 52 ! Each key corresponds to some data calculation 53 ! algorithm. 54 ! To check whether data should be calculated do: 55 ! if ( 'DOS-Gf' .in. save_DATA ) then 56 ! calculate DOS of Gf 57 ! end fi 58 type(dictionary_t) :: save_DATA 59 60 ! Number of eigenchannels to calculate 61 integer :: N_eigen = 0 62 63 ! If the energy-contour is not perfectly divisable by the number of nodes then adjust 64 integer :: BTD_method = 0 ! Optimization method for determining the best tri-diagonal matrix split 65 ! 0 == We optimize for speed 66 ! 1 == We optimize for memory 67 68 ! A quantity describing the accuracy of the coordinates of the 69 ! electrodes. 70 ! * Should only be edited by experienced users * 71 real(dp) :: Elecs_xa_EPS = 1.e-4_dp 72 73 ! Every 5% of the calculation progress it will print an estimation 74 real :: percent_tracker = 5. 75 76#ifdef NCDF_4 77 ! Save file names for data files 78 character(len=250) :: cdf_fname = ' ' 79 character(len=250) :: cdf_fname_sigma = ' ' 80 character(len=250) :: cdf_fname_proj = ' ' 81#endif 82 83 84 ! List of private formats for printing information 85 character(len=*), parameter, private :: f1 ='(''tbt: '',a,t53,''='',tr4,l1)' 86 character(len=*), parameter, private :: f10='(''tbt: '',a,t53,''='',tr4,a)' 87 character(len=*), parameter, private :: f11='(''tbt: '',a)' 88 character(len=*), parameter, private :: f12='(''tbt: '',a,t53,''='',tr2,i0)' 89 character(len=*), parameter, private :: f5 ='(''tbt: '',a,t53,''='',i5,a)' 90 character(len=*), parameter, private :: f20='(''tbt: '',a,t53,''='',i0,'' -- '',i0)' 91 character(len=*), parameter, private :: f6 ='(''tbt: '',a,t53,''='',f10.4,tr1,a)' 92 character(len=*), parameter, private :: f7 ='(''tbt: '',a,t53,''='',f12.6,tr1,a)' 93 character(len=*), parameter, private :: f8 ='(''tbt: '',a,t53,''='',f10.4)' 94 character(len=*), parameter, private :: f9 ='(''tbt: '',a,t53,''='',tr1,e9.3)' 95 character(len=*), parameter, private :: f15='(''tbt: '',a,t53,''='',2(tr1,i0,'' x''),'' '',i0)' 96 97 98contains 99 100 subroutine read_tbt_generic(na_u, lasto) 101 102 use m_region, only: rgn_delete 103 use fdf, only: fdf_defined 104 use m_ts_method, only: ts_init_regions 105 ! This array is never used used, so we delete it 106 use m_ts_method, only: r_pvt 107 108 use m_tbt_diag, only: init_diag 109 110 ! The number of atoms 111 integer, intent(in) :: na_u 112 ! A summated list of last orbitals on atoms. 113 integer, intent(in) :: lasto(0:na_u) 114 115 ! Initialize the buffer regions 116 if ( fdf_defined('TBT.Atoms.Buffer') ) then 117 call ts_init_regions('TBT',na_u,lasto) 118 else 119 call ts_init_regions('TS',na_u,lasto) 120 end if 121 122 call rgn_delete(r_pvt) 123 124 ! Initialize the diagonalization method. 125 call init_diag( ) 126 127 end subroutine read_tbt_generic 128 129 130 ! > Reads the chemical potentials as well as the applied 131 ! Bias. 132 ! The bias is an intricate part of the chemical potential why it 133 ! is read in here. 134 subroutine read_tbt_chem_pot( ) 135 136 use fdf, only : fdf_get 137 use units, only: eV, Kelvin 138 139 use m_ts_chem_pot, only : fdf_nmu, fdffake_mu, fdf_mu, name 140 141 use m_tbt_hs, only: Volt, IsVolt 142 143 implicit none 144 145 ! ******************* 146 ! * LOCAL variables * 147 ! ******************* 148 logical :: err 149 integer :: i 150 151 ! Read in the temperature 152 kT = fdf_get('ElectronicTemperature',1.9e-3_dp,'Ry') 153 kT = fdf_get('TS.ElectronicTemperature',kT,'Ry') 154 kT = fdf_get('TBT.ElectronicTemperature',kT,'Ry') 155 156 ! Read in the chemical potentials 157 N_mu = fdf_nmu('TBT',kT,mus) 158 if ( N_mu < 1 ) then 159 N_mu = fdf_nmu('TS',kT,mus) 160 end if 161 err = .true. 162 if ( N_mu < 1 ) then 163 err = .false. 164 if ( IsVolt ) then 165 ! There is a bias: default 166 ! to two chemical potentials with the 167 ! applied bias 168 N_mu = fdffake_mu(mus,kT,Volt) 169 170 else 171 ! There is no bias. Simply create 172 ! one chemical potential. 173 ! This will make the electrodes 174 ! default to the one chemical potential. 175 N_mu = 1 176 allocate(mus(1)) 177 mus(1)%kT = kT 178 mus(1)%ckT = ' ' 179 mus(1)%name = 'Fermi-level' 180 mus(1)%mu = 0._dp 181 mus(1)%cmu = '0. eV' 182 mus(1)%ID = 1 183 ! These are not used, but we populate 184 ! them anyway 185 mus(1)%N_poles = 1 186 allocate(mus(1)%Eq_seg(1)) 187 mus(1)%Eq_seg(1) = '*NONE' 188 189 end if 190 191 end if 192 do i = 1 , N_mu 193 ! Default things that could be of importance 194 if ( fdf_mu('TBT',mus(i),kT,Volt) ) then 195 ! success 196 else if ( fdf_mu('TS',mus(i),kT,Volt) ) then 197 ! success 198 else if ( err ) then 199 ! only error out if it couldn't be found and forced 200 ! created 201 call die('Could not find chemical potential: ' & 202 //trim(name(mus(i)))) 203 end if 204 end do 205 206#ifdef TBT_PHONON 207 ! Phonon transport cannot define different chemical potentials 208 ! Furthermore, they should be zero 209 do i = 1 , N_mu 210 if ( abs(mus(i)%mu) > 1.e-10 * eV ) then 211 call die('Phonon transport does not define chemical & 212 &potentials. I.e. you cannot lift the frequency spectra.') 213 end if 214 end do 215#endif 216 217 end subroutine read_tbt_chem_pot 218 219 220 ! Reads all information regarding the electrodes, nothing more. 221 subroutine read_tbt_elec( cell, na_u, xa, lasto) 222 223 use fdf, only : fdf_get, fdf_obsolete, fdf_deprecated, leqi 224 use parallel, only : IONode 225 use intrinsic_missing, only : IDX_SPC_PROJ, EYE 226 use intrinsic_missing, only : VEC_PROJ_SCA, VNORM 227 228 use m_os, only : file_exist 229 230 use files, only: slabel 231 use units, only: eV 232 233 use m_tbt_hs, only: spin_idx 234 235 use m_ts_chem_pot, only : copy, chem_pot_add_Elec 236 237 use m_ts_electype, only : fdf_nElec, fdf_Elec 238 use m_ts_electype, only : Name, TotUsedOrbs, TotUsedAtoms 239 use m_ts_electype, only : init_Elec_sim 240 241 use m_ts_method, only : ts_init_electrodes, a_isBuffer 242 243 implicit none 244 245 ! ******************* 246 ! * INPUT variables * 247 ! ******************* 248 real(dp), intent(in) :: cell(3,3) 249 integer, intent(in) :: na_u, lasto(0:na_u) 250 real(dp), intent(in) :: xa(3,na_u) 251 252 ! ******************* 253 ! * LOCAL variables * 254 ! ******************* 255 integer :: i, j 256 real(dp) :: rtmp 257 logical :: err, bool 258 259 if ( N_mu == 0 ) call die('read_tbt_elecs: error in programming') 260 261 ! To determine the same coordinate nature of the electrodes 262 Elecs_xa_EPS = fdf_get('TS.Elecs.Coord.Eps',0.001_dp*Ang, 'Bohr') 263 Elecs_xa_EPS = fdf_get('TBT.Elecs.Coord.Eps',Elecs_xa_EPS,'Bohr') 264 265 ! detect how many electrodes we have 266 N_Elec = fdf_nElec('TBT', Elecs) 267 if ( N_Elec < 1 ) N_Elec = fdf_nElec('TS', Elecs) 268 if ( N_Elec < 1 ) then 269 ! We initialize to 2 electrodes (Left/Right) 270 N_Elec = 2 271 allocate(Elecs(N_Elec)) 272 Elecs(1)%name = 'Left' 273 Elecs(1)%ID = 1 274 Elecs(2)%name = 'Right' 275 Elecs(2)%ID = 2 276 ! if they do-not exist, the user will be told 277 if ( IONode ) then 278 write(*,'(/,''tbt: *** '',a)') 'No electrode names were found, & 279 &default Left/Right are expected' 280 end if 281 end if 282 283 ! Setup default parameters for the electrodes 284 ! first electrode is the "left" 285 ! last electrode is the "right" 286 ! the remaining electrodes have their chemical potential at 0 287 ! Currently the transport direction for all electrodes is the default 288 ! We should probably warn if +2 electrodes are used and t_dir is the 289 ! same for all electrodes... Then the user needs to know what (s)he is doing... 290 Elecs(:)%Bulk = fdf_get('TS.Elecs.Bulk',.true.) ! default everything to bulk electrodes 291 Elecs(:)%Bulk = fdf_get('TBT.Elecs.Bulk',Elecs(1)%Bulk) 292 293 rtmp = fdf_get('TS.Elecs.Eta',0.001_dp*eV,'Ry') 294 rtmp = fdf_get('TBT.Elecs.Eta',rtmp,'Ry') 295#ifdef TBT_PHONON 296 ! eta value needs to be squared as it is phonon spectrum 297 if ( rtmp > 0._dp ) rtmp = rtmp ** 2 298#endif 299 Elecs(:)%Eta = rtmp 300 301 rtmp = fdf_get('TS.Elecs.Accuracy',1.e-13_dp*eV,'Ry') 302 rtmp = fdf_get('TBT.Elecs.Accuracy',rtmp,'Ry') 303 Elecs(:)%accu = rtmp 304 305 ! whether all calculations should be performed 306 ! "out-of-core" i.e. whether the GF files should be created or not 307 ! In tbtrans this is now defaulted to in-core 308 Elecs(:)%out_of_core = fdf_get('TBT.Elecs.Out-of-core',.false.) 309 310 ! Whether we should try and re-use the surface Green function 311 ! files 312 Elecs(:)%ReUseGF = fdf_get('TS.Elecs.GF.ReUse',.true.) 313 Elecs(:)%ReUseGF = fdf_get('TBT.Elecs.GF.ReUse',Elecs(1)%ReUseGF) 314 315 ! Will stop after creating the GF files. 316 stop_after_GS = fdf_get('TBT.Elecs.GF.Only',.false.) 317 318 do i = 1 , N_Elec 319 320 ! If we only have 2 electrodes we take them 321 ! as though the atomic indices are the first and last 322 ! respectively. 323 if ( N_Elec == 2 ) then 324 if ( i == 1 ) then 325 err = fdf_Elec('TBT',slabel,Elecs(i),N_mu,mus,idx_a= 1, & 326 name_prefix = name_prefix) 327 if ( .not. err ) & 328 err = fdf_Elec('TS',slabel,Elecs(i),N_mu,mus,idx_a= 1, & 329 name_prefix = name_prefix) 330 else 331 err = fdf_Elec('TBT',slabel,Elecs(i),N_mu,mus,idx_a=-1, & 332 name_prefix = name_prefix) 333 if ( .not. err ) & 334 err = fdf_Elec('TS',slabel,Elecs(i),N_mu,mus,idx_a=-1, & 335 name_prefix = name_prefix) 336 end if 337 else 338 ! Default things that could be of importance 339 err = fdf_Elec('TBT',slabel,Elecs(i),N_mu,mus) 340 if ( .not. err ) & 341 err = fdf_Elec('TS',slabel,Elecs(i),N_mu,mus, & 342 name_prefix = name_prefix) 343 end if 344 if ( .not. err ) then 345 call die('Could not find electrode: '//trim(name(Elecs(i)))) 346 end if 347 if ( Elecs(i)%idx_a < 0 ) & 348 Elecs(i)%idx_a = na_u + Elecs(i)%idx_a + 1 349 if ( Elecs(i)%idx_a < 1 .or. & 350 na_u < Elecs(i)%idx_a ) then 351 print *,Elecs(i)%idx_a,na_u 352 call die("Electrode position does not exist") 353 end if 354 if ( N_Elec == 2 ) then 355 ! Correct for buffer atoms, first electrode steps "up" 356 ! second electrode steps "down" 357 if ( i == 1 ) then 358 j = Elecs(i)%idx_a 359 do while ( a_isBuffer(j) ) 360 j = j + 1 361 end do 362 Elecs(i)%idx_a = j 363 else 364 j = Elecs(i)%idx_a + TotUsedAtoms(Elecs(i)) - 1 365 do while ( a_isBuffer(j) ) 366 j = j - 1 367 end do 368 Elecs(i)%idx_a = j - TotUsedAtoms(Elecs(i)) + 1 369 end if 370 end if 371 ! set the placement in orbitals 372 Elecs(i)%idx_o = lasto(Elecs(i)%idx_a-1)+1 373 374 ! Initialize electrode parameters 375 call init_Elec_sim(Elecs(i),cell,na_u,xa) 376 377 end do 378 379 ! Initialize the electrode regions 380 call ts_init_electrodes(na_u,lasto,N_Elec,Elecs) 381 382 ! If many electrodes, no transport direction can be specified 383 ! Hence we use this as an error-check (also for N_Elec == 1) 384 if ( any(Elecs(:)%t_dir > 3) ) then 385 ts_tidx = - N_Elec 386 else 387 388 if ( N_Elec /= 2 ) then 389 ! Signals no specific unit-cell direction of transport 390 ts_tidx = - N_Elec 391 else 392 393 ! Retrieve the indices of the unit-cell directions 394 ! according to the electrode transport directions. 395 ! We have already calculated the pivoting table for 396 ! the electrodes 397 i = Elecs(1)%pvt(Elecs(1)%t_dir) 398 j = Elecs(2)%pvt(Elecs(2)%t_dir) 399 400 bool = i == j 401 402 ! For a single transport direction to be true, 403 ! both the projections _has_ to be 1, exactly! 404 rtmp = VEC_PROJ_SCA(cell(:,i), Elecs(1)%cell(:,Elecs(1)%t_dir)) 405 rtmp = rtmp / VNORM(Elecs(1)%cell(:,Elecs(1)%t_dir)) 406 bool = bool .and. abs(abs(rtmp) - 1._dp) < 1.e-5_dp 407 rtmp = VEC_PROJ_SCA(cell(:,j), Elecs(2)%cell(:,Elecs(2)%t_dir)) 408 rtmp = rtmp / VNORM(Elecs(2)%cell(:,Elecs(2)%t_dir)) 409 bool = bool .and. abs(abs(rtmp) - 1._dp) < 1.e-5_dp 410 411 if ( bool ) then 412 413 ! The transport direction for the electrodes are the same... 414 ! And fully encompassed! We have a single transport 415 ! direction. 416 ts_tidx = i 417 418 else 419 420 ! In case we have a skewed transport direction 421 ! we have some restrictions... 422 ts_tidx = -N_Elec 423 424 end if 425 end if 426 end if 427 428 ! Populate the electrodes in the chemical potential type 429 do i = 1 , N_Elec 430 err = .true. 431 do j = 1 , N_mu 432 if ( associated(Elecs(i)%mu,target=mus(j)) ) then 433 call chem_pot_add_Elec(mus(j),i) 434 err = .false. 435 exit 436 end if 437 end do 438 if ( err ) then 439 call die('We could not attribute a chemical potential & 440 &to electrode: '//trim(Elecs(i)%name)) 441 end if 442 end do 443 444 ! check that all electrodes and chemical potentials are paired in 445 ! some way. 446 if ( any(mus(:)%N_El == 0) ) then 447 call die('A/Some chemical potential(s) has/have not been assigned any electrodes. & 448 &All chemical potentials *MUST* be assigned an electrode') 449 end if 450 451 if ( na_u <= sum(TotUsedAtoms(Elecs)) ) then 452 write(*,'(a)') 'Please stop this madness. What where you thinking?' 453 call die('Electrodes occupy the entire device!!!') 454 end if 455 456 end subroutine read_tbt_elec 457 458 459 subroutine read_tbt_after_Elec(nspin, cell, na_u, lasto, xa, no_u, kscell, kdispl) 460 461 use fdf, only: fdf_get, leqi 462 use parallel, only: IONode 463 464 use m_ts_method, only: ts_method, TS_BTD 465 use m_ts_method, only: ts_A_method, TS_BTD_A_COLUMN, TS_BTD_A_PROPAGATION 466 467 use m_tbt_contour, only: read_contour_options 468 469 use m_tbt_sigma_save, only: init_Sigma_options 470 use m_tbt_dH, only: init_dH_options 471 use m_tbt_dSE, only: init_dSE_options 472 473 ! ******************* 474 ! * INPUT variables * 475 ! ******************* 476 integer, intent(in) :: nspin 477 real(dp), intent(in) :: cell(3,3) 478 integer, intent(in) :: na_u, lasto(0:na_u) 479 real(dp), intent(in) :: xa(3,na_u) 480 integer, intent(in) :: no_u 481 integer, intent(in) :: kscell(3,3) 482 real(dp), intent(in) :: kdispl(3) 483 484 integer :: i 485 logical :: ltmp, Gamma3(3), only_T_Gf 486 character(len=150) :: chars 487 488 ! we must have read the electrodes first 489 if ( N_Elec == 0 ) call die('read_tbt_options: Error in programming') 490 491 percent_tracker = fdf_get('TBT.Progress',5.) 492 percent_tracker = max(0., percent_tracker) 493 494 ! Reading the Transiesta solution method 495 chars = fdf_get('TBT.SolutionMethod','BTD') 496 if ( leqi(chars,'BTD').or.leqi(chars,'tri') ) then 497 ts_method = TS_BTD 498 else 499 call die('Unrecognized TBtrans solution method: '//trim(chars)) 500 end if 501 502 chars = fdf_get('TS.BTD.Optimize','speed') 503 chars = fdf_get('TBT.BTD.Optimize',trim(chars)) 504 if ( leqi(chars,'speed') ) then 505 BTD_method = 0 506 else if ( leqi(chars,'memory') ) then 507 BTD_method = 1 508 else 509 call die('Could not determine flag TBT.BTD.Optimize, please & 510 &see manual.') 511 end if 512 513 ! Read spectral calculation method for BTD method 514 if ( N_Elec > 3 ) then 515 chars = fdf_get('TS.BTD.Spectral','column') 516 else 517 chars = fdf_get('TS.BTD.Spectral','propagation') 518 end if 519 chars = fdf_get('TBT.BTD.Spectral',trim(chars)) 520 if ( leqi(chars,'column') ) then 521 ts_A_method = TS_BTD_A_COLUMN 522 else if ( leqi(chars,'propagation') ) then 523 ts_A_method = TS_BTD_A_PROPAGATION 524 else 525 call die('TBT.BTD.Spectral option is not column or propagation. & 526 &Please correct input.') 527 end if 528 529 ! initial 530 only_T_Gf = .true. 531 532 ! Whether we should assert and calculate 533 ! all transmission amplitudes 534 ltmp = fdf_get('TBT.T.Elecs.All',N_Elec == 1) 535 ltmp = fdf_get('TBT.T.All',N_Elec == 1) 536 if ( ltmp ) then 537 save_DATA = save_DATA // ('T-all'.kv.1) 538 end if 539 540 N_eigen = fdf_get('TBT.T.Eig',0) 541 if ( N_eigen > 0 ) then 542 save_DATA = save_DATA // ('T-eig'.kv.N_eigen) 543 only_T_Gf = .false. 544 end if 545 546 ! Should we calculate DOS of electrode bulk Green function 547 ltmp = fdf_get('TBT.DOS.Elecs', .false. ) 548 if ( ltmp ) then 549 save_DATA = save_DATA // ('DOS-Elecs'.kv.1) 550 end if 551 552 ltmp = fdf_get('TBT.T.Bulk', N_Elec == 1) 553 if ( ltmp ) then 554 ! when calculating the DOS for the electrode 555 ! we also get the bulk transmission. 556 save_DATA = save_DATA // ('DOS-Elecs'.kv.1) 557 end if 558 559 ! Should we calculate DOS of Green function 560 ltmp = fdf_get('TBT.DOS.Gf', N_Elec == 1) 561 if ( ltmp ) then 562 save_DATA = save_DATA // ('DOS-Gf'.kv.1) 563 only_T_Gf = .false. 564 end if 565 566 ! Should we calculate DOS of spectral function 567 ltmp = fdf_get('TBT.DOS.A', N_Elec == 1) 568 if ( ltmp ) then 569 save_DATA = save_DATA // ('DOS-A'.kv.1) 570 only_T_Gf = .false. 571 end if 572 573 ! Should we calculate DOS of all spectral functions 574 ltmp = fdf_get('TBT.DOS.A.All', N_Elec == 1) 575 if ( ltmp ) then 576 save_DATA = save_DATA // ('DOS-A'.kv.1) 577 save_DATA = save_DATA // ('DOS-A-all'.kv.1) 578 only_T_Gf = .false. 579 end if 580 581 ! Should we calculate orbital current 582 ltmp = fdf_get('TBT.Current.Orb', .false. ) 583 if ( ltmp ) then 584 save_DATA = save_DATA // ('DOS-A'.kv.1) 585 save_DATA = save_DATA // ('orb-current'.kv.1) 586 only_T_Gf = .false. 587 end if 588 589 ! Options for density-matrix calculations 590 ltmp = fdf_get('TBT.DM.Gf', .false.) 591 if ( ltmp ) then 592 save_DATA = save_DATA // ('DOS-Gf'.kv.1) 593 save_DATA = save_DATA // ('DM-Gf'.kv.1) 594 end if 595 596 ltmp = fdf_get('TBT.DM.A', .false.) 597 if ( ltmp ) then 598 save_DATA = save_DATA // ('DOS-A'.kv.1) 599 save_DATA = save_DATA // ('DM-A'.kv.1) 600 end if 601 602 603 ! Options for COOP and COHP curves. 604 ! These are orbital (energy) populations that can be used to determine the 605 ! bonding nature of the material. 606 ltmp = fdf_get('TBT.COOP.Gf', .false.) 607 if ( ltmp ) then 608 save_DATA = save_DATA // ('DOS-Gf'.kv.1) 609 save_DATA = save_DATA // ('COOP-Gf'.kv.1) 610 end if 611 612 ltmp = fdf_get('TBT.COOP.A', .false. ) 613 if ( ltmp ) then 614 save_DATA = save_DATA // ('DOS-A'.kv.1) 615 save_DATA = save_DATA // ('COOP-A'.kv.1) 616 end if 617 618 ltmp = fdf_get('TBT.COHP.Gf', .false. ) 619 if ( ltmp ) then 620 save_DATA = save_DATA // ('DOS-Gf'.kv.1) 621 save_DATA = save_DATA // ('COHP-Gf'.kv.1) 622 end if 623 624 ltmp = fdf_get('TBT.COHP.A', .false. ) 625 if ( ltmp ) then 626 save_DATA = save_DATA // ('DOS-A'.kv.1) 627 save_DATA = save_DATA // ('COHP-A'.kv.1) 628 end if 629 630 ltmp = fdf_get('TBT.T.Out',N_Elec == 1) 631 if ( ltmp ) then 632 save_DATA = save_DATA // ('T-sum-out'.kv.1) 633 end if 634 635 ! We cannot calculate the transmission for more than 3 636 ! electrodes if using the diagonal 637 if ( only_T_Gf .and. N_Elec > 3 ) then 638 only_T_Gf = .false. 639 end if 640 641 if ( only_T_Gf ) then 642 ! TODO, consider changing this to .true. 643 only_T_Gf = fdf_get('TBT.T.Gf',.false.) 644 end if 645 if ( only_T_Gf ) then 646 save_DATA = save_DATA // ('T-Gf'.kv.1) 647 ts_A_method = TS_BTD_A_COLUMN 648 649 ! We get sum-out for free, so save it 650 save_DATA = save_DATA // ('T-sum-out'.kv.1) 651 652 if ( N_Elec > 2 ) then 653 ! When having more than one electrode 654 ! we *must* calculate all transmissions 655 ! to be able to get the actual transmissions 656 ! using the diagonal Green function 657 save_DATA = save_DATA // ('T-all'.kv.1) 658 end if 659 660 end if 661 662 call init_Sigma_options( save_DATA ) 663 664 if ( IONode ) write(*,*) ! new-line 665 ! The init_dH_options also checks for the 666 ! consecutive sparse patterns. 667 call init_dH_options( no_u ) 668 call init_dSE_options( ) 669 670 671 ! read in contour options 672 call read_contour_options( N_Elec, Elecs, N_mu, mus ) 673 674 ! Check for Gamma in each direction 675 do i = 1 , 3 676 if ( kdispl(i) /= 0._dp ) then 677 Gamma3(i) = .false. 678 else if ( sum(kscell(i,:)) > 1 ) then 679 ! Note it is the off-diagonal for this unit-cell 680 ! direction 681 Gamma3(i) = .false. 682 else 683 Gamma3(i) = .true. 684 end if 685 end do 686 687 do i = 1 , N_Elec 688 ! Initialize the electrode quantities for the stored values 689 call check_Elec_sim(Elecs(i), nspin, cell, na_u, xa, & 690 Elecs_xa_EPS, lasto, Gamma3) 691 end do 692 693 end subroutine read_tbt_after_Elec 694 695 subroutine print_tbt_options(nspin) 696 697 use units, only: Kelvin, eV 698 use parallel, only: IONode 699 use files, only: slabel 700 701 use m_ts_method, only: ts_A_method, TS_BTD_A_COLUMN, TS_BTD_A_PROPAGATION 702 703 use m_tbt_contour, only: print_contour_tbt_options, io_contour_tbt 704 use m_tbt_contour, only: print_contour_tbt_block 705 use m_tbt_save, only: print_save_options 706 use m_tbt_diag, only: print_diag 707 use m_tbt_dH, only: print_dH_options 708 use m_tbt_dSE, only: print_dSE_options 709 use m_tbt_sigma_save, only: print_Sigma_options 710 use m_tbt_proj, only: print_proj_options 711 use m_tbt_hs, only: Volt, IsVolt, spin_idx 712 713 integer, intent(in) :: nspin 714 715 integer :: i 716 717 if ( .not. IONode ) return 718 719 write(*,*) 720 write(*,f11) repeat('*', 62) 721 722 write(*,f6) 'Electronic temperature (reference)',kT/Kelvin,'K' 723 if ( IsVolt ) then 724 write(*,f6) 'Voltage', Volt/eV,'Volts' 725 else 726 write(*,f11) 'No applied bias' 727 end if 728 write(*,f1) 'Calculate transmission only using diag(Gf)',('T-Gf'.in.save_DATA) 729 write(*,f1) 'Saving bulk transmission for electrodes',('DOS-Elecs'.in.save_DATA) 730 write(*,f1) 'Saving DOS from bulk electrodes',('DOS-Elecs'.in.save_DATA) 731 write(*,f1) 'Saving DOS from Green function',('DOS-Gf'.in.save_DATA) 732 if ( 'DOS-A-all' .in. save_DATA ) then 733 write(*,f1) 'Saving DOS from all spectral functions',.true. 734 else 735 write(*,f1) 'Saving DOS from spectral functions',('DOS-A' .in. save_DATA) 736 end if 737 write(*,f1) 'Saving bond currents (orb-orb)',('orb-current'.in.save_DATA) 738 739 ! DM 740 write(*,f1) 'Saving DM from Green function',('DM-Gf'.in.save_DATA) 741 write(*,f1) 'Saving DM from spectral functions',('DM-A'.in.save_DATA) 742 743 ! COOP/COHP curves 744 write(*,f1) 'Saving COOP from Green function',('COOP-Gf'.in.save_DATA) 745 write(*,f1) 'Saving COOP from spectral functions',('COOP-A'.in.save_DATA) 746 write(*,f1) 'Saving COHP from Green function',('COHP-Gf'.in.save_DATA) 747 write(*,f1) 'Saving COHP from spectral functions',('COHP-A'.in.save_DATA) 748 749 write(*,f12) 'Calc. # transmission eigenvalues',N_eigen 750 write(*,f1) 'Calc. T between all electrodes',('T-all'.in.save_DATA) 751 write(*,f1) 'Calc. total T out of electrodes',('T-sum-out'.in.save_DATA) 752 if ( nspin > 1 ) then 753 if ( spin_idx == 0 ) then 754 write(*,f11) 'Calculate spin UP and DOWN' 755 else if ( spin_idx == 1 ) then 756 write(*,f10) 'Calculate spin ','UP' 757 else if ( spin_idx == 2 ) then 758 write(*,f10) 'Calculate spin ','DOWN' 759 else 760 call die('Error in spin_idx') 761 end if 762 else 763 write(*,f11) 'Non-polarized Hamiltonian' 764 end if 765 766 767 ! Algorithm choices... 768 select case ( BTD_method ) 769 case ( 0 ) 770 write(*,f10)'BTD creation algorithm', 'speed' 771 case ( 1 ) 772 write(*,f10)'BTD creation algorithm', 'memory' 773 end select 774 select case ( ts_A_method ) 775 case ( TS_BTD_A_PROPAGATION ) 776 write(*,f10)'BTD spectral function algorithm','propagation' 777 case ( TS_BTD_A_COLUMN ) 778 write(*,f10)'BTD spectral function algorithm','column' 779 case default 780 call die('Error in setup BTD. A calc') 781 end select 782 783 784 call print_diag() 785 call print_Sigma_options( save_DATA ) 786 call print_dH_options( ) 787 call print_dSE_options( ) 788 789 call print_save_options() 790 call print_proj_options( save_DATA ) 791 792 write(*,f11)' >> Electrodes << ' 793 do i = 1 , size(Elecs) 794 call print_settings(Elecs(i),'tbt') 795 end do 796 797 call print_contour_tbt_options( 'TBT' ) 798 799 write(*,f11) repeat('*', 62) 800 write(*,*) 801 802 call io_contour_tbt(slabel) 803 804 write(*,f11) repeat('<', 62) 805 806 call print_contour_tbt_block( 'TBT' ) 807 808 write(*,f11) repeat('<', 62) 809 810 end subroutine print_tbt_options 811 812 subroutine print_tbt_warnings( Gamma ) 813 814 use fdf, only: fdf_get, fdf_defined 815 use units, only: eV 816 use parallel, only: IONode 817 818 use m_tbt_save, only: save_parallel 819 820 use m_tbt_dH, only: print_dH_warnings 821 use m_tbt_hs, only: Volt 822 823 ! Whether the user requests a Gamma calculation 824 logical, intent(in) :: Gamma 825 826 integer :: i 827 logical :: ltmp, rem_DOS_Elecs, rem_T_Gf, has 828 829 ! Removal of keys 830 rem_DOS_Elecs = 'DOS-Elecs' .in. save_DATA 831 if ( rem_DOS_Elecs ) then 832 ! Currently the TBTGF files does not contain the DOS 833 if ( all(Elecs(:)%out_of_core) ) & 834 call delete(save_DATA,'DOS-Elecs') 835 if ( 'Sigma-only' .in. save_DATA ) & 836 call delete(save_DATA,'DOS-Elecs') 837 end if 838 839 rem_T_Gf = 'T-Gf' .in. save_DATA 840 if ( N_Elec > 3 ) then 841 call delete(save_DATA,'T-Gf') 842 end if 843 844 if ( .not. IONode ) return 845 846 write(*,*) ! new-line 847 write(*,'(3a)') repeat('*',24),' Begin: TBT CHECKS AND WARNINGS ',repeat('*',24) 848 849 if ( N_eigen < 0 ) then 850 call die('Number of transmission eigenvalues MUST be & 851 &zero or positive.') 852 end if 853 854 if ( .not. Gamma ) then 855 ltmp = .not. fdf_get('SpinSpiral',.false.) 856 ltmp = fdf_get('TBT.Symmetry.TimeReversal',ltmp) 857 if ( ('orb-current' .in.save_DATA) ) then 858 if ( IONode .and. ltmp ) then 859 write(*,'(a,/,a)') 'WARNING: k-averaging orbital currents with & 860 &time-reversal symmetry will not reproduce','the correct & 861 &orbital current. Set TBT.Symmetry.TimeReversal F' 862 end if 863 end if 864 865 ! Also for COOP analysis 866 has = 'COOP-Gf' .in. save_DATA 867 has = has .or. ('COOP-A' .in. save_DATA) 868 has = has .or. ('COHP-Gf' .in. save_DATA) 869 has = has .or. ('COHP-A' .in. save_DATA) 870 if ( IONode .and. ltmp .and. has ) then 871 write(*,'(a,/,a)') 'WARNING: k-averaging COOP/COHP with & 872 &time-reversal symmetry will not reproduce','the correct & 873 &populations. Set TBT.Symmetry.TimeReversal F' 874 end if 875 end if 876 877 ! CHECK THIS (we could allow it by only checking the difference...) 878 if ( maxval(mus(:)%mu) - minval(mus(:)%mu) - abs(Volt) > 1.e-8_dp * eV ) then 879 if ( IONode ) then 880 write(*,'(a)') 'Chemical potentials [eV]:' 881 do i = 1 , N_Elec 882 write(*,'(a,f10.5,a)') trim(Name(Elecs(i)))//' at ',Elecs(i)%mu%mu/eV,' eV' 883 end do 884 write(*,'(a)') 'The difference must satisfy: "max(ChemPots)-min(ChemPots) - abs(Volt) < 1e-8 eV"' 885 write(*,'(a,f10.5,a)') 'max(ChemPots) at ', maxval(mus(:)%mu)/eV,' eV' 886 write(*,'(a,f10.5,a)') 'min(ChemPots) at ', minval(mus(:)%mu)/eV,' eV' 887 write(*,'(a,f10.5,a)') '|V| at ', abs(Volt)/eV,' eV' 888 end if 889 call die('Chemical potentials are not consistent with the bias applied.') 890 end if 891 892 893 ! Check that the bias does not introduce a gating 894 if ( any(abs(mus(:)%mu) - abs(Volt) > 1.e-9_dp ) ) then 895 write(*,'(a)') 'Chemical potentials must lie in the range [-V;V] with the maximum & 896 &difference being V' 897 call die('Chemical potentials must not introduce consistent Ef shift to the system.') 898 end if 899 900 if ( rem_DOS_Elecs ) then 901 if ( any(Elecs(:)%out_of_core) ) then 902 write(*,'(a)')' Disabling electrode DOS calculation, only & 903 &enabled for in-core self-energy calculations.' 904 end if 905 if ( 'Sigma-only' .in. save_DATA ) then 906 write(*,'(a)')' Disabling electrode DOS calculation, only & 907 &enabled when calculating transmission (not TBT.SelfEnergy.Only).' 908 end if 909 end if 910 911 if ( rem_T_Gf .and. N_Elec > 3 ) then 912 write(*,'(a)')' ** Disabling transport calculation using diagonal, & 913 ¬ possible with N_elec > 3.' 914 end if 915 916 do i = 1, N_Elec 917 if ( Elecs(i)%repeat .and. Elecs(i)%bloch%size() > 1 ) then 918 write(*,'(a)') 'Electrode '//trim(Elecs(i)%name)//' is & 919 &using Bloch unfolding using the repeat scheme! & 920 &Please use the tiling scheme (it is orders of magnitudes faster!).' 921 end if 922 end do 923 924#ifdef MPI 925#ifdef NCDF_PARALLEL 926 if ( .not. save_parallel ) then 927 write(*,'(a)') ' ** Speed up the execution by utilizing parallel I/O' 928 write(*,'(a)') ' > TBT.CDF.MPI true' 929 end if 930#endif 931#endif 932 933 call print_dH_warnings( save_DATA ) 934 935 write(*,'(3a,/)') repeat('*',24),' End: TBT CHECKS AND WARNINGS ',repeat('*',26) 936 937 end subroutine print_tbt_warnings 938 939end module m_tbt_options 940