Line data    Source code

       1             : module m_calc_mpsmat
       2             :    use m_types
       3             :    use m_judft
       4             :    use m_constants
       5             : contains
       6           0 :    subroutine calc_mpsmat(fi, mpdata, smat)
       7             :       implicit none
       8             :       type(t_mat), intent(inout) :: smat
       9             :       type(t_fleurinput), intent(in)    :: fi
      10             :       TYPE(t_mpdata), intent(in)        :: mpdata
      11             : 
      12             :       integer :: igpt2, igpt1
      13             : 
      14           0 :       call timestart("calc smat")
      15           0 :       call smat%alloc(.False., mpdata%num_gpts(),mpdata%num_gpts())
      16             :       !$OMP PARALLEL DO default(shared) schedule(guided)&
      17             :       !$OMP private(igpt2,igpt1)&
      18           0 :       !$OMP shared(mpdata, fi, smat)
      19             :       DO igpt2 = 1, mpdata%num_gpts()
      20             :          DO igpt1 = 1, mpdata%num_gpts() !igpt2
      21             :             smat%data_c(igpt1, igpt2) = calc_smat_elem(fi, mpdata, igpt1, igpt2)
      22             :          END DO
      23             :       END DO
      24             :       !$OMP END PARALLEL DO
      25             : 
      26             :       ! call smat%u2l()
      27           0 :       call timestop("calc smat")
      28           0 :    end subroutine calc_mpsmat
      29             : 
      30       56410 :    function calc_smat_elem(fi, mpdata, igpt1, igpt2) result(elem)
      31             :          ! Calculate the hermitian matrix smat(i,j) = sum(a) integral(MT(a)) exp[i(Gj-Gi)r] dr
      32             :       implicit none
      33             :       type(t_fleurinput), intent(in)  :: fi
      34             :       TYPE(t_mpdata), intent(in)      :: mpdata 
      35             :       integer, intent(in)             :: igpt1, igpt2
      36             :       complex :: elem
      37             : 
      38             :       integer :: g(3), ic, itype
      39             :       real    :: gnorm, rdum 
      40             : 
      41       56410 :       elem = 0.0
      42      225640 :       g = mpdata%g(:, igpt2) - mpdata%g(:, igpt1)
      43       56410 :       gnorm = gptnorm(g, fi%cell%bmat)
      44       56410 :       IF (abs(gnorm) < 1e-12) THEN
      45        2082 :          DO itype = 1, fi%atoms%ntype
      46        2082 :             elem = elem + fi%atoms%neq(itype)*fpi_const*fi%atoms%rmt(itype)**3/3
      47             :          END DO
      48             :       ELSE
      49      165706 :          do ic = 1,fi%atoms%nat
      50      110012 :             itype = fi%atoms%itype(ic)
      51      110012 :             rdum = fi%atoms%rmt(itype)*gnorm
      52      110012 :             rdum = fpi_const*(SIN(rdum) - rdum*COS(rdum))/gnorm**3
      53      495742 :             elem = elem + rdum*EXP(ImagUnit*tpi_const*dot_PRODUCT(fi%atoms%taual(:, ic), g))
      54             :          END DO
      55             :       END IF
      56       56410 :    end function calc_smat_elem
      57             : end module m_calc_mpsmat