diff --git a/src/Dets/H_apply_template.f b/src/Dets/H_apply_template.f index 5e310788..3b6f1e35 100644 --- a/src/Dets/H_apply_template.f +++ b/src/Dets/H_apply_template.f @@ -346,16 +346,16 @@ subroutine $subroutine($params_main) integer :: i_generator, k do i_generator=1,N_det_generators - call $subroutine_monoexc(psi_generators(1,1,i_generator), & - generators_bitmask(1,1,s_hole ,i_bitmask_gen), & - generators_bitmask(1,1,s_part ,i_bitmask_gen), & - i_generator $params_post) call $subroutine_diexc(psi_generators(1,1,i_generator), & generators_bitmask(1,1,d_hole1,i_bitmask_gen), & generators_bitmask(1,1,d_part1,i_bitmask_gen), & generators_bitmask(1,1,d_hole2,i_bitmask_gen), & generators_bitmask(1,1,d_part2,i_bitmask_gen), & i_generator $params_post) + call $subroutine_monoexc(psi_generators(1,1,i_generator), & + generators_bitmask(1,1,s_hole ,i_bitmask_gen), & + generators_bitmask(1,1,s_part ,i_bitmask_gen), & + i_generator $params_post) if (abort_here) then exit endif diff --git a/src/Dets/davidson.irp.f b/src/Dets/davidson.irp.f index 26fa1195..2bc500d2 100644 --- a/src/Dets/davidson.irp.f +++ b/src/Dets/davidson.irp.f @@ -374,8 +374,8 @@ end BEGIN_DOC ! Can be : [ energy | residual | both | wall_time | cpu_time | iterations ] END_DOC - davidson_criterion = 'both' - davidson_threshold = 1.d-8 + davidson_criterion = 'residual' + davidson_threshold = 1.d-6 END_PROVIDER subroutine davidson_converged(energy,residual,wall,iterations,cpu,N_st,converged) diff --git a/src/Dets/determinants.irp.f b/src/Dets/determinants.irp.f index d16f8ce6..b33bf664 100644 --- a/src/Dets/determinants.irp.f +++ b/src/Dets/determinants.irp.f @@ -80,3 +80,54 @@ BEGIN_PROVIDER [ integer, N_det_reference ] N_det_reference = N_det ASSERT (N_det_reference > 0) END_PROVIDER + +BEGIN_PROVIDER [ double precision, psi_average_norm_contrib, (N_det) ] + implicit none + BEGIN_DOC + ! Contribution of determinants to the state-averaged density + END_DOC + integer :: i,j,k + double precision :: f + f = 1.d0/dble(N_states) + do i=1,N_det + psi_average_norm_contrib(i) = psi_coef(i,1)*psi_coef(i,1)*f + enddo + do k=2,N_states + do i=1,N_det + psi_average_norm_contrib(i) = psi_average_norm_contrib(i) + & + psi_coef(i,k)*psi_coef(i,k)*f + enddo + enddo +END_PROVIDER + + BEGIN_PROVIDER [ integer(bit_kind), psi_det_sorted, (N_int,2,N_det) ] +&BEGIN_PROVIDER [ double precision, psi_coef_sorted, (N_det,N_states) ] +&BEGIN_PROVIDER [ double precision, psi_average_norm_contrib_sorted, (N_det) ] + implicit none + BEGIN_DOC + ! Wave function sorted by determinants (state-averaged) + END_DOC + integer :: i,j,k + integer, allocatable :: iorder(:) + allocate ( iorder(N_det) ) + do i=1,N_det + psi_average_norm_contrib_sorted(i) = -psi_average_norm_contrib(i) + iorder(i) = i + enddo + call dsort(psi_average_norm_contrib_sorted,iorder,N_det) + !DIR$ IVDEP + do i=1,N_det + do j=1,N_int + psi_det_sorted(j,1,i) = psi_det(j,1,iorder(i)) + psi_det_sorted(j,2,i) = psi_det(j,2,iorder(i)) + enddo + do k=1,N_states + psi_coef_sorted(i,k) = psi_coef(iorder(i),k) + enddo + psi_average_norm_contrib_sorted(i) = -psi_average_norm_contrib_sorted(i) + enddo + + deallocate(iorder) + +END_PROVIDER + diff --git a/src/Dets/diagonalize_CI.irp.f b/src/Dets/diagonalize_CI.irp.f index 1f1ad1cf..0dc39448 100644 --- a/src/Dets/diagonalize_CI.irp.f +++ b/src/Dets/diagonalize_CI.irp.f @@ -22,8 +22,12 @@ BEGIN_PROVIDER [ double precision, CI_energy, (N_states) ] END_DOC integer :: j + character*(8) :: st + call write_time(output_Dets) do j=1,N_states CI_energy(j) = CI_electronic_energy(j) + nuclear_repulsion + write(st,'(I4)') j + call write_double(output_Dets,CI_energy(j),'Energy of state '//trim(st)) enddo END_PROVIDER diff --git a/src/Full_CI/README.rst b/src/Full_CI/README.rst index f63f5a06..5f0ae69c 100644 --- a/src/Full_CI/README.rst +++ b/src/Full_CI/README.rst @@ -10,6 +10,9 @@ Documentation .. Do not edit this section. It was auto-generated from the .. NEEDED_MODULES file. +`cisd `_ + Undocumented + Needed Modules diff --git a/src/Full_CI/full_ci.irp.f b/src/Full_CI/full_ci.irp.f index 796dbe6a..7d563420 100644 --- a/src/Full_CI/full_ci.irp.f +++ b/src/Full_CI/full_ci.irp.f @@ -3,18 +3,15 @@ program cisd integer :: i,k - double precision, allocatable :: pt2(:), norm_pert(:) - double precision :: H_pert_diag, E_old + double precision, allocatable :: pt2(:), norm_pert(:), H_pert_diag(:) integer :: N_st, degree character*(64) :: perturbation N_st = N_states - allocate (pt2(N_st), norm_pert(N_st)) + allocate (pt2(N_st), norm_pert(N_st),H_pert_diag(N_st)) pt2 = 1.d0 diag_algorithm = "Lapack" do while (maxval(abs(pt2(1:N_st))) > 1.d-6) - print *, '-----' - E_old = CI_energy(1) call H_apply_FCI(pt2, norm_pert, H_pert_diag, N_st) call diagonalize_CI print *, 'N_det = ', N_det @@ -22,6 +19,7 @@ program cisd print *, 'PT2 = ', pt2 print *, 'E = ', CI_energy print *, 'E+PT2 = ', CI_energy+pt2 + print *, '-----' if (abort_all) then exit endif diff --git a/src/Perturbation/selection.irp.f b/src/Perturbation/selection.irp.f index 25f6cf37..400d6976 100644 --- a/src/Perturbation/selection.irp.f +++ b/src/Perturbation/selection.irp.f @@ -31,13 +31,13 @@ subroutine fill_H_apply_buffer_selection(n_selected,det_buffer,e_2_pert_buffer,c l=H_apply_buffer(iproc)%N_det do i=1,n_selected - s = 0.d0 + is_selected = .False. do j=1,N_st - s -= e_2_pert_buffer(j,i) + s = dabs(e_2_pert_buffer(j,i)) + is_selected = s > selection_criterion*selection_criterion_factor .or. is_selected enddo ASSERT (s>=-1.d-8) - is_selected = s > selection_criterion * selection_criterion_factor if (is_selected) then l = l+1 @@ -71,7 +71,7 @@ end BEGIN_DOC ! Threshold to select determinants. Set by selection routines. END_DOC - selection_criterion = .1d0 + selection_criterion = 1.d0 selection_criterion_factor = 0.01d0 selection_criterion_min = selection_criterion @@ -80,35 +80,37 @@ END_PROVIDER subroutine remove_small_contributions implicit none BEGIN_DOC -! Remove determinants with small contributions +! Remove determinants with small contributions. N_states is assumed to be +! provided. END_DOC integer :: i,j,k, N_removed logical keep + double precision :: i_H_psi_array(N_states) + k = 0 N_removed = 0 - do i=N_det,1,-1 + do i=N_det, 50 + call i_H_psi(psi_det_sorted(1,1,i),psi_det_sorted,psi_coef_sorted,N_int,N_det,psi_det_size,N_states,i_H_psi_array) keep = .False. do j=1,N_states - keep = keep .or. (dabs(psi_coef(i,j)) > selection_criterion_min) + keep = keep .or. (dabs(psi_coef_sorted(i,j)*i_H_psi_array(j)) > selection_criterion_min) enddo - if (.not.keep) then - do k=i+1,N_det - do j=1,N_int - psi_det(j,1,k-1) = psi_det(j,1,k) - psi_det(j,2,k-1) = psi_det(j,2,k) - enddo + if (keep) then + k += 1 + do j=1,N_int + psi_det(j,1,k) = psi_det_sorted(j,1,i) + psi_det(j,2,k) = psi_det_sorted(j,2,i) enddo do j=1,N_states - do k=i+1,N_det - psi_coef(k-1,j) = psi_coef(k,j) - enddo + psi_coef(k,j) = psi_coef_sorted(i,j) enddo + else N_removed += 1 endif enddo if (N_removed > 0) then - N_det -= N_removed call write_int(output_dets,N_removed, 'Removed determinants') endif + SOFT_TOUCH N_det psi_det psi_coef end