10
0
mirror of https://github.com/LCPQ/quantum_package synced 2024-12-25 05:43:47 +01:00

corrected bitmask in Full_CI_ZMQ - seemingly working shiftedbk

This commit is contained in:
Yann Garniron 2018-03-15 17:47:32 +01:00
parent b29fb01e6f
commit ade11beb25
8 changed files with 57 additions and 36 deletions

View File

@ -19,6 +19,7 @@ subroutine run
double precision :: E_CI_before, relative_error, absolute_error, eqt double precision :: E_CI_before, relative_error, absolute_error, eqt
allocate (pt2(N_states)) allocate (pt2(N_states))
call diagonalize_CI()
pt2 = 0.d0 pt2 = 0.d0
E_CI_before = pt2_E0_denominator(1) + nuclear_repulsion E_CI_before = pt2_E0_denominator(1) + nuclear_repulsion

View File

@ -342,7 +342,7 @@ subroutine pt2_collector(zmq_socket_pull, E, b, tbc, comb, Ncomb, computed, pt2_
endif endif
end if end if
end do pullLoop end do pullLoop
E0 = sum(pt2_detail(pt2_stoch_istate,:first_det_of_teeth(tooth)-1)) E0 = sum(pt2_detail(pt2_stoch_istate,:first_det_of_teeth(tooth)-1))
prop = ((1d0 - dfloat(comb_teeth - tooth + 1) * comb_step) - pt2_cweight(first_det_of_teeth(tooth)-1)) prop = ((1d0 - dfloat(comb_teeth - tooth + 1) * comb_step) - pt2_cweight(first_det_of_teeth(tooth)-1))
prop = prop * pt2_weight_inv(first_det_of_teeth(tooth)) prop = prop * pt2_weight_inv(first_det_of_teeth(tooth))

View File

@ -41,7 +41,7 @@ subroutine run_pt2_slave(thread,iproc,energy)
buf%N = 0 buf%N = 0
n_tasks = 0 n_tasks = 0
call create_selection_buffer(1, 2, buf) call create_selection_buffer(0, 0, buf)
done = .False. done = .False.
do while (.not.done) do while (.not.done)

View File

@ -253,6 +253,7 @@ subroutine get_m0(gen, phasemask, bannedOrb, vect, mask, h, p, sp, coefs)
deallocate(lbanned) deallocate(lbanned)
end end
subroutine select_singles_and_doubles(i_generator,hole_mask,particle_mask,fock_diag_tmp,E0,pt2,buf,subset) subroutine select_singles_and_doubles(i_generator,hole_mask,particle_mask,fock_diag_tmp,E0,pt2,buf,subset)
use bitmasks use bitmasks
use selection_types use selection_types
@ -292,10 +293,17 @@ subroutine select_singles_and_doubles(i_generator,hole_mask,particle_mask,fock_d
do k=1,N_int do k=1,N_int
hole (k,1) = iand(psi_det_generators(k,1,i_generator), hole_mask(k,1)) !if(buf%N > 0) then
hole (k,2) = iand(psi_det_generators(k,2,i_generator), hole_mask(k,2)) ! hole (k,1) = iand(psi_det_generators(k,1,i_generator), hole_mask(k,1))
particle(k,1) = iand(not(psi_det_generators(k,1,i_generator)), particle_mask(k,1)) ! hole (k,2) = iand(psi_det_generators(k,2,i_generator), hole_mask(k,2))
particle(k,2) = iand(not(psi_det_generators(k,2,i_generator)), particle_mask(k,2)) ! particle(k,1) = iand(not(psi_det_generators(k,1,i_generator)), particle_mask(k,1))
! particle(k,2) = iand(not(psi_det_generators(k,2,i_generator)), particle_mask(k,2))
!else
hole (k,1) = iand(psi_det_generators(k,1,i_generator), full_ijkl_bitmask(k))
hole (k,2) = iand(psi_det_generators(k,2,i_generator), full_ijkl_bitmask(k))
particle(k,1) = iand(not(psi_det_generators(k,1,i_generator)), full_ijkl_bitmask(k))
particle(k,2) = iand(not(psi_det_generators(k,2,i_generator)), full_ijkl_bitmask(k))
!end if
enddo enddo
@ -598,7 +606,6 @@ subroutine fill_buffer_double(i_generator, sp, h1, h2, bannedOrb, banned, fock_d
logical, external :: detEq logical, external :: detEq
if(sp == 3) then if(sp == 3) then
s1 = 1 s1 = 1
s2 = 2 s2 = 2
@ -616,15 +623,13 @@ subroutine fill_buffer_double(i_generator, sp, h1, h2, bannedOrb, banned, fock_d
do p2=ib,mo_tot_num do p2=ib,mo_tot_num
if(bannedOrb(p2, s2)) cycle if(bannedOrb(p2, s2)) cycle
if(banned(p1,p2)) cycle if(banned(p1,p2)) cycle
if( sum(abs(mat(1:N_states, p1, p2))) == 0d0) cycle if( sum(abs(mat(1:N_states, p1, p2))) == 0d0) cycle
call apply_particles(mask, s1, p1, s2, p2, det, ok, N_int) call apply_particles(mask, s1, p1, s2, p2, det, ok, N_int)
Hii = diag_H_mat_elem_fock(psi_det_generators(1,1,i_generator),det,fock_diag_tmp,N_int) Hii = diag_H_mat_elem_fock(psi_det_generators(1,1,i_generator),det,fock_diag_tmp,N_int)
min_e_pert = 0d0 min_e_pert = 0d0
! double precision :: hij
! call i_h_j(psi_det_generators(1,1,i_generator), det, N_int, hij)
do istate=1,N_states do istate=1,N_states
delta_E = E0(istate) - Hii delta_E = E0(istate) - Hii
val = mat(istate, p1, p2) + mat(istate, p1, p2) val = mat(istate, p1, p2) + mat(istate, p1, p2)
@ -635,7 +640,6 @@ subroutine fill_buffer_double(i_generator, sp, h1, h2, bannedOrb, banned, fock_d
e_pert = 0.5d0 * (tmp - delta_E) e_pert = 0.5d0 * (tmp - delta_E)
pt2(istate) = pt2(istate) + e_pert pt2(istate) = pt2(istate) + e_pert
min_e_pert = min(e_pert,min_e_pert) min_e_pert = min(e_pert,min_e_pert)
! ci(istate) = e_pert / hij
end do end do
if(min_e_pert <= buf%mini) then if(min_e_pert <= buf%mini) then

View File

@ -39,7 +39,7 @@ subroutine add_to_selection_buffer(b, det, val)
double precision, intent(in) :: val double precision, intent(in) :: val
integer :: i integer :: i
if(val <= b%mini) then if(b%N > 0 .and. val <= b%mini) then
b%cur += 1 b%cur += 1
b%det(1:N_int,1:2,b%cur) = det(1:N_int,1:2) b%det(1:N_int,1:2,b%cur) = det(1:N_int,1:2)
b%val(b%cur) = val b%val(b%cur) = val
@ -119,7 +119,7 @@ subroutine sort_selection_buffer(b)
integer(bit_kind), pointer :: detmp(:,:,:) integer(bit_kind), pointer :: detmp(:,:,:)
integer :: i, nmwen integer :: i, nmwen
logical, external :: detEq logical, external :: detEq
if (b%cur == 0) return if (b%N == 0 .or. b%cur == 0) return
nmwen = min(b%N, b%cur) nmwen = min(b%N, b%cur)
allocate(iorder(b%cur), detmp(N_int, 2, size(b%det,3))) allocate(iorder(b%cur), detmp(N_int, 2, size(b%det,3)))

View File

@ -286,7 +286,7 @@ subroutine dress_collector(zmq_socket_pull, E, relative_error, delta, delta_s2,
call wall_time(time) call wall_time(time)
if ((dabs(eqt) < relative_error .and. cps_N(cur_cp) >= 5) .or. total_computed == N_det_generators) then if ((dabs(eqt) < relative_error .and. cps_N(cur_cp) >= 30) .or. total_computed == N_det_generators) then
! Termination ! Termination
print '(2X, F16.7, 2X, G16.3, 2X, F16.4, A20)', avg+E(istate)+E0, eqt, time-time0, '' print '(2X, F16.7, 2X, G16.3, 2X, F16.4, A20)', avg+E(istate)+E0, eqt, time-time0, ''
if (zmq_abort(zmq_to_qp_run_socket) == -1) then if (zmq_abort(zmq_to_qp_run_socket) == -1) then

View File

@ -86,7 +86,7 @@ BEGIN_PROVIDER [ double precision, delta_ij, (N_states,N_det,2) ]
! else ! else
! errr = 1d-4 ! errr = 1d-4
! end if ! end if
relative_error = 1.d-3 relative_error = 1.d-4
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))

View File

@ -4,14 +4,19 @@ program mrcc_sto
BEGIN_DOC BEGIN_DOC
! TODO ! TODO
END_DOC END_DOC
call diagonalize_CI()
call dress_zmq() call dress_zmq()
end end
BEGIN_PROVIDER [ double precision, fock_diag_tmp_, (2,mo_tot_num+1,Nproc) ] BEGIN_PROVIDER [ double precision, fock_diag_tmp_, (2,mo_tot_num+1,Nproc) ]
&BEGIN_PROVIDER [ integer, current_generator_, (Nproc) ] &BEGIN_PROVIDER [ integer, current_generator_, (Nproc) ]
&BEGIN_PROVIDER [ double precision, a_h_i, (N_det, Nproc) ]
&BEGIN_PROVIDER [ double precision, a_s2_i, (N_det, Nproc) ]
implicit none implicit none
current_generator_(:) = 0 current_generator_(:) = 0
a_h_i = 0d0
a_s2_i = 0d0
END_PROVIDER END_PROVIDER
@ -31,40 +36,48 @@ subroutine dress_with_alpha_buffer(Nstates,Ndet,Nint,delta_ij_loc, i_gen, minili
integer(bit_kind), intent(in) :: alpha(Nint,2), det_minilist(Nint, 2, n_minilist) integer(bit_kind), intent(in) :: alpha(Nint,2), det_minilist(Nint, 2, n_minilist)
integer,intent(in) :: minilist(n_minilist) integer,intent(in) :: minilist(n_minilist)
double precision, intent(inout) :: delta_ij_loc(Nstates,N_det,2) double precision, intent(inout) :: delta_ij_loc(Nstates,N_det,2)
double precision :: hii, hij, sij, delta_e double precision :: haa, hij, sij
double precision, external :: diag_H_mat_elem_fock double precision, external :: diag_H_mat_elem_fock
integer :: i,j,k,l,m, l_sd integer :: i,j,k,l,m, l_sd
double precision, save :: tot = 0d0 double precision :: hdress, sdress
double precision :: de(N_states), val, tmp double precision :: de, a_h_psi(Nstates), c_alpha
a_h_psi = 0d0
if(current_generator_(iproc) /= i_gen) then if(current_generator_(iproc) /= i_gen) then
current_generator_(iproc) = i_gen current_generator_(iproc) = i_gen
call build_fock_tmp(fock_diag_tmp_(1,1,iproc),psi_det_generators(1,1,i_gen),N_int) call build_fock_tmp(fock_diag_tmp_(1,1,iproc),psi_det_generators(1,1,i_gen),N_int)
end if end if
hii = diag_H_mat_elem_fock(psi_det_generators(1,1,i_gen),alpha,fock_diag_tmp_(1,1,iproc),N_int) haa = diag_H_mat_elem_fock(psi_det_generators(1,1,i_gen),alpha,fock_diag_tmp_(1,1,iproc),N_int)
do i=1,N_states
de(i) = (E0_denominator(i) - hii) do l_sd=1,n_minilist
end do call i_h_j_s2(alpha,det_minilist(1,1,l_sd),N_int,hij, sij)
a_h_i(l_sd, iproc) = hij
do i=1,N_states a_s2_i(l_sd, iproc) = sij
val = 0D0 do i=1,Nstates
do l_sd=1,n_minilist a_h_psi(i) += hij * psi_coef(minilist(l_sd), i)
call i_h_j_s2(alpha,det_minilist(1,1,l_sd),N_int,hij, sij) end do
val += hij end do
do i=1,Nstates
de = E0_denominator(i) - haa
if(DABS(de) < 1D-2) cycle
c_alpha = a_h_psi(i) / de
do l_sd=1,n_minilist
hdress = c_alpha * a_h_i(l_sd, iproc)
sdress = c_alpha * a_s2_i(l_sd, iproc)
delta_ij_loc(i, minilist(l_sd), 1) += hdress
delta_ij_loc(i, minilist(l_sd), 2) += sdress
end do end do
val = 2d0 * val
tmp = dsqrt(de(i)**2 + val**2)
if(de(i) < 0d0) tmp = -tmp
delta_ij_loc(i, minilist(l_sd), 1) += 0.5d0 * (tmp - de(i)) ! * psi_coef(minilist(l_sd), i)
end do end do
end subroutine end subroutine
BEGIN_PROVIDER [ logical, initialize_E0_denominator ] BEGIN_PROVIDER [ logical, initialize_E0_denominator ]
implicit none implicit none
BEGIN_DOC BEGIN_DOC
@ -72,7 +85,8 @@ BEGIN_PROVIDER [ logical, initialize_E0_denominator ]
END_DOC END_DOC
initialize_E0_denominator = .True. initialize_E0_denominator = .True.
END_PROVIDER END_PROVIDER
BEGIN_PROVIDER [ double precision, E0_denominator, (N_states) ] BEGIN_PROVIDER [ double precision, E0_denominator, (N_states) ]
implicit none implicit none
BEGIN_DOC BEGIN_DOC
@ -88,3 +102,5 @@ BEGIN_PROVIDER [ double precision, E0_denominator, (N_states) ]
E0_denominator = -huge(1.d0) E0_denominator = -huge(1.d0)
endif endif
END_PROVIDER END_PROVIDER