@brief Test if an excitation generator generates all and only expected states with the correct pgen.
@details
The pgen_diagnostic is given by
\f[\sum_i \frac{1}{p_i N}\f]
and should be roughly one.
All problematic states with respect to that diagnostic are
printed in the end with their pgen_diagnostic,
excitation type (ic), matrix element, and pgen.
(Does not count determinants with nJ(1) == 0!)
that can be reached from the reference state.
for a given reference an excitation. If a state is never generated,
the pgen cannot be taken from the excitation generator.
Adds an additional column to the output table.
Return true, if an excitation exc from determinant det_I
and a given pgen_diagnostic (sum 1/pgen) is considered
to be problematic.
If it returns true, an entry in the final table is printed.
If there is any problematic excitation, the out parameter
successful will become .false..
By default all states with a pgen_diagnostic that deviate
with more than 5\,\% from 100\,\% and have nonzereo matrix element
are printed.
| Type | Intent | Optional | Attributes | Name | ||
|---|---|---|---|---|---|---|
| procedure(generate_excitation_t) | :: | excit_gen | ||||
| character(len=*), | intent(in) | :: | excit_gen_name | |||
| integer, | intent(in), | optional | :: | opt_nI(nel) | ||
| integer, | intent(in), | optional | :: | opt_n_dets | ||
| procedure(generate_all_excits_t), | optional | :: | gen_all_excits | |||
| procedure(calc_pgen_t), | optional | :: | calc_pgen | |||
| procedure(problem_filter_t), | optional | :: | problem_filter | |||
| integer, | intent(in), | optional | :: | i_unit | ||
| logical, | intent(out), | optional | :: | successful |
subroutine run_excit_gen_tester_function(& excit_gen, excit_gen_name, opt_nI, opt_n_dets, & gen_all_excits, calc_pgen, problem_filter, i_unit, successful) procedure(generate_excitation_t) :: excit_gen character(*), intent(in) :: excit_gen_name integer, intent(in), optional :: opt_nI(nel), opt_n_dets procedure(generate_all_excits_t), optional :: gen_all_excits procedure(calc_pgen_t), optional :: calc_pgen procedure(problem_filter_t), optional :: problem_filter integer, intent(in), optional :: i_unit logical, intent(out), optional :: successful character(*), parameter :: this_routine = "run_excit_gen_tester" integer :: i, nI(nel), n_dets_target integer :: i_unit_ integer, parameter :: default_n_dets_target = 100000 integer(n_int) :: ilut(0:niftot), tgt_ilut(0:niftot) integer :: nJ(nel), n_excits, ex(2, maxExcit), ic, ex_flag, i_unused = 0 type(excit_gen_store_type) :: store logical :: tPar, found_all real(dp) :: pgen, contrib real(dp), allocatable :: pgen_list(:) HElement_t(dp) :: hel integer(n_int), allocatable :: det_list(:, :) real(dp), allocatable :: contrib_list(:) logical, allocatable :: generated_list(:) integer :: n_generated, pos, n_iter procedure(problem_filter_t), pointer :: problem_filter_ ! and also nbasis and stuff.. ASSERT(nbasis > 0) ASSERT(nel <= nbasis) if (present(i_unit)) then i_unit_ = i_unit else i_unit_ = stdout end if if (present(problem_filter)) then problem_filter_ => problem_filter else problem_filter_ => default_predicate end if call init_excit_gen_store(store) ! use some default values if not provided: ! nel must be set! if (present(opt_nI)) then nI = opt_nI else ! use HF-det as default nI = [(i, i=1, nel)] end if if (present(opt_n_dets)) then n_dets_target = opt_n_dets else n_dets_target = default_n_dets_target end if ! Calculate all possible excitations. ! For special systems (hubbard, UEG, GAS, etc.) the calling site ! has to support a function. if (present(gen_all_excits)) then call gen_all_excits(nI, n_excits, det_list) else call gen_all_excits_default(nI, n_excits, det_list) end if allocate(generated_list(n_excits), source=.false.) allocate(contrib_list(n_excits), source=0.0_dp) allocate(pgen_list(n_excits), source=0.0_dp) write(i_unit_, *) "---------------------------------" write(i_unit_, *) "testing: ", excit_gen_name write(i_unit_, *) "for ", size(det_list, 2), " configurations" write(i_unit_, *) " and ", n_dets_target, " determinants to generate" call EncodeBitDet(nI, ilut) write(i_unit_, *) "for starting determinant: ", nI write(i_unit_, *) ! linebreak write(i_unit_, '(A)') 'Progressbar' ! call this below now for the number of specified determinants ! (also use excitations of the first inputted, to be really ! consistent) block integer(int64) :: L type(timer) :: loop_timer loop_timer%timer_name = 'loop '//excit_gen_name L = n_dets_target .div. 100 n_generated = 0 n_iter = 0 contrib = 0.0_dp do while (n_generated < n_dets_target) n_iter = n_iter + 1 call set_timer(loop_timer) call excit_gen(nI, ilut, nJ, tgt_ilut, ex_flag, ic, ex, tpar, pgen, & hel, store) call halt_timer(loop_timer) if (nJ(1) == 0) cycle call EncodeBitDet(nJ, tgt_ilut) pos = binary_search_ilut(det_list, tgt_ilut, nifd + 1) if (pos < 0) then write(i_unit_, *) "nJ: ", nJ write(i_unit_, *) "ilutJ:", tgt_ilut call stop_all(this_routine, 'Unexpected determinant generated') else if (mod(n_generated, L) == 0_int64) then write(i_unit_, '(A)', advance='no') '#' flush (i_unit_) end if generated_list(pos) = .true. n_generated = n_generated + 1 contrib = contrib + 1.0_dp / pgen contrib_list(pos) = contrib_list(pos) + 1.0_dp / pgen pgen_list(pos) = pgen end if end do write(i_unit_, *) ! linebreak write(i_unit_, *) ! linebreak write(i_unit_, *) "==================================================================" write(i_unit_, *) n_generated, " dets generated in ", n_iter, " iterations " write(i_unit_, *) n_generated, (n_iter - n_generated), n_iter, 'valid invalid and total excitations' write(i_unit_, *) 100.0_dp * real(n_iter - n_generated, dp) / real(n_iter, dp), "% abortion rate" write(i_unit_, *) "Averaged contribution: ", contrib / real(n_excits * n_dets_target, dp) write(i_unit_, *) ! linebreak write(i_unit_, *) 'Elapsed Time in seconds:', get_total_time(loop_timer) write(i_unit_, *) 'Elapsed Time in micro seconds per valid excitation:', get_total_time(loop_timer) / dble(n_dets_target) * 10.**6 write(i_unit_, *) 'Elapsed Time in micro seconds per excitation:', get_total_time(loop_timer) / dble(n_iter) * 10.**6 write(i_unit_, *) "==================================================================" write(i_unit_, *) ! linebreak write(i_unit_, *) ! linebreak end block block real(dp) :: pgen_diagnostic integer :: nJ(nEl), exc(2, maxExcit) logical :: tParity write(i_unit_, *) "==================================" write(i_unit_, *) "Problematic contribution List: " write(i_unit_, *) "==================================" write(i_unit_, '("| Determinant | Sum{1 / pgen} / n_iter | ic | <psi_I H psi_J > | pgen |")', advance='no') if (present(calc_pgen)) write(i_unit_, '(" calc_pgen |")', advance='no') write(i_unit_, *) ! linebreak successful = .true. do i = 1, n_excits block logical :: pgens_differ pgen_diagnostic = contrib_list(i) / real(n_iter, dp) ic = findbitexcitlevel(ilut, det_list(:, i)) call decode_bit_det(nJ, det_list(:, i)) call get_excitation(nI, nJ, ic, exc, tParity) if (present(calc_pgen)) then pgens_differ = .not. (calc_pgen(nI, ilut, exc, ic, store%ClassCountOcc, store%ClassCountUnocc) .isclose. pgen_list(i)) else pgens_differ = .false. end if if (problem_filter_(nI, exc, ic, pgen_diagnostic) .or. pgens_differ) then successful = .false. call write_det(i_unit_, nJ, .false.) write(i_unit_, '("|"F10.5"|"I2"|"F10.5"|"F15.10"|")', advance='no') & pgen_diagnostic, ic, get_helement(nI, nJ), pgen_list(i) if (present(calc_pgen)) then write(i_unit_, '(F15.10"|")', advance='no') & calc_pgen(nI, ilut, exc, ic, store%ClassCountOcc, store%ClassCountUnocc) end if write(i_unit_, *) end if end block end do end block call clean_excit_gen_store(store) contains logical function default_predicate(nI, exc, ic, pgen_diagnostic) integer, intent(in) :: nI(nEl), exc(2, maxExcit), ic real(dp), intent(in) :: pgen_diagnostic default_predicate = & (abs(1._dp - pgen_diagnostic) >= 0.05_dp) & .and. .not. near_zero(dyn_sltcnd_excit_old(nI, ic, exc, .true.)) end function end subroutine run_excit_gen_tester_function