LCPQ/quantum_package
10
0
Fork 0
mirror of https://github.com/LCPQ/quantum_package synced 2024-06-26 07:02:14 +02:00
Code Issues Releases Wiki Activity

Merge pull request #195 from scemama/master

Browse Source 
Fixed u0Hu0
This commit is contained in:
Thomas Applencourt 2017-05-16 10:15:42 -05:00 committed by GitHub
commit cebf1594af
26 changed files with 1806 additions and 1547 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

2
config/gfortran_debug.cfg
View File

@ -10,7 +10,7 @@
#
#
[COMMON]
FC : gfortran -g -ffree-line-length-none -I . -static-libgcc
FC : gfortran -g -ffree-line-length-none -I .
LAPACK_LIB : -llapack -lblas
IRPF90 : irpf90
IRPF90_FLAGS : --ninja --align=32

6
config/travis.cfg
View File

@ -13,7 +13,7 @@
FC : gfortran -ffree-line-length-none -I . -g
LAPACK_LIB : -llapack -lblas
IRPF90 : irpf90
IRPF90_FLAGS : --ninja --align=32
IRPF90_FLAGS : --ninja --align=32 --assert
# Global options
################
@ -22,7 +22,7 @@ IRPF90_FLAGS : --ninja --align=32
# 0 : Deactivate
#
[OPTION]
MODE : OPT ; [ OPT | PROFILE | DEBUG ] : Chooses the section below
MODE : DEBUG ; [ OPT | PROFILE | DEBUG ] : Chooses the section below
CACHE : 1 ; Enable cache_compile.py
OPENMP : 1 ; Append OpenMP flags
@ -51,7 +51,7 @@ FCFLAGS : -Ofast
# -g : Extra debugging information
#
[DEBUG]
FCFLAGS : -fcheck=all -g
FCFLAGS : -Ofast -fcheck=all -g -Waliasing -Wampersand -Wconversion -Wsurprising -Wintrinsics-std -Wno-tabs -Wintrinsic-shadow -Wline-truncation -Wreal-q-constant
# OpenMP flags
#################

12
plugins/CAS_SD_ZMQ/run_selection_slave.irp.f
View File

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

2
plugins/CAS_SD_ZMQ/selection.irp.f
View File

@ -112,7 +112,7 @@ double precision function get_phase_bi(phasemask, s1, s2, h1, p1, h2, p2)
if(s1 == s2 .and. max(h1, p1) > min(h2, p2)) np = np + 1_1
get_phase_bi = res(iand(np,1_1))
end subroutine
end function

7
plugins/Full_CI_ZMQ/pt2_slave.irp.f
View File

@ -47,7 +47,12 @@ subroutine run_wf
print *, 'PT2'
call zmq_get_psi(zmq_to_qp_run_socket,1,energy,N_states)
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
!$OMP PARALLEL PRIVATE(i)
i = omp_get_thread_num()
call pt2_slave_tcp(i, energy)

3
plugins/Full_CI_ZMQ/pt2_stoch.irp.f
View File

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

41
plugins/Full_CI_ZMQ/pt2_stoch_routines.irp.f
View File

@ -27,7 +27,7 @@ subroutine ZMQ_pt2(E, pt2,relative_error)
double precision, external :: omp_get_wtime
double precision :: time
allocate(pt2_detail(N_states, N_det_generators), comb(N_det_generators), computed(N_det_generators), tbc(0:size_tbc))
allocate(pt2_detail(N_states, N_det_generators+1), comb(N_det_generators), computed(N_det_generators), tbc(0:size_tbc))
sumabove = 0d0
sum2above = 0d0
Nabove = 0d0
@ -211,7 +211,6 @@ subroutine pt2_collector(E, b, tbc, comb, Ncomb, computed, pt2_detail, sumabove,
call get_first_tooth(actually_computed, tooth)
Nabove_old = Nabove(tooth)
! print *, 'N_deterministic = ', first_det_of_teeth(1)-1
pullLoop : do while (more == 1)
call pull_pt2_results(zmq_socket_pull, Nindex, index, pt2_mwen, task_id, ntask)
@ -256,11 +255,15 @@ subroutine pt2_collector(E, b, tbc, comb, Ncomb, computed, pt2_detail, sumabove,
end do
E0 = sum(pt2_detail(1,: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))
E0 += pt2_detail(1,first_det_of_teeth(tooth)) * prop
avg = E0 + (sumabove(tooth) / Nabove(tooth))
eqt = sqrt(1d0 / (Nabove(tooth)-1) * abs(sum2above(tooth) / Nabove(tooth) - (sumabove(tooth)/Nabove(tooth))**2))
if (tooth <= comb_teeth) then
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))
E0 += pt2_detail(1,first_det_of_teeth(tooth)) * prop
avg = E0 + (sumabove(tooth) / Nabove(tooth))
eqt = sqrt(1d0 / (Nabove(tooth)-1) * abs(sum2above(tooth) / Nabove(tooth) - (sumabove(tooth)/Nabove(tooth))**2))
else
eqt = 0.d0
endif
call wall_time(time)
if (dabs(eqt/avg) < relative_error) then
! Termination
@ -375,7 +378,7 @@ BEGIN_PROVIDER [ integer, size_tbc ]
BEGIN_DOC
! Size of the tbc array
END_DOC
size_tbc = 2*N_det_generators + fragment_count*fragment_first
size_tbc = (comb_teeth+1)*N_det_generators + fragment_count*fragment_first
END_PROVIDER
subroutine get_carlo_workbatch(computed, comb, Ncomb, tbc)
@ -387,16 +390,16 @@ subroutine get_carlo_workbatch(computed, comb, Ncomb, tbc)
integer :: i, j, last_full, dets(comb_teeth)
integer :: icount, n
integer :: k, l
l=1
l=first_det_of_comb
call RANDOM_NUMBER(comb)
do i=1,size(comb)
comb(i) = comb(i) * comb_step
!DIR$ FORCEINLINE
call add_comb(comb(i), computed, tbc, size_tbc, comb_teeth)
Ncomb = i
if (tbc(0) == N_det_generators) return
do while (computed(l))
l=l+1
if (l == N_det_generators+1) return
enddo
k=tbc(0)+1
tbc(k) = l
@ -511,15 +514,25 @@ end subroutine
pt2_weight(1) = psi_coef_generators(1,1)**2
pt2_cweight(1) = psi_coef_generators(1,1)**2
do i=2,N_det_generators
do i=1,N_det_generators
pt2_weight(i) = psi_coef_generators(i,1)**2
pt2_cweight(i) = pt2_cweight(i-1) + psi_coef_generators(i,1)**2
enddo
! Important to loop backwards for numerical precision
pt2_cweight(N_det_generators) = pt2_weight(N_det_generators)
do i=N_det_generators-1,1,-1
pt2_cweight(i) = pt2_weight(i) + pt2_cweight(i+1)
end do
do i=1,N_det_generators
pt2_weight(i) = pt2_weight(i) / pt2_cweight(N_det_generators)
pt2_cweight(i) = pt2_cweight(i) / pt2_cweight(N_det_generators)
pt2_weight(i) = pt2_weight(i) / pt2_cweight(1)
pt2_cweight(i) = pt2_cweight(i) / pt2_cweight(1)
enddo
do i=1,N_det_generators-1
pt2_cweight(i) = 1.d0 - pt2_cweight(i+1)
end do
pt2_cweight(N_det_generators) = 1.d0
norm_left = 1d0

14
plugins/Full_CI_ZMQ/run_pt2_slave.irp.f
View File

@ -113,8 +113,12 @@ subroutine push_pt2_results(zmq_socket_push, N, index, pt2_detail, task_id, ntas
if(rc /= 4*ntask) stop "push"
! Activate is zmq_socket_push is a REQ
! character*(2) :: ok
! rc = f77_zmq_recv( zmq_socket_push, ok, 2, 0)
IRP_IF ZMQ_PUSH
IRP_ELSE
character*(2) :: ok
rc = f77_zmq_recv( zmq_socket_push, ok, 2, 0)
IRP_ENDIF
end subroutine
@ -144,7 +148,11 @@ subroutine pull_pt2_results(zmq_socket_pull, N, index, pt2_detail, task_id, ntas
if(rc /= 4*ntask) stop "pull"
! Activate is zmq_socket_pull is a REP
! rc = f77_zmq_send( zmq_socket_pull, 'ok', 2, 0)
IRP_IF ZMQ_PUSH
IRP_ELSE
rc = f77_zmq_send( zmq_socket_pull, 'ok', 2, 0)
IRP_ENDIF
end subroutine

88
plugins/Full_CI_ZMQ/run_selection_slave.irp.f
View File

@ -21,6 +21,11 @@ subroutine run_selection_slave(thread,iproc,energy)
logical :: done
double precision :: pt2(N_states)
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
zmq_to_qp_run_socket = new_zmq_to_qp_run_socket()
zmq_socket_push = new_zmq_push_socket(thread)
call connect_to_taskserver(zmq_to_qp_run_socket,worker_id,thread)
@ -59,8 +64,8 @@ subroutine run_selection_slave(thread,iproc,energy)
end do
if(ctask > 0) then
call sort_selection_buffer(buf)
call push_selection_results(zmq_socket_push, pt2, buf, task_id(1), ctask)
call merge_selection_buffers(buf,buf2)
call push_selection_results(zmq_socket_push, pt2, buf, task_id(1), ctask)
buf%mini = buf2%mini
pt2 = 0d0
buf%cur = 0
@ -93,24 +98,45 @@ subroutine push_selection_results(zmq_socket_push, pt2, b, task_id, ntask)
integer :: rc
rc = f77_zmq_send( zmq_socket_push, b%cur, 4, ZMQ_SNDMORE)
if(rc /= 4) stop "push"
rc = f77_zmq_send( zmq_socket_push, pt2, 8*N_states, ZMQ_SNDMORE)
if(rc /= 8*N_states) stop "push"
if(rc /= 4) then
print *, 'f77_zmq_send( zmq_socket_push, b%cur, 4, ZMQ_SNDMORE)'
endif
rc = f77_zmq_send( zmq_socket_push, b%val(1), 8*b%cur, ZMQ_SNDMORE)
if(rc /= 8*b%cur) stop "push"
if (b%cur > 0) then
rc = f77_zmq_send( zmq_socket_push, b%det(1,1,1), bit_kind*N_int*2*b%cur, ZMQ_SNDMORE)
if(rc /= bit_kind*N_int*2*b%cur) stop "push"
rc = f77_zmq_send( zmq_socket_push, pt2, 8*N_states, ZMQ_SNDMORE)
if(rc /= 8*N_states) then
print *, 'f77_zmq_send( zmq_socket_push, pt2, 8*N_states, ZMQ_SNDMORE)'
endif
rc = f77_zmq_send( zmq_socket_push, b%val(1), 8*b%cur, ZMQ_SNDMORE)
if(rc /= 8*b%cur) then
print *, 'f77_zmq_send( zmq_socket_push, b%val(1), 8*b%cur, ZMQ_SNDMORE)'
endif
rc = f77_zmq_send( zmq_socket_push, b%det(1,1,1), bit_kind*N_int*2*b%cur, ZMQ_SNDMORE)
if(rc /= bit_kind*N_int*2*b%cur) then
print *, 'f77_zmq_send( zmq_socket_push, b%det(1,1,1), bit_kind*N_int*2*b%cur, ZMQ_SNDMORE)'
endif
endif
rc = f77_zmq_send( zmq_socket_push, ntask, 4, ZMQ_SNDMORE)
if(rc /= 4) stop "push"
if(rc /= 4) then
print *, 'f77_zmq_send( zmq_socket_push, ntask, 4, ZMQ_SNDMORE)'
endif
rc = f77_zmq_send( zmq_socket_push, task_id(1), ntask*4, 0)
if(rc /= 4*ntask) stop "push"
if(rc /= 4*ntask) then
print *, 'f77_zmq_send( zmq_socket_push, task_id(1), ntask*4, 0)'
endif
! Activate is zmq_socket_push is a REQ
! rc = f77_zmq_recv( zmq_socket_push, task_id(1), ntask*4, 0)
IRP_IF ZMQ_PUSH
IRP_ELSE
rc = f77_zmq_recv( zmq_socket_push, task_id(1), ntask*4, 0)
IRP_ENDIF
end subroutine
@ -126,25 +152,45 @@ subroutine pull_selection_results(zmq_socket_pull, pt2, val, det, N, task_id, nt
integer :: rc, rn, i
rc = f77_zmq_recv( zmq_socket_pull, N, 4, 0)
if(rc /= 4) stop "pull"
if(rc /= 4) then
print *, 'f77_zmq_recv( zmq_socket_pull, N, 4, 0)'
endif
rc = f77_zmq_recv( zmq_socket_pull, pt2, N_states*8, 0)
if(rc /= 8*N_states) stop "pull"
if (N>0) then
rc = f77_zmq_recv( zmq_socket_pull, pt2, N_states*8, 0)
if(rc /= 8*N_states) then
print *, 'f77_zmq_recv( zmq_socket_pull, pt2, N_states*8, 0)'
endif
rc = f77_zmq_recv( zmq_socket_pull, val(1), 8*N, 0)
if(rc /= 8*N) stop "pull"
rc = f77_zmq_recv( zmq_socket_pull, val(1), 8*N, 0)
if(rc /= 8*N) then
print *, 'f77_zmq_recv( zmq_socket_pull, val(1), 8*N, 0)'
endif
rc = f77_zmq_recv( zmq_socket_pull, det(1,1,1), bit_kind*N_int*2*N, 0)
if(rc /= bit_kind*N_int*2*N) stop "pull"
rc = f77_zmq_recv( zmq_socket_pull, det(1,1,1), bit_kind*N_int*2*N, 0)
if(rc /= bit_kind*N_int*2*N) then
print *, 'f77_zmq_recv( zmq_socket_pull, det(1,1,1), bit_kind*N_int*2*N, 0)'
endif
else
pt2(:) = 0.d0
endif
rc = f77_zmq_recv( zmq_socket_pull, ntask, 4, 0)
if(rc /= 4) stop "pull"
if(rc /= 4) then
print *, 'f77_zmq_recv( zmq_socket_pull, ntask, 4, 0)'
endif
rc = f77_zmq_recv( zmq_socket_pull, task_id(1), ntask*4, 0)
if(rc /= 4*ntask) stop "pull"
if(rc /= 4*ntask) then
print *, 'f77_zmq_recv( zmq_socket_pull, task_id(1), ntask*4, 0)'
endif
! Activate is zmq_socket_pull is a REP
! rc = f77_zmq_send( zmq_socket_pull, task_id(1), ntask*4, 0)
IRP_IF ZMQ_PUSH
IRP_ELSE
rc = f77_zmq_send( zmq_socket_pull, task_id(1), ntask*4, 0)
IRP_ENDIF
end subroutine

2264
plugins/Full_CI_ZMQ/selection.irp.f
View File

File diff suppressed because it is too large Load Diff

10
plugins/Full_CI_ZMQ/selection_buffer.irp.f
View File

@ -53,11 +53,18 @@ subroutine merge_selection_buffers(b1, b2)
BEGIN_DOC
! Merges the selection buffers b1 and b2 into b2
END_DOC
type(selection_buffer), intent(in) :: b1
type(selection_buffer), intent(inout) :: b1
type(selection_buffer), intent(inout) :: b2
integer(bit_kind), pointer :: detmp(:,:,:)
double precision, pointer :: val(:)
integer :: i, i1, i2, k, nmwen
if (b1%cur == 0) return
do while (b1%val(b1%cur) > b2%mini)
b1%cur = b1%cur-1
if (b1%cur == 0) then
return
endif
enddo
nmwen = min(b1%N, b1%cur+b2%cur)
allocate( val(size(b1%val)), detmp(N_int, 2, size(b1%det,3)) )
i1=1
@ -106,6 +113,7 @@ subroutine sort_selection_buffer(b)
integer(bit_kind), pointer :: detmp(:,:,:)
integer :: i, nmwen
logical, external :: detEq
if (b%cur == 0) return
nmwen = min(b%N, b%cur)
allocate(iorder(b%cur), detmp(N_int, 2, size(b%det,3)))

4
plugins/Full_CI_ZMQ/selection_davidson_slave.irp.f
View File

@ -63,7 +63,9 @@ subroutine run_wf
! --------
print *, 'Davidson'
call davidson_slave_tcp(i)
call omp_set_nested(.True.)
call davidson_slave_tcp(0)
call omp_set_nested(.False.)
print *, 'Davidson done'
else if (trim(zmq_state) == 'pt2') then

16
plugins/Full_CI_ZMQ/zmq_selection.irp.f
View File

@ -17,8 +17,12 @@ subroutine ZMQ_selection(N_in, pt2)
N = max(N_in,1)
if (.True.) then
PROVIDE pt2_e0_denominator
provide nproc
PROVIDE pt2_e0_denominator nproc
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
call new_parallel_job(zmq_to_qp_run_socket,"selection")
call zmq_put_psi(zmq_to_qp_run_socket,1,pt2_e0_denominator,size(pt2_e0_denominator))
call create_selection_buffer(N, N*2, b)
@ -98,7 +102,7 @@ subroutine selection_collector(b, N, pt2)
zmq_to_qp_run_socket = new_zmq_to_qp_run_socket()
zmq_socket_pull = new_zmq_pull_socket()
call create_selection_buffer(N, N*2, b2)
call create_selection_buffer(N, N*8, b2)
allocate(task_id(N_det_generators))
more = 1
pt2(:) = 0d0
@ -107,7 +111,11 @@ subroutine selection_collector(b, N, pt2)
call pull_selection_results(zmq_socket_pull, pt2_mwen, b2%val(1), b2%det(1,1,1), b2%cur, task_id, ntask)
pt2 += pt2_mwen
call merge_selection_buffers(b2,b)
do i=1, b2%cur
call add_to_selection_buffer(b, b2%det(1,1,i), b2%val(i))
if (b2%val(i) > b%mini) exit
end do
do i=1, ntask
if(task_id(i) == 0) then
print *, "Error in collector"

593
plugins/loc_cele/loc_cele.irp.f
View File

@ -1,301 +1,298 @@
program loc_rasorb
implicit none
BEGIN_DOC
! This program performs a localization of the active orbitals
! of a CASSCF wavefunction, reading the orbitals from a RASORB
! file of molcas.
! id1=max is the number of MO in a given symmetry.
END_DOC
integer id1,i_atom,shift,shift_h
parameter (id1=300)
character*1 jobz,uplo
character*64 file1,file2
character*72 string(id1,8),cdum
double precision :: cmo(id1,id1,1),cmoref(id1,id1,1),newcmo(id1,id1,1)
double precision ::s(id1,id1,1),dum,ddum(id1,id1),ovl(id1,id1)
double precision :: w(id1),work(3*id1),t(id1,id1),wi(id1,id1)
integer n,i,j,k,l,nmo(8),isym,nsym,idum,nrot(8),irot(id1,8)
integer ipiv(id1),info,lwork
logical *1 z54
print*,'passed the first copy'
z54=.false.
!Read the name of the RasOrb file
print*,'Entering in the loc program'
! read(5,*) z54
print*,'before = '
accu_norm = 0.d0
do i =1,mo_tot_num
accu_norm += dabs(mo_overlap(i,i))
enddo
print*,'accu_norm = ',accu_norm
nsym = 1
nmo(1) = mo_tot_num
print*,'nmo(1) = ',nmo(1)
cmo = 0.d0
do isym=1,nsym
do i=1,nmo(isym)
do j = 1, ao_num
cmo(j,i,isym) = mo_coef(j,i)
program loc_rasorb
implicit none
BEGIN_DOC
! This program performs a localization of the active orbitals
! of a CASSCF wavefunction, reading the orbitals from a RASORB
! file of molcas.
! id1=max is the number of MO in a given symmetry.
END_DOC
integer id1,i_atom,shift,shift_h
parameter (id1=300)
character*1 jobz,uplo
character*64 file1,file2
character*72 string(id1,8),cdum
double precision :: cmo(id1,id1,1),cmoref(id1,id1,1),newcmo(id1,id1,1)
double precision :: s(id1,id1,1),dum,ddum(id1,id1),ovl(id1,id1)
double precision :: w(id1),work(3*id1),t(id1,id1),wi(id1,id1)
integer n,i,j,k,l,nmo(8),isym,nsym,idum,nrot(8),irot(id1,8)
integer ipiv(id1),info,lwork
logical *1 z54
print*,'passed the first copy'
z54=.false.
!Read the name of the RasOrb file
print*,'Entering in the loc program'
! read(5,*) z54
print*,'before = '
accu_norm = 0.d0
do i =1,mo_tot_num
accu_norm += dabs(mo_overlap(i,i))
enddo
enddo
enddo
print*,'passed the first copy'
do isym=1,nsym
do j=1,mo_tot_num
do i=1,ao_num
newcmo(i,j,isym)=cmo(i,j,isym)
enddo
enddo
enddo
print*,'passed the copy'
nrot(1) = 2 ! number of orbitals to be localized
integer :: index_rot(1000,1)
cmoref = 0.d0
irot = 0
irot(1,1) = 11
irot(2,1) = 12
cmoref(15,1,1) = 1.d0 !
cmoref(14,2,1) = 1.d0 !
! ESATRIENE with 3 bonding and anti bonding orbitals
! First bonding orbital for esa
! cmoref(7,1,1) = 1.d0 !
! cmoref(26,1,1) = 1.d0 !
! Second bonding orbital for esa
! cmoref(45,2,1) = 1.d0 !
! cmoref(64,2,1) = 1.d0 !
! Third bonding orbital for esa
! cmoref(83,3,1) = 1.d0 !
! cmoref(102,3,1) = 1.d0 !
! First anti bonding orbital for esa
! cmoref(7,4,1) = 1.d0 !
! cmoref(26,4,1) = -1.d0 !
! Second anti bonding orbital for esa
! cmoref(45,5,1) = 1.d0 !
! cmoref(64,5,1) = -1.d0 !
! Third anti bonding orbital for esa
! cmoref(83,6,1) = 1.d0 !
! cmoref(102,6,1) = -1.d0 !
! ESATRIENE with 2 bonding and anti bonding orbitals
! AND 2 radical orbitals
! First radical orbital
! cmoref(7,1,1) = 1.d0 !
! First bonding orbital
! cmoref(26,2,1) = 1.d0 !
! cmoref(45,2,1) = 1.d0 !
! Second bonding orbital
! cmoref(64,3,1) = 1.d0 !
! cmoref(83,3,1) = 1.d0 !
! Second radical orbital for esa
! cmoref(102,4,1) = 1.d0 !
! First anti bonding orbital for esa
! cmoref(26,5,1) = 1.d0 !
! cmoref(45,5,1) =-1.d0 !
! Second anti bonding orbital for esa
! cmoref(64,6,1) = 1.d0 !
! cmoref(83,6,1) =-1.d0 !
! ESATRIENE with 1 central bonding and anti bonding orbitals
! AND 4 radical orbitals
! First radical orbital
cmoref(7,1,1) = 1.d0 !
! Second radical orbital
cmoref(26,2,1) = 1.d0 !
! First bonding orbital
cmoref(45,3,1) = 1.d0 !
cmoref(64,3,1) = 1.d0 !
! Third radical orbital for esa
cmoref(83,4,1) = 1.d0 !
! Fourth radical orbital for esa
cmoref(102,5,1) = 1.d0 !
! First anti bonding orbital
cmoref(45,6,1) = 1.d0 !
cmoref(64,6,1) =-1.d0 !
do i = 1, nrot(1)
print*,'irot(i,1) = ',irot(i,1)
enddo
print*,'passed the definition of the referent vectors '
do i = 1, ao_num
do j =1, ao_num
s(i,j,1) = ao_overlap(i,j)
enddo
enddo
!Now big loop over symmetry
do isym=1,nsym
if (nrot(isym).eq.0) cycle
write (6,*)
write (6,*)
write (6,*)
write (6,*) 'WORKING ON SYMMETRY',isym
write (6,*)
!Compute the overlap matrix <ref|vec>
! do i=1,nmo(isym)
do j=1,nrot(isym)
do i=1,ao_num
ddum(i,j)=0.d0
do k=1,ao_num
ddum(i,j)=ddum(i,j)+s(i,k,isym)*cmo(k,irot(j,isym),isym)
enddo
enddo
enddo
do i=1,nrot(isym)
do j=1,nrot(isym)
ovl(i,j)=0.d0
do k=1,ao_num
! do k=1,mo_tot_num
ovl(i,j)=ovl(i,j)+cmoref(k,i,isym)*ddum(k,j)
enddo
enddo
enddo
call maxovl(nrot(isym),nrot(isym),ovl,t,wi)
do i=1,nrot(isym)
do j=1,ao_num
! write (6,*) 'isym,',isym,nrot(isym),nmo(isym)
newcmo(j,irot(i,isym),isym)=0.d0
do k=1,nrot(isym)
newcmo(j,irot(i,isym),isym)=newcmo(j,irot(i,isym),isym) + cmo(j,irot(k,isym),isym)*t(k,i)
print*,'accu_norm = ',accu_norm
nsym = 1
nmo(1) = mo_tot_num
print*,'nmo(1) = ',nmo(1)
cmo = 0.d0
do isym=1,nsym
do i=1,nmo(isym)
do j = 1, ao_num
cmo(j,i,isym) = mo_coef(j,i)
enddo
enddo
enddo
enddo
! if(dabs(newcmo(3,19,1) - mo_coef(3,19)) .gt.1.d-10 )then
! print*,'Something wrong bitch !!'
! print*,'newcmo(3,19,1) = ',newcmo(3,19,1)
! print*,'mo_coef(3,19) = ',mo_coef(3,19)
! stop
! endif
enddo !big loop over symmetry
10 format (4E19.12)
! Now we copyt the newcmo into the mo_coef
mo_coef = 0.d0
do isym=1,nsym
do i=1,nmo(isym)
do j = 1, ao_num
mo_coef(j,i) = newcmo(j,i,isym)
enddo
enddo
enddo
! pause
! we say that it hase been touched, and valid and that everything that
! depends on mo_coef must not be reprovided
double precision :: accu_norm
touch mo_coef
print*,'after = '
accu_norm = 0.d0
do i =1,mo_tot_num
accu_norm += dabs(mo_overlap(i,i))
enddo
print*,'accu_norm = ',accu_norm
! We call the routine that saves mo_coef in the ezfio format
call save_mos
stop
end
print*,'passed the first copy'
do isym=1,nsym
do j=1,mo_tot_num
do i=1,ao_num
newcmo(i,j,isym)=cmo(i,j,isym)
enddo
enddo
enddo
print*,'passed the copy'
nrot(1) = 2 ! number of orbitals to be localized
integer :: index_rot(1000,1)
cmoref = 0.d0
irot = 0
irot(1,1) = 11
irot(2,1) = 12
cmoref(15,1,1) = 1.d0 !
cmoref(14,2,1) = 1.d0 !
! ESATRIENE with 3 bonding and anti bonding orbitals
! First bonding orbital for esa
! cmoref(7,1,1) = 1.d0 !
! cmoref(26,1,1) = 1.d0 !
! Second bonding orbital for esa
! cmoref(45,2,1) = 1.d0 !
! cmoref(64,2,1) = 1.d0 !
! Third bonding orbital for esa
! cmoref(83,3,1) = 1.d0 !
! cmoref(102,3,1) = 1.d0 !
! First anti bonding orbital for esa
! cmoref(7,4,1) = 1.d0 !
! cmoref(26,4,1) = -1.d0 !
! Second anti bonding orbital for esa
! cmoref(45,5,1) = 1.d0 !
! cmoref(64,5,1) = -1.d0 !
! Third anti bonding orbital for esa
! cmoref(83,6,1) = 1.d0 !
! cmoref(102,6,1) = -1.d0 !
! ESATRIENE with 2 bonding and anti bonding orbitals
! AND 2 radical orbitals
! First radical orbital
! cmoref(7,1,1) = 1.d0 !
! First bonding orbital
! cmoref(26,2,1) = 1.d0 !
! cmoref(45,2,1) = 1.d0 !
! Second bonding orbital
! cmoref(64,3,1) = 1.d0 !
! cmoref(83,3,1) = 1.d0 !
! Second radical orbital for esa
! cmoref(102,4,1) = 1.d0 !
! First anti bonding orbital for esa
! cmoref(26,5,1) = 1.d0 !
! cmoref(45,5,1) =-1.d0 !
! Second anti bonding orbital for esa
! cmoref(64,6,1) = 1.d0 !
! cmoref(83,6,1) =-1.d0 !
! ESATRIENE with 1 central bonding and anti bonding orbitals
! AND 4 radical orbitals
! First radical orbital
cmoref(7,1,1) = 1.d0 !
! Second radical orbital
cmoref(26,2,1) = 1.d0 !
! First bonding orbital
cmoref(45,3,1) = 1.d0 !
cmoref(64,3,1) = 1.d0 !
! Third radical orbital for esa
cmoref(83,4,1) = 1.d0 !
! Fourth radical orbital for esa
cmoref(102,5,1) = 1.d0 !
! First anti bonding orbital
cmoref(45,6,1) = 1.d0 !
cmoref(64,6,1) =-1.d0 !
do i = 1, nrot(1)
print*,'irot(i,1) = ',irot(i,1)
enddo
print*,'passed the definition of the referent vectors '
do i = 1, ao_num
do j =1, ao_num
s(i,j,1) = ao_overlap(i,j)
enddo
enddo
!Now big loop over symmetry
do isym=1,nsym
if (nrot(isym).eq.0) cycle
write (6,*)
write (6,*)
write (6,*)
write (6,*) 'WORKING ON SYMMETRY',isym
write (6,*)
!Compute the overlap matrix <ref|vec>
! do i=1,nmo(isym)
do j=1,nrot(isym)
do i=1,ao_num
ddum(i,j)=0.d0
do k=1,ao_num
ddum(i,j)=ddum(i,j)+s(i,k,isym)*cmo(k,irot(j,isym),isym)
enddo
enddo
enddo
do i=1,nrot(isym)
do j=1,nrot(isym)
ovl(i,j)=0.d0
do k=1,ao_num
! do k=1,mo_tot_num
ovl(i,j)=ovl(i,j)+cmoref(k,i,isym)*ddum(k,j)
enddo
enddo
enddo
call maxovl(nrot(isym),nrot(isym),ovl,t,wi)
do i=1,nrot(isym)
do j=1,ao_num
! write (6,*) 'isym,',isym,nrot(isym),nmo(isym)
newcmo(j,irot(i,isym),isym)=0.d0
do k=1,nrot(isym)
newcmo(j,irot(i,isym),isym)=newcmo(j,irot(i,isym),isym) + cmo(j,irot(k,isym),isym)*t(k,i)
enddo
enddo
enddo
! if(dabs(newcmo(3,19,1) - mo_coef(3,19)) .gt.1.d-10 )then
! print*,'Something wrong bitch !!'
! print*,'newcmo(3,19,1) = ',newcmo(3,19,1)
! print*,'mo_coef(3,19) = ',mo_coef(3,19)
! stop
! endif
enddo !big loop over symmetry
10 format (4E19.12)
! Now we copyt the newcmo into the mo_coef
mo_coef = 0.d0
do isym=1,nsym
do i=1,nmo(isym)
do j = 1, ao_num
mo_coef(j,i) = newcmo(j,i,isym)
enddo
enddo
enddo
! pause
! we say that it hase been touched, and valid and that everything that
! depends on mo_coef must not be reprovided
double precision :: accu_norm
touch mo_coef
print*,'after = '
accu_norm = 0.d0
do i =1,mo_tot_num
accu_norm += dabs(mo_overlap(i,i))
enddo
print*,'accu_norm = ',accu_norm
! We call the routine that saves mo_coef in the ezfio format
call save_mos
stop
end

2
plugins/mrcepa0/dressing.irp.f
View File

@ -684,7 +684,7 @@ subroutine getHP(a,h,p,Nint)
end do lh
h = deg
!isInCassd = .true.
end function
end subroutine
BEGIN_PROVIDER [ double precision, delta_mrcepa0_ij, (N_det_ref,N_det_non_ref,N_states) ]

32
plugins/mrcepa0/dressing_slave.irp.f
View File

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

45
src/Davidson/davidson_parallel.irp.f
View File

@ -63,13 +63,14 @@ subroutine davidson_slave_work(zmq_to_qp_run_socket, zmq_socket_push, N_st, sze,
character*(512) :: msg
integer :: imin, imax, ishift, istep
integer, allocatable :: psi_det_read(:,:,:)
double precision, allocatable :: v_0(:,:), s_0(:,:), u_t(:,:)
!DIR$ ATTRIBUTES ALIGN : $IRP_ALIGN :: u_t, v_0, s_0
! Get wave function (u_t)
! -----------------------
integer :: rc
integer :: rc
integer :: N_states_read, N_det_read, psi_det_size_read
integer :: N_det_selectors_read, N_det_generators_read
double precision :: energy(N_st)
@ -107,12 +108,11 @@ subroutine davidson_slave_work(zmq_to_qp_run_socket, zmq_socket_push, N_st, sze,
TOUCH N_det
endif
allocate(v_0(sze,N_st), s_0(sze,N_st),u_t(N_st,N_det_read))
allocate(v_0(sze,N_st), s_0(sze,N_st),u_t(N_st,N_det))
rc = f77_zmq_recv(zmq_to_qp_run_socket,psi_det,N_int*2*N_det*bit_kind,0)
if (rc /= N_int*2*N_det*bit_kind) then
print *, 'f77_zmq_recv(zmq_to_qp_run_socket,psi_det,N_int*2*N_det*bit_kind,0)'
rc = f77_zmq_recv(zmq_to_qp_run_socket,psi_det,N_int*2*N_det_read*bit_kind,0)
if (rc /= N_int*2*N_det_read*bit_kind) then
print *, 'f77_zmq_recv(zmq_to_qp_run_socket,psi_det,N_int*2*N_det_read*bit_kind,0)'
stop 'error'
endif
@ -132,6 +132,7 @@ subroutine davidson_slave_work(zmq_to_qp_run_socket, zmq_socket_push, N_st, sze,
! Run tasks
! ---------
do
v_0 = 0.d0
s_0 = 0.d0
@ -168,12 +169,15 @@ subroutine davidson_push_results(zmq_socket_push, v_0, s_0, task_id)
if(rc /= 4) stop "davidson_push_results failed to push task_id"
! Activate is zmq_socket_push is a REQ
! integer :: idummy
! rc = f77_zmq_recv( zmq_socket_push, idummy, 4, 0)
! if (rc /= 4) then
! print *, irp_here, ': f77_zmq_send( zmq_socket_push, idummy, 4, 0)'
! stop 'error'
! endif
IRP_IF ZMQ_PUSH
IRP_ELSE
integer :: idummy
rc = f77_zmq_recv( zmq_socket_push, idummy, 4, 0)
if (rc /= 4) then
print *, irp_here, ': f77_zmq_send( zmq_socket_push, idummy, 4, 0)'
stop 'error'
endif
IRP_ENDIF
end subroutine
@ -200,11 +204,14 @@ subroutine davidson_pull_results(zmq_socket_pull, v_0, s_0, task_id)
if(rc /= 4) stop "davidson_pull_results failed to pull task_id"
! Activate if zmq_socket_pull is a REP
! rc = f77_zmq_send( zmq_socket_pull, 0, 4, 0)
! if (rc /= 4) then
! print *, irp_here, ' : f77_zmq_send (zmq_socket_pull,...'
! stop 'error'
! endif
IRP_IF ZMQ_PUSH
IRP_ELSE
rc = f77_zmq_send( zmq_socket_pull, 0, 4, 0)
if (rc /= 4) then
print *, irp_here, ' : f77_zmq_send (zmq_socket_pull,...'
stop 'error'
endif
IRP_ENDIF
end subroutine
@ -293,10 +300,6 @@ subroutine H_S2_u_0_nstates_zmq(v_0,s_0,u_0,N_st,sze)
if(N_st /= N_states_diag .or. sze < N_det) stop "assert fail in H_S2_u_0_nstates"
ASSERT (Nint > 0)
ASSERT (Nint == N_int)
ASSERT (n>0)
call new_parallel_job(zmq_to_qp_run_socket,'davidson')
character*(512) :: task

1
src/Davidson/davidson_slave.irp.f
View File

@ -10,6 +10,7 @@ program davidson_slave
call provide_everything
call switch_qp_run_to_master
call omp_set_nested(.True.)
zmq_context = f77_zmq_ctx_new ()
zmq_state = 'davidson'

72
src/Davidson/u0Hu0.irp.f
View File

@ -9,12 +9,13 @@ subroutine u_0_H_u_0(e_0,u_0,n,keys_tmp,Nint,N_st,sze)
END_DOC
integer, intent(in) :: n,Nint, N_st, sze
double precision, intent(out) :: e_0(N_st)
double precision, intent(inout):: u_0(sze,N_st)
double precision, intent(inout) :: u_0(sze,N_st)
integer(bit_kind),intent(in) :: keys_tmp(Nint,2,n)
double precision, allocatable :: v_0(:,:), s_0(:,:)
double precision :: u_dot_u,u_dot_v,diag_H_mat_elem
integer :: i,j
integer :: i,j
allocate (v_0(sze,N_st),s_0(sze,N_st))
call H_S2_u_0_nstates_openmp(v_0,s_0,u_0,N_st,sze)
do i=1,N_st
@ -132,7 +133,7 @@ subroutine H_S2_u_0_nstates_openmp_work_$N_int(v_0,s_0,u_t,N_st,sze,istart,iend,
integer, allocatable :: singles_a(:)
integer, allocatable :: singles_b(:)
integer, allocatable :: idx(:), idx0(:)
integer :: maxab, n_singles_a, n_singles_b, kcol_prev, nmax
integer :: maxab, n_singles_a, n_singles_b, kcol_prev
integer*8 :: k8
double precision, allocatable :: v_t(:,:), s_t(:,:)
!DIR$ ATTRIBUTES ALIGN : $IRP_ALIGN :: v_t, s_t
@ -158,11 +159,11 @@ subroutine H_S2_u_0_nstates_openmp_work_$N_int(v_0,s_0,u_t,N_st,sze,istart,iend,
!$OMP psi_bilinear_matrix_transp_order, N_st, &
!$OMP psi_bilinear_matrix_order_transp_reverse, &
!$OMP psi_bilinear_matrix_columns_loc, &
!$OMP istart, iend, istep, &
!$OMP istart, iend, istep, irp_here, &
!$OMP ishift, idx0, u_t, maxab, v_0, s_0) &
!$OMP PRIVATE(krow, kcol, tmp_det, spindet, k_a, k_b, i, &
!$OMP lcol, lrow, l_a, l_b, nmax, &
!$OMP buffer, doubles, n_doubles, &
!$OMP lcol, lrow, l_a, l_b, &
!$OMP buffer, doubles, n_doubles, &
!$OMP tmp_det2, hij, sij, idx, l, kcol_prev, v_t, &
!$OMP singles_a, n_singles_a, singles_b, &
!$OMP n_singles_b, s_t, k8)
@ -181,12 +182,18 @@ subroutine H_S2_u_0_nstates_openmp_work_$N_int(v_0,s_0,u_t,N_st,sze,istart,iend,
v_t = 0.d0
s_t = 0.d0
ASSERT (iend <= N_det)
ASSERT (istart > 0)
ASSERT (istep > 0)
!$OMP DO SCHEDULE(dynamic,64)
do k_a=istart+ishift,iend,istep
krow = psi_bilinear_matrix_rows(k_a)
ASSERT (krow <= N_det_alpha_unique)
kcol = psi_bilinear_matrix_columns(k_a)
ASSERT (kcol <= N_det_beta_unique)
tmp_det(1:$N_int,1) = psi_det_alpha_unique(1:$N_int, krow)
tmp_det(1:$N_int,2) = psi_det_beta_unique (1:$N_int, kcol)
@ -208,11 +215,15 @@ subroutine H_S2_u_0_nstates_openmp_work_$N_int(v_0,s_0,u_t,N_st,sze,istart,iend,
tmp_det2(1:$N_int,2) = psi_det_beta_unique(1:$N_int, lcol)
l_a = psi_bilinear_matrix_columns_loc(lcol)
ASSERT (l_a <= N_det)
nmax = psi_bilinear_matrix_columns_loc(lcol+1) - l_a
do j=1,nmax
do j=1,psi_bilinear_matrix_columns_loc(lcol+1) - l_a
lrow = psi_bilinear_matrix_rows(l_a)
ASSERT (lrow <= N_det_alpha_unique)
buffer(1:$N_int,j) = psi_det_alpha_unique(1:$N_int, lrow)
ASSERT (l_a <= N_det)
idx(j) = l_a
l_a = l_a+1
enddo
@ -227,7 +238,11 @@ subroutine H_S2_u_0_nstates_openmp_work_$N_int(v_0,s_0,u_t,N_st,sze,istart,iend,
do k = 1,n_singles_a
l_a = singles_a(k)
ASSERT (l_a <= N_det)
lrow = psi_bilinear_matrix_rows(l_a)
ASSERT (lrow <= N_det_alpha_unique)
tmp_det2(1:$N_int,1) = psi_det_alpha_unique(1:$N_int, lrow)
call i_H_j_double_alpha_beta(tmp_det,tmp_det2,$N_int,hij)
call get_s2(tmp_det,tmp_det2,$N_int,sij)
@ -256,7 +271,10 @@ subroutine H_S2_u_0_nstates_openmp_work_$N_int(v_0,s_0,u_t,N_st,sze,istart,iend,
! -----------------------------------------------------------------------
krow = psi_bilinear_matrix_rows(k_a)
ASSERT (krow <= N_det_alpha_unique)
kcol = psi_bilinear_matrix_columns(k_a)
ASSERT (kcol <= N_det_beta_unique)
tmp_det(1:$N_int,1) = psi_det_alpha_unique(1:$N_int, krow)
tmp_det(1:$N_int,2) = psi_det_beta_unique (1:$N_int, kcol)
@ -265,16 +283,19 @@ subroutine H_S2_u_0_nstates_openmp_work_$N_int(v_0,s_0,u_t,N_st,sze,istart,iend,
! ----------------------------------------------------------------------
k_b = psi_bilinear_matrix_order_transp_reverse(k_a)
ASSERT (k_b <= N_det)
spindet(1:$N_int) = tmp_det(1:$N_int,1)
! Loop inside the beta column to gather all the connected alphas
l_a = k_a+1
nmax = min(N_det_alpha_unique, N_det - l_a)
do i=1,nmax
do i=1,N_det_alpha_unique
if (l_a > N_det) exit
lcol = psi_bilinear_matrix_columns(l_a)
if (lcol /= kcol) exit
lrow = psi_bilinear_matrix_rows(l_a)
ASSERT (lrow <= N_det_alpha_unique)
buffer(1:$N_int,i) = psi_det_alpha_unique(1:$N_int, lrow)
idx(i) = l_a
l_a = l_a+1
@ -291,9 +312,14 @@ subroutine H_S2_u_0_nstates_openmp_work_$N_int(v_0,s_0,u_t,N_st,sze,istart,iend,
tmp_det2(1:$N_int,2) = psi_det_beta_unique (1:$N_int, kcol)
do i=1,n_singles_a
l_a = singles_a(i)
ASSERT (l_a <= N_det)
lrow = psi_bilinear_matrix_rows(l_a)
ASSERT (lrow <= N_det_alpha_unique)
tmp_det2(1:$N_int,1) = psi_det_alpha_unique(1:$N_int, lrow)
call i_H_j_mono_spin( tmp_det, tmp_det2, $N_int, 1, hij)
do l=1,N_st
v_t(l,l_a) = v_t(l,l_a) + hij * u_t(l,k_a)
v_t(l,k_a) = v_t(l,k_a) + hij * u_t(l,l_a)
@ -307,7 +333,11 @@ subroutine H_S2_u_0_nstates_openmp_work_$N_int(v_0,s_0,u_t,N_st,sze,istart,iend,
do i=1,n_doubles
l_a = doubles(i)
ASSERT (l_a <= N_det)
lrow = psi_bilinear_matrix_rows(l_a)
ASSERT (lrow <= N_det_alpha_unique)
call i_H_j_double_spin( tmp_det(1,1), psi_det_alpha_unique(1, lrow), $N_int, hij)
do l=1,N_st
v_t(l,l_a) = v_t(l,l_a) + hij * u_t(l,k_a)
@ -317,7 +347,6 @@ subroutine H_S2_u_0_nstates_openmp_work_$N_int(v_0,s_0,u_t,N_st,sze,istart,iend,
enddo
! Single and double beta excitations
! ==================================
@ -337,14 +366,17 @@ subroutine H_S2_u_0_nstates_openmp_work_$N_int(v_0,s_0,u_t,N_st,sze,istart,iend,
! -----------------------------------------------------------------------
k_b = psi_bilinear_matrix_order_transp_reverse(k_a)
ASSERT (k_b <= N_det)
! Loop inside the alpha row to gather all the connected betas
l_b = k_b+1
nmax = min(N_det_beta_unique, N_det - l_b)
do i=1,nmax
do i=1,N_det_beta_unique
if (l_b > N_det) exit
lrow = psi_bilinear_matrix_transp_rows(l_b)
if (lrow /= krow) exit
lcol = psi_bilinear_matrix_transp_columns(l_b)
ASSERT (lcol <= N_det_beta_unique)
buffer(1:$N_int,i) = psi_det_beta_unique(1:$N_int, lcol)
idx(i) = l_b
l_b = l_b+1
@ -361,10 +393,15 @@ subroutine H_S2_u_0_nstates_openmp_work_$N_int(v_0,s_0,u_t,N_st,sze,istart,iend,
tmp_det2(1:$N_int,1) = psi_det_alpha_unique(1:$N_int, krow)
do i=1,n_singles_b
l_b = singles_b(i)
ASSERT (l_b <= N_det)
lcol = psi_bilinear_matrix_transp_columns(l_b)
ASSERT (lcol <= N_det_beta_unique)
tmp_det2(1:$N_int,2) = psi_det_beta_unique (1:$N_int, lcol)
call i_H_j_mono_spin( tmp_det, tmp_det2, $N_int, 2, hij)
l_a = psi_bilinear_matrix_transp_order(l_b)
ASSERT (l_a <= N_det)
do l=1,N_st
v_t(l,l_a) = v_t(l,l_a) + hij * u_t(l,k_a)
v_t(l,k_a) = v_t(l,k_a) + hij * u_t(l,l_a)
@ -377,9 +414,15 @@ subroutine H_S2_u_0_nstates_openmp_work_$N_int(v_0,s_0,u_t,N_st,sze,istart,iend,
do i=1,n_doubles
l_b = doubles(i)
ASSERT (l_b <= N_det)
lcol = psi_bilinear_matrix_transp_columns(l_b)
ASSERT (lcol <= N_det_beta_unique)
call i_H_j_double_spin( tmp_det(1,2), psi_det_beta_unique(1, lcol), $N_int, hij)
l_a = psi_bilinear_matrix_transp_order(l_b)
ASSERT (l_a <= N_det)
do l=1,N_st
v_t(l,l_a) = v_t(l,l_a) + hij * u_t(l,k_a)
v_t(l,k_a) = v_t(l,k_a) + hij * u_t(l,l_a)
@ -396,7 +439,10 @@ subroutine H_S2_u_0_nstates_openmp_work_$N_int(v_0,s_0,u_t,N_st,sze,istart,iend,
! -----------------------------------------------------------------------
krow = psi_bilinear_matrix_rows(k_a)
ASSERT (krow <= N_det_alpha_unique)
kcol = psi_bilinear_matrix_columns(k_a)
ASSERT (kcol <= N_det_beta_unique)
tmp_det(1:$N_int,1) = psi_det_alpha_unique(1:$N_int, krow)
tmp_det(1:$N_int,2) = psi_det_beta_unique (1:$N_int, kcol)

3
src/Davidson/u0Hu0_old.irp.f
View File

@ -496,10 +496,11 @@ subroutine H_S2_u_0_nstates_test(v_0,s_0,u_0,H_jj,S2_jj,n,keys_tmp,Nint,N_st,sze
! endif
! if ((degree == 2).and.(exc(0,1,1)==1)) cycle
! if ((degree > 1)) cycle
! if ((degree == 1)) cycle
! if (exc(0,1,2) /= 0) cycle
! if (exc(0,1,1) == 2) cycle
! if (exc(0,1,2) == 2) cycle
! if ((degree==1).and.(exc(0,1,2) == 1)) cycle
if ((degree==1).and.(exc(0,1,1) == 1)) cycle
call i_H_j(keys_tmp(1,1,j),keys_tmp(1,1,i),Nint,hij)
do l=1,N_st
!$OMP ATOMIC

29
src/Determinants/H_apply.irp.f
View File

@ -362,12 +362,16 @@ subroutine push_pt2(zmq_socket_push,pt2,norm_pert,H_pert_diag,i_generator,N_st,t
endif
! Activate if zmq_socket_push is a REQ
! integer :: idummy
! rc = f77_zmq_recv( zmq_socket_push, idummy, 4, 0)
! if (rc /= 4) then
! print *, irp_here, 'f77_zmq_send( zmq_socket_push, idummy, 4, 0)'
! stop 'error'
! endif
IRP_IF ZMQ_PUSH
IRP_ELSE
integer :: idummy
rc = f77_zmq_recv( zmq_socket_push, idummy, 4, 0)
if (rc /= 4) then
print *, irp_here, 'f77_zmq_send( zmq_socket_push, idummy, 4, 0)'
stop 'error'
endif
IRP_ENDIF
end
subroutine pull_pt2(zmq_socket_pull,pt2,norm_pert,H_pert_diag,i_generator,N_st,n,task_id)
@ -433,11 +437,14 @@ subroutine pull_pt2(zmq_socket_pull,pt2,norm_pert,H_pert_diag,i_generator,N_st,n
endif
! Activate if zmq_socket_pull is a REP
! rc = f77_zmq_send( zmq_socket_pull, 0, 4, 0)
! if (rc /= 4) then
! print *, irp_here, 'f77_zmq_send( zmq_socket_pull, 0, 4, 0)'
! stop 'error'
! endif
IRP_IF ZMQ_PUSH
IRP_ELSE
rc = f77_zmq_send( zmq_socket_pull, 0, 4, 0)
if (rc /= 4) then
print *, irp_here, 'f77_zmq_send( zmq_socket_pull, 0, 4, 0)'
stop 'error'
endif
IRP_ENDIF
end

7
src/Determinants/determinants.irp.f
View File

@ -324,9 +324,12 @@ END_PROVIDER
BEGIN_PROVIDER [ integer(bit_kind), psi_det_sorted, (N_int,2,psi_det_size) ]
&BEGIN_PROVIDER [ double precision, psi_coef_sorted, (psi_det_size,N_states) ]
&BEGIN_PROVIDER [ double precision, psi_average_norm_contrib_sorted, (psi_det_size) ]
&BEGIN_PROVIDER [ integer, psi_det_sorted_order, (psi_det_size) ]
implicit none
BEGIN_DOC
! Wave function sorted by determinants contribution to the norm (state-averaged)
!
! psi_det_sorted_order(i) -> k : index in psi_det
END_DOC
integer :: i,j,k
integer, allocatable :: iorder(:)
@ -346,6 +349,10 @@ END_PROVIDER
enddo
psi_average_norm_contrib_sorted(i) = -psi_average_norm_contrib_sorted(i)
enddo
do i=1,N_det
psi_det_sorted_order(iorder(i)) = i
enddo
deallocate(iorder)

7
src/Determinants/slater_rules.irp.f
View File

@ -1929,7 +1929,7 @@ subroutine a_operator(iorb,ispin,key,hjj,Nint,na,nb)
ASSERT (Nint > 0)
k = ishft(iorb-1,-bit_kind_shift)+1
ASSERT (k > 0)
ASSERT (k>0)
l = iorb - ishft(k-1,bit_kind_shift)-1
key(k,ispin) = ibclr(key(k,ispin),l)
other_spin = iand(ispin,1)+1
@ -1977,11 +1977,12 @@ subroutine ac_operator(iorb,ispin,key,hjj,Nint,na,nb)
!DIR$ FORCEINLINE
call bitstring_to_list_ab(key, occ, tmp, Nint)
ASSERT (tmp(1) == elec_alpha_num)
ASSERT (tmp(2) == elec_beta_num)
ASSERT (tmp(2) == elec_alpha_num)
k = ishft(iorb-1,-bit_kind_shift)+1
ASSERT (k > 0)
ASSERT (k >0)
l = iorb - ishft(k-1,bit_kind_shift)-1
ASSERT (l >= 0)
key(k,ispin) = ibset(key(k,ispin),l)
other_spin = iand(ispin,1)+1

55
src/Determinants/spindeterminants.irp.f
View File

@ -20,7 +20,7 @@ integer*8 function spin_det_search_key(det,Nint)
do i=2,Nint
spin_det_search_key = ieor(spin_det_search_key,det(i))
enddo
spin_det_search_key = spin_det_search_key+1_bit_kind-unsigned_shift
spin_det_search_key = spin_det_search_key+unsigned_shift
end
@ -189,9 +189,7 @@ integer function get_index_in_psi_det_alpha_unique(key,Nint)
enddo
i += 1
if (i > N_det_alpha_unique) then
return
endif
ASSERT (i <= N_det_alpha_unique)
!DIR$ FORCEINLINE
do while (spin_det_search_key(psi_det_alpha_unique(1,i),Nint) == det_ref)
@ -213,6 +211,7 @@ integer function get_index_in_psi_det_alpha_unique(key,Nint)
endif
i += 1
if (i > N_det_alpha_unique) then
ASSERT (get_index_in_psi_det_alpha_unique > 0)
return
endif
@ -270,9 +269,7 @@ integer function get_index_in_psi_det_beta_unique(key,Nint)
enddo
i += 1
if (i > N_det_beta_unique) then
return
endif
ASSERT (i <= N_det_beta_unique)
!DIR$ FORCEINLINE
do while (spin_det_search_key(psi_det_beta_unique(1,i),Nint) == det_ref)
@ -294,6 +291,7 @@ integer function get_index_in_psi_det_beta_unique(key,Nint)
endif
i += 1
if (i > N_det_beta_unique) then
ASSERT (get_index_in_psi_det_beta_unique > 0)
return
endif
@ -413,7 +411,12 @@ BEGIN_PROVIDER [ double precision, psi_bilinear_matrix_values, (N_det,N_states)
!$OMP PARALLEL DO DEFAULT(SHARED) PRIVATE(i,j,k,l)
do k=1,N_det
i = get_index_in_psi_det_alpha_unique(psi_det(1,1,k),N_int)
ASSERT (i>0)
ASSERT (i<=N_det_alpha_unique)
j = get_index_in_psi_det_beta_unique (psi_det(1,2,k),N_int)
ASSERT (j>0)
ASSERT (j<=N_det_alpha_unique)
do l=1,N_states
psi_bilinear_matrix_values(k,l) = psi_coef(k,l)
@ -421,6 +424,7 @@ BEGIN_PROVIDER [ double precision, psi_bilinear_matrix_values, (N_det,N_states)
psi_bilinear_matrix_rows(k) = i
psi_bilinear_matrix_columns(k) = j
to_sort(k) = int(N_det_alpha_unique,8) * int(j-1,8) + int(i,8)
ASSERT (to_sort(k) > 0_8)
psi_bilinear_matrix_order(k) = k
enddo
!$OMP END PARALLEL DO
@ -431,6 +435,12 @@ BEGIN_PROVIDER [ double precision, psi_bilinear_matrix_values, (N_det,N_states)
call dset_order(psi_bilinear_matrix_values(1,l),psi_bilinear_matrix_order,N_det)
enddo
deallocate(to_sort)
ASSERT (minval(psi_bilinear_matrix_rows) == 1)
ASSERT (minval(psi_bilinear_matrix_columns) == 1)
ASSERT (minval(psi_bilinear_matrix_order) == 1)
ASSERT (maxval(psi_bilinear_matrix_rows) == N_det_alpha_unique)
ASSERT (maxval(psi_bilinear_matrix_columns) == N_det_beta_unique)
ASSERT (maxval(psi_bilinear_matrix_order) == N_det)
END_PROVIDER
@ -438,7 +448,7 @@ BEGIN_PROVIDER [ integer, psi_bilinear_matrix_order_reverse , (N_det) ]
use bitmasks
implicit none
BEGIN_DOC
! Order which allors to go from psi_bilinear_matrix to psi_det
! Order which allows to go from psi_bilinear_matrix to psi_det
END_DOC
integer :: k
!$OMP PARALLEL DO DEFAULT(SHARED) PRIVATE(k)
@ -446,6 +456,8 @@ BEGIN_PROVIDER [ integer, psi_bilinear_matrix_order_reverse , (N_det) ]
psi_bilinear_matrix_order_reverse(psi_bilinear_matrix_order(k)) = k
enddo
!$OMP END PARALLEL DO
ASSERT (minval(psi_bilinear_matrix_order) == 1)
ASSERT (maxval(psi_bilinear_matrix_order) == N_det)
END_PROVIDER
@ -471,8 +483,13 @@ BEGIN_PROVIDER [ integer, psi_bilinear_matrix_columns_loc, (N_det_beta_unique+1)
l = psi_bilinear_matrix_columns(k)
psi_bilinear_matrix_columns_loc(l) = k
endif
if (psi_bilinear_matrix_columns(k) < 1) then
stop '(psi_bilinear_matrix_columns(k) < 1)'
endif
enddo
psi_bilinear_matrix_columns_loc(N_det_beta_unique+1) = N_det+1
ASSERT (minval(psi_bilinear_matrix_columns_loc) == 1)
ASSERT (maxval(psi_bilinear_matrix_columns_loc) == N_det+1)
END_PROVIDER
BEGIN_PROVIDER [ double precision, psi_bilinear_matrix_transp_values, (N_det,N_states) ]
@ -496,20 +513,27 @@ BEGIN_PROVIDER [ double precision, psi_bilinear_matrix_transp_values, (N_det,N_
integer*8, allocatable :: to_sort(:)
allocate(to_sort(N_det))
!$OMP PARALLEL DEFAULT(SHARED) PRIVATE(i,j,k,l)
!$OMP DO COLLAPSE(2)
do l=1,N_states
!$OMP DO
do k=1,N_det
psi_bilinear_matrix_transp_values (k,l) = psi_bilinear_matrix_values (k,l)
enddo
!$OMP ENDDO
enddo
!$OMP ENDDO
!$OMP DO
do k=1,N_det
psi_bilinear_matrix_transp_columns(k) = psi_bilinear_matrix_columns(k)
ASSERT (psi_bilinear_matrix_transp_columns(k) > 0)
ASSERT (psi_bilinear_matrix_transp_columns(k) <= N_det)
psi_bilinear_matrix_transp_rows (k) = psi_bilinear_matrix_rows (k)
ASSERT (psi_bilinear_matrix_transp_rows(k) > 0)
ASSERT (psi_bilinear_matrix_transp_rows(k) <= N_det)
i = psi_bilinear_matrix_transp_columns(k)
j = psi_bilinear_matrix_transp_rows (k)
to_sort(k) = int(N_det_beta_unique,8) * int(j-1,8) + int(i,8)
ASSERT (to_sort(k) > 0)
psi_bilinear_matrix_transp_order(k) = k
enddo
!$OMP ENDDO
@ -521,6 +545,12 @@ BEGIN_PROVIDER [ double precision, psi_bilinear_matrix_transp_values, (N_det,N_
call dset_order(psi_bilinear_matrix_transp_values(1,l),psi_bilinear_matrix_transp_order,N_det)
enddo
deallocate(to_sort)
ASSERT (minval(psi_bilinear_matrix_transp_columns) == 1)
ASSERT (minval(psi_bilinear_matrix_transp_rows) == 1)
ASSERT (minval(psi_bilinear_matrix_transp_order) == 1)
ASSERT (maxval(psi_bilinear_matrix_transp_columns) == N_det_beta_unique)
ASSERT (maxval(psi_bilinear_matrix_transp_rows) == N_det_alpha_unique)
ASSERT (maxval(psi_bilinear_matrix_transp_order) == N_det)
END_PROVIDER
BEGIN_PROVIDER [ integer, psi_bilinear_matrix_transp_rows_loc, (N_det_alpha_unique+1) ]
@ -542,6 +572,8 @@ BEGIN_PROVIDER [ integer, psi_bilinear_matrix_transp_rows_loc, (N_det_alpha_uniq
endif
enddo
psi_bilinear_matrix_transp_rows_loc(N_det_alpha_unique+1) = N_det+1
ASSERT (minval(psi_bilinear_matrix_transp_rows_loc) == 1)
ASSERT (maxval(psi_bilinear_matrix_transp_rows_loc) == N_det+1)
END_PROVIDER
BEGIN_PROVIDER [ integer, psi_bilinear_matrix_order_transp_reverse , (N_det) ]
@ -552,11 +584,14 @@ BEGIN_PROVIDER [ integer, psi_bilinear_matrix_order_transp_reverse , (N_det) ]
END_DOC
integer :: k
psi_bilinear_matrix_order_transp_reverse = -1
!$OMP PARALLEL DO DEFAULT(SHARED) PRIVATE(k)
do k=1,N_det
psi_bilinear_matrix_order_transp_reverse(psi_bilinear_matrix_transp_order(k)) = k
enddo
!$OMP END PARALLEL DO
ASSERT (minval(psi_bilinear_matrix_order_transp_reverse) == 1)
ASSERT (maxval(psi_bilinear_matrix_order_transp_reverse) == N_det)
END_PROVIDER

28
src/Integrals_Bielec/ao_bielec_integrals_in_map_slave.irp.f
View File

@ -57,12 +57,15 @@ subroutine push_integrals(zmq_socket_push, n_integrals, buffer_i, buffer_value,
endif
! Activate is zmq_socket_push is a REQ
! integer :: idummy
! rc = f77_zmq_recv( zmq_socket_push, idummy, 4, 0)
! if (rc /= 4) then
! print *, irp_here, ': f77_zmq_send( zmq_socket_push, idummy, 4, 0)'
! stop 'error'
! endif
IRP_IF ZMQ_PUSH
IRP_ELSE
integer :: idummy
rc = f77_zmq_recv( zmq_socket_push, idummy, 4, 0)
if (rc /= 4) then
print *, irp_here, ': f77_zmq_send( zmq_socket_push, idummy, 4, 0)'
stop 'error'
endif
IRP_ENDIF
end
@ -187,11 +190,14 @@ subroutine ao_bielec_integrals_in_map_collector
rc = f77_zmq_recv( zmq_socket_pull, task_id, 4, 0)
! Activate if zmq_socket_pull is a REP
! rc = f77_zmq_send( zmq_socket_pull, 0, 4, 0)
! if (rc /= 4) then
! print *, irp_here, ' : f77_zmq_send (zmq_socket_pull,...'
! stop 'error'
! endif
IRP_IF ZMQ_PUSH
IRP_ELSE
rc = f77_zmq_send( zmq_socket_pull, 0, 4, 0)
if (rc /= 4) then
print *, irp_here, ' : f77_zmq_send (zmq_socket_pull,...'
stop 'error'
endif
IRP_ENDIF
call insert_into_ao_integrals_map(n_integrals,buffer_i,buffer_value)

10
src/ZMQ/utils.irp.f
View File

@ -232,8 +232,11 @@ function new_zmq_pull_socket()
if (zmq_context == 0_ZMQ_PTR) then
stop 'zmq_context is uninitialized'
endif
IRP_IF ZMQ_PUSH
new_zmq_pull_socket = f77_zmq_socket(zmq_context, ZMQ_PULL)
! new_zmq_pull_socket = f77_zmq_socket(zmq_context, ZMQ_REP)
IRP_ELSE
new_zmq_pull_socket = f77_zmq_socket(zmq_context, ZMQ_REP)
IRP_ENDIF
call omp_unset_lock(zmq_lock)
if (new_zmq_pull_socket == 0_ZMQ_PTR) then
stop 'Unable to create zmq pull socket'
@ -309,8 +312,11 @@ function new_zmq_push_socket(thread)
if (zmq_context == 0_ZMQ_PTR) then
stop 'zmq_context is uninitialized'
endif
IRP_IF ZMQ_PUSH
new_zmq_push_socket = f77_zmq_socket(zmq_context, ZMQ_PUSH)
! new_zmq_push_socket = f77_zmq_socket(zmq_context, ZMQ_REQ)
IRP_ELSE
new_zmq_push_socket = f77_zmq_socket(zmq_context, ZMQ_REQ)
IRP_ENDIF
call omp_unset_lock(zmq_lock)
if (new_zmq_push_socket == 0_ZMQ_PTR) then
stop 'Unable to create zmq push socket'