10
0
mirror of https://github.com/LCPQ/quantum_package synced 2024-11-04 05:03:54 +01:00

Simplified PT2 stoch

This commit is contained in:
Anthony Scemama 2017-05-09 23:11:45 +02:00
parent cbafcb5f55
commit ff20894479
4 changed files with 80 additions and 106 deletions

View File

@ -111,7 +111,7 @@ subroutine push_selection_results(zmq_socket_push, pt2, b, task_id, ntask)
if(rc /= 4*ntask) stop "push" if(rc /= 4*ntask) stop "push"
! Activate is zmq_socket_push is a REQ ! Activate is zmq_socket_push is a REQ
rc = f77_zmq_recv( zmq_socket_push, task_id(1), ntask*4, 0) ! rc = f77_zmq_recv( zmq_socket_push, task_id(1), ntask*4, 0)
end subroutine end subroutine
@ -145,7 +145,7 @@ subroutine pull_selection_results(zmq_socket_pull, pt2, val, det, N, task_id, nt
if(rc /= 4*ntask) stop "pull" if(rc /= 4*ntask) stop "pull"
! Activate is zmq_socket_pull is a REP ! Activate is zmq_socket_pull is a REP
rc = f77_zmq_send( zmq_socket_pull, task_id(1), ntask*4, 0) ! rc = f77_zmq_send( zmq_socket_pull, task_id(1), ntask*4, 0)
end subroutine end subroutine

View File

@ -24,8 +24,8 @@ subroutine run
E_CI_before = pt2_E0_denominator(1) + nuclear_repulsion E_CI_before = pt2_E0_denominator(1) + nuclear_repulsion
threshold_selectors = 1.d0 threshold_selectors = 1.d0
threshold_generators = 1d0 threshold_generators = 1d0
! relative_error = 1.d-3 relative_error = 1.d-3
relative_error = 1.d-8 ! relative_error = 1.d-8
call ZMQ_pt2(E_CI_before, pt2, relative_error) call ZMQ_pt2(E_CI_before, pt2, relative_error)
print *, 'Final step' print *, 'Final step'
print *, 'N_det = ', N_det print *, 'N_det = ', N_det

View File

@ -47,31 +47,16 @@ subroutine ZMQ_pt2(E, pt2,relative_error)
print *, '========== ================= ================= =================' print *, '========== ================= ================= ================='
print *, ' Samples Energy Stat. Error Seconds ' print *, ' Samples Energy Stat. Error Seconds '
print *, '========== ================= ================= =================' print *, '========== ================= ================= ================='
do while(.true.)
call new_parallel_job(zmq_to_qp_run_socket,'pt2') call new_parallel_job(zmq_to_qp_run_socket,'pt2')
call zmq_put_psi(zmq_to_qp_run_socket,1,pt2_e0_denominator,size(pt2_e0_denominator)) call zmq_put_psi(zmq_to_qp_run_socket,1,pt2_e0_denominator,size(pt2_e0_denominator))
call create_selection_buffer(1, 1*2, b) call create_selection_buffer(1, 1*2, b)
Ncomb=size(comb) Ncomb=size(comb)
! i=N_det_generators
! do while (tbc(0) < i)
call get_carlo_workbatch(computed, comb, Ncomb, tbc) call get_carlo_workbatch(computed, comb, Ncomb, tbc)
! i=0
! do j=1,N_det_generators
! if (.not.computed(j)) then
! i = i+1
! endif
! enddo
! i = i/2
! enddo
integer(ZMQ_PTR), external :: new_zmq_to_qp_run_socket integer(ZMQ_PTR), external :: new_zmq_to_qp_run_socket
integer :: ipos integer :: ipos
logical :: tasks
tasks = .False.
ipos=1 ipos=1
do i=1,tbc(0) do i=1,tbc(0)
@ -81,7 +66,6 @@ subroutine ZMQ_pt2(E, pt2,relative_error)
if (ipos > 63980) then if (ipos > 63980) then
call add_task_to_taskserver(zmq_to_qp_run_socket,trim(task(1:ipos))) call add_task_to_taskserver(zmq_to_qp_run_socket,trim(task(1:ipos)))
ipos=1 ipos=1
tasks = .True.
endif endif
else else
do j=1,fragment_count do j=1,fragment_count
@ -90,17 +74,14 @@ subroutine ZMQ_pt2(E, pt2,relative_error)
if (ipos > 63980) then if (ipos > 63980) then
call add_task_to_taskserver(zmq_to_qp_run_socket,trim(task(1:ipos))) call add_task_to_taskserver(zmq_to_qp_run_socket,trim(task(1:ipos)))
ipos=1 ipos=1
tasks = .True.
endif endif
end do end do
end if end if
end do end do
if (ipos > 1) then if (ipos > 1) then
call add_task_to_taskserver(zmq_to_qp_run_socket,trim(task(1:ipos))) call add_task_to_taskserver(zmq_to_qp_run_socket,trim(task(1:ipos)))
tasks = .True.
endif endif
if (tasks) then
call zmq_set_running(zmq_to_qp_run_socket) call zmq_set_running(zmq_to_qp_run_socket)
!$OMP PARALLEL DEFAULT(shared) NUM_THREADS(nproc+1) & !$OMP PARALLEL DEFAULT(shared) NUM_THREADS(nproc+1) &
@ -115,20 +96,6 @@ subroutine ZMQ_pt2(E, pt2,relative_error)
call delete_selection_buffer(b) call delete_selection_buffer(b)
call end_parallel_job(zmq_to_qp_run_socket, 'pt2') call end_parallel_job(zmq_to_qp_run_socket, 'pt2')
else
pt2 = 0.d0
do i=1,N_det_generators
do k=1,N_states
pt2(k) = pt2(k) + pt2_detail(k,i)
enddo
enddo
endif
tbc(0) = 0
if (pt2(1) /= 0.d0) then
exit
endif
end do
print *, '========== ================= ================= =================' print *, '========== ================= ================= ================='
deallocate(pt2_detail, comb, computed, tbc) deallocate(pt2_detail, comb, computed, tbc)
@ -208,6 +175,7 @@ subroutine pt2_collector(E, b, tbc, comb, Ncomb, computed, pt2_detail, sumabove,
integer :: tooth, firstTBDcomb, orgTBDcomb integer :: tooth, firstTBDcomb, orgTBDcomb
integer, allocatable :: parts_to_get(:) integer, allocatable :: parts_to_get(:)
logical, allocatable :: actually_computed(:) logical, allocatable :: actually_computed(:)
character*(512) :: task
Nabove_old = -1.d0 Nabove_old = -1.d0
allocate(actually_computed(N_det_generators), parts_to_get(N_det_generators), & allocate(actually_computed(N_det_generators), parts_to_get(N_det_generators), &
@ -295,7 +263,19 @@ subroutine pt2_collector(E, b, tbc, comb, Ncomb, computed, pt2_detail, sumabove,
eqt = sqrt(1d0 / (Nabove(tooth)-1) * abs(sum2above(tooth) / Nabove(tooth) - (sumabove(tooth)/Nabove(tooth))**2)) eqt = sqrt(1d0 / (Nabove(tooth)-1) * abs(sum2above(tooth) / Nabove(tooth) - (sumabove(tooth)/Nabove(tooth))**2))
call wall_time(time) call wall_time(time)
if (dabs(eqt/avg) < relative_error) then if (dabs(eqt/avg) < relative_error) then
! Termination
pt2(1) = avg pt2(1) = avg
print '(G10.3, 2X, F16.10, 2X, G16.3, 2X, F16.4, A20)', Nabove(tooth), avg+E, eqt, time-time0, ''
integer :: worker_id
call connect_to_taskserver(zmq_to_qp_run_socket,worker_id,0)
if(worker_id /= -1) then
do
call get_task_from_taskserver(zmq_to_qp_run_socket,worker_id, task_id(1), task)
if (task_id(1) == 0) exit
call task_done_to_taskserver(zmq_to_qp_run_socket,worker_id,task_id(1))
call zmq_delete_task(zmq_to_qp_run_socket,zmq_socket_pull,task_id(1),more)
enddo
end if
else else
if (Nabove(tooth) > Nabove_old) then if (Nabove(tooth) > Nabove_old) then
print '(G10.3, 2X, F16.10, 2X, G16.3, 2X, F16.4, A20)', Nabove(tooth), avg+E, eqt, time-time0, '' print '(G10.3, 2X, F16.10, 2X, G16.3, 2X, F16.4, A20)', Nabove(tooth), avg+E, eqt, time-time0, ''
@ -397,7 +377,7 @@ BEGIN_PROVIDER [ integer, size_tbc ]
BEGIN_DOC BEGIN_DOC
! Size of the tbc array ! Size of the tbc array
END_DOC END_DOC
size_tbc = N_det_generators + fragment_count*fragment_first size_tbc = 2*N_det_generators + fragment_count*fragment_first
END_PROVIDER END_PROVIDER
subroutine get_carlo_workbatch(computed, comb, Ncomb, tbc) subroutine get_carlo_workbatch(computed, comb, Ncomb, tbc)
@ -406,30 +386,24 @@ subroutine get_carlo_workbatch(computed, comb, Ncomb, tbc)
double precision, intent(out) :: comb(Ncomb) double precision, intent(out) :: comb(Ncomb)
integer, intent(inout) :: tbc(0:size_tbc) integer, intent(inout) :: tbc(0:size_tbc)
logical, intent(inout) :: computed(N_det_generators) logical, intent(inout) :: computed(N_det_generators)
integer :: i, j, last_full, dets(comb_teeth), tbc_save integer :: i, j, last_full, dets(comb_teeth)
integer :: icount, n integer :: icount, n
integer :: k, l integer :: k, l
l=1 l=1
call RANDOM_NUMBER(comb) call RANDOM_NUMBER(comb)
do i=1,size(comb) do i=1,size(comb)
comb(i) = comb(i) * comb_step comb(i) = comb(i) * comb_step
tbc_save = tbc(0)
!DIR$ FORCEINLINE !DIR$ FORCEINLINE
call add_comb(comb(i), computed, tbc, size_tbc, comb_teeth) call add_comb(comb(i), computed, tbc, size_tbc, comb_teeth)
if ( (tbc(0) < size(tbc)-1).and.(l < first_det_of_teeth(comb_teeth)) ) then
Ncomb = i Ncomb = i
do while (computed(l)) do while (computed(l))
l=l+1 l=l+1
if (l == size(computed)) exit if (l == N_det_generators+1) return
enddo enddo
k=tbc(0)+1 k=tbc(0)+1
tbc(k) = l tbc(k) = l
computed(l) = .True. computed(l) = .True.
tbc(0) = k tbc(0) = k
else
tbc(0) = tbc_save
return
endif
enddo enddo
end subroutine end subroutine

View File

@ -316,12 +316,12 @@ subroutine push_mrsc2_results(zmq_socket_push, I_i, J, delta, delta_s2, task_id)
endif endif
! Activate is zmq_socket_push is a REQ ! Activate is zmq_socket_push is a REQ
integer :: idummy ! integer :: idummy
rc = f77_zmq_recv( zmq_socket_push, idummy, 4, 0) ! rc = f77_zmq_recv( zmq_socket_push, idummy, 4, 0)
if (rc /= 4) then ! if (rc /= 4) then
print *, irp_here, 'f77_zmq_send( zmq_socket_push, idummy, 4, 0)' ! print *, irp_here, 'f77_zmq_send( zmq_socket_push, idummy, 4, 0)'
stop 'error' ! stop 'error'
endif ! endif
end end
@ -390,12 +390,12 @@ subroutine pull_mrsc2_results(zmq_socket_pull, I_i, J, n, idx, delta, delta_s2,
! Activate is zmq_socket_pull is a REP ! Activate is zmq_socket_pull is a REP
integer :: idummy ! integer :: idummy
rc = f77_zmq_send( zmq_socket_pull, idummy, 4, 0) ! rc = f77_zmq_send( zmq_socket_pull, idummy, 4, 0)
if (rc /= 4) then ! if (rc /= 4) then
print *, irp_here, 'f77_zmq_send( zmq_socket_pull, idummy, 4, 0)' ! print *, irp_here, 'f77_zmq_send( zmq_socket_pull, idummy, 4, 0)'
stop 'error' ! stop 'error'
endif ! endif
end end