1!
2! Copyright (C) 2001 PWSCF 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!
9!-----------------------------------------------------------------------
10subroutine hdiag( npw, max_iter, avg_iter, et_ )
11  !
12  ! Diagonalizes the unperturbed Hamiltonian in a non-selfconsistent way
13  ! by Conjugate Gradient (band-by-band)
14  !
15  USE kinds,     ONLY : DP
16  USE gvect,     ONLY: g, gstart
17  USE wvfct,     ONLY: g2kin, nbnd, npwx
18  USE uspp,      ONLY: vkb, okvan
19  USE noncollin_module,    ONLY: npol
20  USE wavefunctions,ONLY: evc
21  USE ramanm,    ONLY: eth_ns
22  implicit none
23  !
24  !     I/O variables:
25  !
26  integer :: npw, max_iter
27  ! maximum number of iterations
28  real(DP) :: avg_iter, et_(nbnd)
29  ! iteration number in the diagonalization
30  ! eigenvalues of the diagonalization
31  !
32  !     Local variables:
33  !
34  REAL(DP) :: cg_iter
35  ! number of iteration in CG
36  INTEGER  :: ig, ntry, notconv
37  ! counter on G vectors
38  ! number or repeated call to diagonalization in case of non convergence
39  ! number of notconverged elements
40  INTEGER, ALLOCATABLE :: btype(:)
41    ! type of band: valence (1) or conduction (0)
42  REAl(DP), ALLOCATABLE :: h_prec(:)
43    ! preconditioning matrix (diagonal)
44! CG diagonalization uses these external routines on a single band
45   external hs_1psi, s_1psi
46!  subroutine hs_1psi(npwx,npw,psi,hpsi,spsi)  computes H*psi and S*psi
47!  subroutine s_1psi(npwx,npw,psi,spsi)  computes S*psi (if needed)
48
49  call start_clock ('hdiag')
50
51  allocate (h_prec( npwx), btype(nbnd))
52  !
53  !   various initializations
54  !
55  btype(:) = 1
56  !
57  ! Conjugate-Gradient diagonalization
58  !
59  h_prec=1.0_DP
60  do ig = 1, npw
61     h_prec (ig) = max (1.d0, g2kin (ig) )
62  enddo
63  ntry = 0
6410 continue
65  if (ntry > 0) then
66     call rotate_wfc &
67       ( npwx, npw, nbnd, gstart, nbnd, evc, npol, okvan, evc, et_ )
68     avg_iter = avg_iter + 1.d0
69  endif
70  CALL ccgdiagg( hs_1psi, s_1psi, h_prec, &
71       npwx, npw, nbnd, npol, evc, et_, btype, eth_ns, &
72       max_iter, .true., notconv, cg_iter)
73  avg_iter = avg_iter + cg_iter
74  ntry = ntry + 1
75
76  if (ntry.le.5.and.notconv.gt.0) goto 10
77
78  deallocate (btype, h_prec)
79  call stop_clock ('hdiag')
80
81  return
82end subroutine hdiag
83
84