9
1
mirror of https://github.com/QuantumPackage/qp2.git synced 2024-11-20 12:22:21 +01:00

Introduced pt2_type

This commit is contained in:
Anthony Scemama 2020-08-28 00:10:46 +02:00
parent 4ac514dd63
commit e0f17d516b
7 changed files with 126 additions and 85 deletions

View File

@ -1,18 +1,20 @@
subroutine run_cipsi
implicit none
use selection_types
BEGIN_DOC
! Selected Full Configuration Interaction with deterministic selection and
! stochastic PT2.
END_DOC
integer :: i,j,k
double precision, allocatable :: pt2(:), variance(:), norm(:), rpt2(:), zeros(:)
type(pt2_type) :: pt2_data
double precision, allocatable :: rpt2(:), zeros(:)
integer :: to_select
logical, external :: qp_stop
double precision :: threshold_generators_save
double precision :: rss
double precision, external :: memory_of_double
PROVIDE H_apply_buffer_allocated
PROVIDE H_apply_buffer_allocated
N_iter = 1
threshold_generators = 1.d0
@ -21,7 +23,8 @@ subroutine run_cipsi
rss = memory_of_double(N_states)*4.d0
call check_mem(rss,irp_here)
allocate (pt2(N_states), zeros(N_states), rpt2(N_states), norm(N_states), variance(N_states))
allocate (zeros(N_states), rpt2(N_states))
allocate (pt2_data % pt2(N_states), pt2_data % norm2(N_states), pt2_data % variance(N_states))
double precision :: hf_energy_ref
logical :: has
@ -30,10 +33,10 @@ subroutine run_cipsi
relative_error=PT2_relative_error
zeros = 0.d0
pt2 = -huge(1.e0)
rpt2 = -huge(1.e0)
norm = 0.d0
variance = huge(1.e0)
pt2_data % pt2 = -huge(1.e0)
pt2_data % norm2 = 0.d0
pt2_data % variance = huge(1.e0)
if (s2_eig) then
call make_s2_eigenfunction
@ -67,8 +70,8 @@ subroutine run_cipsi
do while ( &
(N_det < N_det_max) .and. &
(maxval(abs(rpt2(1:N_states))) > pt2_max) .and. &
(maxval(abs(variance(1:N_states))) > variance_max) .and. &
(maxval(abs(pt2_data % pt2(1:N_states))) > pt2_max) .and. &
(maxval(abs(pt2_data % variance(1:N_states))) > variance_max) .and. &
(correlation_energy_ratio <= correlation_energy_ratio_max) &
)
write(*,'(A)') '--------------------------------------------------------------------------------'
@ -77,22 +80,21 @@ subroutine run_cipsi
to_select = int(sqrt(dble(N_states))*dble(N_det)*selection_factor)
to_select = max(N_states_diag, to_select)
if (do_pt2) then
pt2 = 0.d0
variance = 0.d0
norm = 0.d0
pt2_data % pt2 = 0.d0
pt2_data % variance = 0.d0
pt2_data % norm2 = 0.d0
threshold_generators_save = threshold_generators
threshold_generators = 1.d0
SOFT_TOUCH threshold_generators
call ZMQ_pt2(psi_energy_with_nucl_rep,pt2,relative_error,error, variance, &
norm, 0) ! Stochastic PT2
call ZMQ_pt2(psi_energy_with_nucl_rep,pt2_data,relative_error,error, 0) ! Stochastic PT2
threshold_generators = threshold_generators_save
SOFT_TOUCH threshold_generators
else
call ZMQ_selection(to_select, pt2, variance, norm)
call ZMQ_selection(to_select, pt2_data)
endif
do k=1,N_states
rpt2(k) = pt2(k)/(1.d0 + norm(k))
rpt2(k) = pt2_data % pt2(k)/(1.d0 + pt2_data % norm2(k))
enddo
correlation_energy_ratio = (psi_energy_with_nucl_rep(1) - hf_energy_ref) / &
@ -100,7 +102,11 @@ subroutine run_cipsi
correlation_energy_ratio = min(1.d0,correlation_energy_ratio)
call write_double(6,correlation_energy_ratio, 'Correlation ratio')
call print_summary(psi_energy_with_nucl_rep,pt2,error,variance,norm,N_det,N_occ_pattern,N_states,psi_s2)
call print_summary(psi_energy_with_nucl_rep, &
pt2_data % pt2, error, &
pt2_data % variance, &
pt2_data % norm2, &
N_det,N_occ_pattern,N_states,psi_s2)
call save_energy(psi_energy_with_nucl_rep, rpt2)
@ -121,7 +127,7 @@ subroutine run_cipsi
call diagonalize_CI
call save_wavefunction
call save_energy(psi_energy_with_nucl_rep, zeros)
if (qp_stop()) exit
if (qp_stop()) exit
enddo
if (.not.qp_stop()) then
@ -132,13 +138,12 @@ subroutine run_cipsi
endif
if (do_pt2) then
pt2(:) = 0.d0
variance(:) = 0.d0
norm(:) = 0.d0
pt2_data % pt2(:) = 0.d0
pt2_data % variance(:) = 0.d0
pt2_data % norm2(:) = 0.d0
threshold_generators = 1d0
SOFT_TOUCH threshold_generators
call ZMQ_pt2(psi_energy_with_nucl_rep, pt2,relative_error,error,variance, &
norm,0) ! Stochastic PT2
call ZMQ_pt2(psi_energy_with_nucl_rep, pt2_data, relative_error, error, 0) ! Stochastic PT2
SOFT_TOUCH threshold_generators
endif
print *, 'N_det = ', N_det
@ -148,11 +153,15 @@ subroutine run_cipsi
do k=1,N_states
rpt2(k) = pt2(k)/(1.d0 + norm(k))
rpt2(k) = pt2_data % pt2(k)/(1.d0 + pt2_data % norm2(k))
enddo
call save_energy(psi_energy_with_nucl_rep, rpt2)
call print_summary(psi_energy_with_nucl_rep(1:N_states),pt2,error,variance,norm,N_det,N_occ_pattern,N_states,psi_s2)
call print_summary(psi_energy_with_nucl_rep(1:N_states), &
pt2_data % pt2, error, &
pt2_data % variance, &
pt2_data % norm2, &
N_det,N_occ_pattern,N_states,psi_s2)
call save_iterations(psi_energy_with_nucl_rep(1:N_states),rpt2,N_det)
call print_extrapolated_energy()
endif

View File

@ -107,7 +107,7 @@ end function
subroutine ZMQ_pt2(E, pt2,relative_error, error, variance, norm2, N_in)
subroutine ZMQ_pt2(E, pt2_data, relative_error, error, N_in)
use f77_zmq
use selection_types
@ -117,10 +117,9 @@ subroutine ZMQ_pt2(E, pt2,relative_error, error, variance, norm2, N_in)
integer, intent(in) :: N_in
! integer, intent(inout) :: N_in
double precision, intent(in) :: relative_error, E(N_states)
double precision, intent(out) :: pt2(N_states),error(N_states)
double precision, intent(out) :: variance(N_states),norm2(N_states)
double precision, intent(out) :: error(N_states)
type(pt2_type), intent(inout) :: pt2_data
!
integer :: i, N
double precision :: state_average_weight_save(N_states), w(N_states,4)
@ -138,10 +137,10 @@ subroutine ZMQ_pt2(E, pt2,relative_error, error, variance, norm2, N_in)
endif
if (N_det <= max(4,N_states) .or. pt2_N_teeth < 2) then
pt2=0.d0
variance=0.d0
norm2=0.d0
call ZMQ_selection(N_in, pt2, variance, norm2)
pt2_data % pt2=0.d0
pt2_data % variance=0.d0
pt2_data % norm2=0.d0
call ZMQ_selection(N_in, pt2_data % pt2, pt2_data % variance, pt2_data % norm2)
error(:) = 0.d0
else
@ -304,10 +303,11 @@ subroutine ZMQ_pt2(E, pt2,relative_error, error, variance, norm2, N_in)
if (i==0) then
call pt2_collector(zmq_socket_pull, E(pt2_stoch_istate),relative_error, w(1,1), w(1,2), w(1,3), w(1,4), b, N)
pt2(pt2_stoch_istate) = w(pt2_stoch_istate,1)
pt2_data % pt2(pt2_stoch_istate) = w(pt2_stoch_istate,1)
error(pt2_stoch_istate) = w(pt2_stoch_istate,2)
variance(pt2_stoch_istate) = w(pt2_stoch_istate,3)
norm2(pt2_stoch_istate) = w(pt2_stoch_istate,4)
pt2_data % variance(pt2_stoch_istate) = w(pt2_stoch_istate,3)
pt2_data % norm2(pt2_stoch_istate) = w(pt2_stoch_istate,4)
!TODO SEGV
else
call pt2_slave_inproc(i)
@ -335,10 +335,10 @@ subroutine ZMQ_pt2(E, pt2,relative_error, error, variance, norm2, N_in)
TOUCH state_average_weight
endif
do k=N_det+1,N_states
pt2(k) = 0.d0
pt2_data % pt2(k) = 0.d0
enddo
call update_pt2_and_variance_weights(pt2, variance, norm2, N_states)
call update_pt2_and_variance_weights(pt2_data, N_states)
end subroutine

View File

@ -19,15 +19,17 @@ BEGIN_PROVIDER [ double precision, variance_match_weight, (N_states) ]
variance_match_weight(:) = 1.d0
END_PROVIDER
subroutine update_pt2_and_variance_weights(pt2, variance, norm2, N_st)
subroutine update_pt2_and_variance_weights(pt2_data, N_st)
implicit none
use selection_types
BEGIN_DOC
! Updates the PT2- and Variance- matching weights.
END_DOC
integer, intent(in) :: N_st
double precision, intent(in) :: pt2(N_st)
double precision, intent(in) :: variance(N_st)
double precision, intent(in) :: norm2(N_st)
type(pt2_type), intent(in) :: pt2_data
double precision :: pt2(N_st)
double precision :: variance(N_st)
double precision :: norm2(N_st)
double precision :: avg, rpt2(N_st), element, dt, x
integer :: k
@ -35,6 +37,10 @@ subroutine update_pt2_and_variance_weights(pt2, variance, norm2, N_st)
integer, parameter :: i_itermax = 1
double precision, allocatable, save :: memo_variance(:,:), memo_pt2(:,:)
pt2(:) = pt2_data % pt2(:)
variance(:) = pt2_data % variance(:)
norm2(:) = pt2_data % norm2(:)
if (i_iter == 0) then
allocate(memo_variance(N_st,i_itermax), memo_pt2(N_st,i_itermax))
memo_pt2(:,:) = 1.d0

View File

@ -5,5 +5,11 @@ module selection_types
double precision, pointer :: val(:)
double precision :: mini
endtype
type pt2_type
double precision, allocatable :: pt2(:)
double precision, allocatable :: variance(:)
double precision, allocatable :: norm2(:)
endtype
end module

View File

@ -1,12 +1,15 @@
subroutine run_stochastic_cipsi
use selection_types
implicit none
BEGIN_DOC
! Selected Full Configuration Interaction with Stochastic selection and PT2.
END_DOC
integer :: i,j,k
double precision, allocatable :: pt2(:), variance(:), norm2(:), rpt2(:), zeros(:)
double precision, allocatable :: rpt2(:), zeros(:)
integer :: to_select
logical, external :: qp_stop
type(pt2_type) :: pt2_data
logical, external :: qp_stop
double precision :: rss
double precision, external :: memory_of_double
@ -19,7 +22,11 @@ subroutine run_stochastic_cipsi
rss = memory_of_double(N_states)*4.d0
call check_mem(rss,irp_here)
allocate (pt2(N_states), zeros(N_states), rpt2(N_states), norm2(N_states), variance(N_states))
allocate (zeros(N_states), rpt2(N_states))
allocate( pt2_data % pt2(N_states) )
allocate( pt2_data % variance(N_states) )
allocate( pt2_data % norm2(N_states) )
double precision :: hf_energy_ref
logical :: has
@ -28,10 +35,10 @@ subroutine run_stochastic_cipsi
relative_error=PT2_relative_error
zeros = 0.d0
pt2 = -huge(1.e0)
pt2_data % pt2 = -huge(1.e0)
rpt2 = -huge(1.e0)
norm2 = 0.d0
variance = huge(1.e0)
pt2_data % norm2 = 0.d0
pt2_data % variance = huge(1.e0)
if (s2_eig) then
call make_s2_eigenfunction
@ -65,8 +72,8 @@ subroutine run_stochastic_cipsi
do while ( &
(N_det < N_det_max) .and. &
(maxval(abs(rpt2(1:N_states))) > pt2_max) .and. &
(maxval(abs(variance(1:N_states))) > variance_max) .and. &
(maxval(abs(pt2_data % pt2(1:N_states))) > pt2_max) .and. &
(maxval(abs(pt2_data % variance(1:N_states))) > variance_max) .and. &
(correlation_energy_ratio <= correlation_energy_ratio_max) &
)
write(*,'(A)') '--------------------------------------------------------------------------------'
@ -75,14 +82,15 @@ subroutine run_stochastic_cipsi
to_select = int(sqrt(dble(N_states))*dble(N_det)*selection_factor)
to_select = max(N_states_diag, to_select)
pt2 = 0.d0
variance = 0.d0
norm2 = 0.d0
call ZMQ_pt2(psi_energy_with_nucl_rep,pt2,relative_error,error, variance, &
norm2, to_select) ! Stochastic PT2 and selection
pt2_data % pt2 = 0.d0
pt2_data % variance = 0.d0
pt2_data % norm2 = 0.d0
call ZMQ_pt2(psi_energy_with_nucl_rep,pt2_data,relative_error,error, &
to_select) ! Stochastic PT2 and selection
do k=1,N_states
rpt2(k) = pt2(k)/(1.d0 + norm2(k))
rpt2(k) = pt2_data % pt2(k)/(1.d0 + pt2_data % norm2(k))
enddo
correlation_energy_ratio = (psi_energy_with_nucl_rep(1) - hf_energy_ref) / &
@ -90,7 +98,11 @@ subroutine run_stochastic_cipsi
correlation_energy_ratio = min(1.d0,correlation_energy_ratio)
call write_double(6,correlation_energy_ratio, 'Correlation ratio')
call print_summary(psi_energy_with_nucl_rep,pt2,error,variance,norm2,N_det,N_occ_pattern,N_states,psi_s2)
call print_summary(psi_energy_with_nucl_rep, &
pt2_data % pt2, error, &
pt2_data % variance, &
pt2_data % norm2, &
N_det,N_occ_pattern,N_states,psi_s2)
call save_energy(psi_energy_with_nucl_rep, rpt2)
@ -121,18 +133,21 @@ subroutine run_stochastic_cipsi
call save_energy(psi_energy_with_nucl_rep, zeros)
endif
pt2(:) = 0.d0
variance(:) = 0.d0
norm2(:) = 0.d0
call ZMQ_pt2(psi_energy_with_nucl_rep, pt2,relative_error,error,variance, &
norm2,0) ! Stochastic PT2
pt2_data % pt2(:) = 0.d0
pt2_data % variance(:) = 0.d0
pt2_data % norm2(:) = 0.d0
call ZMQ_pt2(psi_energy_with_nucl_rep, pt2_data, relative_error, error, 0) ! Stochastic PT2
do k=1,N_states
rpt2(k) = pt2(k)/(1.d0 + norm2(k))
rpt2(k) = pt2_data % pt2(k)/(1.d0 + pt2_data % norm2(k))
enddo
call save_energy(psi_energy_with_nucl_rep, rpt2)
call print_summary(psi_energy_with_nucl_rep(1:N_states),pt2,error,variance,norm2,N_det,N_occ_pattern,N_states,psi_s2)
call print_summary(psi_energy_with_nucl_rep, &
pt2_data % pt2, error, &
pt2_data % variance, &
pt2_data % norm2, &
N_det,N_occ_pattern,N_states,psi_s2)
call save_iterations(psi_energy_with_nucl_rep(1:N_states),rpt2,N_det)
call print_extrapolated_energy()
endif

View File

@ -1,4 +1,4 @@
subroutine ZMQ_selection(N_in, pt2, variance, norm2)
subroutine ZMQ_selection(N_in, pt2_data)
use f77_zmq
use selection_types
@ -9,9 +9,7 @@ subroutine ZMQ_selection(N_in, pt2, variance, norm2)
type(selection_buffer) :: b
integer :: i, N
integer, external :: omp_get_thread_num
double precision, intent(out) :: pt2(N_states)
double precision, intent(out) :: variance(N_states)
double precision, intent(out) :: norm2(N_states)
type(pt2_type), intent(inout) :: pt2_data
! PROVIDE psi_det psi_coef N_det qp_max_mem N_states pt2_F s2_eig N_det_generators
@ -112,19 +110,20 @@ subroutine ZMQ_selection(N_in, pt2, variance, norm2)
enddo
endif
!$OMP PARALLEL DEFAULT(shared) SHARED(b, pt2, variance, norm2) PRIVATE(i) NUM_THREADS(nproc_target+1)
!$OMP PARALLEL DEFAULT(shared) SHARED(b, pt2_data) PRIVATE(i) NUM_THREADS(nproc_target+1)
i = omp_get_thread_num()
if (i==0) then
call selection_collector(zmq_socket_pull, b, N, pt2, variance, norm2)
call selection_collector(zmq_socket_pull, b, N, &
pt2_data % pt2, pt2_data % variance, pt2_data % norm2)
else
call selection_slave_inproc(i)
endif
!$OMP END PARALLEL
call end_parallel_job(zmq_to_qp_run_socket, zmq_socket_pull, 'selection')
do i=N_det+1,N_states
pt2(i) = 0.d0
variance(i) = 0.d0
norm2(i) = 0.d0
pt2_data % pt2(i) = 0.d0
pt2_data % variance(i) = 0.d0
pt2_data % norm2(i) = 0.d0
enddo
if (N_in > 0) then
if (s2_eig) then
@ -134,12 +133,12 @@ subroutine ZMQ_selection(N_in, pt2, variance, norm2)
endif
call delete_selection_buffer(b)
do k=1,N_states
pt2(k) = pt2(k) * f(k)
variance(k) = variance(k) * f(k)
norm2(k) = norm2(k) * f(k)
pt2_data % pt2(k) = pt2_data % pt2(k) * f(k)
pt2_data % variance(k) = pt2_data % variance(k) * f(k)
pt2_data % norm2(k) = pt2_data % norm2(k) * f(k)
enddo
call update_pt2_and_variance_weights(pt2, variance, norm2, N_states)
call update_pt2_and_variance_weights(pt2_data, N_states)
end subroutine

View File

@ -28,35 +28,41 @@ end
subroutine run
implicit none
use selection_types
integer :: i,j,k
logical, external :: detEq
double precision :: pt2(N_states)
type(pt2_type) :: pt2_data
integer :: degree
integer :: n_det_before, to_select
double precision :: threshold_davidson_in
double precision :: E_CI_before(N_states), relative_error, error(N_states), variance(N_states), norm2(N_states), rpt2(N_states)
pt2(:) = 0.d0
allocate( pt2_data % pt2(N_states) )
allocate( pt2_data % variance(N_states) )
allocate( pt2_data % norm2(N_states) )
E_CI_before(:) = psi_energy(:) + nuclear_repulsion
relative_error=PT2_relative_error
if (do_pt2) then
call ZMQ_pt2(psi_energy_with_nucl_rep,pt2,relative_error,error, variance, &
norm2,0) ! Stochastic PT2
call ZMQ_pt2(psi_energy_with_nucl_rep,pt2_data,relative_error,error,0) ! Stochastic PT2
else
call ZMQ_selection(0, pt2, variance, norm2)
call ZMQ_selection(0, pt2_data)
endif
do k=1,N_states
rpt2(k) = pt2(k)/(1.d0 + norm2(k))
rpt2(k) = pt2_data % pt2(k)/(1.d0 + pt2_data % norm2(k))
enddo
call print_summary(psi_energy_with_nucl_rep(1:N_states),pt2,error,variance,norm2,N_det,N_occ_pattern,N_states,psi_s2)
call print_summary(psi_energy_with_nucl_rep(1:N_states), &
pt2_data % pt2, error, &
pt2_data % variance, &
pt2_data % norm2, &
N_det,N_occ_pattern,N_states,psi_s2)
call save_energy(E_CI_before,pt2)
call save_energy(E_CI_before,pt2_data % pt2)
end