10
0
mirror of https://github.com/LCPQ/quantum_package synced 2024-11-04 21:24:02 +01:00

multistate shifted_bk

This commit is contained in:
Yann Garniron 2018-05-17 16:02:51 +02:00
parent 7cc33f1ab3
commit f45fc46b4f
5 changed files with 66 additions and 33 deletions

View File

@ -28,6 +28,11 @@ subroutine run_wf
double precision :: energy(N_states_diag) double precision :: energy(N_states_diag)
character*(64) :: states(1) character*(64) :: states(1)
integer :: rc, i integer :: rc, i
integer, external :: zmq_get_dvector, zmq_get_N_det_generators
integer, external :: zmq_get_psi, zmq_get_N_det_selectors
integer, external :: zmq_get_N_states_diag
double precision :: tmp
call provide_everything call provide_everything
@ -43,10 +48,22 @@ subroutine run_wf
exit exit
else if (zmq_state(:5) == 'dress') then else if (zmq_state(:5) == 'dress') then
! Dress
! Selection
! --------- ! ---------
call zmq_get_psi(zmq_to_qp_run_socket,1,energy,N_states) !call zmq_get_psi(zmq_to_qp_run_socket,1,energy,N_states)
if (zmq_get_psi(zmq_to_qp_run_socket,1) == -1) cycle
!TOUCH psi_det
if (zmq_get_N_det_generators (zmq_to_qp_run_socket, 1) == -1) cycle
if (zmq_get_N_det_selectors(zmq_to_qp_run_socket, 1) == -1) cycle
if (zmq_get_dvector(zmq_to_qp_run_socket,1,'state_average_weight',state_average_weight,N_states) == -1) cycle
if (zmq_get_dvector(zmq_to_qp_run_socket,1,'energy',energy,N_states) == -1) cycle
if (zmq_get_dvector(zmq_to_qp_run_socket,1,'dress_stoch_istate',tmp,1) == -1) cycle
dress_stoch_istate = int(tmp)
TOUCH dress_stoch_istate
TOUCH state_average_weight
PROVIDE psi_bilinear_matrix_columns_loc psi_det_alpha_unique psi_det_beta_unique 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_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_rows_loc psi_bilinear_matrix_transp_columns
@ -57,7 +74,6 @@ subroutine run_wf
call dress_slave_tcp(0, energy) call dress_slave_tcp(0, energy)
!!$OMP END PARALLEL !!$OMP END PARALLEL
endif endif
end do end do
end end

View File

@ -4,11 +4,12 @@ BEGIN_PROVIDER [ integer, fragment_first ]
END_PROVIDER END_PROVIDER
subroutine ZMQ_dress(E, dress, delta_out, delta_s2_out, relative_error) subroutine ZMQ_dress(E, dress, delta_out, delta_s2_out, relative_error, lndet)
use f77_zmq use f77_zmq
implicit none implicit none
integer, intent(in) :: lndet
character(len=64000) :: task character(len=64000) :: task
character(len=3200) :: temp character(len=3200) :: temp
integer(ZMQ_PTR) :: zmq_to_qp_run_socket, zmq_socket_pull integer(ZMQ_PTR) :: zmq_to_qp_run_socket, zmq_socket_pull
@ -27,11 +28,9 @@ subroutine ZMQ_dress(E, dress, delta_out, delta_s2_out, relative_error)
double precision :: time double precision :: time
integer, external :: add_task_to_taskserver integer, external :: add_task_to_taskserver
double precision :: state_average_weight_save(N_states) double precision :: state_average_weight_save(N_states)
task(:) = CHAR(0) task(:) = CHAR(0)
temp(:) = CHAR(0) temp(:) = CHAR(0)
allocate(delta(N_states,N_det), delta_s2(N_det,N_states)) allocate(delta(N_states,N_det), delta_s2(N_states, N_det))
state_average_weight_save(:) = state_average_weight(:) state_average_weight_save(:) = state_average_weight(:)
do dress_stoch_istate=1,N_states do dress_stoch_istate=1,N_states
SOFT_TOUCH dress_stoch_istate SOFT_TOUCH dress_stoch_istate
@ -39,14 +38,14 @@ subroutine ZMQ_dress(E, dress, delta_out, delta_s2_out, relative_error)
state_average_weight(dress_stoch_istate) = 1.d0 state_average_weight(dress_stoch_istate) = 1.d0
TOUCH state_average_weight TOUCH state_average_weight
!provide psi_coef_generators
provide nproc fragment_first fragment_count mo_bielec_integrals_in_map mo_mono_elec_integral dress_weight psi_selectors provide nproc fragment_first fragment_count mo_bielec_integrals_in_map mo_mono_elec_integral dress_weight psi_selectors
!print *, dress_e0_denominator
print *, '========== ================= ================= =================' print *, '========== ================= ================= ================='
print *, ' Samples Energy Stat. Error Seconds ' print *, ' Samples Energy Stat. Error Seconds '
print *, '========== ================= ================= =================' print *, '========== ================= ================= ================='
call new_parallel_job(zmq_to_qp_run_socket,zmq_socket_pull, 'dress') call new_parallel_job(zmq_to_qp_run_socket,zmq_socket_pull, 'dress')
integer, external :: zmq_put_psi integer, external :: zmq_put_psi
@ -54,6 +53,7 @@ subroutine ZMQ_dress(E, dress, delta_out, delta_s2_out, relative_error)
integer, external :: zmq_put_N_det_selectors integer, external :: zmq_put_N_det_selectors
integer, external :: zmq_put_dvector integer, external :: zmq_put_dvector
integer, external :: zmq_set_running integer, external :: zmq_set_running
if (zmq_put_psi(zmq_to_qp_run_socket,1) == -1) then if (zmq_put_psi(zmq_to_qp_run_socket,1) == -1) then
stop 'Unable to put psi on ZMQ server' stop 'Unable to put psi on ZMQ server'
endif endif
@ -66,6 +66,13 @@ subroutine ZMQ_dress(E, dress, delta_out, delta_s2_out, relative_error)
if (zmq_put_dvector(zmq_to_qp_run_socket,1,'energy',dress_e0_denominator,size(dress_e0_denominator)) == -1) then if (zmq_put_dvector(zmq_to_qp_run_socket,1,'energy',dress_e0_denominator,size(dress_e0_denominator)) == -1) then
stop 'Unable to put energy on ZMQ server' stop 'Unable to put energy on ZMQ server'
endif endif
if (zmq_put_dvector(zmq_to_qp_run_socket,1,"state_average_weight",state_average_weight,N_states) == -1) then
stop 'Unable to put state_average_weight on ZMQ server'
endif
if (zmq_put_dvector(zmq_to_qp_run_socket,1,"dress_stoch_istate",real(dress_stoch_istate,8),1) == -1) then
stop 'Unable to put dress_stoch_istate on ZMQ server'
endif
integer(ZMQ_PTR), external :: new_zmq_to_qp_run_socket integer(ZMQ_PTR), external :: new_zmq_to_qp_run_socket
integer :: ipos, sz integer :: ipos, sz
@ -131,13 +138,13 @@ subroutine ZMQ_dress(E, dress, delta_out, delta_s2_out, relative_error)
!i = omp_get_thread_num() !i = omp_get_thread_num()
!if (i==0) then !if (i==0) then
call dress_collector(zmq_socket_pull,E, relative_error, delta, delta_s2, dress,& call dress_collector(zmq_socket_pull,E, relative_error, delta, delta_s2, dress,&
dress_stoch_istate) dress_stoch_istate)
!else !else
! call dress_slave_inproc(i) ! call dress_slave_inproc(i)
!endif !endif
!!$OMP END PARALLEL !!$OMP END PARALLEL
delta_out(dress_stoch_istate,1:N_det) = delta(dress_stoch_istate,1:N_det) delta_out(dress_stoch_istate,1:N_det) = delta(dress_stoch_istate,1:N_det)
delta_s2_out(dress_stoch_istate,1:N_det) = delta_s2_out(dress_stoch_istate,1:N_det) delta_s2_out(dress_stoch_istate,1:N_det) = delta_s2(dress_stoch_istate,1:N_det)
call end_parallel_job(zmq_to_qp_run_socket, zmq_socket_pull, 'dress') call end_parallel_job(zmq_to_qp_run_socket, zmq_socket_pull, 'dress')
print *, '========== ================= ================= =================' print *, '========== ================= ================= ================='
@ -194,7 +201,7 @@ subroutine dress_collector(zmq_socket_pull, E, relative_error, delta, delta_s2,
double precision, allocatable :: double_buf(:) double precision, allocatable :: double_buf(:)
integer(bit_kind), allocatable :: det_buf(:,:,:) integer(bit_kind), allocatable :: det_buf(:,:,:)
integer, external :: zmq_delete_tasks integer, external :: zmq_delete_tasks
last_cp = 10000000
allocate(agreg_for_cp(N_cp)) allocate(agreg_for_cp(N_cp))
agreg_for_cp = 0 agreg_for_cp = 0
allocate(int_buf(N_dress_int_buffer), double_buf(N_dress_double_buffer), det_buf(N_int,2,N_dress_det_buffer)) allocate(int_buf(N_dress_int_buffer), double_buf(N_dress_double_buffer), det_buf(N_int,2,N_dress_det_buffer))
@ -222,6 +229,7 @@ subroutine dress_collector(zmq_socket_pull, E, relative_error, delta, delta_s2,
call pull_dress_results(zmq_socket_pull, ind, cur_cp, delta_loc, int_buf, double_buf, det_buf, N_buf, task_id, dress_mwen) call pull_dress_results(zmq_socket_pull, ind, cur_cp, delta_loc, int_buf, double_buf, det_buf, N_buf, task_id, dress_mwen)
if(floop) then if(floop) then
call wall_time(time) call wall_time(time)
print *, "first_pull", time-time0
time0 = time time0 = time
floop = .false. floop = .false.
end if end if
@ -237,6 +245,7 @@ subroutine dress_collector(zmq_socket_pull, E, relative_error, delta, delta_s2,
endif endif
if(more == 0) loop = .false. !stop 'loop = .false.' !!!!!!!!!!!!!!!! if(more == 0) loop = .false. !stop 'loop = .false.' !!!!!!!!!!!!!!!!
dress_detail(:, ind) = dress_mwen(:) dress_detail(:, ind) = dress_mwen(:)
!print *, "DETAIL", ind, dress_mwen
else if(cur_cp > 0) then else if(cur_cp > 0) then
if(ind == 0) cycle if(ind == 0) cycle
!$OMP PARALLEL DO DEFAULT(SHARED) PRIVATE(i) !$OMP PARALLEL DO DEFAULT(SHARED) PRIVATE(i)
@ -248,8 +257,9 @@ subroutine dress_collector(zmq_socket_pull, E, relative_error, delta, delta_s2,
do i=1,N_det do i=1,N_det
cp(:,i,cur_cp,2) += delta_loc(:,i,2) cp(:,i,cur_cp,2) += delta_loc(:,i,2)
end do end do
!$OMP END PARALLEL DO
agreg_for_cp(cur_cp) += ind agreg_for_cp(cur_cp) += ind
!print *, agreg_for_cp(cur_cp), ind, needed_by_cp(cur_cp), cur_cp
if(agreg_for_cp(cur_cp) > needed_by_cp(cur_cp)) then if(agreg_for_cp(cur_cp) > needed_by_cp(cur_cp)) then
stop "too much results..." stop "too much results..."
end if end if
@ -270,6 +280,7 @@ subroutine dress_collector(zmq_socket_pull, E, relative_error, delta, delta_s2,
su += val su += val
su2 += val*val su2 += val*val
end do end do
!$OMP END PARALLEL DO
avg = su / cps_N(cur_cp) avg = su / cps_N(cur_cp)
eqt = dsqrt( ((su2 / cps_N(cur_cp)) - avg*avg) / cps_N(cur_cp) ) eqt = dsqrt( ((su2 / cps_N(cur_cp)) - avg*avg) / cps_N(cur_cp) )
@ -278,8 +289,9 @@ subroutine dress_collector(zmq_socket_pull, E, relative_error, delta, delta_s2,
E0 = E0 + dress_detail(istate, first_det_of_teeth(cp_first_tooth(cur_cp))) * (1d0-fractage(cp_first_tooth(cur_cp))) E0 = E0 + dress_detail(istate, first_det_of_teeth(cp_first_tooth(cur_cp))) * (1d0-fractage(cp_first_tooth(cur_cp)))
end if end if
print '(2X, F16.7, 2X, G16.3, 2X, F16.4, A20)', avg+E(istate)+E0, eqt, time-time0, '' !print '(I2X, F16.7, 2X, G16.3, 2X, F16.4, A20)', avg+E(istate)+E0, eqt, time-time0, ''
if ((dabs(eqt) < relative_error .and. cps_N(cur_cp) >= 30) .or. cur_cp == cur_cp-2) then print '(G10.3, 2X, F16.10, 2X, G16.3, 2X, F16.4, A20)', cps_N(cur_cp), avg+E0+E(istate), eqt, time-time0, ''
if ((dabs(eqt) < relative_error .and. cps_N(cur_cp) >= 30)) then
! Termination ! Termination
print *, "TERMINATE" print *, "TERMINATE"
if (zmq_abort(zmq_to_qp_run_socket) == -1) then if (zmq_abort(zmq_to_qp_run_socket) == -1) then
@ -365,6 +377,7 @@ END_PROVIDER
logical :: fracted logical :: fracted
integer :: first_suspect integer :: first_suspect
provide psi_coef_generators
first_suspect = 1 first_suspect = 1
allocate(filler(n_det_generators)) allocate(filler(n_det_generators))
@ -397,7 +410,7 @@ END_PROVIDER
end do end do
l=first_det_of_comb l=first_det_of_comb
call random_seed(put=(/321,654,65,321,65,321,654,65,321,65321,654,65,321,65321,654,65,321,65321,654,65,321,65/)) call random_seed(put=(/321,654,65,321,65,321,654,65,321,6321,654,65,321,621,654,65,321,65,654,65,321,65/))
call RANDOM_NUMBER(comb) call RANDOM_NUMBER(comb)
lfiller = 1 lfiller = 1
nfiller = 1 nfiller = 1

View File

@ -100,14 +100,15 @@ BEGIN_PROVIDER [ double precision, delta_ij_tmp, (N_states,N_det_delta_ij,2) ]
! else ! else
! errr = 1d-4 ! errr = 1d-4
! end if ! end if
relative_error = 1.d-5 relative_error = 5.d-5
call write_double(6,relative_error,"Convergence of the stochastic algorithm") call write_double(6,relative_error,"Convergence of the stochastic algorithm")
call ZMQ_dress(E_CI_before, dress, del, del_s2, abs(relative_error)) call ZMQ_dress(E_CI_before, dress, del, del_s2, abs(relative_error), N_det_delta_ij)
delta_ij_tmp(:,:,1) = del(:,:) delta_ij_tmp(:,:,1) = del(:,:)
delta_ij_tmp(:,:,2) = del_s2(:,:) delta_ij_tmp(:,:,2) = del_s2(:,:)
deallocate(dress, del, del_s2) deallocate(dress, del, del_s2)
end if end if
END_PROVIDER END_PROVIDER

View File

@ -116,7 +116,7 @@ subroutine run_dress_slave(thread,iproce,energy)
!$OMP ATOMIC !$OMP ATOMIC
done_for(done_cp_at_det(i_generator)) += 1 done_for(done_cp_at_det(i_generator)) += 1
! print *, "IGEN", i_generator, done_cp_at_det(i_generator)
delta_ij_loc(:,:,:) = 0d0 delta_ij_loc(:,:,:) = 0d0
call generator_start(i_generator, iproc) call generator_start(i_generator, iproc)
call alpha_callback(delta_ij_loc, i_generator, subset, iproc) call alpha_callback(delta_ij_loc, i_generator, subset, iproc)
@ -175,12 +175,14 @@ subroutine run_dress_slave(thread,iproce,energy)
delta_ij_loc = 0d0 delta_ij_loc = 0d0
if(cur_cp < 1) stop "cur_cp < 1" if(cur_cp < 1) stop "cur_cp < 1"
do i=1,cur_cp do i=1,cur_cp
delta_ij_loc(:,:,:) += cp(:,:,i,:) delta_ij_loc(:,:,1) += cp(:,:,i,1)
delta_ij_loc(:,:,2) += cp(:,:,i,2)
end do end do
delta_ij_loc(:,:,:) = delta_ij_loc(:,:,:) / cps_N(cur_cp) delta_ij_loc(:,:,:) = delta_ij_loc(:,:,:) / cps_N(cur_cp)
do i=cp_first_tooth(cur_cp)-1,0,-1 do i=cp_first_tooth(cur_cp)-1,0,-1
delta_ij_loc(:,:,:) = delta_ij_loc(:,:,:) +delta_det(:,:,i,:) delta_ij_loc(:,:,1) = delta_ij_loc(:,:,1) +delta_det(:,:,i,1)
delta_ij_loc(:,:,2) = delta_ij_loc(:,:,2) +delta_det(:,:,i,2)
end do end do
call push_dress_results(zmq_socket_push, done_for(cur_cp), cur_cp, delta_ij_loc, int_buf, double_buf, det_buf, N_buf, -1) call push_dress_results(zmq_socket_push, done_for(cur_cp), cur_cp, delta_ij_loc, int_buf, double_buf, det_buf, N_buf, -1)
end if end if

View File

@ -196,16 +196,17 @@ subroutine undress_with_alpha(old_generators, old_det_gen, alpha, n_alpha)
c_alpha(:,1) += c_alpha(:,i) c_alpha(:,1) += c_alpha(:,i)
end do end do
delta_ij_tmp(:,:,1) -= delta_ij_loc(:,:,1,1)
delta_ij_tmp(:,:,2) -= delta_ij_loc(:,:,2,1) delta_ij_tmp(:,:,1) -= delta_ij_loc(:,:,1,1)
delta_ij_tmp(:,:,2) -= delta_ij_loc(:,:,2,1)
!print *, "SUM ALPHA2 PRE", global_sum_alpha2 !print *, "SUM ALPHA2 PRE", global_sum_alpha2
!global_sum_alpha2(:) -= c_alpha(:,1) !global_sum_alpha2(:) -= c_alpha(:,1)
print *, "SUM C_ALPHA^2 ", global_sum_alpha2(:) print *, "SUM C_ALPHA^2 =", global_sum_alpha2(:)
print *, "*** DRESSINS DIVIDED BY 1+SUM C_ALPHA^2 ***" !print *, "*** DRESSINS DIVIDED BY 1+SUM C_ALPHA^2 ***"
do i=1,N_states !do i=1,N_states
delta_ij_tmp(i,:,:) = delta_ij_tmp(i,:,:) / (1d0 + global_sum_alpha2(i)) ! delta_ij_tmp(i,:,:) = delta_ij_tmp(i,:,:) / (1d0 + global_sum_alpha2(i))
end do !end do
global_sum_alpha2 = 0d0 global_sum_alpha2 = 0d0
end subroutine end subroutine