LCPQ
/
qp2
mirror of https://github.com/QuantumPackage/qp2.git
9
1
Fork 0
Code Releases Activity

Merge branch 'develop' of https://github.com/QuantumPackage/qp2 into develop

This commit is contained in:
eginer 2019-02-06 15:27:59 +01:00
parent 6ffdccb26d 2016edd734
commit 4759309d11
27 changed files with 532 additions and 313 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

13
INSTALL.rst
View File

@ -210,6 +210,12 @@ Zlib
make
make install
With Debian or Ubuntu, you can use
.. code:: bash
sudo apt install zlib1g-dev
OCaml
@ -217,6 +223,13 @@ OCaml
*OCaml* is a general purpose programming language with an emphasis on expressiveness and safety.
* The following packages are required (Debian or Ubuntu):
.. code:: bash
sudo apt install libncurses5-dev pkg-config libgmp3-dev m4
* Download the installer of the OPAM package manager here :
`<https://raw.githubusercontent.com/ocaml/opam/master/shell/install.sh>`_
and move it in the :file:`${QP_ROOT}/external` directory

14
REPLACE
View File

@ -183,3 +183,17 @@ qp_name save_one_body_dm -r save_one_e_dm
qp_name ezfio_set_aux_quantities_data_one_e_alpha_dm_mo -r ezfio_set_aux_quantities_data_one_e_dm_alpha_mo
qp_name ezfio_set_aux_quantities_data_one_e_beta_dm_mo -r ezfio_set_aux_quantities_data_one_e_dm_beta_mo
qp_name two_electron_energy -r two_e_energy
qp_name do_mono_excitation -r do_single_excitation
qp_name get_mono_excitation -r get_single_excitation
qp_name get_mono_excitation_from_fock -r get_single_excitation_from_fock
qp_name is_connected_to_by_mono -r is_connected_to_by_single
qp_name connected_to_ref_by_mono -r connected_to_ref_by_single
qp_name mono_excitation_wee -r single_excitation_wee
qp_name get_mono_excitation_spin
qp_name get_mono_excitation_spin -r get_single_excitation_spin
qp_name get_excitation_degree_vector_mono -r get_excitation_degree_vector_single
qp_name get_excitation_degree_vector_mono_or_exchange -r get_excitation_degree_vector_single_or_exchange_or_exchange
qp_name get_excitation_degree_vector_single_or_exchange_or_exchange -r get_excitation_degree_vector_single_or_exchange
qp_name get_excitation_degree_vector_mono_or_exchange_verbose -r get_excitation_degree_vector_single_or_exchange_verbose
qp_name i_h_j_mono_spin -r i_h_j_single_spin
qp_name i_Wee_j_mono -r i_Wee_j_single

4
TODO
View File

@ -49,8 +49,6 @@ Refaire les benchmarks
# Documentation de qpsh
# Documentation de /etc
# Extrapolation qui prend aussi en compe la variance? a tester
Parler dans le papier de rPT2
# Toto
Re-design de qp command
@ -58,3 +56,5 @@ Re-design de qp command
Doc: plugins et qp_plugins
Ajouter les symetries dans devel
Compiler ezfio avec openmp

19
configure vendored
View File

@ -361,23 +361,24 @@ EOF
done
source quantum_package.rc
NINJA=$(find_exe ninja)
if [[ ${NINJA} = $(not_found) ]] ; then
error "Ninja is not installed."
error "Ninja (ninja) is not installed."
fail
fi
IRPF90=$(find_exe irpf90)
if [[ ${IRPF90} = $(not_found) ]] ; then
error "IRPf90 is not installed."
error "IRPf90 (irpf90) is not installed."
fail
fi
ZEROMQ=$(find_lib -lzmq)
if [[ ${ZEROMQ} = $(not_found) ]] ; then
error "ZeroMQ is not installed."
error "ZeroMQ (zeromq) is not installed."
fail
fi
@ -395,31 +396,31 @@ fi
OCAML=$(find_exe ocaml)
if [[ ${OCAML} = $(not_found) ]] ; then
error "OCaml compiler is not installed."
error "OCaml (ocaml) compiler is not installed."
fail
fi
EZFIO=$(find_dir "${QP_ROOT}"/external/ezfio)
if [[ ${EZFIO} = $(not_found) ]] ; then
error "EZFIO is not installed."
error "EZFIO (ezfio) is not installed."
fail
fi
ZLIB=$(find_lib -lz)
if [[ ${ZLIB} = $(not_found) ]] ; then
error "Zlib is not installed."
error "Zlib (zlib) is not installed."
fail
fi
DOCOPT=$(find_python_lib docopt)
if [[ ${DOCOPT} = $(not_found) ]] ; then
error "docopt is not installed."
error "docopt (docopt) is not installed."
fail
fi
RESULTSFILE=$(find_python_lib resultsFile)
if [[ ${RESULTSFILE} = $(not_found) ]] ; then
error "resultsFile is not installed."
error "resultsFile (resultsFile) is not installed."
fail
fi

69
src/cipsi/pt2_stoch_routines.irp.f
View File

@ -10,13 +10,20 @@ END_PROVIDER
&BEGIN_PROVIDER [ integer, pt2_n_tasks_max ]
implicit none
logical, external :: testTeethBuilding
integer :: i
integer :: e
e = elec_num - n_core_orb * 2
pt2_n_tasks_max = 1+min((e*(e-1)), int(dsqrt(dble(N_det_selectors)))/10)
do i=1,N_det_generators
pt2_F(i) = 1 + int(dble(pt2_n_tasks_max)*dabs(maxval(psi_coef_sorted_gen(i,:)))**(0.75d0))
integer :: i,j
pt2_n_tasks_max = elec_beta_num*elec_beta_num + elec_alpha_num*elec_beta_num - n_core_orb*2
pt2_n_tasks_max = min(pt2_n_tasks_max,1+N_det_generators/10000)
call write_int(6,pt2_n_tasks_max,'pt2_n_tasks_max')
pt2_F(:) = int(sqrt(float(pt2_n_tasks_max)))
do i=1,pt2_n_0(1+pt2_N_teeth/4)
pt2_F(i) = pt2_n_tasks_max
enddo
do i=1+pt2_n_0(pt2_N_teeth-pt2_N_teeth/4), N_det_generators
pt2_F(i) = 1
enddo
END_PROVIDER
BEGIN_PROVIDER [ integer, pt2_N_teeth ]
@ -54,17 +61,16 @@ logical function testTeethBuilding(minF, N)
allocate(tilde_w(N_det_generators), tilde_cW(0:N_det_generators))
do i=1,N_det_generators
tilde_w(i) = psi_coef_sorted_gen(i,pt2_stoch_istate)**2 !+ 1.d-20
enddo
double precision :: norm
norm = 0.d0
double precision :: norm
do i=N_det_generators,1,-1
norm += tilde_w(i)
tilde_w(i) = psi_coef_sorted_gen(i,pt2_stoch_istate) * &
psi_coef_sorted_gen(i,pt2_stoch_istate)
norm = norm + tilde_w(i)
enddo
tilde_w(:) = tilde_w(:) / norm
f = 1.d0/norm
tilde_w(:) = tilde_w(:) * f
tilde_cW(0) = -1.d0
do i=1,N_det_generators
@ -74,10 +80,14 @@ logical function testTeethBuilding(minF, N)
n0 = 0
testTeethBuilding = .false.
double precision :: f
integer :: minFN
minFN = N_det_generators - minF * N
f = 1.d0/dble(N)
do
u0 = tilde_cW(n0)
r = tilde_cW(n0 + minF)
Wt = (1d0 - u0) / dble(N)
Wt = (1d0 - u0) * f
if (dabs(Wt) <= 1.d-3) then
return
endif
@ -86,7 +96,8 @@ logical function testTeethBuilding(minF, N)
return
end if
n0 += 1
if(N_det_generators - n0 < minF * N) then
! if(N_det_generators - n0 < minF * N) then
if(n0 > minFN) then
return
end if
end do
@ -103,7 +114,6 @@ subroutine ZMQ_pt2(E, pt2,relative_error, error, variance, norm, N_in)
integer(ZMQ_PTR) :: zmq_to_qp_run_socket, zmq_socket_pull
integer, intent(in) :: N_in
integer, external :: omp_get_thread_num
double precision, intent(in) :: relative_error, E(N_states)
double precision, intent(out) :: pt2(N_states),error(N_states)
double precision, intent(out) :: variance(N_states),norm(N_states)
@ -111,7 +121,6 @@ subroutine ZMQ_pt2(E, pt2,relative_error, error, variance, norm, N_in)
integer :: i, N
double precision, external :: omp_get_wtime
double precision :: state_average_weight_save(N_states), w(N_states,4)
integer(ZMQ_PTR), external :: new_zmq_to_qp_run_socket
type(selection_buffer) :: b
@ -120,9 +129,9 @@ subroutine ZMQ_pt2(E, pt2,relative_error, error, variance, norm, N_in)
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 psi_selectors_coef_transp psi_det_sorted
PROVIDE psi_det_hii
PROVIDE psi_det_hii N_generators_bitmask
if (s2_eig) then
if (h0_type == 'SOP') then
PROVIDE psi_occ_pattern_hii det_to_occ_pattern
endif
@ -136,6 +145,10 @@ subroutine ZMQ_pt2(E, pt2,relative_error, error, variance, norm, N_in)
N = max(N_in,1) * N_states
state_average_weight_save(:) = state_average_weight(:)
if (int(N,8)*2_8 > huge(1)) then
print *, irp_here, ': integer too large'
stop -1
endif
call create_selection_buffer(N, N*2, b)
ASSERT (associated(b%det))
ASSERT (associated(b%val))
@ -279,6 +292,7 @@ subroutine ZMQ_pt2(E, pt2,relative_error, error, variance, norm, N_in)
print '(A)', ' Samples Energy Stat. Err Variance Norm Seconds '
print '(A)', '========== ================= =========== =============== =============== ================='
PROVIDE global_selection_buffer
!$OMP PARALLEL DEFAULT(shared) NUM_THREADS(nproc_target+1) &
!$OMP PRIVATE(i)
i = omp_get_thread_num()
@ -326,6 +340,7 @@ subroutine pt2_slave_inproc(i)
implicit none
integer, intent(in) :: i
PROVIDE global_selection_buffer
call run_pt2_slave(1,i,pt2_e0_denominator)
end
@ -359,7 +374,6 @@ subroutine pt2_collector(zmq_socket_pull, E, relative_error, pt2, error, varianc
integer, allocatable :: task_id(:)
integer, allocatable :: index(:)
double precision, external :: omp_get_wtime
double precision :: v, x, x2, x3, avg, avg2, avg3, eqt, E0, v0, n0
double precision :: time, time1, time0
@ -425,7 +439,6 @@ subroutine pt2_collector(zmq_socket_pull, E, relative_error, pt2, error, varianc
stop_now = .false.
do while (n <= N_det_generators)
if(f(pt2_J(n)) == 0) then
!print *, 'f(pt2_J(n)) == 0'
d(pt2_J(n)) = .true.
do while(d(U+1))
U += 1
@ -478,6 +491,7 @@ subroutine pt2_collector(zmq_socket_pull, E, relative_error, pt2, error, varianc
pt2(pt2_stoch_istate) = avg
variance(pt2_stoch_istate) = avg2
norm(pt2_stoch_istate) = avg3
call wall_time(time)
! 1/(N-1.5) : see Brugger, The American Statistician (23) 4 p. 32 (1969)
if(c > 2) then
eqt = dabs((S2(t) / c) - (S(t)/c)**2) ! dabs for numerical stability
@ -498,7 +512,6 @@ subroutine pt2_collector(zmq_socket_pull, E, relative_error, pt2, error, varianc
endif
endif
endif
call wall_time(time)
end if
n += 1
else if(more == 0) then
@ -506,21 +519,21 @@ subroutine pt2_collector(zmq_socket_pull, E, relative_error, pt2, error, varianc
else
call pull_pt2_results(zmq_socket_pull, index, eI_task, vI_task, nI_task, task_id, n_tasks, b2)
if (zmq_delete_tasks_async_send(zmq_to_qp_run_socket,task_id,n_tasks,sending) == -1) then
stop 'Unable to delete tasks'
stop 'PT2: Unable to delete tasks (send)'
endif
do i=1,n_tasks
eI(:, index(i)) += eI_task(:,i)
vI(:, index(i)) += vI_task(:,i)
nI(:, index(i)) += nI_task(:,i)
eI(1:N_states, index(i)) += eI_task(1:N_states,i)
vI(1:N_states, index(i)) += vI_task(1:N_states,i)
nI(1:N_states, index(i)) += nI_task(1:N_states,i)
f(index(i)) -= 1
end do
do i=1, b2%cur
call add_to_selection_buffer(b, b2%det(1,1,i), b2%val(i))
! We assume the pulled buffer is sorted
if (b2%val(i) > b%mini) exit
call add_to_selection_buffer(b, b2%det(1,1,i), b2%val(i))
end do
if (zmq_delete_tasks_async_recv(zmq_to_qp_run_socket,more,sending) == -1) then
stop 'Unable to delete tasks'
stop 'PT2: Unable to delete tasks (recv)'
endif
end if
end do

414
src/cipsi/run_pt2_slave.irp.f
View File

@ -1,6 +1,46 @@
use omp_lib
use selection_types
use f77_zmq
BEGIN_PROVIDER [ integer(omp_lock_kind), global_selection_buffer_lock ]
use omp_lib
implicit none
BEGIN_DOC
! Global buffer for the OpenMP selection
END_DOC
call omp_init_lock(global_selection_buffer_lock)
END_PROVIDER
BEGIN_PROVIDER [ type(selection_buffer), global_selection_buffer ]
use omp_lib
implicit none
BEGIN_DOC
! Global buffer for the OpenMP selection
END_DOC
call omp_set_lock(global_selection_buffer_lock)
call delete_selection_buffer(global_selection_buffer)
call create_selection_buffer(N_det_generators, 2*N_det_generators, &
global_selection_buffer)
call omp_unset_lock(global_selection_buffer_lock)
END_PROVIDER
subroutine run_pt2_slave(thread,iproc,energy)
use f77_zmq
use selection_types
use selection_types
use f77_zmq
implicit none
double precision, intent(in) :: energy(N_states_diag)
integer, intent(in) :: thread, iproc
if (N_det > nproc*(elec_alpha_num * (mo_num-elec_alpha_num))**2) then
call run_pt2_slave_large(thread,iproc,energy)
else
call run_pt2_slave_small(thread,iproc,energy)
endif
end
subroutine run_pt2_slave_small(thread,iproc,energy)
use selection_types
use f77_zmq
implicit none
double precision, intent(in) :: energy(N_states_diag)
@ -24,14 +64,9 @@ subroutine run_pt2_slave(thread,iproc,energy)
integer :: n_tasks, k, N
integer, allocatable :: i_generator(:), subset(:)
double precision :: rss
double precision, external :: memory_of_double, memory_of_int
integer :: bsize ! Size of selection buffers
logical :: sending
rss = memory_of_int(pt2_n_tasks_max)*67.d0
rss += memory_of_double(pt2_n_tasks_max)*(N_states*3)
call check_mem(rss,irp_here)
! logical :: sending
allocate(task_id(pt2_n_tasks_max), task(pt2_n_tasks_max))
allocate(pt2(N_states,pt2_n_tasks_max), i_generator(pt2_n_tasks_max), subset(pt2_n_tasks_max))
@ -52,9 +87,8 @@ subroutine run_pt2_slave(thread,iproc,energy)
buffer_ready = .False.
n_tasks = 1
sending = .False.
! sending = .False.
done = .False.
n_tasks = 1
do while (.not.done)
n_tasks = max(1,n_tasks)
@ -79,7 +113,130 @@ subroutine run_pt2_slave(thread,iproc,energy)
call create_selection_buffer(bsize, bsize*2, b)
buffer_ready = .True.
else
ASSERT (N == b%N)
ASSERT (b%N == bsize)
endif
double precision :: time0, time1
call wall_time(time0)
do k=1,n_tasks
pt2(:,k) = 0.d0
variance(:,k) = 0.d0
norm(:,k) = 0.d0
b%cur = 0
!double precision :: time2
!call wall_time(time2)
call select_connected(i_generator(k),energy,pt2(1,k),variance(1,k),norm(1,k),b,subset(k),pt2_F(i_generator(k)))
!call wall_time(time1)
!print *, i_generator(1), time1-time2, n_tasks, pt2_F(i_generator(1))
enddo
call wall_time(time1)
!print *, '-->', i_generator(1), time1-time0, n_tasks
integer, external :: tasks_done_to_taskserver
if (tasks_done_to_taskserver(zmq_to_qp_run_socket,worker_id,task_id,n_tasks) == -1) then
done = .true.
endif
call sort_selection_buffer(b)
call push_pt2_results(zmq_socket_push, i_generator, pt2, variance, norm, b, task_id, n_tasks)
b%cur=0
! ! Try to adjust n_tasks around nproc/2 seconds per job
! n_tasks = min(2*n_tasks,int( dble(n_tasks * nproc/2) / (time1 - time0 + 1.d0)))
n_tasks = 1
end do
integer, external :: disconnect_from_taskserver
do i=1,300
if (disconnect_from_taskserver(zmq_to_qp_run_socket,worker_id) /= -2) exit
call sleep(1)
print *, 'Retry disconnect...'
end do
call end_zmq_push_socket(zmq_socket_push,thread)
call end_zmq_to_qp_run_socket(zmq_to_qp_run_socket)
if (buffer_ready) then
call delete_selection_buffer(b)
endif
end subroutine
subroutine run_pt2_slave_large(thread,iproc,energy)
use selection_types
use f77_zmq
implicit none
double precision, intent(in) :: energy(N_states_diag)
integer, intent(in) :: thread, iproc
integer :: rc, i
integer :: worker_id, ctask, ltask
character*(512), allocatable :: task(:)
integer, allocatable :: task_id(:)
integer(ZMQ_PTR),external :: new_zmq_to_qp_run_socket
integer(ZMQ_PTR) :: zmq_to_qp_run_socket
integer(ZMQ_PTR), external :: new_zmq_push_socket
integer(ZMQ_PTR) :: zmq_socket_push
type(selection_buffer) :: b
logical :: done, buffer_ready
double precision,allocatable :: pt2(:,:), variance(:,:), norm(:,:)
integer :: n_tasks, k, N
integer, allocatable :: i_generator(:), subset(:)
integer :: bsize ! Size of selection buffers
logical :: sending
PROVIDE global_selection_buffer global_selection_buffer_lock
allocate(task_id(pt2_n_tasks_max), task(pt2_n_tasks_max))
allocate(pt2(N_states,pt2_n_tasks_max), i_generator(pt2_n_tasks_max), subset(pt2_n_tasks_max))
allocate(variance(N_states,pt2_n_tasks_max))
allocate(norm(N_states,pt2_n_tasks_max))
zmq_to_qp_run_socket = new_zmq_to_qp_run_socket()
integer, external :: connect_to_taskserver
if (connect_to_taskserver(zmq_to_qp_run_socket,worker_id,thread) == -1) then
call end_zmq_to_qp_run_socket(zmq_to_qp_run_socket)
return
endif
zmq_socket_push = new_zmq_push_socket(thread)
b%N = 0
buffer_ready = .False.
n_tasks = 1
sending = .False.
done = .False.
do while (.not.done)
n_tasks = max(1,n_tasks)
n_tasks = min(pt2_n_tasks_max,n_tasks)
integer, external :: get_tasks_from_taskserver
if (get_tasks_from_taskserver(zmq_to_qp_run_socket,worker_id, task_id, task, n_tasks) == -1) then
exit
endif
done = task_id(n_tasks) == 0
if (done) then
n_tasks = n_tasks-1
endif
if (n_tasks == 0) exit
do k=1,n_tasks
read (task(k),*) subset(k), i_generator(k), N
enddo
if (b%N == 0) then
! Only first time
bsize = min(N, (elec_alpha_num * (mo_num-elec_alpha_num))**2)
call create_selection_buffer(bsize, bsize*2, b)
buffer_ready = .True.
else
ASSERT (b%N == bsize)
endif
double precision :: time0, time1
@ -104,11 +261,23 @@ subroutine run_pt2_slave(thread,iproc,energy)
endif
call sort_selection_buffer(b)
call push_pt2_results_async_recv(zmq_socket_push,b%mini,sending)
call push_pt2_results_async_send(zmq_socket_push, i_generator, pt2, variance, norm, b, task_id, n_tasks,sending)
call omp_set_lock(global_selection_buffer_lock)
global_selection_buffer%mini = b%mini
call merge_selection_buffers(b,global_selection_buffer)
b%cur=0
call omp_unset_lock(global_selection_buffer_lock)
if ( iproc == 1 ) then
call omp_set_lock(global_selection_buffer_lock)
call push_pt2_results_async_send(zmq_socket_push, i_generator, pt2, variance, norm, global_selection_buffer, task_id, n_tasks,sending)
global_selection_buffer%cur = 0
call omp_unset_lock(global_selection_buffer_lock)
else
call push_pt2_results_async_send(zmq_socket_push, i_generator, pt2, variance, norm, b, task_id, n_tasks,sending)
endif
! Try to adjust n_tasks around nproc/2 seconds per job
n_tasks = min(2*n_tasks,int( dble(n_tasks * nproc/2) / (time1 - time0 + 1.d0)))
! ! Try to adjust n_tasks around nproc/2 seconds per job
! n_tasks = min(2*n_tasks,int( dble(n_tasks * nproc/2) / (time1 - time0 + 1.d0)))
n_tasks = 1
end do
call push_pt2_results_async_recv(zmq_socket_push,b%mini,sending)
@ -124,12 +293,13 @@ subroutine run_pt2_slave(thread,iproc,energy)
if (buffer_ready) then
call delete_selection_buffer(b)
endif
FREE global_selection_buffer
end subroutine
subroutine push_pt2_results(zmq_socket_push, index, pt2, variance, norm, b, task_id, n_tasks)
use f77_zmq
use selection_types
use selection_types
use f77_zmq
implicit none
integer(ZMQ_PTR), intent(in) :: zmq_socket_push
@ -138,99 +308,17 @@ subroutine push_pt2_results(zmq_socket_push, index, pt2, variance, norm, b, task
double precision, intent(in) :: norm(N_states,n_tasks)
integer, intent(in) :: n_tasks, index(n_tasks), task_id(n_tasks)
type(selection_buffer), intent(inout) :: b
integer :: rc
rc = f77_zmq_send( zmq_socket_push, n_tasks, 4, ZMQ_SNDMORE)
if (rc == -1) then
return
else if(rc /= 4) then
stop 'push'
endif
rc = f77_zmq_send( zmq_socket_push, index, 4*n_tasks, ZMQ_SNDMORE)
if (rc == -1) then
return
else if(rc /= 4*n_tasks) then
stop 'push'
endif
rc = f77_zmq_send( zmq_socket_push, pt2, 8*N_states*n_tasks, ZMQ_SNDMORE)
if (rc == -1) then
return
else if(rc /= 8*N_states*n_tasks) then
stop 'push'
endif
rc = f77_zmq_send( zmq_socket_push, variance, 8*N_states*n_tasks, ZMQ_SNDMORE)
if (rc == -1) then
return
else if(rc /= 8*N_states*n_tasks) then
stop 'push'
endif
rc = f77_zmq_send( zmq_socket_push, norm, 8*N_states*n_tasks, ZMQ_SNDMORE)
if (rc == -1) then
return
else if(rc /= 8*N_states*n_tasks) then
stop 'push'
endif
rc = f77_zmq_send( zmq_socket_push, task_id, n_tasks*4, ZMQ_SNDMORE)
if (rc == -1) then
return
else if(rc /= 4*n_tasks) then
stop 'push'
endif
rc = f77_zmq_send( zmq_socket_push, b%cur, 4, ZMQ_SNDMORE)
if (rc == -1) then
return
else if(rc /= 4) then
stop 'push'
endif
rc = f77_zmq_send( zmq_socket_push, b%val, 8*b%cur, ZMQ_SNDMORE)
if (rc == -1) then
return
else if(rc /= 8*b%cur) then
stop 'push'
endif
rc = f77_zmq_send( zmq_socket_push, b%det, bit_kind*N_int*2*b%cur, 0)
if (rc == -1) then
return
else if(rc /= N_int*2*8*b%cur) then
stop 'push'
endif
! Activate is zmq_socket_push is a REQ
IRP_IF ZMQ_PUSH
IRP_ELSE
character*(2) :: ok
rc = f77_zmq_recv( zmq_socket_push, ok, 2, 0)
if (rc == -1) then
return
else if ((rc /= 2).and.(ok(1:2) /= 'ok')) then
print *, irp_here//': error in receiving ok'
stop -1
endif
IRP_ENDIF
logical :: sending
sending = .False.
call push_pt2_results_async_send(zmq_socket_push, index, pt2, variance, norm, b, task_id, n_tasks, sending)
call push_pt2_results_async_recv(zmq_socket_push, b%mini, sending)
end subroutine
subroutine push_pt2_results_async_send(zmq_socket_push, index, pt2, variance, norm, b, task_id, n_tasks, sending)
use f77_zmq
use selection_types
use selection_types
use f77_zmq
implicit none
integer(ZMQ_PTR), intent(in) :: zmq_socket_push
@ -241,6 +329,7 @@ subroutine push_pt2_results_async_send(zmq_socket_push, index, pt2, variance, no
type(selection_buffer), intent(inout) :: b
logical, intent(inout) :: sending
integer :: rc
integer*8 :: rc8
if (sending) then
print *, irp_here, ': sending is true'
@ -250,6 +339,8 @@ subroutine push_pt2_results_async_send(zmq_socket_push, index, pt2, variance, no
rc = f77_zmq_send( zmq_socket_push, n_tasks, 4, ZMQ_SNDMORE)
if (rc == -1) then
print *, irp_here, ': error sending result'
stop 1
return
else if(rc /= 4) then
stop 'push'
@ -258,6 +349,8 @@ subroutine push_pt2_results_async_send(zmq_socket_push, index, pt2, variance, no
rc = f77_zmq_send( zmq_socket_push, index, 4*n_tasks, ZMQ_SNDMORE)
if (rc == -1) then
print *, irp_here, ': error sending result'
stop 2
return
else if(rc /= 4*n_tasks) then
stop 'push'
@ -266,6 +359,8 @@ subroutine push_pt2_results_async_send(zmq_socket_push, index, pt2, variance, no
rc = f77_zmq_send( zmq_socket_push, pt2, 8*N_states*n_tasks, ZMQ_SNDMORE)
if (rc == -1) then
print *, irp_here, ': error sending result'
stop 3
return
else if(rc /= 8*N_states*n_tasks) then
stop 'push'
@ -274,6 +369,8 @@ subroutine push_pt2_results_async_send(zmq_socket_push, index, pt2, variance, no
rc = f77_zmq_send( zmq_socket_push, variance, 8*N_states*n_tasks, ZMQ_SNDMORE)
if (rc == -1) then
print *, irp_here, ': error sending result'
stop 4
return
else if(rc /= 8*N_states*n_tasks) then
stop 'push'
@ -282,6 +379,8 @@ subroutine push_pt2_results_async_send(zmq_socket_push, index, pt2, variance, no
rc = f77_zmq_send( zmq_socket_push, norm, 8*N_states*n_tasks, ZMQ_SNDMORE)
if (rc == -1) then
print *, irp_here, ': error sending result'
stop 5
return
else if(rc /= 8*N_states*n_tasks) then
stop 'push'
@ -290,40 +389,63 @@ subroutine push_pt2_results_async_send(zmq_socket_push, index, pt2, variance, no
rc = f77_zmq_send( zmq_socket_push, task_id, n_tasks*4, ZMQ_SNDMORE)
if (rc == -1) then
print *, irp_here, ': error sending result'
stop 6
return
else if(rc /= 4*n_tasks) then
stop 'push'
endif
rc = f77_zmq_send( zmq_socket_push, b%cur, 4, ZMQ_SNDMORE)
if (rc == -1) then
return
else if(rc /= 4) then
stop 'push'
endif
if (b%cur == 0) then
rc = f77_zmq_send( zmq_socket_push, b%cur, 4, 0)
if (rc == -1) then
print *, irp_here, ': error sending result'
stop 7
return
else if(rc /= 4) then
stop 'push'
endif
else
rc = f77_zmq_send( zmq_socket_push, b%cur, 4, ZMQ_SNDMORE)
if (rc == -1) then
print *, irp_here, ': error sending result'
stop 7
return
else if(rc /= 4) then
stop 'push'
endif
rc = f77_zmq_send( zmq_socket_push, b%val, 8*b%cur, ZMQ_SNDMORE)
if (rc == -1) then
return
else if(rc /= 8*b%cur) then
stop 'push'
endif
rc8 = f77_zmq_send8( zmq_socket_push, b%val, 8_8*int(b%cur,8), ZMQ_SNDMORE)
if (rc8 == -1_8) then
print *, irp_here, ': error sending result'
stop 8
return
else if(rc8 /= 8_8*int(b%cur,8)) then
stop 'push'
endif
rc = f77_zmq_send( zmq_socket_push, b%det, bit_kind*N_int*2*b%cur, 0)
if (rc == -1) then
return
else if(rc /= N_int*2*8*b%cur) then
stop 'push'
rc8 = f77_zmq_send8( zmq_socket_push, b%det, int(bit_kind*N_int*2,8)*int(b%cur,8), 0)
if (rc8 == -1_8) then
print *, irp_here, ': error sending result'
stop 9
return
else if(rc8 /= int(N_int*2*8,8)*int(b%cur,8)) then
stop 'push'
endif
endif
end subroutine
subroutine push_pt2_results_async_recv(zmq_socket_push,mini,sending)
use f77_zmq
use selection_types
use selection_types
use f77_zmq
implicit none
integer(ZMQ_PTR), intent(in) :: zmq_socket_push
@ -339,12 +461,22 @@ IRP_ELSE
character*(2) :: ok
rc = f77_zmq_recv( zmq_socket_push, ok, 2, 0)
if (rc == -1) then
print *, irp_here, ': error sending result'
stop 10
return
else if ((rc /= 2).and.(ok(1:2) /= 'ok')) then
print *, irp_here//': error in receiving ok'
stop -1
endif
rc = f77_zmq_recv( zmq_socket_push, mini, 8, 0)
if (rc == -1) then
print *, irp_here, ': error sending result'
stop 11
return
else if (rc /= 8) then
print *, irp_here//': error in receiving mini'
stop 12
endif
IRP_ENDIF
sending = .False.
end subroutine
@ -352,8 +484,8 @@ end subroutine
subroutine pull_pt2_results(zmq_socket_pull, index, pt2, variance, norm, task_id, n_tasks, b)
use f77_zmq
use selection_types
use selection_types
use f77_zmq
implicit none
integer(ZMQ_PTR), intent(in) :: zmq_socket_pull
double precision, intent(inout) :: pt2(N_states,*)
@ -363,6 +495,7 @@ subroutine pull_pt2_results(zmq_socket_pull, index, pt2, variance, norm, task_id
integer, intent(out) :: index(*)
integer, intent(out) :: n_tasks, task_id(*)
integer :: rc, rn, i
integer*8 :: rc8
rc = f77_zmq_recv( zmq_socket_pull, n_tasks, 4, 0)
if (rc == -1) then
@ -420,22 +553,25 @@ subroutine pull_pt2_results(zmq_socket_pull, index, pt2, variance, norm, task_id
stop 'pull'
endif
rc = f77_zmq_recv( zmq_socket_pull, b%val, 8*b%cur, 0)
if (rc == -1) then
n_tasks = 1
task_id(1) = 0
else if(rc /= 8*b%cur) then
stop 'pull'
endif
if (b%cur > 0) then
rc = f77_zmq_recv( zmq_socket_pull, b%det, bit_kind*N_int*2*b%cur, 0)
if (rc == -1) then
n_tasks = 1
task_id(1) = 0
else if(rc /= N_int*2*8*b%cur) then
stop 'pull'
endif
rc8 = f77_zmq_recv8( zmq_socket_pull, b%val, 8_8*int(b%cur,8), 0)
if (rc8 == -1_8) then
n_tasks = 1
task_id(1) = 0
else if(rc8 /= 8_8*int(b%cur,8)) then
stop 'pull'
endif
rc8 = f77_zmq_recv8( zmq_socket_pull, b%det, int(bit_kind*N_int*2,8)*int(b%cur,8), 0)
if (rc8 == -1_8) then
n_tasks = 1
task_id(1) = 0
else if(rc8 /= int(N_int*2*8,8)*int(b%cur,8)) then
stop 'pull'
endif
endif
! Activate is zmq_socket_pull is a REP
IRP_IF ZMQ_PUSH

12
src/cipsi/selection_buffer.irp.f
View File

@ -37,6 +37,11 @@ subroutine delete_selection_buffer(b)
if (associated(b%val)) then
deallocate(b%val)
endif
NULLIFY(b%det)
NULLIFY(b%val)
b%cur = 0
b%mini = 0.d0
b%N = 0
end
@ -69,7 +74,7 @@ subroutine merge_selection_buffers(b1, b2)
type(selection_buffer), intent(inout) :: b2
integer(bit_kind), pointer :: detmp(:,:,:)
double precision, pointer :: val(:)
integer :: i, i1, i2, k, nmwen
integer :: i, i1, i2, k, nmwen, sze
if (b1%cur == 0) return
do while (b1%val(b1%cur) > b2%mini)
b1%cur = b1%cur-1
@ -80,9 +85,10 @@ subroutine merge_selection_buffers(b1, b2)
nmwen = min(b1%N, b1%cur+b2%cur)
double precision :: rss
double precision, external :: memory_of_double
rss = memory_of_double(size(b1%val)) + 2*N_int*memory_of_double(size(b1%det,3))
sze = max(size(b1%val), size(b2%val))
rss = memory_of_double(sze) + 2*N_int*memory_of_double(sze)
call check_mem(rss,irp_here)
allocate( val(size(b1%val)), detmp(N_int, 2, size(b1%det,3)) )
allocate(val(sze), detmp(N_int, 2, sze))
i1=1
i2=1
do i=1,nmwen

49
src/cipsi/slave_cipsi.irp.f
View File

@ -232,6 +232,7 @@ subroutine run_slave_main
endif
IRP_ENDIF
IRP_IF MPI_DEBUG
call mpi_print('Entering OpenMP section')
IRP_ENDIF
@ -251,7 +252,7 @@ subroutine run_slave_main
+ 64.d0*pt2_n_tasks_max & ! task
+ 3.d0*pt2_n_tasks_max*N_states & ! pt2, variance, norm
+ 1.d0*pt2_n_tasks_max & ! i_generator, subset
+ 2.d0*(N_int*2.d0*ii+ ii) & ! selection buffer
+ 3.d0*(N_int*2.d0*ii+ ii) & ! selection buffer
+ 1.d0*(N_int*2.d0*ii+ ii) & ! sort selection buffer
+ 2.0d0*(ii) & ! preinteresting, interesting,
! prefullinteresting, fullinteresting
@ -272,29 +273,37 @@ subroutine run_slave_main
nproc_target = nproc_target - 1
enddo
if (N_det > 100000) then
if (mpi_master) then
print *, 'N_det', N_det
print *, 'N_det_generators', N_det_generators
print *, 'N_det_selectors', N_det_selectors
print *, 'pt2_e0_denominator', pt2_e0_denominator
print *, 'pt2_stoch_istate', pt2_stoch_istate
print *, 'state_average_weight', state_average_weight
print *, 'Number of threads', nproc_target
PROVIDE psi_det_hii
if (s2_eig) then
PROVIDE psi_occ_pattern_hii det_to_occ_pattern
if (mpi_master) then
print *, 'N_det', N_det
print *, 'N_det_generators', N_det_generators
print *, 'N_det_selectors', N_det_selectors
print *, 'pt2_e0_denominator', pt2_e0_denominator
print *, 'pt2_stoch_istate', pt2_stoch_istate
print *, 'state_average_weight', state_average_weight
print *, 'Number of threads', nproc_target
endif
endif
!$OMP PARALLEL PRIVATE(i) NUM_THREADS(nproc_target+1)
i = omp_get_thread_num()
call run_pt2_slave(0,i,pt2_e0_denominator)
!$OMP END PARALLEL
if (h0_type == 'SOP') then
PROVIDE det_to_occ_pattern
endif
PROVIDE global_selection_buffer
if (mpi_master) then
print *, 'Running PT2'
endif
!$OMP PARALLEL PRIVATE(i) NUM_THREADS(nproc_target+1)
i = omp_get_thread_num()
call run_pt2_slave(0,i,pt2_e0_denominator)
!$OMP END PARALLEL
FREE state_average_weight
print *, mpi_rank, ': PT2 done'
print *, '-------'
endif
endif
FREE state_average_weight
print *, mpi_rank, ': PT2 done'
IRP_IF MPI
call MPI_BARRIER(MPI_COMM_WORLD, ierr)

4
src/davidson/davidson_parallel.irp.f
View File

@ -359,7 +359,7 @@ subroutine davidson_collector(zmq_to_qp_run_socket, zmq_socket_pull, v0, s0, sze
do while (more == 1)
call davidson_pull_results(zmq_socket_pull, v_t, s_t, imin, imax, task_id)
if (zmq_delete_task_async_send(zmq_to_qp_run_socket,task_id,sending) == -1) then
stop 'Unable to delete task'
stop 'davidson: Unable to delete task (send)'
endif
do j=1,N_st
do i=imin,imax
@ -368,7 +368,7 @@ subroutine davidson_collector(zmq_to_qp_run_socket, zmq_socket_pull, v0, s0, sze
enddo
enddo
if (zmq_delete_task_async_recv(zmq_to_qp_run_socket,more,sending) == -1) then
stop 'Unable to delete task'
stop 'davidson: Unable to delete task (recv)'
endif
end do
deallocate(v_t,s_t)

4
src/davidson/diagonalize_ci.irp.f
View File

@ -202,6 +202,8 @@ subroutine diagonalize_CI
psi_coef(i,j) = CI_eigenvectors(i,j)
enddo
enddo
psi_energy(1:N_states) = CI_electronic_energy(1:N_states)
psi_s2(1:N_states) = CI_s2(1:N_states)
SOFT_TOUCH psi_coef CI_electronic_energy CI_energy CI_eigenvectors CI_s2
SOFT_TOUCH psi_coef CI_electronic_energy CI_energy CI_eigenvectors CI_s2 psi_energy psi_s2
end

20
src/davidson/u0_h_u0.irp.f
View File

@ -41,13 +41,27 @@ subroutine u_0_H_u_0(e_0,s_0,u_0,n,keys_tmp,Nint,N_st,sze)
double precision, allocatable :: v_0(:,:), s_vec(:,:), u_1(:,:)
double precision :: u_dot_u,u_dot_v,diag_H_mat_elem
integer :: i,j
integer :: i,j, istate
if ((n > 100000).and.distributed_davidson) then
allocate (v_0(n,N_states_diag),s_vec(n,N_states_diag), u_1(n,N_states_diag))
u_1(:,:) = 0.d0
u_1(1:n,1:N_st) = u_0(1:n,1:N_st)
call H_S2_u_0_nstates_zmq(v_0,s_vec,u_1,N_states_diag,n)
else if (n < n_det_max_full) then
allocate (v_0(n,N_st),s_vec(n,N_st), u_1(n,N_st))
v_0(:,:) = 0.d0
u_1(:,:) = 0.d0
s_vec(:,:) = 0.d0
u_1(1:n,1:N_st) = u_0(1:n,1:N_st)
do istate = 1,N_st
do j=1,n
do i=1,n
v_0(i,istate) = v_0(i,istate) + h_matrix_all_dets(i,j) * u_0(j,istate)
s_vec(i,istate) = s_vec(i,istate) + S2_matrix_all_dets(i,j) * u_0(j,istate)
enddo
enddo
enddo
else
allocate (v_0(n,N_st),s_vec(n,N_st),u_1(n,N_st))
u_1(:,:) = 0.d0
@ -469,7 +483,7 @@ compute_singles=.True.
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)
call i_h_j_single_spin( tmp_det, tmp_det2, $N_int, 1, hij)
!DIR$ LOOP COUNT AVG(4)
do l=1,N_st
@ -554,7 +568,7 @@ compute_singles=.True.
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)
call i_h_j_single_spin( tmp_det, tmp_det2, $N_int, 2, hij)
l_a = psi_bilinear_matrix_transp_order(l_b)
ASSERT (l_a <= N_det)
!DIR$ LOOP COUNT AVG(4)

4
src/davidson/u0_wee_u0.irp.f
View File

@ -304,7 +304,7 @@ subroutine H_S2_u_0_two_e_nstates_openmp_work_$N_int(v_t,s_t,u_t,N_st,sze,istart
ASSERT (lrow <= N_det_alpha_unique)
tmp_det2(1:$N_int,1) = psi_det_alpha_unique(1:$N_int, lrow)
call i_Wee_j_mono( tmp_det, tmp_det2, $N_int, 1, hij)
call i_Wee_j_single( tmp_det, tmp_det2, $N_int, 1, hij)
do l=1,N_st
v_t(l,k_a) = v_t(l,k_a) + hij * u_t(l,l_a)
@ -384,7 +384,7 @@ subroutine H_S2_u_0_two_e_nstates_openmp_work_$N_int(v_t,s_t,u_t,N_st,sze,istart
ASSERT (lcol <= N_det_beta_unique)
tmp_det2(1:$N_int,2) = psi_det_beta_unique (1:$N_int, lcol)
call i_Wee_j_mono( tmp_det, tmp_det2, $N_int, 2, hij)
call i_Wee_j_single( 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

18
src/determinants/connected_to_ref.irp.f
View File

@ -185,7 +185,7 @@ end
logical function is_connected_to_by_mono(key,keys,Nint,Ndet)
logical function is_connected_to_by_single(key,keys,Nint,Ndet)
use bitmasks
implicit none
BEGIN_DOC
@ -202,7 +202,7 @@ logical function is_connected_to_by_mono(key,keys,Nint,Ndet)
ASSERT (Nint > 0)
ASSERT (Nint == N_int)
is_connected_to_by_mono = .false.
is_connected_to_by_single = .false.
do i=1,Ndet
degree_x2 = popcnt(xor( key(1,1), keys(1,1,i))) + &
@ -214,7 +214,7 @@ logical function is_connected_to_by_mono(key,keys,Nint,Ndet)
if (degree_x2 > 2) then
cycle
else
is_connected_to_by_mono = .true.
is_connected_to_by_single = .true.
return
endif
enddo
@ -333,7 +333,7 @@ end
integer function connected_to_ref_by_mono(key,keys,Nint,N_past_in,Ndet)
integer function connected_to_ref_by_single(key,keys,Nint,N_past_in,Ndet)
use bitmasks
implicit none
BEGIN_DOC
@ -368,7 +368,7 @@ integer function connected_to_ref_by_mono(key,keys,Nint,N_past_in,Ndet)
ASSERT (Nint > 0)
ASSERT (Nint == N_int)
connected_to_ref_by_mono = 0
connected_to_ref_by_single = 0
N_past = max(1,N_past_in)
if (Nint == 1) then
@ -380,7 +380,7 @@ integer function connected_to_ref_by_mono(key,keys,Nint,N_past_in,Ndet)
else if (degree_x2 == 4)then
cycle
else if(degree_x2 == 2)then
connected_to_ref_by_mono = i
connected_to_ref_by_single = i
return
endif
enddo
@ -400,7 +400,7 @@ integer function connected_to_ref_by_mono(key,keys,Nint,N_past_in,Ndet)
else if (degree_x2 == 4)then
cycle
else if(degree_x2 == 2)then
connected_to_ref_by_mono = i
connected_to_ref_by_single = i
return
endif
enddo
@ -421,7 +421,7 @@ integer function connected_to_ref_by_mono(key,keys,Nint,N_past_in,Ndet)
else if (degree_x2 == 4)then
cycle
else if(degree_x2 == 2)then
connected_to_ref_by_mono = i
connected_to_ref_by_single = i
return
endif
enddo
@ -442,7 +442,7 @@ integer function connected_to_ref_by_mono(key,keys,Nint,N_past_in,Ndet)
else if (degree_x2 == 4)then
cycle
else if(degree_x2 == 2)then
connected_to_ref_by_mono = i
connected_to_ref_by_single = i
return
endif
enddo

2
src/determinants/create_excitations.irp.f
View File

@ -1,4 +1,4 @@
subroutine do_mono_excitation(key_in,i_hole,i_particle,ispin,i_ok)
subroutine do_single_excitation(key_in,i_hole,i_particle,ispin,i_ok)
implicit none
BEGIN_DOC
! Apply the single excitation operator : a^{dager}_(i_particle) a_(i_hole) of spin = ispin

4
src/determinants/density_matrix.irp.f
View File

@ -150,7 +150,7 @@ END_PROVIDER
call get_excitation_degree_spin(tmp_det(1,1),tmp_det2,degree,N_int)
if (degree == 1) then
exc = 0
call get_mono_excitation_spin(tmp_det(1,1),tmp_det2,exc,phase,N_int)
call get_single_excitation_spin(tmp_det(1,1),tmp_det2,exc,phase,N_int)
call decode_exc_spin(exc,h1,p1,h2,p2)
do m=1,N_states
ckl = psi_bilinear_matrix_values(k_a,m)*psi_bilinear_matrix_values(l,m) * phase
@ -206,7 +206,7 @@ END_PROVIDER
call get_excitation_degree_spin(tmp_det(1,2),tmp_det2,degree,N_int)
if (degree == 1) then
exc = 0
call get_mono_excitation_spin(tmp_det(1,2),tmp_det2,exc,phase,N_int)
call get_single_excitation_spin(tmp_det(1,2),tmp_det2,exc,phase,N_int)
call decode_exc_spin(exc,h1,p1,h2,p2)
do m=1,N_states
ckl = psi_bilinear_matrix_transp_values(k_b,m)*psi_bilinear_matrix_transp_values(l,m) * phase

4
src/determinants/example.irp.f
View File

@ -29,7 +29,7 @@ subroutine example_determinants
print*,'h1 --> p1 of spin s1'
print*,'i_ok == +1 : excitation is possible '
print*,'i_ok == -1 : excitation is NOT possible '
call do_mono_excitation(det_i,h1,p1,s1,i_ok)
call do_single_excitation(det_i,h1,p1,s1,i_ok)
print*,'h1,p1,s1,i_ok'
print*, h1,p1,s1,i_ok
if(i_ok == -1)then
@ -54,7 +54,7 @@ subroutine example_determinants
h1 = elec_alpha_num
p1 = elec_alpha_num + 1
s1 = 2
call do_mono_excitation(det_i,h1,p1,s1,i_ok)
call do_single_excitation(det_i,h1,p1,s1,i_ok)
print*,'h1,p1,s1,i_ok'
print*, h1,p1,s1,i_ok
call i_H_j(det_i,det_i,N_int,h0i)

4
src/determinants/h_apply.template.f
View File

@ -236,7 +236,7 @@ subroutine $subroutine_diexcOrg(key_in,key_mask,hole_1,particl_1,hole_2, particl
call bitstring_to_list_ab(particle_tmp,occ_particle_tmp,N_elec_in_key_part_2,N_int)
call bitstring_to_list_ab(hole_tmp,occ_hole_tmp,N_elec_in_key_hole_2,N_int)
! hole = a^(+)_j_a(ispin) a_i_a(ispin)|key_in> : mono exc :: orb(i_a,ispin) --> orb(j_a,ispin)
! hole = a^(+)_j_a(ispin) a_i_a(ispin)|key_in> : single exc :: orb(i_a,ispin) --> orb(j_a,ispin)
hole_save = hole
! Build array of the non-zero integrals of second excitation
@ -297,7 +297,7 @@ subroutine $subroutine_diexcOrg(key_in,key_mask,hole_1,particl_1,hole_2, particl
enddo
endif
! does all the mono excitations of the same spin
! does all the single excitations of the same spin
i=0
do kk = 1,N_elec_in_key_hole_2(ispin)
i_b = occ_hole_tmp(kk,ispin)

8
src/determinants/mono_excitations_bielec.irp.f → src/determinants/single_excitation_two_e.irp.f
View File

@ -1,5 +1,5 @@
use bitmasks
subroutine mono_excitation_wee(det_1,det_2,h,p,spin,phase,hij)
subroutine single_excitation_wee(det_1,det_2,h,p,spin,phase,hij)
use bitmasks
implicit none
integer,intent(in) :: h,p,spin
@ -79,7 +79,7 @@ BEGIN_PROVIDER [double precision, fock_wee_closed_shell, (mo_num, mo_num) ]
enddo
double precision :: array_coulomb(mo_num),array_exchange(mo_num)
call bitstring_to_list_ab(key_virt, occ_virt, n_occ_ab_virt, N_int)
! docc ---> virt mono excitations
! docc ---> virt single excitations
do i0 = 1, n_occ_ab(1)
i=occ(i0,1)
do j0 = 1, n_occ_ab_virt(1)
@ -97,7 +97,7 @@ BEGIN_PROVIDER [double precision, fock_wee_closed_shell, (mo_num, mo_num) ]
enddo
enddo
! virt ---> virt mono excitations
! virt ---> virt single excitations
do i0 = 1, n_occ_ab_virt(1)
i=occ_virt(i0,1)
do j0 = 1, n_occ_ab_virt(1)
@ -114,7 +114,7 @@ BEGIN_PROVIDER [double precision, fock_wee_closed_shell, (mo_num, mo_num) ]
enddo
enddo
! docc ---> docc mono excitations
! docc ---> docc single excitations
do i0 = 1, n_occ_ab(1)
i=occ(i0,1)
do j0 = 1, n_occ_ab(1)

8
src/determinants/single_excitations.irp.f
View File

@ -39,7 +39,7 @@ BEGIN_PROVIDER [double precision, fock_operator_closed_shell_ref_bitmask, (mo_nu
double precision, allocatable :: array_coulomb(:),array_exchange(:)
allocate (array_coulomb(mo_num),array_exchange(mo_num))
call bitstring_to_list_ab(key_virt, occ_virt, n_occ_ab_virt, N_int)
! docc ---> virt mono excitations
! docc ---> virt single excitations
do i0 = 1, n_occ_ab(1)
i=occ(i0,1)
do j0 = 1, n_occ_ab_virt(1)
@ -57,7 +57,7 @@ BEGIN_PROVIDER [double precision, fock_operator_closed_shell_ref_bitmask, (mo_nu
enddo
enddo
! virt ---> virt mono excitations
! virt ---> virt single excitations
do i0 = 1, n_occ_ab_virt(1)
i=occ_virt(i0,1)
do j0 = 1, n_occ_ab_virt(1)
@ -74,7 +74,7 @@ BEGIN_PROVIDER [double precision, fock_operator_closed_shell_ref_bitmask, (mo_nu
enddo
enddo
! docc ---> docc mono excitations
! docc ---> docc single excitations
do i0 = 1, n_occ_ab(1)
i=occ(i0,1)
do j0 = 1, n_occ_ab(1)
@ -94,7 +94,7 @@ BEGIN_PROVIDER [double precision, fock_operator_closed_shell_ref_bitmask, (mo_nu
END_PROVIDER
subroutine get_mono_excitation_from_fock(det_1,det_2,h,p,spin,phase,hij)
subroutine get_single_excitation_from_fock(det_1,det_2,h,p,spin,phase,hij)
use bitmasks
implicit none
integer,intent(in) :: h,p,spin

73
src/determinants/slater_rules.irp.f
View File

@ -93,7 +93,7 @@ subroutine get_excitation(det1,det2,exc,degree,phase,Nint)
return
case (1)
call get_mono_excitation(det1,det2,exc,phase,Nint)
call get_single_excitation(det1,det2,exc,phase,Nint)
return
case(0)
@ -336,7 +336,7 @@ end
subroutine get_mono_excitation(det1,det2,exc,phase,Nint)
subroutine get_single_excitation(det1,det2,exc,phase,Nint)
use bitmasks
implicit none
BEGIN_DOC
@ -499,7 +499,7 @@ subroutine i_H_j_s2(key_i,key_j,Nint,hij,s2)
select case (degree)
case (2)
call get_double_excitation(key_i,key_j,exc,phase,Nint)
! Mono alpha, mono beta
! Single alpha, single beta
if (exc(0,1,1) == 1) then
if ( (exc(1,1,1) == exc(1,2,2)).and.(exc(1,1,2) == exc(1,2,1)) ) then
s2 = -phase
@ -541,21 +541,21 @@ subroutine i_H_j_s2(key_i,key_j,Nint,hij,s2)
exc(1,2,2) ,mo_integrals_map) )
endif
case (1)
call get_mono_excitation(key_i,key_j,exc,phase,Nint)
call get_single_excitation(key_i,key_j,exc,phase,Nint)
!DIR$ FORCEINLINE
call bitstring_to_list_ab(key_i, occ, n_occ_ab, Nint)
! Mono alpha
! Single alpha
if (exc(0,1,1) == 1) then
m = exc(1,1,1)
p = exc(1,2,1)
spin = 1
! Mono beta
! Single beta
else
m = exc(1,1,2)
p = exc(1,2,2)
spin = 2
endif
call get_mono_excitation_from_fock(key_i,key_j,p,m,spin,phase,hij)
call get_single_excitation_from_fock(key_i,key_j,p,m,spin,phase,hij)
case (0)
double precision, external :: diag_S_mat_elem
@ -602,7 +602,7 @@ subroutine i_H_j(key_i,key_j,Nint,hij)
case (2)
call get_double_excitation(key_i,key_j,exc,phase,Nint)
if (exc(0,1,1) == 1) then
! Mono alpha, mono beta
! Single alpha, single beta
if(exc(1,1,1) == exc(1,2,2) )then
hij = phase * big_array_exchange_integrals(exc(1,1,1),exc(1,1,2),exc(1,2,1))
else if (exc(1,2,1) ==exc(1,1,2))then
@ -640,21 +640,21 @@ subroutine i_H_j(key_i,key_j,Nint,hij)
exc(1,2,2) ,mo_integrals_map) )
endif
case (1)
call get_mono_excitation(key_i,key_j,exc,phase,Nint)
call get_single_excitation(key_i,key_j,exc,phase,Nint)
!DIR$ FORCEINLINE
call bitstring_to_list_ab(key_i, occ, n_occ_ab, Nint)
if (exc(0,1,1) == 1) then
! Mono alpha
! Single alpha
m = exc(1,1,1)
p = exc(1,2,1)
spin = 1
else
! Mono beta
! Single beta
m = exc(1,1,2)
p = exc(1,2,2)
spin = 2
endif
call get_mono_excitation_from_fock(key_i,key_j,p,m,spin,phase,hij)
call get_single_excitation_from_fock(key_i,key_j,p,m,spin,phase,hij)
case (0)
hij = diag_H_mat_elem(key_i,Nint)
@ -703,7 +703,7 @@ subroutine i_H_j_verbose(key_i,key_j,Nint,hij,hmono,hdouble,phase)
case (2)
call get_double_excitation(key_i,key_j,exc,phase,Nint)
if (exc(0,1,1) == 1) then
! Mono alpha, mono beta
! Single alpha, single beta
hij = phase*get_two_e_integral( &
exc(1,1,1), &
exc(1,1,2), &
@ -736,12 +736,12 @@ subroutine i_H_j_verbose(key_i,key_j,Nint,hij,hmono,hdouble,phase)
exc(1,2,2) ,mo_integrals_map) )
endif
case (1)
call get_mono_excitation(key_i,key_j,exc,phase,Nint)
call get_single_excitation(key_i,key_j,exc,phase,Nint)
!DIR$ FORCEINLINE
call bitstring_to_list_ab(key_i, occ, n_occ_ab, Nint)
has_mipi = .False.
if (exc(0,1,1) == 1) then
! Mono alpha
! Single alpha
m = exc(1,1,1)
p = exc(1,2,1)
do k = 1, elec_alpha_num
@ -768,7 +768,7 @@ subroutine i_H_j_verbose(key_i,key_j,Nint,hij,hmono,hdouble,phase)
enddo
else
! Mono beta
! Single beta
m = exc(1,1,2)
p = exc(1,2,2)
do k = 1, elec_beta_num
@ -1060,7 +1060,7 @@ end
subroutine get_excitation_degree_vector_mono(key1,key2,degree,Nint,sze,idx)
subroutine get_excitation_degree_vector_single(key1,key2,degree,Nint,sze,idx)
use bitmasks
implicit none
BEGIN_DOC
@ -1154,7 +1154,7 @@ subroutine get_excitation_degree_vector_mono(key1,key2,degree,Nint,sze,idx)
end
subroutine get_excitation_degree_vector_mono_or_exchange(key1,key2,degree,Nint,sze,idx)
subroutine get_excitation_degree_vector_single_or_exchange(key1,key2,degree,Nint,sze,idx)
use bitmasks
implicit none
BEGIN_DOC
@ -1202,7 +1202,7 @@ subroutine get_excitation_degree_vector_mono_or_exchange(key1,key2,degree,Nint,s
enddo
else
print*, 'get_excitation_degree_vector_mono_or_exchange not yet implemented for N_int > 1 ...'
print*, 'get_excitation_degree_vector_single_or_exchange not yet implemented for N_int > 1 ...'
stop
endif
@ -1322,7 +1322,7 @@ subroutine get_excitation_degree_vector_double_alpha_beta(key1,key2,degree,Nint,
end
subroutine get_excitation_degree_vector_mono_or_exchange_verbose(key1,key2,degree,Nint,sze,idx)
subroutine get_excitation_degree_vector_single_or_exchange_verbose(key1,key2,degree,Nint,sze,idx)
use bitmasks
implicit none
BEGIN_DOC
@ -1635,7 +1635,7 @@ double precision function diag_H_mat_elem_fock(det_ref,det_pert,fock_diag_tmp,Ni
endif
else if (degree == 1) then
call get_mono_excitation(det_ref,det_pert,exc,phase,Nint)
call get_single_excitation(det_ref,det_pert,exc,phase,Nint)
call decode_exc(exc,1,h1,p1,h2,p2,s1,s2)
if (s1 == 1) then
diag_H_mat_elem_fock = E0 - fock_diag_tmp(1,h1) &
@ -1775,7 +1775,17 @@ subroutine ac_operator(iorb,ispin,key,hjj,Nint,na,nb)
integer :: other_spin
integer :: k,l,i
ASSERT (iorb > 0)
if (iorb < 1) then
print *, irp_here, 'iorb < 1'
print *, iorb, mo_num
stop -1
endif
if (iorb > mo_num) then
print *, irp_here, 'iorb > mo_num'
print *, iorb, mo_num
stop -1
endif
ASSERT (ispin > 0)
ASSERT (ispin < 3)
ASSERT (Nint > 0)
@ -1793,11 +1803,6 @@ subroutine ac_operator(iorb,ispin,key,hjj,Nint,na,nb)
key(k,ispin) = ibset(key(k,ispin),l)
other_spin = iand(ispin,1)+1
! if (iorb > mo_num) then
! print *, irp_here, 'iorb > mo_num'
! print *, iorb, mo_num
! stop -1
! endif
hjj = hjj + mo_one_e_integrals(iorb,iorb)
! Same spin
@ -1921,7 +1926,7 @@ subroutine get_excitation_spin(det1,det2,exc,degree,phase,Nint)
return
case (1)
call get_mono_excitation_spin(det1,det2,exc,phase,Nint)
call get_single_excitation_spin(det1,det2,exc,phase,Nint)
return
case(0)
@ -2093,7 +2098,7 @@ subroutine get_double_excitation_spin(det1,det2,exc,phase,Nint)
end
subroutine get_mono_excitation_spin(det1,det2,exc,phase,Nint)
subroutine get_single_excitation_spin(det1,det2,exc,phase,Nint)
use bitmasks
implicit none
BEGIN_DOC
@ -2169,7 +2174,7 @@ subroutine get_mono_excitation_spin(det1,det2,exc,phase,Nint)
enddo
end
subroutine i_H_j_mono_spin(key_i,key_j,Nint,spin,hij)
subroutine i_H_j_single_spin(key_i,key_j,Nint,spin,hij)
use bitmasks
implicit none
BEGIN_DOC
@ -2185,8 +2190,8 @@ subroutine i_H_j_mono_spin(key_i,key_j,Nint,spin,hij)
PROVIDE big_array_exchange_integrals mo_two_e_integrals_in_map
call get_mono_excitation_spin(key_i(1,spin),key_j(1,spin),exc,phase,Nint)
call get_mono_excitation_from_fock(key_i,key_j,exc(1,1),exc(1,2),spin,phase,hij)
call get_single_excitation_spin(key_i(1,spin),key_j(1,spin),exc,phase,Nint)
call get_single_excitation_from_fock(key_i,key_j,exc(1,1),exc(1,2),spin,phase,hij)
end
subroutine i_H_j_double_spin(key_i,key_j,Nint,hij)
@ -2235,8 +2240,8 @@ subroutine i_H_j_double_alpha_beta(key_i,key_j,Nint,hij)
PROVIDE big_array_exchange_integrals mo_two_e_integrals_in_map
call get_mono_excitation_spin(key_i(1,1),key_j(1,1),exc(0,1,1),phase,Nint)
call get_mono_excitation_spin(key_i(1,2),key_j(1,2),exc(0,1,2),phase2,Nint)
call get_single_excitation_spin(key_i(1,1),key_j(1,1),exc(0,1,1),phase,Nint)
call get_single_excitation_spin(key_i(1,2),key_j(1,2),exc(0,1,2),phase2,Nint)
phase = phase*phase2
if (exc(1,1,1) == exc(1,2,2)) then
hij = phase * big_array_exchange_integrals(exc(1,1,1),exc(1,1,2),exc(1,2,1))

14
src/determinants/slater_rules_wee_mono.irp.f
View File

@ -1,5 +1,5 @@
subroutine i_Wee_j_mono(key_i,key_j,Nint,spin,hij)
subroutine i_Wee_j_single(key_i,key_j,Nint,spin,hij)
use bitmasks
implicit none
BEGIN_DOC
@ -15,8 +15,8 @@ subroutine i_Wee_j_mono(key_i,key_j,Nint,spin,hij)
PROVIDE big_array_exchange_integrals mo_two_e_integrals_in_map
call get_mono_excitation_spin(key_i(1,spin),key_j(1,spin),exc,phase,Nint)
call mono_excitation_wee(key_i,key_j,exc(1,1),exc(1,2),spin,phase,hij)
call get_single_excitation_spin(key_i(1,spin),key_j(1,spin),exc,phase,Nint)
call single_excitation_wee(key_i,key_j,exc(1,1),exc(1,2),spin,phase,hij)
end
@ -188,7 +188,7 @@ subroutine i_H_j_mono_spin_one_e(key_i,key_j,Nint,spin,hij)
integer :: exc(0:2,2)
double precision :: phase
call get_mono_excitation_spin(key_i(1,spin),key_j(1,spin),exc,phase,Nint)
call get_single_excitation_spin(key_i(1,spin),key_j(1,spin),exc,phase,Nint)
integer :: m,p
m = exc(1,1)
p = exc(1,2)
@ -252,7 +252,7 @@ subroutine i_H_j_one_e(key_i,key_j,Nint,hij)
if(degree==0)then
hij = diag_H_mat_elem_one_e(key_i,N_int)
else
call get_mono_excitation(key_i,key_j,exc,phase,Nint)
call get_single_excitation(key_i,key_j,exc,phase,Nint)
if (exc(0,1,1) == 1) then
! Mono alpha
m = exc(1,1,1)
@ -340,7 +340,7 @@ subroutine i_H_j_two_e(key_i,key_j,Nint,hij)
exc(1,2,2) ,mo_integrals_map) )
endif
case (1)
call get_mono_excitation(key_i,key_j,exc,phase,Nint)
call get_single_excitation(key_i,key_j,exc,phase,Nint)
!DIR$ FORCEINLINE
call bitstring_to_list_ab(key_i, occ, n_occ_ab, Nint)
if (exc(0,1,1) == 1) then
@ -354,7 +354,7 @@ subroutine i_H_j_two_e(key_i,key_j,Nint,hij)
p = exc(1,2,2)
spin = 2
endif
call mono_excitation_wee(key_i,key_j,p,m,spin,phase,hij)
call single_excitation_wee(key_i,key_j,p,m,spin,phase,hij)
case (0)
double precision :: diag_wee_mat_elem
hij = diag_wee_mat_elem(key_i,Nint)

51
src/fci/40.fci.bats
View File

@ -18,8 +18,27 @@ function run() {
}
function run_stoch() {
thresh=$2
test_exe fci || skip
qp set perturbation do_pt2 True
qp set determinants n_det_max 100000
qp set determinants n_states 1
qp set davidson threshold_davidson 1.e-10
qp set davidson n_states_diag 1
qp run fci
energy1="$(ezfio get fci energy_pt2 | tr '[]' ' ' | cut -d ',' -f 1)"
eq $energy1 $1 $thresh
}
@test "F2" { # 4.07m
[[ -n $TRAVIS ]] && skip
qp set_file f2.ezfio
qp set_frozen_core
run_stoch -199.30496 1.e-4
}
@test "NH3" { # 10.6657s
qp set_file nh3.ezfio
qp set_mo_class --core="[1-4]" --act="[5-72]"
@ -34,13 +53,13 @@ function run() {
@test "HCO" { # 12.2868s
qp set_file hco.ezfio
run -113.296806579881 1.e-05
run -113.296794171915 2.e-05
}
@test "H2O2" { # 12.9214s
qp set_file h2o2.ezfio
qp set_mo_class --core="[1-2]" --act="[3-24]" --del="[25-38]"
run -151.004935161155 1.e-5
run -151.004888189874 2.e-5
}
@test "HBO" { # 13.3144s
@ -52,19 +71,19 @@ function run() {
@test "H2O" { # 11.3727s
[[ -n $TRAVIS ]] && skip
qp set_file h2o.ezfio
run -76.2359268957699 1.e-5
run -76.2359268957699 2.e-5
}
@test "ClO" { # 13.3755s
[[ -n $TRAVIS ]] && skip
qp set_file clo.ezfio
run -534.546053053143 1.e-5
run -534.546005867797 5.e-5
}
@test "SO" { # 13.4952s
[[ -n $TRAVIS ]] && skip
qp set_file so.ezfio
run -26.0126370436611 1.e-5
run -26.0144622194831 1.e-5
}
@test "H2S" { # 13.6745s
@ -94,20 +113,20 @@ function run() {
@test "SiH3" { # 15.99s
[[ -n $TRAVIS ]] && skip
qp set_file sih3.ezfio
run -5.57267383364177 1.e-05
run -5.57269434557089 2.e-05
}
@test "CH4" { # 16.1612s
[[ -n $TRAVIS ]] && skip
qp set_file ch4.ezfio
qp set_mo_class --core="[1]" --act="[2-30]" --del="[31-59]"
run -40.2409672510721 1.e-5
run -40.2409858175829 2.e-5
}
@test "ClF" { # 16.8864s
[[ -n $TRAVIS ]] && skip
qp set_file clf.ezfio
run -559.168731496312 1.e-5
run -559.170116079903 1.e-5
}
@test "SO2" { # 17.5645s
@ -121,36 +140,30 @@ function run() {
[[ -n $TRAVIS ]] && skip
qp set_file c2h2.ezfio
qp set_mo_class --act="[1-30]" --del="[31-36]"
run -12.3671467643742 1.e-5
run -12.3678973551285 2.e-5
}
@test "N2" { # 18.0198s
[[ -n $TRAVIS ]] && skip
qp set_file n2.ezfio
qp set_mo_class --core="[1,2]" --act="[3-40]" --del="[41-60]"
run -109.291382745482 1.e-5
run -109.291310557766 1.e-4
}
@test "N2H4" { # 18.5006s
[[ -n $TRAVIS ]] && skip
qp set_file n2h4.ezfio
qp set_mo_class --core="[1-2]" --act="[3-24]" --del="[25-48]"
run -111.367234092521 1.e-5
run -111.367234092521 2.e-5
}
@test "CO2" { # 21.1748s
[[ -n $TRAVIS ]] && skip
qp set_file co2.ezfio
qp set_mo_class --core="[1,2]" --act="[3-30]" --del="[31-42]"
run -187.969721107603 1.e-5
run -187.969556614801 1.e-5
}
@test "F2" { # 21.331s
[[ -n $TRAVIS ]] && skip
qp set_file f2.ezfio
qp set_mo_class --core="[1,2]" --act="[3-30]" --del="[31-62]"
run -199.068518708366 1.e-5
}
@test "[Cu(NH3)4]2+" { # 25.0417s
[[ -n $TRAVIS ]] && skip
@ -163,6 +176,6 @@ function run() {
[[ -n $TRAVIS ]] && skip
qp set_file hcn.ezfio
qp set_mo_class --core="[1,2]" --act="[3-40]" --del="[41-55]"
run -93.0779744802522 1.e-5
run -93.0794109423741 2.e-5
}

6
src/mo_one_e_ints/mo_one_e_ints.irp.f
View File

@ -2,10 +2,10 @@ BEGIN_PROVIDER [ double precision, mo_one_e_integrals,(mo_num,mo_num)]
implicit none
integer :: i,j,n,l
BEGIN_DOC
! array of the mono electronic hamiltonian on the MOs basis :
! sum of the kinetic and nuclear electronic potential (and pseudo potential if needed)
! array of the one-electron Hamiltonian on the |MO| basis :
! sum of the kinetic and nuclear electronic potentials (and pseudo potential if needed)
END_DOC
print*,'Providing the mono electronic integrals'
print*,'Providing the one-electron integrals'
IF (read_mo_one_e_integrals) THEN
call ezfio_get_mo_one_e_ints_mo_one_e_integrals(mo_one_e_integrals)

6
src/perturbation/EZFIO.cfg
View File

@ -22,4 +22,8 @@ doc: The selection process stops at a fixed correlation ratio (useful for gettin
interface: ezfio,provider,ocaml
default: 1.00
[h0_type]
type: character*(32)
doc: Type of denominator in PT2. [EN | SOP | HF]
interface: ezfio,provider,ocaml
default: EN

13
src/perturbation/h0_type.irp.f
View File

@ -1,13 +0,0 @@
BEGIN_PROVIDER [ character*32,h0_type ]
implicit none
BEGIN_DOC
! Type of zeroth-order Hamiltonian
END_DOC
if (s2_eig) then
h0_type = 'SOP'
else
h0_type = 'EN'
endif
! h0_type = 'HF'
END_PROVIDER

4
src/perturbation/perturbation.template.f
View File

@ -198,7 +198,7 @@ subroutine perturb_buffer_by_mono_$PERT(i_generator,buffer,buffer_size,e_2_pert_
double precision, intent(inout) :: coef_pert_buffer(N_st,buffer_size),e_2_pert_buffer(N_st,buffer_size),sum_H_pert_diag(N_st)
double precision :: c_pert(N_st), e_2_pert(N_st), H_pert_diag(N_st)
integer :: i,k, c_ref, ni, ex
integer, external :: connected_to_ref_by_mono
integer, external :: connected_to_ref_by_single
logical, external :: is_in_wavefunction
integer(bit_kind), allocatable :: minilist(:,:,:)
@ -232,7 +232,7 @@ subroutine perturb_buffer_by_mono_$PERT(i_generator,buffer,buffer_size,e_2_pert_
do i=1,buffer_size
c_ref = connected_to_ref_by_mono(buffer(1,1,i),psi_det_generators,Nint,i_generator,N_det)
c_ref = connected_to_ref_by_single(buffer(1,1,i),psi_det_generators,Nint,i_generator,N_det)
if (c_ref /= 0) then
cycle

4
src/zmq/utils.irp.f
View File

@ -604,11 +604,11 @@ subroutine end_parallel_job(zmq_to_qp_run_socket,zmq_socket_pull,name_in)
rc = f77_zmq_send(zmq_to_qp_run_socket, 'end_job '//trim(zmq_state),8+len(trim(zmq_state)),0)
rc = f77_zmq_recv(zmq_to_qp_run_socket, message, 512, 0)
if (trim(message(1:13)) == 'error waiting') then
call sleep(1)
cycle
else if (message(1:2) == 'ok') then
exit
endif
call sleep(1)
end do
if (i==0) then
print *, '.. Forcing kill ..'
@ -1127,6 +1127,7 @@ integer function zmq_delete_task_async_recv(zmq_to_qp_run_socket,more,sending)
integer :: rc
character*(512) :: message
character*(64) :: reply
zmq_delete_task_async_recv = 0
if (.not.sending) return
sending = .False.
reply = ''
@ -1136,6 +1137,7 @@ integer function zmq_delete_task_async_recv(zmq_to_qp_run_socket,more,sending)
else if (reply(16:19) == 'done') then
more = 0
else
print *, reply(1:rc)
zmq_delete_task_async_recv = -1
return
endif