subroutine init_spectral_lanczos()
use bit_rep_data, only: nifd, extract_sign
use DetBitOps, only: DetBitEq, EncodeBitDet, ilut_lt, ilut_gt
use gndts_mod, only: gndts
use load_balance_calcnodes, only: DetermineDetNode
use sort_mod, only: sort
use SystemData, only: nbasis, nel, BRR, nBasisMax, G1, tSpn, LMS, tParity, SymRestrict
integer :: ndets_this_proc, ndets_tot
integer(MPIArg) :: mpi_temp
integer, allocatable :: nI_list(:, :)
integer(n_int) :: ilut(0:NIfTot)
character(len=*), parameter :: t_r = 'init_spectral_lanczos'
integer :: i, j, ndets, proc, ierr, temp(1, 1), hf_ind
real(dp) :: real_sign(lenof_sign)
real(dp) :: norm_pert
write(stdout, '(/,1x,a39,/)') "Beginning spectral Lanczos calculation."
call neci_flush(stdout)
allocate(ndets_sl(0:nProcessors - 1))
allocate(disps_sl(0:nProcessors - 1))
write(stdout, '(1x,a56)', advance='yes') "Enumerating and storing all determinants in the space..."
call neci_flush(stdout)
! Determine the total number of determinants.
call gndts(nel, nbasis, BRR, nBasisMax, temp, .true., G1, tSpn, lms, tParity, SymRestrict, ndets, hf_ind)
allocate(nI_list(nel, ndets))
! Now generate them again and store them this time.
call gndts(nel, nbasis, BRR, nBasisMax, nI_list, .false., G1, tSpn, lms, tParity, SymRestrict, ndets, hf_ind)
! Count the number of determinants belonging to this processor.
ndets_this_proc = 0
do i = 1, ndets
proc = DetermineDetNode(nel, nI_list(:, i), 0)
if (proc == iProcIndex) ndets_this_proc = ndets_this_proc + 1
end do
allocate(sl_ilut_list(0:NIfTot, ndets_this_proc))
sl_ilut_list = 0_n_int
j = 0
do i = 1, ndets
proc = DetermineDetNode(nel, nI_list(:, i), 0)
if (proc == iProcIndex) then
call EncodeBitDet(nI_list(:, i), ilut)
j = j + 1
sl_ilut_list(0:nifd, j) = ilut(0:nifd)
end if
end do
call sort(sl_ilut_list, ilut_lt, ilut_gt)
write(stdout, '(1x,a9)') "Complete."
call neci_flush(stdout)
mpi_temp = int(ndets_this_proc, MPIArg)
call MPIAllGather(mpi_temp, ndets_sl, ierr)
disps_sl(0) = 0
do i = 1, nProcessors - 1
disps_sl(i) = disps_sl(i - 1) + ndets_sl(i - 1)
end do
ndets_tot = int(sum(ndets_sl))
write(stdout, '(1x,a44)', advance='no') "Allocating arrays to hold Lanczos vectors..."
call neci_flush(stdout)
allocate(sl_vecs(ndets_this_proc, n_lanc_vecs_sl))
sl_vecs = 0.0_dp
allocate(full_vec_sl(ndets_tot))
write(stdout, '(1x,a9)') "Complete."
call neci_flush(stdout)
allocate(pert_ground_left(ndets_this_proc), stat=ierr)
if (ierr /= 0) then
write(stdout, '(1x,a11,1x,i5)') "Error code:", ierr
call stop_all(t_r, "Error allocating array to hold left perturbed ground state components.")
end if
allocate(sl_hamil(n_lanc_vecs_sl, n_lanc_vecs_sl))
allocate(sl_h_eigv(n_lanc_vecs_sl))
allocate(sl_overlaps(n_lanc_vecs_sl))
allocate(trans_amps_left(n_lanc_vecs_sl))
allocate(trans_amps_right(n_lanc_vecs_sl))
sl_hamil = 0.0_dp
sl_h_eigv = 0.0_dp
sl_overlaps = 0.0_dp
trans_amps_left = 0.0_dp
trans_amps_right = 0.0_dp
write(stdout, '(1x,a48)') "Allocating and calculating Hamiltonian matrix..."
call neci_flush(stdout)
call calculate_sparse_ham_par(ndets_sl, sl_ilut_list, .true.)
write(stdout, '(1x,a48,/)') "Hamiltonian allocation and calculation complete."
call neci_flush(stdout)
call return_perturbed_ground_spec(left_perturb_spectral, sl_ilut_list, pert_ground_left, left_pert_norm)
call return_perturbed_ground_spec(right_perturb_spectral, sl_ilut_list, sl_vecs(:, 1), right_pert_norm)
! Normalise the perturbed wave function, currently stored in sl_vecs(:,1),
! so that it can be used as the first lanczos vector.
sl_vecs(:, 1) = sl_vecs(:, 1) * sqrt(pops_norm) / right_pert_norm
end subroutine init_spectral_lanczos