10
0
mirror of https://github.com/LCPQ/quantum_package synced 2024-07-03 09:55:59 +02:00

minilisted dress_zmq - mrcc_sto as unittest

This commit is contained in:
Yann Garniron 2018-02-14 10:33:11 +01:00
parent 66f7019ad1
commit eacc63624c
9 changed files with 285 additions and 99 deletions

View File

@ -43,7 +43,7 @@ subroutine alpha_callback(delta_ij_loc, i_generator, subset)
use bitmasks
implicit none
integer, intent(in) :: i_generator, subset
double precision,intent(inout) :: delta_ij_loc(N_states,N_det_non_ref,2)
double precision,intent(inout) :: delta_ij_loc(N_states,N_det,2)
integer :: k,l
@ -63,7 +63,7 @@ subroutine generate_singles_and_doubles(delta_ij_loc, i_generator, bitmask_index
! WARNING /!\ : It is assumed that the generators and selectors are psi_det_sorted
END_DOC
double precision,intent(inout) :: delta_ij_loc(N_states,N_det_non_ref,2)
double precision,intent(inout) :: delta_ij_loc(N_states,N_det,2)
integer, intent(in) :: i_generator, subset, bitmask_index
integer :: h1,h2,s1,s2,s3,i1,i2,ib,sp,k,i,j,nt,ii,n
@ -587,6 +587,115 @@ end
subroutine count_pq(mask, sp, det, i_gen, N_sel, bannedOrb, banned, mat, interesting)
use bitmasks
implicit none
integer, intent(in) :: sp, i_gen, N_sel
integer, intent(in) :: interesting(0:N_sel)
integer(bit_kind),intent(in) :: mask(N_int, 2), det(N_int, 2, N_sel)
logical, intent(inout) :: bannedOrb(mo_tot_num, 2), banned(mo_tot_num, mo_tot_num, 2)
integer, intent(inout) :: mat(mo_tot_num, mo_tot_num)
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 :: phasemask(2,N_int*bit_kind_size)
PROVIDE psi_selectors_coef_transp
mat = 0
do i=1,N_int
negMask(i,1) = not(mask(i,1))
negMask(i,2) = not(mask(i,2))
end do
do i=1, N_sel ! interesting(0)
!i = interesting(ii)
if (interesting(i) < 0) then
stop 'prefetch interesting(i)'
endif
mobMask(1,1) = iand(negMask(1,1), det(1,1,i))
mobMask(1,2) = iand(negMask(1,2), det(1,2,i))
nt = popcnt(mobMask(1, 1)) + popcnt(mobMask(1, 2))
if(nt > 4) cycle
do j=2,N_int
mobMask(j,1) = iand(negMask(j,1), det(j,1,i))
mobMask(j,2) = iand(negMask(j,2), det(j,2,i))
nt = nt + popcnt(mobMask(j, 1)) + popcnt(mobMask(j, 2))
end do
if(nt > 4) cycle
if (interesting(i) == i_gen) then
if(sp == 3) then
do j=1,mo_tot_num
do k=1,mo_tot_num
banned(j,k,2) = banned(k,j,1)
enddo
enddo
else
do k=1,mo_tot_num
do l=k+1,mo_tot_num
banned(l,k,1) = banned(k,l,1)
end do
end do
end if
end if
call bitstring_to_list_in_selection(mobMask(1,1), p(1,1), p(0,1), N_int)
call bitstring_to_list_in_selection(mobMask(1,2), p(1,2), p(0,2), N_int)
perMask(1,1) = iand(mask(1,1), not(det(1,1,i)))
perMask(1,2) = iand(mask(1,2), not(det(1,2,i)))
do j=2,N_int
perMask(j,1) = iand(mask(j,1), not(det(j,1,i)))
perMask(j,2) = iand(mask(j,2), not(det(j,2,i)))
end do
call bitstring_to_list_in_selection(perMask(1,1), h(1,1), h(0,1), N_int)
call bitstring_to_list_in_selection(perMask(1,2), h(1,2), h(0,2), N_int)
if (interesting(i) >= i_gen) then
call get_mask_phase(psi_det_sorted(1,1,interesting(i)), phasemask)
if(nt == 4) then
call count_d2(mat, p, sp)
else if(nt == 3) then
call count_d1(mat, p, sp)
else
mat(:,:) = mat(:,:) + 1
end if
else
if(nt == 4) call past_d2(banned, p, sp)
if(nt == 3) call past_d1(bannedOrb, p)
end if
end do
do i=1,mo_tot_num
do j=1,mo_tot_num
if(banned(i,j,1)) mat(i,j) = 0
end do
end do
if(sp == 3) then
do i=1,mo_tot_num
if(bannedOrb(i, 1)) mat(i, :) = 0
if(bannedOrb(i, 2)) mat(:, i) = 0
end do
else
do i=1,mo_tot_num
if(bannedOrb(i, sp)) then
mat(:,i) = 0
mat(i,:) = 0
end if
end do
end if
end
subroutine splash_pq(mask, sp, det, i_gen, N_sel, bannedOrb, banned, mat, interesting)
use bitmasks
@ -1040,6 +1149,65 @@ subroutine past_d1(bannedOrb, p)
end
subroutine count_d1(mat, p, sp)
use bitmasks
implicit none
integer, intent(inout) :: mat(mo_tot_num, mo_tot_num)
integer, intent(in) :: p(0:4, 2), sp
integer :: i,s,j
if(sp == 3) then
do i=1,p(0,1)
mat(p(i,1), :) += 1
end do
do i=1,p(0,2)
mat(:, p(i,2)) += 1
end do
do i=1,p(0,1)
do j=1,p(0,2)
mat(p(i,1), p(j,2)) -= 1
end do
end do
else
if(sp == 1 .and. p(0,2) /= 0) stop "count_d1 bug"
if(sp == 2 .and. p(0,1) /= 0) stop "count_d1 bug"
do i=1,p(0,sp)
mat(:p(i,sp), p(i,sp)) += 1
mat(p(i,sp), p(i,sp):) += 1
mat(p(i,sp), p(i,sp)) -= 1
end do
end if
end
subroutine count_d2(mat, p, sp)
use bitmasks
implicit none
integer, intent(inout) :: mat(mo_tot_num, mo_tot_num)
integer, intent(in) :: p(0:4, 2), sp
integer :: i,j
if(sp == 3) then
do i=1,p(0,1)
do j=1,p(0,2)
mat(p(i,1), p(j,2)) += 1
end do
end do
else
do i=1,p(0, sp)
do j=1,i-1
mat(p(j,sp), p(i,sp)) += 1
mat(p(i,sp), p(j,sp)) += 1
end do
end do
end if
end
subroutine past_d2(banned, p, sp)
use bitmasks
implicit none

View File

@ -15,8 +15,8 @@ subroutine ZMQ_dress(E, dress, delta, delta_s2, relative_error)
integer, external :: omp_get_thread_num
double precision, intent(in) :: relative_error, E
double precision, intent(out) :: dress(N_states)
double precision, intent(out) :: delta(N_states, N_det_non_ref)
double precision, intent(out) :: delta_s2(N_states, N_det_non_ref)
double precision, intent(out) :: delta(N_states, N_det)
double precision, intent(out) :: delta_s2(N_states, N_det)
integer :: i, j, k, Ncp
@ -32,7 +32,7 @@ subroutine ZMQ_dress(E, dress, delta, delta_s2, relative_error)
!!!!!!!!!!!!!!! demander a TOTO !!!!!!!
w(:) = 0.d0
w(dress_stoch_istate) = 1.d0
call update_psi_average_norm_contrib(w)
!call update_psi_average_norm_contrib(w)
@ -135,8 +135,8 @@ subroutine dress_collector(zmq_socket_pull, E, relative_error, delta, delta_s2,
double precision, intent(out) :: dress(N_states)
double precision, allocatable :: cp(:,:,:,:)
double precision, intent(out) :: delta(N_states, N_det_non_ref)
double precision, intent(out) :: delta_s2(N_states, N_det_non_ref)
double precision, intent(out) :: delta(N_states, N_det)
double precision, intent(out) :: delta_s2(N_states, N_det)
double precision, allocatable :: delta_loc(:,:,:), delta_det(:,:,:,:)
double precision, allocatable :: dress_detail(:,:)
double precision :: dress_mwen(N_states)
@ -158,9 +158,9 @@ subroutine dress_collector(zmq_socket_pull, E, relative_error, delta, delta_s2,
delta = 0d0
delta_s2 = 0d0
allocate(delta_det(N_states, N_det_non_ref, 0:comb_teeth+1, 2))
allocate(cp(N_states, N_det_non_ref, N_cp, 2), dress_detail(N_states, N_det_generators))
allocate(delta_loc(N_states, N_det_non_ref, 2))
allocate(delta_det(N_states, N_det, 0:comb_teeth+1, 2))
allocate(cp(N_states, N_det, N_cp, 2), dress_detail(N_states, N_det))
allocate(delta_loc(N_states, N_det, 2))
dress_detail = 0d0
delta_det = 0d0
cp = 0d0
@ -192,8 +192,8 @@ subroutine dress_collector(zmq_socket_pull, E, relative_error, delta, delta_s2,
loop = .true.
pullLoop : do while (loop)
call pull_dress_results(zmq_socket_pull, ind, dress_mwen, delta_loc, task_id)
call pull_dress_results(zmq_socket_pull, ind, delta_loc, task_id)
dress_mwen(:) = 0d0 !!!!!!!! A CALCULER ICI
dress_detail(:, ind) += dress_mwen(:)
do j=1,N_cp !! optimizable
if(cps(ind, j) > 0d0) then
@ -202,7 +202,7 @@ subroutine dress_collector(zmq_socket_pull, E, relative_error, delta, delta_s2,
integer :: toothMwen
logical :: fracted
fac = cps(ind, j) / cps_N(j) * dress_weight_inv(ind) * comb_step
do k=1,N_det_non_ref
do k=1,N_det
do i_state=1,N_states
cp(i_state,k,j,1) += delta_loc(i_state,k,1) * fac
cp(i_state,k,j,2) += delta_loc(i_state,k,2) * fac
@ -499,7 +499,7 @@ subroutine add_comb(com, computed, cp, N, tbc)
implicit none
double precision, intent(in) :: com
integer, intent(inout) :: N
double precision, intent(inout) :: cp(N_det_non_ref)
double precision, intent(inout) :: cp(N_det)
logical, intent(inout) :: computed(N_det_generators)
integer, intent(inout) :: tbc(N_det_generators)
integer :: i, k, l, dets(comb_teeth)

View File

@ -64,10 +64,7 @@ END_PROVIDER
BEGIN_PROVIDER [ double precision, delta_ij, (N_states,N_det_non_ref, N_det_ref) ]
&BEGIN_PROVIDER [ double precision, delta_ii, (N_states, N_det_ref) ]
&BEGIN_PROVIDER [ double precision, delta_ij_s2, (N_states,N_det_non_ref, N_det_ref) ]
&BEGIN_PROVIDER [ double precision, delta_ii_s2, (N_states, N_det_ref) ]
BEGIN_PROVIDER [ double precision, delta_ij, (N_states,N_det,2) ]
use bitmasks
implicit none
@ -77,18 +74,15 @@ END_PROVIDER
double precision :: E_CI_before, relative_error
double precision, save :: errr = 0d0
allocate(dress(N_states), del(N_states, N_det_non_ref), del_s2(N_states, N_det_non_ref))
allocate(dress(N_states), del(N_states, N_det), del_s2(N_states, N_det))
delta_ij = 0d0
delta_ii = 0d0
delta_ij_s2 = 0d0
delta_ii_s2 = 0d0
E_CI_before = dress_E0_denominator(1) + nuclear_repulsion
threshold_selectors = 1.d0
threshold_generators = 1d0
if(errr /= 0d0) then
errr = errr / 2d0 !
errr = errr / 2d0
else
errr = 1d-4
end if
@ -97,11 +91,12 @@ END_PROVIDER
call ZMQ_dress(E_CI_before, dress, del, del_s2, abs(relative_error))
delta_ij(:,:,1) = del(:,:)
delta_ij_s2(:,:,1) = del_s2(:,:)
do i=N_det_non_ref,1,-1
delta_ii(dress_stoch_istate,1) -= delta_ij(dress_stoch_istate, i, 1) / psi_ref_coef(1,dress_stoch_istate) * psi_non_ref_coef(i, dress_stoch_istate)
delta_ii_s2(dress_stoch_istate,1) -= delta_ij_s2(dress_stoch_istate, i, 1) / psi_ref_coef(1,dress_stoch_istate) * psi_non_ref_coef(i, dress_stoch_istate)
end do
!delta_ij_s2(:,:,1) = del_s2(:,:)
delta_ij(:,:,2) = del_s2(:,:)
!do i=N_det,1,-1
! delta_ii(dress_stoch_istate,1) -= delta_ij(dress_stoch_istate, i, 1) / psi_ref_coef(1,dress_stoch_istate) * psi_coef(i, dress_stoch_istate)
! delta_ii_s2(dress_stoch_istate,1) -= delta_ij_s2(dress_stoch_istate, i, 1) / psi_ref_coef(1,dress_stoch_istate) * psi_coef(i, dress_stoch_istate)
!end do
END_PROVIDER

View File

@ -0,0 +1,29 @@
BEGIN_PROVIDER [ double precision, dressing_column_h, (N_det,N_states) ]
&BEGIN_PROVIDER [ double precision, dressing_column_s, (N_det,N_states) ]
implicit none
BEGIN_DOC
! Null dressing vectors
END_DOC
dressing_column_h(:,:) = 0.d0
dressing_column_s(:,:) = 0.d0
integer :: i,ii,k,j,jj, l
double precision :: f, tmp
double precision, external :: u_dot_v
do k=1,N_states
l = dressed_column_idx(k)
f = 1.d0/psi_coef(l,k)
do jj = 1, n_det
j = jj !idx_non_ref(jj)
dressing_column_h(j,k) = delta_ij (k,jj,1)
dressing_column_s(j,k) = delta_ij (k,jj,2)!delta_ij_s2(k,jj)
enddo
tmp = u_dot_v(dressing_column_h(1,k), psi_coef(1,k), N_det)
dressing_column_h(l,k) -= tmp * f
tmp = u_dot_v(dressing_column_s(1,k), psi_coef(1,k), N_det)
dressing_column_s(l,k) -= tmp * f
enddo
END_PROVIDER

View File

@ -29,19 +29,11 @@ subroutine run_dress_slave(thread,iproc,energy)
double precision,allocatable :: dress_detail(:)
integer :: ind
integer(bit_kind),allocatable :: abuf(:,:,:)
integer(bit_kind) :: mask(N_int,2), omask(N_int,2)
double precision,allocatable :: delta_ij_loc(:,:,:)
double precision,allocatable :: delta_ii_loc(:,:)
integer :: h,p,n
logical :: ok
double precision :: contrib(N_states)
allocate(delta_ij_loc(N_states,N_det_non_ref,2) &
,delta_ii_loc(N_states,2))
allocate(abuf(N_int, 2, N_det_non_ref))
allocate(dress_detail(N_states))
allocate(delta_ij_loc(N_states,N_det,2))
zmq_to_qp_run_socket = new_zmq_to_qp_run_socket()
zmq_socket_push = new_zmq_push_socket(thread)
@ -52,47 +44,15 @@ subroutine run_dress_slave(thread,iproc,energy)
call end_zmq_push_socket(zmq_socket_push,thread)
return
end if
dress_detail = 0d0
do
call get_task_from_taskserver(zmq_to_qp_run_socket,worker_id, task_id, task)
if(task_id /= 0) then
read (task,*) subset, i_generator
contrib = 0d0
delta_ij_loc = 0d0
delta_ii_loc = 0d0
if(do_dress_with_alpha_buffer .or. do_dress_with_alpha) then
do h=1, hh_shortcut(0)
call apply_hole_local(psi_det_generators(1,1,i_generator), hh_exists(1, h), mask, ok, N_int)
if(.not. ok) cycle
omask = 0_bit_kind
if(hh_exists(1, h) /= 0) omask = mask
n = 1
do p=hh_shortcut(h), hh_shortcut(h+1)-1
call apply_particle_local(mask, pp_exists(1, p), abuf(1,1,n), ok, N_int)
if(ok) n = n + 1
if(n > N_det_non_ref) stop "Buffer too small in dress..."
end do
n = n - 1
if(n /= 0) then
if(do_dress_with_alpha_buffer) then
call dress_with_alpha_buffer(delta_ij_loc(1,1,1), delta_ii_loc(1,1), delta_ij_loc(1,1,2), delta_ii_loc(1,2), &
i_generator,n,abuf,N_int,omask,contrib)
else
stop 'dress_with_alpha not implemented yet'
end if
endif
end do
else if(do_dress_with_generator) then
stop 'dress_with_generator not implemented yet'
else
stop 'no dressing level defined'
end if
dress_detail(:) = contrib
call alpha_callback(delta_ij_loc, i_generator, subset)
call task_done_to_taskserver(zmq_to_qp_run_socket,worker_id,task_id)
call push_dress_results(zmq_socket_push, i_generator, dress_detail, delta_ij_loc(1,1,1), task_id)
dress_detail(:) = 0d0
call push_dress_results(zmq_socket_push, i_generator, delta_ij_loc, task_id)
else
exit
end if
@ -103,13 +63,12 @@ subroutine run_dress_slave(thread,iproc,energy)
end subroutine
subroutine push_dress_results(zmq_socket_push, ind, dress_detail, delta_loc, task_id)
subroutine push_dress_results(zmq_socket_push, ind, delta_loc, task_id)
use f77_zmq
implicit none
integer(ZMQ_PTR), intent(in) :: zmq_socket_push
double precision, intent(in) :: dress_detail(N_states, N_det_generators)
double precision, intent(in) :: delta_loc(N_states, N_det_non_ref, 2)
double precision, intent(in) :: delta_loc(N_states, N_det, 2)
integer, intent(in) :: ind, task_id
integer :: rc, i
@ -118,16 +77,9 @@ subroutine push_dress_results(zmq_socket_push, ind, dress_detail, delta_loc, tas
if(rc /= 4) stop "push"
rc = f77_zmq_send( zmq_socket_push, dress_detail, 8*N_states, ZMQ_SNDMORE)
if(rc /= 8*N_states) stop "push"
rc = f77_zmq_send( zmq_socket_push, delta_loc(1,1,1), 8*N_states*N_det_non_ref, ZMQ_SNDMORE)
if(rc /= 8*N_states*N_det_non_ref) stop "push"
rc = f77_zmq_send( zmq_socket_push, delta_loc, 8*N_states*N_det*2, ZMQ_SNDMORE)
if(rc /= 8*N_states*N_det*2) stop "push"
rc = f77_zmq_send( zmq_socket_push, delta_loc(1,1,2), 8*N_states*N_det_non_ref, ZMQ_SNDMORE)
if(rc /= 8*N_states*N_det_non_ref) stop "push"
rc = f77_zmq_send( zmq_socket_push, task_id, 4, 0)
if(rc /= 4) stop "push"
@ -141,30 +93,21 @@ IRP_ENDIF
end subroutine
subroutine pull_dress_results(zmq_socket_pull, ind, dress_detail, delta_loc, task_id)
subroutine pull_dress_results(zmq_socket_pull, ind, delta_loc, task_id)
use f77_zmq
implicit none
integer(ZMQ_PTR), intent(in) :: zmq_socket_pull
double precision, intent(inout) :: dress_detail(N_states)
double precision, intent(inout) :: delta_loc(N_states, N_det_non_ref, 2)
double precision, allocatable :: dress_dress_mwen(:,:)
double precision, intent(inout) :: delta_loc(N_states, N_det, 2)
integer, intent(out) :: ind
integer, intent(out) :: task_id
integer :: rc, rn, i
integer :: rc, i
allocate(dress_dress_mwen(N_states, N_det_non_ref))
rc = f77_zmq_recv( zmq_socket_pull, ind, 4, 0)
if(rc /= 4) stop "pull"
rc = f77_zmq_recv( zmq_socket_pull, dress_detail, N_states*8, 0)
if(rc /= 8*N_states) stop "pull"
rc = f77_zmq_recv( zmq_socket_pull, delta_loc(1,1,1), N_states*8*N_det_non_ref, 0)
if(rc /= 8*N_states*N_det_non_ref) stop "pull"
rc = f77_zmq_recv( zmq_socket_pull, delta_loc(1,1,2), N_states*8*N_det_non_ref, 0)
if(rc /= 8*N_states*N_det_non_ref) stop "pull"
rc = f77_zmq_recv( zmq_socket_pull, delta_loc, N_states*8*N_det*2, 0)
if(rc /= 8*N_states*N_det*2) stop "pull"
rc = f77_zmq_recv( zmq_socket_pull, task_id, 4, 0)
if(rc /= 4) stop "pull"

View File

@ -0,0 +1,37 @@
! BEGIN_PROVIDER [ logical, do_dress_with_alpha ]
!&BEGIN_PROVIDER [ logical, do_dress_with_alpha_buffer ]
!&BEGIN_PROVIDER [ logical, do_dress_with_generator ]
! implicit none
! do_dress_with_alpha = .false.
! do_dress_with_alpha_buffer = .true.
! do_dress_with_generator = .false.
!END_PROVIDER
subroutine dress_with_alpha_buffer(delta_ij_loc, minilist, n_minilist, abuf, n_abuf)
use bitmasks
implicit none
double precision,intent(inout) :: delta_ij_loc(N_states,N_det,2)
integer, intent(in) :: n_minilist, n_abuf
integer(bit_kind),intent(in) :: abuf(N_int, 2, n_abuf), minilist(N_int, 2, n_minilist)
integer :: a, i, nref, nobt, deg
do a=1,n_abuf
nref=0
do i=1,N_det
call get_excitation_degree(psi_det(1,1,i), abuf(1,1,a), deg, N_int)
if(deg <= 2) nref = nref + 1
end do
nobt=0
do i=1,n_minilist
call get_excitation_degree(minilist(1,1,i), abuf(1,1,a), deg, N_int)
if(deg <= 2) nobt = nobt + 1
end do
if(nref /= nobt) stop "foireous minilist"
end do
delta_ij_loc = 1d0
end subroutine

View File

@ -0,0 +1,14 @@
program mrcc_sto
implicit none
BEGIN_DOC
! TODO
END_DOC
print *, "!!!!!!========================!!!!!!"
print *, "!!!!!!========================!!!!!!"
print *, "!!!!!!========================!!!!!!"
print *, "MRCC_STO not implemented - acts as a unittest for dress_zmq"
print *, "!!!!!!========================!!!!!!"
print *, "!!!!!!========================!!!!!!"
print *, "!!!!!!========================!!!!!!"
call dress_zmq()
end

View File

@ -1 +1 @@
Perturbation Selectors_full Generators_full Psiref_CAS MRCC_Utils ZMQ DavidsonDressed
DavidsonDressed Perturbation Selectors_full Generators_full Psiref_CAS MRCC_Utils ZMQ

View File

@ -44,7 +44,7 @@ subroutine ZMQ_mrcc(E, mrcc, delta, delta_s2, relative_error)
w = 0.d0
w = 1.d0
call update_psi_average_norm_contrib(w)
!call update_psi_average_norm_contrib(w)