Commit 3bf3af90 authored by Stefano de Gironcoli's avatar Stefano de Gironcoli
Browse files

update of the different solvers to report similar quantities: avg_iter, total nhpsi, notcnv, ethr

parent 4969f4bc
......@@ -25,9 +25,7 @@ PARA_POSTFIX=" -nk 1 -nd 4 -nb 4 -nt 1 "
#
BIN_DIR="$EXAMPLE_DIR/../../bin"
BIN_LIST="cb_cg cb_cg_gamma cb_davidson cb_davidson_gamma cb_ppcg cb_ppcg_gamma cb_paro cb_paro_gamma"
BIN_LIST="cb_davidson_gamma cb_paro_gamma cb_paro_new_gamma"
BIN_LIST="cb_paro_new"
BIN_LIST="cb_cg cb_cg_gamma cb_davidson cb_davidson_gamma cb_ppcg cb_ppcg_gamma cb_paro cb_paro_gamma cb_paro_new cb_paro_gamma_new"
for j in $BIN_LIST
do
......
......@@ -14,7 +14,8 @@
complex(DP), allocatable :: evc(:,:)
real(dp), allocatable :: eig(:)
integer, parameter :: npol=1
integer :: notconv, cg_iter, ig
integer :: notconv, nhpsi, ig
real(dp):: cg_iter
integer :: max_cg_iter = 100 !default of QE
real(dp), allocatable :: h_diag(:,:) !in case of CG - the preconditioner
logical :: overlap = .false. , lrot =.false. , lscf = .true. ! lscf is true for scf calc
......@@ -75,6 +76,7 @@ call set_mpi_comm_4_solvers( world_comm, intra_bgrp_comm, inter_bgrp_comm )
npwx, npw, nbnd, nbnd, evc, evc, eig )
endif
#endif
nhpsi = nbnd
write(stdout,*) ' then cg diagonalization '
h_diag = 1.D0
......@@ -86,6 +88,7 @@ call set_mpi_comm_4_solvers( world_comm, intra_bgrp_comm, inter_bgrp_comm )
npwx, npw, nbnd, evc, eig, btype, &
ethr, max_cg_iter, .NOT. lscf, notconv, cg_iter )
nhpsi = nhpsi + cg_iter*nbnd
!--------------------------------------------------------------------------------------------------------------!
call stop_clock('cg')
......@@ -93,7 +96,7 @@ call set_mpi_comm_4_solvers( world_comm, intra_bgrp_comm, inter_bgrp_comm )
if (energy_shift .and. current_k==1) ref=eig(4*ncell**3)
call write_bands(eig,ref)
write (6,*) 'cg_iter, notconv, ethr ',cg_iter,notconv, ethr
write (6,*) 'cg_iter, nhpsi, notconv, ethr ',int(cg_iter), nhpsi, notconv, ethr
deallocate( eig )
end do
......@@ -109,6 +112,12 @@ call set_mpi_comm_4_solvers( world_comm, intra_bgrp_comm, inter_bgrp_comm )
call print_clock('s_psi')
call print_clock('hs_1psi')
call print_clock('s_1psi')
!
write (6,*)
write (6,*) ' general FFT routines'
call print_clock('fftw')
call print_clock('ffts')
#if defined(__MPI)
call unset_mpi_comm_4_solvers
call mp_global_end( )
......
......@@ -14,7 +14,8 @@
complex(DP), allocatable :: evc(:,:)
real(dp), allocatable :: eig(:)
integer, parameter :: npol=1
integer :: notconv, cg_iter, ig
integer :: notconv, nhpsi, ig
real(dp):: cg_iter
integer :: max_cg_iter = 100 !default of QE
real(dp), allocatable :: h_diag(:,:) !in case of CG - the preconditioner
logical :: overlap = .false. , lrot =.false. , lscf = .true. ! lscf is true for scf calc
......@@ -75,6 +76,7 @@ call set_mpi_comm_4_solvers( world_comm, intra_bgrp_comm, inter_bgrp_comm )
npwx, npw, nbnd, nbnd, npol, evc, evc, eig )
endif
#endif
nhpsi = nbnd
write(stdout,*) ' then cg diagonalization '
h_diag = 1.D0
......@@ -85,7 +87,7 @@ call set_mpi_comm_4_solvers( world_comm, intra_bgrp_comm, inter_bgrp_comm )
CALL ccgdiagg( cb_hs_1psi, cb_s_1psi, h_diag, &
npwx, npw, nbnd, npol, evc, eig, btype, &
ethr, max_cg_iter, .NOT. lscf, notconv, cg_iter )
nhpsi = nhpsi + cg_iter*nbnd
!--------------------------------------------------------------------------------------------------------------!
call stop_clock('cg')
......@@ -93,7 +95,7 @@ call set_mpi_comm_4_solvers( world_comm, intra_bgrp_comm, inter_bgrp_comm )
if (energy_shift .and. current_k==1) ref=eig(4*ncell**3)
call write_bands(eig,ref)
write (6,*) 'cg_iter, notconv, ethr ',cg_iter,notconv, ethr
write (6,*) 'cg_iter, nhpsi, notconv, ethr ',int(cg_iter), nhpsi, notconv, ethr
deallocate( eig )
end do
......@@ -109,6 +111,12 @@ call set_mpi_comm_4_solvers( world_comm, intra_bgrp_comm, inter_bgrp_comm )
call print_clock('s_psi')
call print_clock('hs_1psi')
call print_clock('s_1psi')
!
write (6,*)
write (6,*) ' general FFT routines'
call print_clock('fftw')
call print_clock('ffts')
#if defined(__MPI)
call unset_mpi_comm_4_solvers
call mp_global_end( )
......
......@@ -97,6 +97,12 @@
call print_clock('h_psi')
call print_clock('s_psi')
call print_clock('g_psi')
!
write (6,*)
write (6,*) ' general FFT routines'
call print_clock('fftw')
call print_clock('ffts')
#if defined(__MPI)
call unset_mpi_comm_4_solvers
call mp_global_end( )
......
......@@ -14,7 +14,7 @@
complex(DP), allocatable :: evc(:,:)
real(dp), allocatable :: eig(:)
integer, parameter :: npol=1
integer :: notcnv, dav_iter
integer :: notcnv, dav_iter, nhpsi
logical :: overlap = .false. , lrot =.false.
! additional local variables
real(dp) :: ref=0.d0
......@@ -64,12 +64,12 @@
#endif
call cegterg( cb_h_psi, cb_s_psi, overlap, cb_g_psi, &
npw, npwx, nbnd, nbndx, npol, evc, ethr, &
eig, btype, notcnv, lrot, dav_iter )
eig, btype, notcnv, lrot, dav_iter, nhpsi )
#if defined(__MPI)
else
call pcegterg(cb_h_psi, cb_s_psi, overlap, cb_g_psi, &
npw, npwx, nbnd, nbndx, npol, evc, ethr, &
eig, btype, notcnv, lrot, dav_iter )
eig, btype, notcnv, lrot, dav_iter, nhpsi )
end if
#endif
!--------------------------------------------------------------------------------------------------------------!
......@@ -79,7 +79,7 @@
if (energy_shift .and. current_k==1) ref=eig(4*ncell**3)
call write_bands(eig,ref)
write (stdout,*) 'dav_iter, notcnv, ethr ',dav_iter,notcnv, ethr
write (stdout,*) 'dav_iter, nhpsi, notcnv, ethr ', dav_iter, nhpsi, notcnv, ethr
deallocate( eig )
end do
......@@ -96,6 +96,11 @@
call print_clock('h_psi')
call print_clock('s_psi')
call print_clock('g_psi')
!
write (6,*)
write (6,*) ' general FFT routines'
call print_clock('fftw')
call print_clock('ffts')
#if defined(__MPI)
call unset_mpi_comm_4_solvers
......
......@@ -90,13 +90,10 @@ PROGRAM cb_paro_gamma_main
#endif
WRITE(stdout,*) ' then paro diagonalization '
!in QE iter and ntry would be used here to determine the lrot, which would call rotate_wfc.
CALL paro_gamma(cb_h_psi, cb_s_psi, cb_hs_1psi, cb_g_1psi, overlap, &
npwx, npw, nbnd, evc, eig, btype, ethr, notconv, nhpsi)
! CALL rotate_wfc_gamma(cb_h_psi, cb_s_psi, &
! npwx, npw, nbnd, nbnd, evc, overlap, evc, eig )
avg_iter = avg_iter + nhpsi/float(nbnd) ; !if (lrot) avg_iter = avg_iter +1.d0
avg_iter = avg_iter + nhpsi/float(nbnd) ; nhpsi = nhpsi + nbnd
CALL stop_clock('paro')
DEALLOCATE ( evc )
......
......@@ -82,22 +82,18 @@ PROGRAM cb_paro_main
ENDIF
#endif
WRITE(stdout,*) ' then paro diagonalization '
!support only serial for now
!in QE iter and ntry would be used here to determine the lrot, which would call rotate_wfc.
CALL paro_k(cb_h_psi, cb_s_psi, cb_hs_1psi, cb_g_1psi, overlap, &
npwx, npw, nbnd, npol, evc, eig, btype, ethr, notconv, nhpsi)
! CALL rotate_wfc_k(cb_h_psi, cb_s_psi, &
! npwx, npw, nbnd, nbnd, npol, evc, overlap, evc, eig )
avg_iter = avg_iter + nhpsi/float(nbnd) ; !if (lrot) avg_iter = avg_iter +1.d0
avg_iter = avg_iter + nhpsi/float(nbnd) ; nhpsi = nhpsi + nbnd
CALL stop_clock('paro')
DEALLOCATE ( evc )
IF (energy_shift .AND. current_k==1) ref=eig(4*ncell**3)
CALL write_bands(eig,ref)
WRITE (6,*) 'avg_iter, notconv, ethr ',avg_iter,notconv, ethr
WRITE (6,*) 'avg_iter, nhpsi, notconv, ethr ', avg_iter ,nhpsi, notconv, ethr
DEALLOCATE( eig )
END DO
......@@ -111,6 +107,12 @@ PROGRAM cb_paro_main
CALL print_clock('pcg:hs_1psi')
CALL print_clock('pcg:ortho')
CALL print_clock('s_psi')
!
write (6,*)
write (6,*) ' general FFT routines'
call print_clock('fftw')
call print_clock('ffts')
#if defined(__MPI)
CALL unset_mpi_comm_4_solvers
CALL mp_global_end( )
......
......@@ -14,7 +14,8 @@
complex(DP), allocatable :: evc(:,:)
real(dp), allocatable :: eig(:)
integer, parameter :: npol=1
integer :: notconv, ppcg_iter, ig
integer :: notconv, nhpsi, ig
real(dp):: ppcg_iter
integer :: max_ppcg_iter = 100 !default of QE
real(dp), allocatable :: h_diag(:,:) !in case of CG - the preconditioner
logical :: overlap = .false. , lrot =.false. , lscf = .true. ! lscf is true for scf calc
......@@ -81,6 +82,7 @@ call set_mpi_comm_4_solvers( world_comm, intra_bgrp_comm, inter_bgrp_comm )
npwx, npw, nbnd, nbnd, evc, evc, eig )
endif
#endif
nhpsi = nbnd
write(stdout,*) ' then ppcg diagonalization '
!support only serial for now
h_diag = 1.D0
......@@ -91,7 +93,7 @@ call set_mpi_comm_4_solvers( world_comm, intra_bgrp_comm, inter_bgrp_comm )
CALL ppcg_gamma( cb_h_psi, cb_s_psi, overlap, h_diag, &
npwx, npw, nbnd, evc, eig, btype, &
ethr, max_ppcg_iter, notconv, ppcg_iter, sbsize , rrstep, iter )
nhpsi = nhpsi + ppcg_iter * nbnd
!--------------------------------------------------------------------------------------------------------------!
call stop_clock('ppcg')
......@@ -99,7 +101,7 @@ call set_mpi_comm_4_solvers( world_comm, intra_bgrp_comm, inter_bgrp_comm )
if (energy_shift .and. current_k==1) ref=eig(4*ncell**3)
call write_bands(eig,ref)
write (6,*) 'ppcg_iter, notconv, ethr ',ppcg_iter,notconv, ethr
write (6,*) 'ppcg_iter, nhpsi, notconv, ethr ', int(ppcg_iter), nhpsi, notconv, ethr
deallocate( eig )
end do
......@@ -118,6 +120,12 @@ call set_mpi_comm_4_solvers( world_comm, intra_bgrp_comm, inter_bgrp_comm )
call print_clock('ppcg:RR')
call print_clock('ppcg:DTRSM')
call print_clock('ppcg:lock')
!
write (6,*)
write (6,*) ' general FFT routines'
call print_clock('fftw')
call print_clock('ffts')
#if defined(__MPI)
call unset_mpi_comm_4_solvers
call mp_global_end( )
......
......@@ -14,7 +14,8 @@
complex(DP), allocatable :: evc(:,:)
real(dp), allocatable :: eig(:)
integer, parameter :: npol=1
integer :: notconv, ppcg_iter, ig
integer :: notconv, nhpsi, ig
real(dp):: ppcg_iter
integer :: max_ppcg_iter = 100 !default of QE
real(dp), allocatable :: h_diag(:,:) !in case of CG - the preconditioner
logical :: overlap = .false. , lrot =.false. , lscf = .true. ! lscf is true for scf calc
......@@ -81,6 +82,7 @@ call set_mpi_comm_4_solvers( world_comm, intra_bgrp_comm, inter_bgrp_comm )
npwx, npw, nbnd, nbnd, npol, evc, evc, eig )
endif
#endif
nhpsi = nbnd
write(stdout,*) ' then ppcg diagonalization '
!support only serial for now
h_diag = 1.D0
......@@ -91,7 +93,7 @@ call set_mpi_comm_4_solvers( world_comm, intra_bgrp_comm, inter_bgrp_comm )
CALL ppcg_k( cb_h_psi, cb_s_psi, overlap, h_diag, &
npwx, npw, nbnd, npol, evc, eig, btype, &
ethr, max_ppcg_iter, notconv, ppcg_iter, sbsize , rrstep, iter )
nhpsi = nhpsi + ppcg_iter*nbnd
!--------------------------------------------------------------------------------------------------------------!
call stop_clock('ppcg')
......@@ -99,7 +101,7 @@ call set_mpi_comm_4_solvers( world_comm, intra_bgrp_comm, inter_bgrp_comm )
if (energy_shift .and. current_k==1) ref=eig(4*ncell**3)
call write_bands(eig,ref)
write (6,*) 'ppcg_iter, notconv, ethr ',ppcg_iter,notconv, ethr
write (6,*) 'ppcg_iter, nhpsi, notconv, ethr ', int(ppcg_iter), nhpsi, notconv, ethr
deallocate( eig )
end do
......@@ -118,6 +120,13 @@ call set_mpi_comm_4_solvers( world_comm, intra_bgrp_comm, inter_bgrp_comm )
call print_clock('ppcg:RR')
call print_clock('ppcg:ZTRSM')
call print_clock('ppcg:lock')
!
write (6,*)
write (6,*) ' general FFT routines'
call print_clock('fftw')
call print_clock('ffts')
#if defined(__MPI)
call unset_mpi_comm_4_solvers
call mp_global_end( )
......
......@@ -11,7 +11,7 @@
!----------------------------------------------------------------------------
SUBROUTINE cegterg( h_psi, s_psi, uspp, g_psi, &
npw, npwx, nvec, nvecx, npol, evc, ethr, &
e, btype, notcnv, lrot, dav_iter )
e, btype, notcnv, lrot, dav_iter, nhpsi)
!----------------------------------------------------------------------------
!
! ... iterative solution of the eigenvalue problem:
......@@ -56,6 +56,8 @@ SUBROUTINE cegterg( h_psi, s_psi, uspp, g_psi, &
INTEGER, INTENT(OUT) :: dav_iter, notcnv
! integer number of iterations performed
! number of unconverged roots
INTEGER, INTENT(OUT), OPTIONAL :: nhpsi
! total number of indivitual hpsi
!
! ... LOCAL variables
!
......@@ -101,6 +103,7 @@ SUBROUTINE cegterg( h_psi, s_psi, uspp, g_psi, &
! calculates (diag(h)-e)^-1 * psi, diagonal approx. to (h-e)^-1*psi
! the first nvec columns contain the trial eigenvectors
!
nhpsi = 0
CALL start_clock( 'cegterg' ); !write(*,*) 'start cegterg' ; FLUSH(6)
!
IF ( nvec > nvecx / 2 ) CALL errore( 'cegterg', 'nvecx is too small', 1 )
......@@ -162,7 +165,7 @@ SUBROUTINE cegterg( h_psi, s_psi, uspp, g_psi, &
!
! ... hpsi contains h times the basis vectors
!
CALL h_psi( npwx, npw, nvec, psi, hpsi )
CALL h_psi( npwx, npw, nvec, psi, hpsi ) ; if (present(nhpsi)) nhpsi = nhpsi + nvec
!
! ... spsi contains s times the basis vectors
!
......@@ -256,7 +259,7 @@ SUBROUTINE cegterg( h_psi, s_psi, uspp, g_psi, &
!
iterate: DO kter = 1, maxter
!
dav_iter = kter
dav_iter = kter ; !write(*,*) kter, notcnv, conv
!
CALL start_clock( 'cegterg:update' )
!
......@@ -368,7 +371,7 @@ SUBROUTINE cegterg( h_psi, s_psi, uspp, g_psi, &
!
! ... here compute the hpsi and spsi of the new functions
!
CALL h_psi( npwx, npw, notcnv, psi(1,1,nb1), hpsi(1,1,nb1) )
CALL h_psi( npwx, npw, notcnv, psi(1,1,nb1), hpsi(1,1,nb1) ) ; if (present(nhpsi)) nhpsi = nhpsi + notcnv
!
IF ( uspp ) CALL s_psi( npwx, npw, notcnv, psi(1,1,nb1), spsi(1,1,nb1) )
!
......@@ -568,7 +571,7 @@ END SUBROUTINE cegterg
!----------------------------------------------------------------------------
SUBROUTINE pcegterg(h_psi, s_psi, uspp, g_psi, &
npw, npwx, nvec, nvecx, npol, evc, ethr, &
e, btype, notcnv, lrot, dav_iter )
e, btype, notcnv, lrot, dav_iter, nhpsi )
!----------------------------------------------------------------------------
!
! ... iterative solution of the eigenvalue problem:
......@@ -615,6 +618,8 @@ SUBROUTINE pcegterg(h_psi, s_psi, uspp, g_psi, &
INTEGER, INTENT(OUT) :: dav_iter, notcnv
! integer number of iterations performed
! number of unconverged roots
INTEGER, INTENT(OUT), OPTIONAL :: nhpsi
! total number of indivitual hpsi
!
! ... LOCAL variables
!
......@@ -665,7 +670,7 @@ SUBROUTINE pcegterg(h_psi, s_psi, uspp, g_psi, &
! calculates (diag(h)-e)^-1 * psi, diagonal approx. to (h-e)^-1*psi
! the first nvec columns contain the trial eigenvectors
!
!
nhpsi = 0
CALL start_clock( 'cegterg' )
!
IF ( nvec > nvecx / 2 ) CALL errore( 'pcegterg', 'nvecx is too small', 1 )
......@@ -776,7 +781,7 @@ SUBROUTINE pcegterg(h_psi, s_psi, uspp, g_psi, &
!
! ... hpsi contains h times the basis vectors
!
CALL h_psi( npwx, npw, nvec, psi, hpsi )
CALL h_psi( npwx, npw, nvec, psi, hpsi ) ; if (present(nhpsi)) nhpsi = nhpsi + nvec
!
IF ( uspp ) CALL s_psi( npwx, npw, nvec, psi, spsi )
!
......@@ -832,7 +837,7 @@ SUBROUTINE pcegterg(h_psi, s_psi, uspp, g_psi, &
!
iterate: DO kter = 1, maxter
!
dav_iter = kter
dav_iter = kter ; !write(*,*) kter, notcnv, conv
!
CALL start_clock( 'cegterg:update' )
!
......@@ -888,7 +893,7 @@ SUBROUTINE pcegterg(h_psi, s_psi, uspp, g_psi, &
!
! ... here compute the hpsi and spsi of the new functions
!
CALL h_psi( npwx, npw, notcnv, psi(1,1,nb1), hpsi(1,1,nb1) )
CALL h_psi( npwx, npw, notcnv, psi(1,1,nb1), hpsi(1,1,nb1) ) ; if(present(nhpsi)) nhpsi = nhpsi + notcnv
!
IF ( uspp ) CALL s_psi( npwx, npw, notcnv, psi(1,1,nb1), spsi(1,1,nb1) )
!
......
......@@ -224,7 +224,7 @@ SUBROUTINE regterg( h_psi, s_psi, uspp, g_psi, &
!
iterate: DO kter = 1, maxter
!
dav_iter = kter
dav_iter = kter ; !write(*,*) kter, notcnv, conv
!
CALL start_clock( 'regterg:update' )
!
......@@ -767,7 +767,7 @@ SUBROUTINE pregterg(h_psi, s_psi, uspp, g_psi, &
!
iterate: DO kter = 1, maxter
!
dav_iter = kter
dav_iter = kter ; !write(*,*) kter, notcnv, conv
!
CALL start_clock( 'regterg:update' )
!
......
......@@ -98,7 +98,7 @@ SUBROUTINE paro_gamma( h_psi, s_psi, hs_1psi, g_1psi, overlap, &
ParO_loop : &
DO itry = 1,paro_ntr
write (6,*) ' paro_itry =', itry, ethr
!write (6,*) ' paro_itry =', itry, ethr
nactive = nbnd - (nconv+1)/2 ! number of correction vectors to be computed (<nbnd)
notconv = nbnd - nconv ! number of needed roots
......@@ -106,8 +106,8 @@ SUBROUTINE paro_gamma( h_psi, s_psi, hs_1psi, g_1psi, overlap, &
nbase = nconv + nactive ! number of orbitals the correction should be orthogonal to (<2*nbnd)
ndiag = nbase + nactive ! dimension of the matrix to be diagonalized at this iteration (<2*nbnd)
write (*,*) itry, conv
write (6,*) ' nbnd, nconv, notconv, nextra, nactive, nbase, ndiag =', nbnd, nconv, notconv, nextra, nactive, nbase, ndiag
!write (*,*) itry, notconv, conv
!write (6,*) ' nbnd, nconv, notconv, nextra, nactive, nbase, ndiag =', nbnd, nconv, notconv, nextra, nactive, nbase, ndiag
call s_psi (npwx,npw,nbnd,psi2,evc) ! computes S*psi needed to ortogonalize to nbase
lbnd = nbase
......
......@@ -85,9 +85,7 @@ SUBROUTINE paro_gamma_new( h_psi, s_psi, hs_1psi, g_1psi, overlap, &
iter = 0
paro_ntr = 20
!
write (6,*) ' enter paro diag'
! write (6,*) ' paro_flag = ', paro_flag
! if (paro_flag /= 1) WRITE(stdout,*) 'wrong setting of paro_flag!! '
!write (6,*) ' enter paro diag'
ALLOCATE ( psi(npwx,2*nbnd), hpsi(npwx,2*nbnd), spsi(npwx,2*nbnd), ew(2*nbnd), conv(nbnd) )
......@@ -112,7 +110,7 @@ SUBROUTINE paro_gamma_new( h_psi, s_psi, hs_1psi, g_1psi, overlap, &
ParO_loop : &
DO itry = 1,paro_ntr
write (6,*) ' paro_itry =', itry, ethr
!write (6,*) ' paro_itry =', itry, ethr
nactive = nbnd - (nconv+1)/2 ! number of correction vectors to be computed (<nbnd)
notconv = nbnd - nconv ! number of needed roots
......@@ -120,8 +118,8 @@ SUBROUTINE paro_gamma_new( h_psi, s_psi, hs_1psi, g_1psi, overlap, &
nbase = nconv + nactive ! number of orbitals the correction should be orthogonal to (<2*nbnd)
ndiag = nbase + nactive ! dimension of the matrix to be diagonalized at this iteration (<2*nbnd)
write (*,*) itry, conv
write (6,*) ' nbnd, nconv, notconv, nextra, nactive, nbase, ndiag =', nbnd, nconv, notconv, nextra, nactive, nbase, ndiag
!write (*,*) itry, notconv, conv
!write (6,*) ' nbnd, nconv, notconv, nextra, nactive, nbase, ndiag =', nbnd, nconv, notconv, nextra, nactive, nbase, ndiag
lbnd = nbase
DO ibnd = 1, nbnd ! pack unconverged roots
......
......@@ -105,6 +105,7 @@ SUBROUTINE paro_k( h_psi, s_psi, hs_1psi, g_1psi, overlap, &
nbase = nconv + nactive ! number of orbitals the correction should be orthogonal to (<2*nbnd)
ndiag = nbase + nactive ! dimension of the matrix to be diagonalized at this iteration (<2*nbnd)
!write (*,*) itry, notconv, conv
!write (6,*) ' nbnd, nconv, notconv, nextra, nactive, nbase, ndiag =', nbnd, nconv, notconv, nextra, nactive, nbase, ndiag
call s_psi (npwx,npw,nbnd,psi2,evc) ! computes S*psi needed to ortogonalize to nbase
......@@ -122,7 +123,6 @@ SUBROUTINE paro_k( h_psi, s_psi, hs_1psi, g_1psi, overlap, &
eig(lbnd-nbase) = eig(lbnd-nbase-1)
END DO
!write (*,*) itry, conv
!write (6,*) ' check nactive = ', lbnd-nbase, nactive
if (lbnd .ne. nbase+nactive ) stop ' nactive check FAILED '
......
......@@ -85,7 +85,7 @@ SUBROUTINE paro_k_new( h_psi, s_psi, hs_1psi, g_1psi, overlap, &
iter = 0
paro_ntr = 20
!
write (6,*) ' enter paro diag'
!write (6,*) ' enter paro diag'
ALLOCATE ( psi(npwx*npol,2*nbnd), hpsi(npwx*npol, 2*nbnd), spsi(npwx*npol,2*nbnd), ew(2*nbnd), conv(nbnd) )
......@@ -111,7 +111,7 @@ SUBROUTINE paro_k_new( h_psi, s_psi, hs_1psi, g_1psi, overlap, &
ParO_loop : &
DO itry = 1,paro_ntr
write (6,*) ' paro_itry =', itry, ethr
!write (6,*) ' paro_itry =', itry, ethr
nactive = nbnd - (nconv+1)/2 ! number of correction vectors to be computed (<nbnd)
notconv = nbnd - nconv ! number of needed roots
......@@ -119,8 +119,8 @@ SUBROUTINE paro_k_new( h_psi, s_psi, hs_1psi, g_1psi, overlap, &
nbase = nconv + nactive ! number of orbitals the correction should be orthogonal to (<2*nbnd)
ndiag = nbase + nactive ! dimension of the matrix to be diagonalized at this iteration (<2*nbnd)
write (*,*) itry, conv
write (6,*) ' nbnd, nconv, notconv, nextra, nactive, nbase, ndiag =', nbnd, nconv, notconv, nextra, nactive, nbase, ndiag
!write (*,*) itry, notconv, conv
!write (6,*) ' nbnd, nconv, notconv, nextra, nactive, nbase, ndiag =', nbnd, nconv, notconv, nextra, nactive, nbase, ndiag
lbnd = nbase
DO ibnd = 1, nbnd ! pack unconverged roots
......
Supports Markdown
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment