10
0
mirror of https://github.com/LCPQ/quantum_package synced 2024-06-02 11:25:26 +02:00

Using rPT2 for extrapolation

This commit is contained in:
Anthony Scemama 2018-11-21 12:33:17 +01:00
parent de90c8cec6
commit 57f2a5597a
10 changed files with 47 additions and 50 deletions

View File

@ -2,7 +2,7 @@
type: Threshold
doc: Thresholds of Davidson's algorithm
interface: ezfio,provider,ocaml
default: 1.e-12
default: 1.e-10
[n_states_diag]
type: States_number

View File

@ -280,7 +280,7 @@ subroutine remove_duplicates_in_psi_det(found_duplicates)
N_det = k
psi_det_sorted_bit(:,:,1:N_det) = psi_det(:,:,1:N_det)
psi_coef_sorted_bit(1:N_det,:) = psi_coef(1:N_det,:)
TOUCH N_det psi_det psi_coef psi_det_sorted_bit psi_coef_sorted_bit
TOUCH N_det psi_det psi_coef psi_det_sorted_bit psi_coef_sorted_bit c0_weight
endif
psi_det = psi_det_sorted
psi_coef = psi_coef_sorted

View File

@ -393,21 +393,20 @@ BEGIN_PROVIDER [ double precision, c0_weight, (N_states) ]
BEGIN_DOC
! Weight of the states in the selection : 1/c_0^2
END_DOC
integer :: i,k
double precision :: c
do i=1,N_states
c0_weight(i) = 1.d-31
c = maxval(psi_coef(:,i) * psi_coef(:,i))
c0_weight(i) = 1.d0/(c+1.d-20)
c0_weight(i) = min(c0_weight(i), 100.d0)
enddo
if (mpi_master) then
print *, ''
print *, 'c0 weights'
print *, '----------'
print *, ''
print *, c0_weight(1:N_states)
print *, ''
if (N_states > 1) then
integer :: i
double precision :: c
do i=1,N_states
c0_weight(i) = 1.d-31
c = maxval(psi_coef(:,i) * psi_coef(:,i))
c0_weight(i) = 1.d0/(c+1.d-20)
enddo
c = 1.d0/minval(c0_weight(:))
do i=1,N_states
c0_weight(i) = c0_weight(i) * c
enddo
else
c0_weight = 1.d0
endif
END_PROVIDER

View File

@ -497,6 +497,10 @@ subroutine save_wavefunction
BEGIN_DOC
! Save the wave function into the EZFIO file
END_DOC
! Trick to avoid re-reading the wave function every time N_det changes
read_wf = .False.
if (N_det < N_states) then
return
endif
@ -555,17 +559,10 @@ subroutine save_wavefunction_general(ndet,nstates,psidet,dim_psicoef,psicoef)
allocate (psi_coef_save(ndet,nstates))
double precision :: accu_norm
do k=1,nstates
accu_norm = 0.d0
do i=1,ndet
accu_norm = accu_norm + psicoef(i,k) * psicoef(i,k)
enddo
if (accu_norm == 0.d0) then
accu_norm = 1.e-12
endif
accu_norm = 1.d0/dsqrt(accu_norm)
do i=1,ndet
psi_coef_save(i,k) = psicoef(i,k) * accu_norm
psi_coef_save(i,k) = psicoef(i,k)
enddo
call normalize(psi_coef_save(1,k),ndet)
enddo
call ezfio_set_determinants_psi_coef(psi_coef_save)

View File

@ -1,11 +1,10 @@
program fci_zmq
implicit none
integer :: i,j,k
double precision, allocatable :: pt2(:), variance(:), norm(:)
double precision, allocatable :: pt2(:), variance(:), norm(:), rpt2(:)
integer :: n_det_before, to_select
double precision :: threshold_davidson_in
allocate (pt2(N_states), norm(N_states), variance(N_states))
allocate (pt2(N_states), rpt2(N_states), norm(N_states), variance(N_states))
double precision :: hf_energy_ref
logical :: has
@ -14,11 +13,9 @@ program fci_zmq
relative_error=PT2_relative_error
pt2 = -huge(1.e0)
rpt2 = -huge(1.e0)
norm = 0.d0
variance = 0.d0
threshold_davidson_in = threshold_davidson
threshold_davidson = threshold_davidson_in * 100.d0
SOFT_TOUCH threshold_davidson
call diagonalize_CI
call save_wavefunction
@ -78,8 +75,13 @@ program fci_zmq
call ezfio_set_fci_energy_pt2(psi_energy_with_nucl_rep+pt2)
call write_double(6,correlation_energy_ratio, 'Correlation ratio')
call print_summary(psi_energy_with_nucl_rep(1:N_states),pt2,error,variance,norm)
call save_iterations(psi_energy_with_nucl_rep(1:N_states),pt2,N_det)
call print_extrapolated_energy(psi_energy_with_nucl_rep(1:N_states),pt2)
do k=1,N_states
rpt2(:) = pt2(:)/(1.d0 + norm(k))
enddo
call save_iterations(psi_energy_with_nucl_rep(1:N_states),rpt2,N_det)
call print_extrapolated_energy(psi_energy_with_nucl_rep(1:N_states),rpt2)
N_iter += 1
n_det_before = N_det
@ -98,16 +100,12 @@ program fci_zmq
PROVIDE psi_det
PROVIDE psi_det_sorted
if (N_det >= N_det_max) then
threshold_davidson = threshold_davidson_in
end if
call diagonalize_CI
call save_wavefunction
call ezfio_set_fci_energy(psi_energy_with_nucl_rep(1:N_states))
enddo
if (N_det < N_det_max) then
threshold_davidson = threshold_davidson_in
call diagonalize_CI
call save_wavefunction
call ezfio_set_fci_energy(psi_energy_with_nucl_rep(1:N_states))
@ -133,8 +131,8 @@ program fci_zmq
print*, 'correlation_ratio = ', correlation_energy_ratio
call save_iterations(psi_energy_with_nucl_rep(1:N_states),pt2,N_det)
call write_double(6,correlation_energy_ratio, 'Correlation ratio')
call print_summary(psi_energy_with_nucl_rep(1:N_states),pt2,error,variance,norm)
call save_iterations(psi_energy_with_nucl_rep(1:N_states),rpt2,N_det)
call print_extrapolated_energy(psi_energy_with_nucl_rep(1:N_states),rpt2)
end

View File

@ -6,10 +6,8 @@ BEGIN_PROVIDER [ integer, pt2_stoch_istate ]
pt2_stoch_istate = 1
END_PROVIDER
BEGIN_PROVIDER [ integer, pt2_N_teeth ]
&BEGIN_PROVIDER [ integer, pt2_minDetInFirstTeeth ]
BEGIN_PROVIDER [ integer, pt2_F, (N_det_generators) ]
&BEGIN_PROVIDER [ integer, pt2_n_tasks_max ]
&BEGIN_PROVIDER [ integer, pt2_F, (N_det_generators) ]
implicit none
logical, external :: testTeethBuilding
integer :: i
@ -19,6 +17,12 @@ END_PROVIDER
do i=1,N_det_generators
pt2_F(i) = 1 + int(dble(pt2_n_tasks_max)*maxval(dsqrt(dabs(psi_coef_sorted_gen(i,1:N_states)))))
enddo
END_PROVIDER
BEGIN_PROVIDER [ integer, pt2_N_teeth ]
&BEGIN_PROVIDER [ integer, pt2_minDetInFirstTeeth ]
implicit none
logical, external :: testTeethBuilding
if(N_det_generators < 1024) then
pt2_minDetInFirstTeeth = 1

View File

@ -25,7 +25,7 @@ subroutine run_selection_slave(thread,iproc,energy)
PROVIDE psi_bilinear_matrix_columns_loc psi_det_alpha_unique psi_det_beta_unique
PROVIDE psi_bilinear_matrix_rows psi_det_sorted_order psi_bilinear_matrix_order
PROVIDE psi_bilinear_matrix_transp_rows_loc psi_bilinear_matrix_transp_columns
PROVIDE psi_bilinear_matrix_transp_order
PROVIDE psi_bilinear_matrix_transp_order N_int pt2_F
zmq_to_qp_run_socket = new_zmq_to_qp_run_socket()

View File

@ -727,9 +727,7 @@ subroutine splash_pq(mask, sp, det, i_gen, N_sel, bannedOrb, banned, mat, intere
integer :: i, ii, j, k, l, h(0:2,2), p(0:4,2), nt
integer(bit_kind) :: perMask(N_int, 2), mobMask(N_int, 2), negMask(N_int, 2)
integer(bit_kind) :: phasemask(N_int,2)
! logical :: bandon
!
! bandon = .false.
PROVIDE psi_selectors_coef_transp
mat = 0d0

View File

@ -13,6 +13,7 @@ subroutine ZMQ_selection(N_in, pt2, variance, norm)
double precision, intent(out) :: variance(N_states)
double precision, intent(out) :: norm(N_states)
PROVIDE psi_det psi_coef N_det qp_max_mem N_states pt2_F s2_eig N_det_generators
N = max(N_in,1)
if (.True.) then
@ -60,8 +61,8 @@ subroutine ZMQ_selection(N_in, pt2, variance, norm)
task = ' '
do i= 1, N_det_generators
do j=1,pt2_F(pt2_J(i))
write(task(ipos:ipos+30),'(I9,1X,I9,1X,I9,''|'')') j, pt2_J(i), N
do j=1,pt2_F(i)
write(task(ipos:ipos+30),'(I9,1X,I9,1X,I9,''|'')') j, i, N
ipos += 30
if (ipos > 100000-30) then
if (add_task_to_taskserver(zmq_to_qp_run_socket,trim(task(1:ipos))) == -1) then

View File

@ -20,7 +20,7 @@ subroutine print_summary(e_,pt2_,error_,variance_,norm_)
N_states_p = min(N_det,N_states)
do i=1,N_states_p
f(i) = 1.d0/(1+norm_(i))
f(i) = 1.d0/(1.d0+norm_(i))
enddo
print *, ''