From f75ce67a8749f92cce038324e74502151f840cc5 Mon Sep 17 00:00:00 2001 From: Anthony Scemama Date: Sun, 25 Sep 2016 22:14:17 +0200 Subject: [PATCH] FIxed Davidson --- plugins/Full_CI_ZMQ/fci_zmq.irp.f | 84 ++++++++++++----------- plugins/Full_CI_ZMQ/selection_slave.irp.f | 2 +- plugins/Generators_full/generators.irp.f | 2 + src/Davidson/diagonalize_CI.irp.f | 20 +++--- src/Davidson/parameters.irp.f | 2 +- 5 files changed, 56 insertions(+), 54 deletions(-) diff --git a/plugins/Full_CI_ZMQ/fci_zmq.irp.f b/plugins/Full_CI_ZMQ/fci_zmq.irp.f index dbf8414b..5c4f15a4 100644 --- a/plugins/Full_CI_ZMQ/fci_zmq.irp.f +++ b/plugins/Full_CI_ZMQ/fci_zmq.irp.f @@ -1,11 +1,11 @@ program fci_zmq implicit none integer :: i,j,k - logical, external :: detEq + logical, external :: detEq double precision, allocatable :: pt2(:) integer :: degree - + allocate (pt2(N_states)) pt2 = 1.d0 @@ -30,10 +30,10 @@ program fci_zmq print *, '-----' enddo endif - double precision :: E_CI_before(N_states) + double precision :: E_CI_before(N_states) - - integer :: n_det_before + + integer :: n_det_before print*,'Beginning the selection ...' E_CI_before(1:N_states) = CI_energy(1:N_states) @@ -44,16 +44,16 @@ program fci_zmq PROVIDE psi_coef PROVIDE psi_det PROVIDE psi_det_sorted - + if (N_det > N_det_max) then - psi_det = psi_det_sorted - psi_coef = psi_coef_sorted - N_det = N_det_max - soft_touch N_det psi_det psi_coef + psi_det = psi_det_sorted + psi_coef = psi_coef_sorted + N_det = N_det_max + soft_touch N_det psi_det psi_coef endif call diagonalize_CI call save_wavefunction - + print *, 'N_det = ', N_det print *, 'N_states = ', N_states do k=1, N_states @@ -64,42 +64,44 @@ program fci_zmq enddo print *, '-----' if(N_states.gt.1)then - print*,'Variational Energy difference' - do i = 2, N_states - print*,'Delta E = ',CI_energy(i) - CI_energy(1) - enddo + print*,'Variational Energy difference' + do i = 2, N_states + print*,'Delta E = ',CI_energy(i) - CI_energy(1) + enddo endif if(N_states.gt.1)then - print*,'Variational + perturbative Energy difference' - do i = 2, N_states - print*,'Delta E = ',E_CI_before(i)+ pt2(i) - (E_CI_before(1) + pt2(1)) - enddo + print*,'Variational + perturbative Energy difference' + do i = 2, N_states + print*,'Delta E = ',E_CI_before(i)+ pt2(i) - (E_CI_before(1) + pt2(1)) + enddo endif E_CI_before(1:N_states) = CI_energy(1:N_states) call ezfio_set_full_ci_energy(CI_energy) enddo - N_det = min(N_det_max,N_det) - touch N_det psi_det psi_coef - call diagonalize_CI - if(do_pt2_end)then - print*,'Last iteration only to compute the PT2' - threshold_selectors = 1.d0 - threshold_generators = 0.9999d0 - E_CI_before(1:N_states) = CI_energy(1:N_states) - call ZMQ_selection(1, pt2) - print *, 'Final step' - print *, 'N_det = ', N_det - print *, 'N_states = ', N_states - do k=1,N_states - print *, 'State', k - print *, 'PT2 = ', pt2 - print *, 'E = ', E_CI_before - print *, 'E+PT2 = ', E_CI_before+pt2 - print *, '-----' - enddo - call ezfio_set_full_ci_energy_pt2(E_CI_before+pt2) - endif - call save_wavefunction + if (N_det > N_det_max) then + N_det = N_det_max + touch N_det psi_det psi_coef + call diagonalize_CI + endif + if(do_pt2_end)then + print*,'Last iteration only to compute the PT2' + threshold_selectors = 1.d0 + threshold_generators = 0.9999d0 + E_CI_before(1:N_states) = CI_energy(1:N_states) + call ZMQ_selection(1, pt2) + print *, 'Final step' + print *, 'N_det = ', N_det + print *, 'N_states = ', N_states + do k=1,N_states + print *, 'State', k + print *, 'PT2 = ', pt2 + print *, 'E = ', E_CI_before + print *, 'E+PT2 = ', E_CI_before+pt2 + print *, '-----' + enddo + call ezfio_set_full_ci_energy_pt2(E_CI_before+pt2) + endif + call save_wavefunction end diff --git a/plugins/Full_CI_ZMQ/selection_slave.irp.f b/plugins/Full_CI_ZMQ/selection_slave.irp.f index 438892d3..4c365238 100644 --- a/plugins/Full_CI_ZMQ/selection_slave.irp.f +++ b/plugins/Full_CI_ZMQ/selection_slave.irp.f @@ -58,7 +58,7 @@ subroutine update_energy(energy) ! Update energy when it is received from ZMQ END_DOC integer :: j,k - do j=1,N_states_diag + do j=1,N_states do k=1,N_det CI_eigenvectors(k,j) = psi_coef(k,j) enddo diff --git a/plugins/Generators_full/generators.irp.f b/plugins/Generators_full/generators.irp.f index 14bc18d4..eea5821b 100644 --- a/plugins/Generators_full/generators.irp.f +++ b/plugins/Generators_full/generators.irp.f @@ -30,6 +30,8 @@ END_PROVIDER ! Hartree-Fock determinant END_DOC integer :: i, k + psi_coef_generators = 0.d0 + psi_det_generators = 0_bit_kind do i=1,N_det_generators do k=1,N_int psi_det_generators(k,1,i) = psi_det_sorted(k,1,i) diff --git a/src/Davidson/diagonalize_CI.irp.f b/src/Davidson/diagonalize_CI.irp.f index b6e57cb9..e6b230b5 100644 --- a/src/Davidson/diagonalize_CI.irp.f +++ b/src/Davidson/diagonalize_CI.irp.f @@ -8,8 +8,10 @@ BEGIN_PROVIDER [ double precision, CI_energy, (N_states_diag) ] integer :: j character*(8) :: st call write_time(output_determinants) - do j=1,min(N_det,N_states) + do j=1,min(N_det,N_states_diag) CI_energy(j) = CI_electronic_energy(j) + nuclear_repulsion + enddo + do j=1,min(N_det,N_states) write(st,'(I4)') j call write_double(output_determinants,CI_energy(j),'Energy of state '//trim(st)) call write_double(output_determinants,CI_eigenvectors_s2(j),'S^2 of state '//trim(st)) @@ -38,14 +40,14 @@ END_PROVIDER double precision, allocatable :: e_array(:) integer, allocatable :: iorder(:) - ! Guess values for the "N_states_diag" states of the CI_eigenvectors + ! Guess values for the "N_states" states of the CI_eigenvectors do j=1,min(N_states,N_det) do i=1,N_det CI_eigenvectors(i,j) = psi_coef(i,j) enddo enddo - do j=N_det+1,N_states_diag + do j=min(N_states,N_det)+1,N_states_diag do i=1,N_det CI_eigenvectors(i,j) = 0.d0 enddo @@ -143,14 +145,15 @@ END_PROVIDER endif - if( s2_eig.and.(n_states_diag > 1).and.(n_det >= n_states_diag) )then + if( s2_eig.and.(N_states_diag > 1).and.(N_det >= N_states_diag) )then ! Diagonalizing S^2 within the "n_states_diag" states found allocate(s2_eigvalues(N_states_diag), e_array(N_states_diag)) - call diagonalize_s2_betweenstates(psi_det,CI_eigenvectors,n_det,size(psi_det,3),size(CI_eigenvectors,1),min(n_states_diag,n_det),s2_eigvalues) + call diagonalize_s2_betweenstates(psi_det,CI_eigenvectors,N_det,size(psi_det,3), & + size(CI_eigenvectors,1),min(n_states_diag,n_det),s2_eigvalues) double precision, allocatable :: psi_coef_tmp(:,:) allocate(psi_coef_tmp(psi_det_size,N_states_diag)) - do j = 1, N_states + do j = 1, N_states_diag do i = 1, N_det psi_coef_tmp(i,j) = CI_eigenvectors(i,j) enddo @@ -200,11 +203,6 @@ END_PROVIDER CI_eigenvectors_s2(i_state + i_other_state) = s2_eigvalues(j) CI_electronic_energy(i_state + i_other_state) = e_array(i_state + i_other_state) enddo - do j=1,N_states - do i=1,N_det - psi_coef(i,j) = psi_coef_tmp(i,j) - enddo - enddo deallocate(iorder,e_array,index_good_state_array,good_state_array,psi_coef_tmp) deallocate(s2_eigvalues) diff --git a/src/Davidson/parameters.irp.f b/src/Davidson/parameters.irp.f index d9c82e3c..82315495 100644 --- a/src/Davidson/parameters.irp.f +++ b/src/Davidson/parameters.irp.f @@ -12,7 +12,7 @@ BEGIN_PROVIDER [ integer, davidson_sze_max ] ! Max number of Davidson sizes END_DOC ASSERT (davidson_sze_max <= davidson_iter_max) - davidson_sze_max = 8*N_states + davidson_sze_max = N_states+7 END_PROVIDER