1!
2! Copyright (C) 2001-2007 Quantum ESPRESSO group
3! This file is distributed under the terms of the
4! GNU General Public License. See the file `License'
5! in the root directory of the present distribution,
6! or http://www.gnu.org/copyleft/gpl.txt .
7!
8!---------------------------------------------------------------------
9subroutine compute_becalp (becq, alpq)
10  !---------------------------------------------------------------------
11  !
12  !     This routine is used only at finite q and in this case
13  !     computes the scalar product of vkb and psi_{k+q}, and of
14  !     the derivative of vkb and psi_{k+q}. Eq. B8 and B10 (at k+q)
15  !     of PRB 64 235118 (2001).
16  !
17
18  USE kinds, only : DP
19  USE cell_base, ONLY : tpiba
20  USE klist,     ONLY : xk, ngk, igk_k
21  USE gvect,     ONLY : g
22  USE becmod, ONLY: calbec, bec_type, becscal
23  USE buffers, ONLY: get_buffer
24  USE uspp, ONLY: nkb, vkb
25  USE noncollin_module, ONLY : noncolin, npol
26  USE wvfct,    ONLY : nbnd, npwx
27  USE paw_variables, ONLY : okpaw
28
29  USE units_lr, ONLY : lrwfc, iuwfc
30  USE control_ph, ONLY : rec_code_read
31  USE control_lr, ONLY : lgamma
32  USE eqv, ONLY : evq
33  USE qpoint, ONLY : nksq, ikqs
34
35  implicit none
36
37  type (bec_type) :: becq(nksq), alpq(3,nksq)
38  ! the becp with psi_{k+q}
39  ! the alphap with psi_{k+q}
40
41  integer :: ik, ikq, ipol, ibnd, ig, ios, npwq
42  ! counter on k points
43  ! counter on polarizations, bands and
44  ! used for i/o control
45
46  complex(DP) :: fact
47  complex(DP), allocatable :: aux (:,:)
48  !
49  if (lgamma) return
50  IF (rec_code_read >= -20.AND..NOT.okpaw) RETURN
51
52  allocate (aux ( npwx*npol , nbnd))
53  do ik = 1, nksq
54     ikq = ikqs(ik)
55     npwq = ngk(ikq)
56     call init_us_2 (npwq, igk_k(1,ikq), xk (1, ikq), vkb)
57     call get_buffer (evq, lrwfc, iuwfc, ikq)
58     call calbec ( npwq, vkb, evq, becq(ik) )
59     do ipol = 1, 3
60        aux=(0.d0,0.d0)
61        do ibnd = 1, nbnd
62           do ig = 1, npwq
63              aux (ig, ibnd) = evq (ig, ibnd) * &
64                   (xk (ipol, ikq) + g (ipol, igk_k(ig,ikq) ) )
65           enddo
66           IF (noncolin) THEN
67              do ig = 1, npwq
68                 aux (ig+npwx, ibnd) = evq (ig+npwx, ibnd) * &
69                   (xk (ipol, ikq) + g (ipol, igk_k(ig,ikq) ) )
70              enddo
71           ENDIF
72        enddo
73
74        call calbec ( npwq, vkb, aux, alpq(ipol,ik) )
75     enddo
76  enddo
77  fact = CMPLX(0.d0, tpiba,kind=DP)
78
79  DO ik=1,nksq
80     DO ipol=1,3
81        CALL becscal(fact,alpq(ipol,ik),nkb,nbnd)
82     ENDDO
83  ENDDO
84
85  deallocate (aux)
86  return
87end subroutine compute_becalp
88