tran_bulk Subroutine

private subroutine tran_bulk()

Uses

  • proc~~tran_bulk~~UsesGraph proc~tran_bulk tran_bulk module~w90_constants w90_constants proc~tran_bulk->module~w90_constants module~w90_parameters w90_parameters proc~tran_bulk->module~w90_parameters module~w90_io w90_io proc~tran_bulk->module~w90_io module~w90_parameters->module~w90_constants module~w90_parameters->module~w90_io module~w90_io->module~w90_constants

Arguments

None

Calls

proc~~tran_bulk~~CallsGraph proc~tran_bulk tran_bulk proc~io_date io_date proc~tran_bulk->proc~io_date zgemm zgemm proc~tran_bulk->zgemm proc~tran_transfer tran_transfer proc~tran_bulk->proc~tran_transfer proc~tran_green tran_green proc~tran_bulk->proc~tran_green proc~tran_read_htx tran_read_htX proc~tran_bulk->proc~tran_read_htx proc~io_file_unit io_file_unit proc~tran_bulk->proc~io_file_unit zaxpy zaxpy proc~tran_transfer->zaxpy zcopy zcopy proc~tran_transfer->zcopy zgesv zgesv proc~tran_transfer->zgesv proc~io_error io_error proc~tran_transfer->proc~io_error proc~tran_green->zgemm proc~tran_green->zgesv proc~tran_green->proc~io_error proc~tran_read_htx->proc~io_file_unit

Called by

proc~~tran_bulk~~CalledByGraph proc~tran_bulk tran_bulk proc~tran_main tran_main proc~tran_main->proc~tran_bulk program~wannier wannier program~wannier->proc~tran_main proc~wannier_run wannier_run proc~wannier_run->proc~tran_main

Contents

Source Code


Source Code

  subroutine tran_bulk()
    !==================================================================!

    use w90_constants, only: dp, cmplx_0, cmplx_1, cmplx_i, pi
    use w90_io, only: io_error, io_stopwatch, seedname, io_date, &
      io_file_unit, stdout
    use w90_parameters, only: tran_num_bb, tran_read_ht, &
      tran_win_min, tran_win_max, tran_energy_step, &
      timing_level

    implicit none

    integer :: qc_unit, dos_unit
    integer :: ierr
    integer :: n_e, n, i
    real(kind=dp) ::  qc, dos
    real(kind=dp) ::  e_scan
    complex(kind=dp) :: e_scan_cmp
    complex(kind=dp), allocatable, dimension(:, :) :: tot, tott
    complex(kind=dp), allocatable, dimension(:, :) :: g_B, gR, gL
    complex(kind=dp), allocatable, dimension(:, :) :: sLr, sRr
    complex(kind=dp), allocatable, dimension(:, :) :: s1, s2, c1
    character(len=50) :: filename
    character(len=9)  :: cdate, ctime

    if (timing_level > 1) call io_stopwatch('tran: bulk', 1)

    allocate (tot(tran_num_bb, tran_num_bb), stat=ierr)
    if (ierr /= 0) call io_error('Error in allocating tot in tran_bulk')
    allocate (tott(tran_num_bb, tran_num_bb), stat=ierr)
    if (ierr /= 0) call io_error('Error in allocating tott in tran_bulk')
    allocate (g_B(tran_num_bb, tran_num_bb), stat=ierr)
    if (ierr /= 0) call io_error('Error in allocating g_B in tran_bulk')
    allocate (gL(tran_num_bb, tran_num_bb), stat=ierr)
    if (ierr /= 0) call io_error('Error in allocating gL in tran_bulk')
    allocate (gR(tran_num_bb, tran_num_bb), stat=ierr)
    if (ierr /= 0) call io_error('Error in allocating gR in tran_bulk')
    allocate (sLr(tran_num_bb, tran_num_bb), stat=ierr)
    if (ierr /= 0) call io_error('Error in allocating sLr in tran_bulk')
    allocate (sRr(tran_num_bb, tran_num_bb), stat=ierr)
    if (ierr /= 0) call io_error('Error in allocating sRr in tran_bulk')
    allocate (s1(tran_num_bb, tran_num_bb), stat=ierr)
    if (ierr /= 0) call io_error('Error in allocating s1 in tran_bulk')
    allocate (s2(tran_num_bb, tran_num_bb), stat=ierr)
    if (ierr /= 0) call io_error('Error in allocating s2 in tran_bulk')
    allocate (c1(tran_num_bb, tran_num_bb), stat=ierr)
    if (ierr /= 0) call io_error('Error in allocating c1 in tran_bulk')

    call io_date(cdate, ctime)

    qc_unit = io_file_unit()
    open (qc_unit, file=trim(seedname)//'_qc.dat', status='unknown', &
          form='formatted', action='write')
    write (qc_unit, *) '## written on '//cdate//' at '//ctime ! Date and time

    dos_unit = io_file_unit()
    open (dos_unit, file=trim(seedname)//'_dos.dat', status='unknown', &
          form='formatted', action='write')
    write (dos_unit, *) '## written on '//cdate//' at '//ctime ! Date and time

    !   set up the layer hamiltonians

    if (tran_read_ht) then
      allocate (hB0(tran_num_bb, tran_num_bb), stat=ierr)
      if (ierr /= 0) call io_error('Error in allocating hB0 in tran_bulk')
      allocate (hB1(tran_num_bb, tran_num_bb), stat=ierr)
      if (ierr /= 0) call io_error('Error in allocating hB1 in tran_bulk')
      filename = trim(seedname)//'_htB.dat'
      call tran_read_htX(tran_num_bb, hB0, hB1, filename)
    end if

    !   loop over the energies

    n_e = floor((tran_win_max - tran_win_min)/tran_energy_step) + 1

    write (stdout, '(/1x,a)', advance='no') 'Calculating quantum&
        & conductance and density of states...'

    do n = 1, n_e
      e_scan = tran_win_min + real(n - 1, dp)*tran_energy_step

!       if (mod(n,nint(0.1*n_e)).eq.0) write(stdout,'(a)',advance='no') '.'

      ! compute conductance according to Fisher and Lee
      ! retarded Green

      e_scan_cmp = e_scan + eta
      call tran_transfer(tot, tott, hB0, hB1, e_scan_cmp, tran_num_bb)
      call tran_green(tot, tott, hB0, hB1, e_scan, g_B, 0, 1, tran_num_bb)

      ! compute S_Lr and S_Rr

      c1(:, :) = cmplx(hB1(:, :), kind=dp)

      ! Self-energy (Sigma_L^r) : sLr = (hB1)^+ * tott
      ! Self-energy (Sigma_R^r) : sRr = (hB1)   * tot
      sLr = cmplx_0
      sRr = cmplx_0
      call ZGEMM('C', 'N', tran_num_bb, tran_num_bb, tran_num_bb, cmplx_1, c1, &
                 tran_num_bb, tott, tran_num_bb, cmplx_0, sLr, tran_num_bb)
      call ZGEMM('N', 'N', tran_num_bb, tran_num_bb, tran_num_bb, cmplx_1, c1, &
                 tran_num_bb, tot, tran_num_bb, cmplx_0, sRr, tran_num_bb)

      ! Gamma_L = i(Sigma_L^r-Sigma_L^a)
      gL = cmplx_i*(sLr - conjg(transpose(sLr)))
      ! Gamma_R = i(Sigma_R^r-Sigma_R^a)
      gR = cmplx_i*(sRr - conjg(transpose(sRr)))

      s1 = cmplx_0
      s2 = cmplx_0
      c1 = cmplx_0
      ! s1 = Gamma_L * g_B^r
      call ZGEMM('N', 'N', tran_num_bb, tran_num_bb, tran_num_bb, cmplx_1, gL, &
                 tran_num_bb, g_B, tran_num_bb, cmplx_0, s1, tran_num_bb)
      ! s2 = Gamma_L * g_B^r * Gamma_R
      call ZGEMM('N', 'N', tran_num_bb, tran_num_bb, tran_num_bb, cmplx_1, s1, &
                 tran_num_bb, gR, tran_num_bb, cmplx_0, s2, tran_num_bb)
      ! c1 = Gamma_L * g_B^r * Gamma_R * g_B^a
      call ZGEMM('N', 'C', tran_num_bb, tran_num_bb, tran_num_bb, cmplx_1, s2, &
                 tran_num_bb, g_B, tran_num_bb, cmplx_0, c1, tran_num_bb)

      qc = 0.0_dp
      do i = 1, tran_num_bb
        qc = qc + real(c1(i, i), dp)
      end do
!       write(qc_unit,'(f12.6,f15.6)') e_scan, qc
      write (qc_unit, '(f15.9,f18.9)') e_scan, qc

      dos = 0.0_dp
      do i = 1, tran_num_bb
        dos = dos - aimag(g_B(i, i))
      end do
      dos = dos/pi
      write (dos_unit, '(f15.9,f18.9)') e_scan, dos

    end do

    write (stdout, '(a/)') ' done'

    close (qc_unit)
    close (dos_unit)

    deallocate (c1, stat=ierr)
    if (ierr /= 0) call io_error('Error in deallocating c1 in tran_bulk')
    deallocate (s2, stat=ierr)
    if (ierr /= 0) call io_error('Error in deallocating s2 in tran_bulk')
    deallocate (s1, stat=ierr)
    if (ierr /= 0) call io_error('Error in deallocating s1 in tran_bulk')
    deallocate (sRr, stat=ierr)
    if (ierr /= 0) call io_error('Error in deallocating sRr in tran_bulk')
    deallocate (sLr, stat=ierr)
    if (ierr /= 0) call io_error('Error in deallocating sLr in tran_bulk')
    deallocate (gR, stat=ierr)
    if (ierr /= 0) call io_error('Error in deallocating gR in tran_bulk')
    deallocate (gL, stat=ierr)
    if (ierr /= 0) call io_error('Error in deallocating gL in tran_bulk')
    deallocate (g_B, stat=ierr)
    if (ierr /= 0) call io_error('Error in deallocating g_B in tran_bulk')
    deallocate (tott, stat=ierr)
    if (ierr /= 0) call io_error('Error in deallocating tott in tran_bulk')
    deallocate (tot, stat=ierr)
    if (ierr /= 0) call io_error('Error in deallocating tot in tran_bulk')

    if (timing_level > 1) call io_stopwatch('tran: bulk', 2)

    return

  end subroutine tran_bulk