subroutine create_direct_ci_arrays(ras, classes, ras_strings, ras_iluts, ras_excit)
! This routine should be used before perform_multiplication is called. It will
! create some arrays which will make this multiplication quicker.
type(ras_parameters), intent(in) :: ras
type(ras_class_data), intent(in) :: classes(ras%num_classes)
integer, intent(out) :: ras_strings(-1:tot_nelec, ras%num_strings)
integer(n_int), intent(out) :: ras_iluts(0:NIfD, ras%num_strings)
type(direct_ci_excit), intent(out) :: ras_excit(ras%num_strings)
integer(n_int) :: ilut_i(0:NIfD)
integer :: class_i, ind, new_ind, counter
integer :: nras1, nras3, par, k
integer :: ex(2)
integer :: string_i(tot_nelec)
logical :: none_left, tgen
type(simple_excit_store), target :: gen_store_1
ilut_i = 0_n_int
! Loop over all classes.
do class_i = 1, ras%num_classes
nras1 = classes(class_i)%nelec_1
nras3 = classes(class_i)%nelec_3
call generate_first_full_string(string_i, ras, classes(class_i))
! Loop over all strings.
do
ind = classes(class_i)%address_map(get_address(classes(class_i), string_i))
do k = 1, class_i - 1
ind = ind + classes(k)%class_size
end do
! Store the string, along with the symmetry, RAS class and ilut.
ras_strings(1:tot_nelec, ind) = string_i
ras_strings(-1, ind) = get_abelian_sym(string_i)
ras_strings(0, ind) = class_i
call encode_string(string_i, ilut_i)
ras_iluts(:, ind) = ilut_i
! Now count the number of excitations.
counter = 0
! Initialise the excitation generator.
tgen = .true.
do
call gen_next_single_ex(string_i, ilut_i, nras1, nras3, new_ind, par, &
ex, ras, classes, gen_store_1, tgen, .true.)
if (tgen) exit
counter = counter + 1
end do
! Store the number of excitations.
ras_excit(ind)%nexcit = counter
! Now we know how many excitations there are, allocate the excitation array and store them.
allocate(ras_excit(ind)%excit_ind(ras_excit(ind)%nexcit))
allocate(ras_excit(ind)%par(ras_excit(ind)%nexcit))
allocate(ras_excit(ind)%orbs(2, ras_excit(ind)%nexcit))
counter = 0
tgen = .true.
do
call gen_next_single_ex(string_i, ilut_i, nras1, nras3, new_ind, par, &
ex, ras, classes, gen_store_1, tgen, .false.)
if (tgen) exit
counter = counter + 1
! Store the index...
ras_excit(ind)%excit_ind(counter) = new_ind
! ...and the parity of the excitation...
ras_excit(ind)%par(counter) = par
! ...and the two orbitals involved in the excitation.
ras_excit(ind)%orbs(:, counter) = ex
end do
call generate_next_string(string_i, ras, classes(class_i), none_left)
! If no strings left in this class, then go to the next class.
if (none_left) exit
end do
end do
end subroutine create_direct_ci_arrays