From 5e7978304d5a258aa49905813c3205d08e2f2d88 Mon Sep 17 00:00:00 2001 From: Anthony Scemama Date: Tue, 11 Jan 2022 11:25:48 +0100 Subject: [PATCH 01/17] Documentation --- src/two_body_rdm/two_e_dm_mo.irp.f | 7 +++---- 1 file changed, 3 insertions(+), 4 deletions(-) diff --git a/src/two_body_rdm/two_e_dm_mo.irp.f b/src/two_body_rdm/two_e_dm_mo.irp.f index 4dadd2e6..a4dea15f 100644 --- a/src/two_body_rdm/two_e_dm_mo.irp.f +++ b/src/two_body_rdm/two_e_dm_mo.irp.f @@ -1,9 +1,8 @@ BEGIN_PROVIDER [double precision, two_e_dm_mo, (mo_num,mo_num,mo_num,mo_num)] implicit none BEGIN_DOC - ! two_e_dm_bb_mo(i,j,k,l,istate) = STATE SPECIFIC physicist notation for 2RDM of beta/beta electrons - ! - ! + ! \sum_{\sigma \sigma'} + ! ! ! where the indices (i,j,k,l) belong to all MOs. ! @@ -12,7 +11,7 @@ BEGIN_PROVIDER [double precision, two_e_dm_mo, (mo_num,mo_num,mo_num,mo_num)] ! !!!!! WARNING !!!!! IF "no_core_density" then all elements involving at least one CORE MO are set to zero ! The state-averaged two-electron energy : ! - ! \sum_{i,j,k,l = 1, mo_num} two_e_dm_mo(i,j,k,l) * < ii jj | kk ll > + ! \sum_{i,j,k,l = 1, mo_num} two_e_dm_mo(i,j,k,l) * < kk ll | ii jj > END_DOC two_e_dm_mo = 0.d0 integer :: i,j,k,l,iorb,jorb,korb,lorb,istate From 188ccf0d06579319346ca7f20cf2f7b969fe36a7 Mon Sep 17 00:00:00 2001 From: Emmanuel Giner Date: Thu, 20 Jan 2022 20:10:50 +0100 Subject: [PATCH 02/17] removed bug in subroutine of dav_general_mat --- src/csf/configurations.irp.f | 24 ++++++++++++++ .../dav_diag_dressed_ext_rout.irp.f | 32 +------------------ src/determinants/s2.irp.f | 6 +++- src/mo_guess/h_core_guess_routine.irp.f | 2 +- 4 files changed, 31 insertions(+), 33 deletions(-) diff --git a/src/csf/configurations.irp.f b/src/csf/configurations.irp.f index 8e2a513c..ce5d48ab 100644 --- a/src/csf/configurations.irp.f +++ b/src/csf/configurations.irp.f @@ -779,6 +779,7 @@ subroutine binary_search_cfg(cfgInp,addcfg) end subroutine BEGIN_PROVIDER [ integer, psi_configuration_to_psi_det, (2,N_configuration) ] +&BEGIN_PROVIDER [ integer, psi_configuration_n_det, (N_configuration) ] &BEGIN_PROVIDER [ integer, psi_configuration_to_psi_det_data, (N_det) ] implicit none @@ -867,6 +868,29 @@ end subroutine enddo deallocate(dets, old_order) + integer :: ndet_conf + do i = 1, N_configuration + ndet_conf = psi_configuration_to_psi_det(2,i) - psi_configuration_to_psi_det(1,i) + 1 + psi_configuration_n_det(i) = ndet_conf + enddo END_PROVIDER + +BEGIN_PROVIDER [ integer, n_elec_alpha_for_psi_configuration, (N_configuration)] + implicit none + integer :: i,j,k,l + integer(bit_kind) :: det_tmp(N_int,2),det_alpha(N_int) + n_elec_alpha_for_psi_configuration = 0 + do i = 1, N_configuration + j = psi_configuration_to_psi_det(2,i) + det_tmp(:,:) = psi_det(:,:,j) + k = 0 + do l = 1, N_int + det_alpha(N_int) = iand(det_tmp(l,1),psi_configuration(l,1,i)) + k += popcnt(det_alpha(l)) + enddo + n_elec_alpha_for_psi_configuration(i) = k + enddo + +END_PROVIDER diff --git a/src/dav_general_mat/dav_diag_dressed_ext_rout.irp.f b/src/dav_general_mat/dav_diag_dressed_ext_rout.irp.f index 2f3d7f80..243e9995 100644 --- a/src/dav_general_mat/dav_diag_dressed_ext_rout.irp.f +++ b/src/dav_general_mat/dav_diag_dressed_ext_rout.irp.f @@ -1,5 +1,5 @@ -subroutine davidson_general_ext_rout(u_in,H_jj,Dress_jj,energies,sze,N_st,N_st_diag_in,converged,hcalc) +subroutine davidson_general_ext_rout_diag_dressed(u_in,H_jj,Dress_jj,energies,sze,N_st,N_st_diag_in,converged,hcalc) use mmap_module implicit none BEGIN_DOC @@ -412,36 +412,6 @@ subroutine davidson_general_ext_rout(u_in,H_jj,Dress_jj,energies,sze,N_st,N_st_d FREE nthreads_davidson end -subroutine hcalc_template(v,u,N_st,sze) - use bitmasks - implicit none - BEGIN_DOC - ! Template of routine for the application of H - ! - ! Here, it is done with the Hamiltonian matrix - ! - ! on the set of determinants of psi_det - ! - ! Computes $v = H | u \rangle$ - ! - END_DOC - integer, intent(in) :: N_st,sze - double precision, intent(in) :: u(sze,N_st) - double precision, intent(inout) :: v(sze,N_st) - integer :: i,j,istate - v = 0.d0 - do istate = 1, N_st - do i = 1, sze - do j = 1, sze - v(i,istate) += H_matrix_all_dets(j,i) * u(j,istate) - enddo - enddo - do i = 1, sze - v(i,istate) += u(i,istate) * nuclear_repulsion - enddo - enddo -end - subroutine dressing_diag_uv(v,u,dress_diag,N_st,sze) implicit none BEGIN_DOC diff --git a/src/determinants/s2.irp.f b/src/determinants/s2.irp.f index d73b2dbf..2c1a8757 100644 --- a/src/determinants/s2.irp.f +++ b/src/determinants/s2.irp.f @@ -103,13 +103,17 @@ BEGIN_PROVIDER [ double precision, expected_s2] END_PROVIDER -BEGIN_PROVIDER [ double precision, s2_values, (N_states) ] + BEGIN_PROVIDER [ double precision, s2_values, (N_states) ] +&BEGIN_PROVIDER [ double precision, s_values, (N_states) ] implicit none BEGIN_DOC ! array of the averaged values of the S^2 operator on the various states END_DOC integer :: i call u_0_S2_u_0(s2_values,psi_coef,n_det,psi_det,N_int,N_states,psi_det_size) + do i = 1, N_states + s_values(i) = 0.5d0 *(-1.d0 + dsqrt(1.d0 + 4 * s2_values(i))) + enddo END_PROVIDER diff --git a/src/mo_guess/h_core_guess_routine.irp.f b/src/mo_guess/h_core_guess_routine.irp.f index cbf23a9a..fcbdde49 100644 --- a/src/mo_guess/h_core_guess_routine.irp.f +++ b/src/mo_guess/h_core_guess_routine.irp.f @@ -7,7 +7,7 @@ subroutine hcore_guess label = 'Guess' call mo_as_eigvectors_of_mo_matrix(mo_one_e_integrals, & size(mo_one_e_integrals,1), & - size(mo_one_e_integrals,2),label,1,.false.) + size(mo_one_e_integrals,2),label,1,.true.) call nullify_small_elements(ao_num, mo_num, mo_coef, size(mo_coef,1), 1.d-12 ) call save_mos TOUCH mo_coef mo_label From a65902aa33da727fa48ce69c81a0fd4caad7d75d Mon Sep 17 00:00:00 2001 From: Emmanuel Giner Date: Tue, 25 Jan 2022 18:20:03 +0100 Subject: [PATCH 03/17] fixed bug in Laplacians --- src/dft_utils_in_r/ao_in_r.irp.f | 34 +++++++++++++++++++++----------- src/dft_utils_in_r/mo_in_r.irp.f | 2 +- 2 files changed, 24 insertions(+), 12 deletions(-) diff --git a/src/dft_utils_in_r/ao_in_r.irp.f b/src/dft_utils_in_r/ao_in_r.irp.f index 6fa6a4c7..38478d21 100644 --- a/src/dft_utils_in_r/ao_in_r.irp.f +++ b/src/dft_utils_in_r/ao_in_r.irp.f @@ -91,7 +91,19 @@ enddo END_PROVIDER - BEGIN_PROVIDER[double precision, aos_lapl_in_r_array, (ao_num,n_points_final_grid,3)] + BEGIN_PROVIDER [double precision, aos_lapl_in_r_array_transp, (ao_num, n_points_final_grid,3)] + implicit none + integer :: i,j,m + do i = 1, n_points_final_grid + do j = 1, ao_num + do m = 1, 3 + aos_lapl_in_r_array_transp(j,i,m) = aos_lapl_in_r_array(m,j,i) + enddo + enddo + enddo + END_PROVIDER + + BEGIN_PROVIDER [double precision, aos_lapl_in_r_array, (3,ao_num,n_points_final_grid)] implicit none BEGIN_DOC ! aos_lapl_in_r_array(i,j,k) = value of the kth component of the laplacian of jth ao on the ith grid point @@ -100,20 +112,20 @@ END_DOC integer :: i,j,m double precision :: aos_array(ao_num), r(3) - double precision :: aos_grad_array(ao_num,3) - double precision :: aos_lapl_array(ao_num,3) + double precision :: aos_grad_array(3,ao_num) + double precision :: aos_lapl_array(3,ao_num) !$OMP PARALLEL DO & !$OMP DEFAULT (NONE) & !$OMP PRIVATE (i,r,aos_array,aos_grad_array,aos_lapl_array,j,m) & !$OMP SHARED(aos_lapl_in_r_array,n_points_final_grid,ao_num,final_grid_points) - do m = 1, 3 - do i = 1, n_points_final_grid - r(1) = final_grid_points(1,i) - r(2) = final_grid_points(2,i) - r(3) = final_grid_points(3,i) - call give_all_aos_and_grad_and_lapl_at_r(r,aos_array,aos_grad_array,aos_lapl_array) - do j = 1, ao_num - aos_lapl_in_r_array(j,i,m) = aos_lapl_array(j,m) + do i = 1, n_points_final_grid + r(1) = final_grid_points(1,i) + r(2) = final_grid_points(2,i) + r(3) = final_grid_points(3,i) + call give_all_aos_and_grad_and_lapl_at_r(r,aos_array,aos_grad_array,aos_lapl_array) + do j = 1, ao_num + do m = 1, 3 + aos_lapl_in_r_array(m,j,i) = aos_lapl_array(m,j) enddo enddo enddo diff --git a/src/dft_utils_in_r/mo_in_r.irp.f b/src/dft_utils_in_r/mo_in_r.irp.f index 0a8b4d52..192cb25a 100644 --- a/src/dft_utils_in_r/mo_in_r.irp.f +++ b/src/dft_utils_in_r/mo_in_r.irp.f @@ -138,7 +138,7 @@ integer :: m mos_lapl_in_r_array = 0.d0 do m=1,3 - call dgemm('N','N',mo_num,n_points_final_grid,ao_num,1.d0,mo_coef_transp,mo_num,aos_lapl_in_r_array(1,1,m),ao_num,0.d0,mos_lapl_in_r_array(1,1,m),mo_num) + call dgemm('N','N',mo_num,n_points_final_grid,ao_num,1.d0,mo_coef_transp,mo_num,aos_lapl_in_r_array_transp(1,1,m),ao_num,0.d0,mos_lapl_in_r_array(1,1,m),mo_num) enddo END_PROVIDER From 173ea9eab1c461b38ac50b6efc1ff715645b805a Mon Sep 17 00:00:00 2001 From: Emmanuel Giner Date: Tue, 15 Feb 2022 11:02:39 +0100 Subject: [PATCH 04/17] minor modifs in print_wf.irp.f --- src/tools/print_wf.irp.f | 6 +++++- 1 file changed, 5 insertions(+), 1 deletion(-) diff --git a/src/tools/print_wf.irp.f b/src/tools/print_wf.irp.f index 7e51caaf..64eb1a1f 100644 --- a/src/tools/print_wf.irp.f +++ b/src/tools/print_wf.irp.f @@ -32,8 +32,9 @@ subroutine routine double precision :: norm_mono_a,norm_mono_b double precision :: norm_mono_a_2,norm_mono_b_2 double precision :: norm_mono_a_pert_2,norm_mono_b_pert_2 - double precision :: norm_mono_a_pert,norm_mono_b_pert + double precision :: norm_mono_a_pert,norm_mono_b_pert,norm_double_1 double precision :: delta_e,coef_2_2 + norm_mono_a = 0.d0 norm_mono_b = 0.d0 norm_mono_a_2 = 0.d0 @@ -42,6 +43,7 @@ subroutine routine norm_mono_b_pert = 0.d0 norm_mono_a_pert_2 = 0.d0 norm_mono_b_pert_2 = 0.d0 + norm_double_1 = 0.d0 do i = 1, min(N_det_print_wf,N_det) print*,'' print*,'i = ',i @@ -93,6 +95,7 @@ subroutine routine print*,'h1,p1 = ',h1,p1 print*,'s2',s2 print*,'h2,p2 = ',h2,p2 + norm_double_1 += dabs(psi_coef_sorted(i,1)/psi_coef_sorted(1,1)) endif print*,' = ',hij @@ -109,6 +112,7 @@ subroutine routine print*,'' print*,'L1 norm of mono alpha = ',norm_mono_a print*,'L1 norm of mono beta = ',norm_mono_b + print*,'L1 norm of double exc = ',norm_double_1 print*, '---' print*,'L2 norm of mono alpha = ',norm_mono_a_2 print*,'L2 norm of mono beta = ',norm_mono_b_2 From 39c00dd990b359dd3b17d7048a1b03633bd87d14 Mon Sep 17 00:00:00 2001 From: Anthony Scemama Date: Sat, 5 Mar 2022 15:31:16 +0100 Subject: [PATCH 05/17] Better sorting in spindeterminants --- external/qp2-dependencies | 2 +- src/determinants/spindeterminants.irp.f | 2 +- 2 files changed, 2 insertions(+), 2 deletions(-) diff --git a/external/qp2-dependencies b/external/qp2-dependencies index bc856147..90ee61f5 160000 --- a/external/qp2-dependencies +++ b/external/qp2-dependencies @@ -1 +1 @@ -Subproject commit bc856147f6e626a6616b20344e5b8e3f30f44a92 +Subproject commit 90ee61f5041c7c94a0c605625a264860292813a0 diff --git a/src/determinants/spindeterminants.irp.f b/src/determinants/spindeterminants.irp.f index dea4a566..dd55e112 100644 --- a/src/determinants/spindeterminants.irp.f +++ b/src/determinants/spindeterminants.irp.f @@ -585,7 +585,7 @@ END_PROVIDER enddo !$OMP ENDDO !$OMP END PARALLEL - call i8radix_sort(to_sort, psi_bilinear_matrix_transp_order, N_det,-1) + call i8sort(to_sort, psi_bilinear_matrix_transp_order, N_det) call iset_order(psi_bilinear_matrix_transp_rows,psi_bilinear_matrix_transp_order,N_det) call iset_order(psi_bilinear_matrix_transp_columns,psi_bilinear_matrix_transp_order,N_det) !$OMP PARALLEL DO DEFAULT(SHARED) PRIVATE(l) From 980b48d9061c9ac74cf11e95df166f3740124d65 Mon Sep 17 00:00:00 2001 From: eginer Date: Mon, 7 Mar 2022 16:36:35 +0100 Subject: [PATCH 06/17] added the possibility to unset frozen core orbitals --- bin/qp_set_frozen_core | 5 ++++- 1 file changed, 4 insertions(+), 1 deletion(-) diff --git a/bin/qp_set_frozen_core b/bin/qp_set_frozen_core index 879c71de..bc6f6834 100755 --- a/bin/qp_set_frozen_core +++ b/bin/qp_set_frozen_core @@ -7,12 +7,13 @@ setting all MOs as Active, except the n/2 first ones which are set as Core. If pseudo-potentials are used, all the MOs are set as Active. Usage: - qp_set_frozen_core [-q|--query] [(-l|-s|--large|--small)] EZFIO_DIR + qp_set_frozen_core [-q|--query] [(-l|-s|-u|--large|--small|--unset)] EZFIO_DIR Options: -q --query Prints in the standard output the number of frozen MOs -l --large Use a small core -s --small Use a large core + -u --unset Unset frozen core Default numbers of frozen electrons: @@ -88,7 +89,9 @@ def main(arguments): elif charge <= 54: n_frozen += 9 elif charge <= 86: n_frozen += 18 elif charge <= 118: n_frozen += 27 + elif arguments["--unset"]: + n_frozen = 0 else: # default for charge in ezfio.nuclei_nucl_charge: if charge <= 4: pass From cef2ab8a91effadc642f08f80807802a9a90152f Mon Sep 17 00:00:00 2001 From: Anthony Scemama Date: Tue, 8 Mar 2022 11:24:17 +0100 Subject: [PATCH 07/17] Accelerated PT2 --- external/qp2-dependencies | 2 +- src/cipsi/run_pt2_slave.irp.f | 8 ++++---- src/davidson/u0_hs2_u0.irp.f | 16 +++++++++------- 3 files changed, 14 insertions(+), 12 deletions(-) diff --git a/external/qp2-dependencies b/external/qp2-dependencies index 90ee61f5..bc856147 160000 --- a/external/qp2-dependencies +++ b/external/qp2-dependencies @@ -1 +1 @@ -Subproject commit 90ee61f5041c7c94a0c605625a264860292813a0 +Subproject commit bc856147f6e626a6616b20344e5b8e3f30f44a92 diff --git a/src/cipsi/run_pt2_slave.irp.f b/src/cipsi/run_pt2_slave.irp.f index f1001f89..731e40ac 100644 --- a/src/cipsi/run_pt2_slave.irp.f +++ b/src/cipsi/run_pt2_slave.irp.f @@ -31,11 +31,11 @@ subroutine run_pt2_slave(thread,iproc,energy) double precision, intent(in) :: energy(N_states_diag) integer, intent(in) :: thread, iproc - if (N_det > 100000 ) then - call run_pt2_slave_large(thread,iproc,energy) - else +! if (N_det > 100000 ) then +! call run_pt2_slave_large(thread,iproc,energy) +! else call run_pt2_slave_small(thread,iproc,energy) - endif +! endif end subroutine run_pt2_slave_small(thread,iproc,energy) diff --git a/src/davidson/u0_hs2_u0.irp.f b/src/davidson/u0_hs2_u0.irp.f index 8f7bf06b..38fb56bd 100644 --- a/src/davidson/u0_hs2_u0.irp.f +++ b/src/davidson/u0_hs2_u0.irp.f @@ -203,7 +203,7 @@ subroutine H_S2_u_0_nstates_openmp_work_$N_int(v_t,s_t,u_t,N_st,sze,istart,iend, integer, allocatable :: doubles(:) integer, allocatable :: singles_a(:) integer, allocatable :: singles_b(:) - integer, allocatable :: idx(:), idx0(:) + integer, allocatable :: idx(:), buffer_lrow(:), idx0(:) integer :: maxab, n_singles_a, n_singles_b, kcol_prev integer*8 :: k8 logical :: compute_singles @@ -253,7 +253,7 @@ compute_singles=.True. !$OMP PRIVATE(krow, kcol, tmp_det, spindet, k_a, k_b, i, & !$OMP lcol, lrow, l_a, l_b, utl, kk, u_is_sparse, & !$OMP buffer, doubles, n_doubles, umax, & - !$OMP tmp_det2, hij, sij, idx, l, kcol_prev, & + !$OMP tmp_det2, hij, sij, idx, buffer_lrow, l, kcol_prev, & !$OMP singles_a, n_singles_a, singles_b, ratio, & !$OMP n_singles_b, k8, last_found,left,right,right_max) @@ -264,7 +264,7 @@ compute_singles=.True. singles_a(maxab), & singles_b(maxab), & doubles(maxab), & - idx(maxab), utl(N_st,block_size)) + idx(maxab), buffer_lrow(maxab), utl(N_st,block_size)) kcol_prev=-1 @@ -332,18 +332,20 @@ compute_singles=.True. l_a = psi_bilinear_matrix_columns_loc(lcol) ASSERT (l_a <= N_det) - !DIR$ UNROLL(8) - !DIR$ LOOP COUNT avg(50000) do j=1,psi_bilinear_matrix_columns_loc(lcol+1) - psi_bilinear_matrix_columns_loc(lcol) 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) ! hot spot + buffer_lrow(j) = lrow ASSERT (l_a <= N_det) idx(j) = l_a l_a = l_a+1 enddo + + do j=1,psi_bilinear_matrix_columns_loc(lcol+1) - psi_bilinear_matrix_columns_loc(lcol) + buffer(1:$N_int,j) = psi_det_alpha_unique(1:$N_int, buffer_lrow(j)) ! hot spot + enddo j = j-1 call get_all_spin_singles_$N_int( & @@ -789,7 +791,7 @@ compute_singles=.True. end do !$OMP END DO - deallocate(buffer, singles_a, singles_b, doubles, idx, utl) + deallocate(buffer, singles_a, singles_b, doubles, idx, buffer_lrow, utl) !$OMP END PARALLEL end From 6b118362df31a7ed647274d2864e737b171bcef2 Mon Sep 17 00:00:00 2001 From: Anthony Scemama Date: Tue, 8 Mar 2022 23:30:13 +0100 Subject: [PATCH 08/17] Accelerating PT2 again --- src/cipsi/run_pt2_slave.irp.f | 57 +++++++++++++++++++++-------------- 1 file changed, 35 insertions(+), 22 deletions(-) diff --git a/src/cipsi/run_pt2_slave.irp.f b/src/cipsi/run_pt2_slave.irp.f index 731e40ac..f9171a42 100644 --- a/src/cipsi/run_pt2_slave.irp.f +++ b/src/cipsi/run_pt2_slave.irp.f @@ -31,11 +31,11 @@ subroutine run_pt2_slave(thread,iproc,energy) double precision, intent(in) :: energy(N_states_diag) integer, intent(in) :: thread, iproc -! if (N_det > 100000 ) then -! call run_pt2_slave_large(thread,iproc,energy) -! else + if (N_det > 100000 ) then + call run_pt2_slave_large(thread,iproc,energy) + else call run_pt2_slave_small(thread,iproc,energy) -! endif + endif end subroutine run_pt2_slave_small(thread,iproc,energy) @@ -116,10 +116,10 @@ subroutine run_pt2_slave_small(thread,iproc,energy) do k=1,n_tasks call pt2_alloc(pt2_data(k),N_states) b%cur = 0 - double precision :: time2 - call wall_time(time2) +! double precision :: time2 +! call wall_time(time2) call select_connected(i_generator(k),energy,pt2_data(k),b,subset(k),pt2_F(i_generator(k))) - call wall_time(time1) +! call wall_time(time1) ! print *, i_generator(1), time1-time2, n_tasks, pt2_F(i_generator(1)) enddo call wall_time(time1) @@ -172,8 +172,8 @@ subroutine run_pt2_slave_large(thread,iproc,energy) integer :: rc, i integer :: worker_id, ctask, ltask - character*(512) :: task - integer :: task_id(1) + 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 @@ -184,15 +184,16 @@ subroutine run_pt2_slave_large(thread,iproc,energy) type(selection_buffer) :: b logical :: done, buffer_ready - type(pt2_type) :: pt2_data + type(pt2_type), allocatable :: pt2_data(:) integer :: n_tasks, k, N - integer :: i_generator, subset - + 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_data(pt2_n_tasks_max), i_generator(pt2_n_tasks_max), subset(pt2_n_tasks_max)) zmq_to_qp_run_socket = new_zmq_to_qp_run_socket() integer, external :: connect_to_taskserver @@ -211,6 +212,9 @@ subroutine run_pt2_slave_large(thread,iproc,energy) 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 @@ -221,11 +225,9 @@ subroutine run_pt2_slave_large(thread,iproc,energy) endif if (n_tasks == 0) exit - call sscanf_ddd(task, subset, i_generator, N) - if( pt2_F(i_generator) <= 0 .or. pt2_F(i_generator) > N_det ) then - print *, irp_here - stop 'bug in selection' - endif + do k=1,n_tasks + call sscanf_ddd(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) @@ -235,9 +237,14 @@ subroutine run_pt2_slave_large(thread,iproc,energy) ASSERT (b%N == bsize) endif - call pt2_alloc(pt2_data,N_states) - b%cur = 0 - call select_connected(i_generator,energy,pt2_data,b,subset,pt2_F(i_generator)) + double precision :: time0, time1 + call wall_time(time0) + do k=1,n_tasks + call pt2_alloc(pt2_data(k),N_states) + b%cur = 0 + call select_connected(i_generator(k),energy,pt2_data(k),b,subset(k),pt2_F(i_generator(k))) + enddo + call wall_time(time1) integer, external :: tasks_done_to_taskserver if (tasks_done_to_taskserver(zmq_to_qp_run_socket,worker_id,task_id,n_tasks) == -1) then @@ -249,7 +256,7 @@ subroutine run_pt2_slave_large(thread,iproc,energy) call merge_selection_buffers(b,global_selection_buffer) b%cur=0 call omp_unset_lock(global_selection_buffer_lock) - if ( iproc == 1 .or. i_generator < 100 .or. done) then + if ( iproc == 1 .or. i_generator(1) < 100 .or. done) then call omp_set_lock(global_selection_buffer_lock) call push_pt2_results_async_recv(zmq_socket_push,b%mini,sending) call push_pt2_results_async_send(zmq_socket_push, (/i_generator/), (/pt2_data/), global_selection_buffer, (/task_id/), 1,sending) @@ -260,7 +267,13 @@ subroutine run_pt2_slave_large(thread,iproc,energy) call push_pt2_results_async_send(zmq_socket_push, (/i_generator/), (/pt2_data/), b, (/task_id/), 1,sending) endif - call pt2_dealloc(pt2_data) + do k=1,n_tasks + call pt2_dealloc(pt2_data(k)) + enddo + b%cur=0 +! ! Try to adjust n_tasks at least 5 seconds per task + n_tasks = min(2*n_tasks,int( dble(5*n_tasks) / (time1 - time0 + 1.d0))) + n_tasks = min(n_tasks, pt2_n_tasks_max) end do call push_pt2_results_async_recv(zmq_socket_push,b%mini,sending) From 0c8f5e5f0b8537a3aa01e692ba21cbbd07509a42 Mon Sep 17 00:00:00 2001 From: Anthony Scemama Date: Tue, 8 Mar 2022 23:43:29 +0100 Subject: [PATCH 09/17] Accelerating PT2 again --- src/cipsi/run_pt2_slave.irp.f | 43 ++++++++++++----------------------- 1 file changed, 15 insertions(+), 28 deletions(-) diff --git a/src/cipsi/run_pt2_slave.irp.f b/src/cipsi/run_pt2_slave.irp.f index f9171a42..083865f3 100644 --- a/src/cipsi/run_pt2_slave.irp.f +++ b/src/cipsi/run_pt2_slave.irp.f @@ -172,8 +172,8 @@ subroutine run_pt2_slave_large(thread,iproc,energy) integer :: rc, i integer :: worker_id, ctask, ltask - character*(512), allocatable :: task(:) - integer, allocatable :: task_id(:) + character*(512) :: task + integer :: task_id(1) integer(ZMQ_PTR),external :: new_zmq_to_qp_run_socket integer(ZMQ_PTR) :: zmq_to_qp_run_socket @@ -184,16 +184,15 @@ subroutine run_pt2_slave_large(thread,iproc,energy) type(selection_buffer) :: b logical :: done, buffer_ready - type(pt2_type), allocatable :: pt2_data(:) + type(pt2_type) :: pt2_data integer :: n_tasks, k, N - integer, allocatable :: i_generator(:), subset(:) + integer :: 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_data(pt2_n_tasks_max), i_generator(pt2_n_tasks_max), subset(pt2_n_tasks_max)) zmq_to_qp_run_socket = new_zmq_to_qp_run_socket() integer, external :: connect_to_taskserver @@ -212,9 +211,6 @@ subroutine run_pt2_slave_large(thread,iproc,energy) 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 @@ -225,9 +221,11 @@ subroutine run_pt2_slave_large(thread,iproc,energy) endif if (n_tasks == 0) exit - do k=1,n_tasks - call sscanf_ddd(task(k), subset(k), i_generator(k), N) - enddo + call sscanf_ddd(task, subset, i_generator, N) + if( pt2_F(i_generator) <= 0 .or. pt2_F(i_generator) > N_det ) then + print *, irp_here + stop 'bug in selection' + endif if (b%N == 0) then ! Only first time bsize = min(N, (elec_alpha_num * (mo_num-elec_alpha_num))**2) @@ -237,14 +235,9 @@ subroutine run_pt2_slave_large(thread,iproc,energy) ASSERT (b%N == bsize) endif - double precision :: time0, time1 - call wall_time(time0) - do k=1,n_tasks - call pt2_alloc(pt2_data(k),N_states) - b%cur = 0 - call select_connected(i_generator(k),energy,pt2_data(k),b,subset(k),pt2_F(i_generator(k))) - enddo - call wall_time(time1) + call pt2_alloc(pt2_data,N_states) + b%cur = 0 + call select_connected(i_generator,energy,pt2_data,b,subset,pt2_F(i_generator)) integer, external :: tasks_done_to_taskserver if (tasks_done_to_taskserver(zmq_to_qp_run_socket,worker_id,task_id,n_tasks) == -1) then @@ -256,7 +249,7 @@ subroutine run_pt2_slave_large(thread,iproc,energy) call merge_selection_buffers(b,global_selection_buffer) b%cur=0 call omp_unset_lock(global_selection_buffer_lock) - if ( iproc == 1 .or. i_generator(1) < 100 .or. done) then + if ( iproc == 1 .or. i_generator < 100 .or. done) then call omp_set_lock(global_selection_buffer_lock) call push_pt2_results_async_recv(zmq_socket_push,b%mini,sending) call push_pt2_results_async_send(zmq_socket_push, (/i_generator/), (/pt2_data/), global_selection_buffer, (/task_id/), 1,sending) @@ -267,13 +260,7 @@ subroutine run_pt2_slave_large(thread,iproc,energy) call push_pt2_results_async_send(zmq_socket_push, (/i_generator/), (/pt2_data/), b, (/task_id/), 1,sending) endif - do k=1,n_tasks - call pt2_dealloc(pt2_data(k)) - enddo - b%cur=0 -! ! Try to adjust n_tasks at least 5 seconds per task - n_tasks = min(2*n_tasks,int( dble(5*n_tasks) / (time1 - time0 + 1.d0))) - n_tasks = min(n_tasks, pt2_n_tasks_max) + call pt2_dealloc(pt2_data) end do call push_pt2_results_async_recv(zmq_socket_push,b%mini,sending) From 25c50aab59aebe83ed2ad2b110eb1e943526ba16 Mon Sep 17 00:00:00 2001 From: Anthony Scemama Date: Wed, 9 Mar 2022 10:23:27 +0100 Subject: [PATCH 10/17] Improving PT2 --- src/cipsi/pt2_stoch_routines.irp.f | 4 ++-- src/cipsi/run_pt2_slave.irp.f | 36 ++++++++++++++++++++---------- src/cipsi/selection.irp.f | 10 ++------- src/ezfio_files/output.irp.f | 2 +- src/utils/memory.irp.f | 2 +- 5 files changed, 30 insertions(+), 24 deletions(-) diff --git a/src/cipsi/pt2_stoch_routines.irp.f b/src/cipsi/pt2_stoch_routines.irp.f index 14b1d060..5019c957 100644 --- a/src/cipsi/pt2_stoch_routines.irp.f +++ b/src/cipsi/pt2_stoch_routines.irp.f @@ -291,7 +291,7 @@ subroutine ZMQ_pt2(E, pt2_data, pt2_data_err, relative_error, N_in) print '(A)', '========== ======================= ===================== ===================== ===========' - print '(A)', ' Samples Energy Variance Norm^2 Seconds' + print '(A)', ' Samples Energy Variance Norm^2 Seconds' print '(A)', '========== ======================= ===================== ===================== ===========' PROVIDE global_selection_buffer @@ -537,7 +537,7 @@ subroutine pt2_collector(zmq_socket_pull, E, relative_error, pt2_data, pt2_data_ if ((time - time1 > 1.d0) .or. (n==N_det_generators)) then time1 = time - print '(I10, X, F12.6, X, G10.3, X, F10.6, X, G10.3, X, F10.6, X, G10.3, X, F10.4)', c, & + print '(I10, X, F12.6, X, G10.3, X, F10.6, X, G10.3, X, F10.6, X, G10.3, X, F10.1)', c, & pt2_data % pt2(pt2_stoch_istate) +E, & pt2_data_err % pt2(pt2_stoch_istate), & pt2_data % variance(pt2_stoch_istate), & diff --git a/src/cipsi/run_pt2_slave.irp.f b/src/cipsi/run_pt2_slave.irp.f index 083865f3..a8215a98 100644 --- a/src/cipsi/run_pt2_slave.irp.f +++ b/src/cipsi/run_pt2_slave.irp.f @@ -190,8 +190,12 @@ subroutine run_pt2_slave_large(thread,iproc,energy) integer :: bsize ! Size of selection buffers logical :: sending + double precision :: time_shift + PROVIDE global_selection_buffer global_selection_buffer_lock + call random_number(time_shift) + time_shift = time_shift*15.d0 zmq_to_qp_run_socket = new_zmq_to_qp_run_socket() @@ -209,6 +213,9 @@ subroutine run_pt2_slave_large(thread,iproc,energy) sending = .False. done = .False. + double precision :: time0, time1 + call wall_time(time0) + time0 = time0+time_shift do while (.not.done) integer, external :: get_tasks_from_taskserver @@ -244,20 +251,25 @@ subroutine run_pt2_slave_large(thread,iproc,energy) done = .true. endif call sort_selection_buffer(b) - 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 .or. i_generator < 100 .or. done) then + + call wall_time(time1) + if (time1-time0 > 15.d0) then call omp_set_lock(global_selection_buffer_lock) - call push_pt2_results_async_recv(zmq_socket_push,b%mini,sending) - call push_pt2_results_async_send(zmq_socket_push, (/i_generator/), (/pt2_data/), global_selection_buffer, (/task_id/), 1,sending) - global_selection_buffer%cur = 0 + 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) - else - call push_pt2_results_async_recv(zmq_socket_push,b%mini,sending) - call push_pt2_results_async_send(zmq_socket_push, (/i_generator/), (/pt2_data/), b, (/task_id/), 1,sending) + call wall_time(time0) + if ( iproc == 1 .or. i_generator < 100 .or. done) then + call omp_set_lock(global_selection_buffer_lock) + call push_pt2_results_async_recv(zmq_socket_push,b%mini,sending) + call push_pt2_results_async_send(zmq_socket_push, (/i_generator/), (/pt2_data/), global_selection_buffer, (/task_id/), 1,sending) + global_selection_buffer%cur = 0 + call omp_unset_lock(global_selection_buffer_lock) + else + call push_pt2_results_async_recv(zmq_socket_push,b%mini,sending) + call push_pt2_results_async_send(zmq_socket_push, (/i_generator/), (/pt2_data/), b, (/task_id/), 1,sending) + endif endif call pt2_dealloc(pt2_data) diff --git a/src/cipsi/selection.irp.f b/src/cipsi/selection.irp.f index f1ec6ff6..acb91fb5 100644 --- a/src/cipsi/selection.irp.f +++ b/src/cipsi/selection.irp.f @@ -474,17 +474,11 @@ subroutine select_singles_and_doubles(i_generator,hole_mask,particle_mask,fock_d ! endif do i=1,fullinteresting(0) - do k=1,N_int - fullminilist(k,1,i) = psi_det_sorted(k,1,fullinteresting(i)) - fullminilist(k,2,i) = psi_det_sorted(k,2,fullinteresting(i)) - enddo + fullminilist(:,:,i) = psi_det_sorted(:,:,fullinteresting(i)) enddo do i=1,interesting(0) - do k=1,N_int - minilist(k,1,i) = psi_det_sorted(k,1,interesting(i)) - minilist(k,2,i) = psi_det_sorted(k,2,interesting(i)) - enddo + minilist(:,:,i) = psi_det_sorted(:,:,interesting(i)) enddo do s2=s1,2 diff --git a/src/ezfio_files/output.irp.f b/src/ezfio_files/output.irp.f index 48512f92..7b2663a0 100644 --- a/src/ezfio_files/output.irp.f +++ b/src/ezfio_files/output.irp.f @@ -25,7 +25,7 @@ subroutine write_time(iunit) ct = ct - output_cpu_time_0 call wall_time(wt) wt = wt - output_wall_time_0 - write(6,'(A,F14.6,A,F14.6,A)') & + write(6,'(A,F14.2,A,F14.2,A)') & '.. >>>>> [ WALL TIME: ', wt, ' s ] [ CPU TIME: ', ct, ' s ] <<<<< ..' write(6,*) end diff --git a/src/utils/memory.irp.f b/src/utils/memory.irp.f index 3ea242b0..d5a066a1 100644 --- a/src/utils/memory.irp.f +++ b/src/utils/memory.irp.f @@ -114,7 +114,7 @@ subroutine print_memory_usage() call resident_memory(rss) call total_memory(mem) - write(*,'(A,F14.6,A,F14.6,A)') & + write(*,'(A,F14.3,A,F14.3,A)') & '.. >>>>> [ RES MEM : ', rss , & ' GB ] [ VIRT MEM : ', mem, ' GB ] <<<<< ..' end From 416ff24ff44cbb7dd5fabea1c3415ea720432146 Mon Sep 17 00:00:00 2001 From: Anthony Scemama Date: Wed, 9 Mar 2022 10:40:30 +0100 Subject: [PATCH 11/17] Fixed previous commit --- src/cipsi/run_pt2_slave.irp.f | 6 +++--- 1 file changed, 3 insertions(+), 3 deletions(-) diff --git a/src/cipsi/run_pt2_slave.irp.f b/src/cipsi/run_pt2_slave.irp.f index a8215a98..9e046877 100644 --- a/src/cipsi/run_pt2_slave.irp.f +++ b/src/cipsi/run_pt2_slave.irp.f @@ -266,10 +266,10 @@ subroutine run_pt2_slave_large(thread,iproc,energy) call push_pt2_results_async_send(zmq_socket_push, (/i_generator/), (/pt2_data/), global_selection_buffer, (/task_id/), 1,sending) global_selection_buffer%cur = 0 call omp_unset_lock(global_selection_buffer_lock) - else - call push_pt2_results_async_recv(zmq_socket_push,b%mini,sending) - call push_pt2_results_async_send(zmq_socket_push, (/i_generator/), (/pt2_data/), b, (/task_id/), 1,sending) endif + else + call push_pt2_results_async_recv(zmq_socket_push,b%mini,sending) + call push_pt2_results_async_send(zmq_socket_push, (/i_generator/), (/pt2_data/), b, (/task_id/), 1,sending) endif call pt2_dealloc(pt2_data) From 27dc0c06aea0ae63ed55621106df40353201e3a2 Mon Sep 17 00:00:00 2001 From: Anthony Scemama Date: Thu, 10 Mar 2022 00:55:23 +0100 Subject: [PATCH 12/17] Optimize PT2 --- src/cipsi/run_pt2_slave.irp.f | 18 ++++---- src/cipsi/selection_buffer.irp.f | 71 +++++++++++++++++++------------- 2 files changed, 51 insertions(+), 38 deletions(-) diff --git a/src/cipsi/run_pt2_slave.irp.f b/src/cipsi/run_pt2_slave.irp.f index 9e046877..30fc7ce0 100644 --- a/src/cipsi/run_pt2_slave.irp.f +++ b/src/cipsi/run_pt2_slave.irp.f @@ -253,22 +253,22 @@ subroutine run_pt2_slave_large(thread,iproc,energy) call sort_selection_buffer(b) call wall_time(time1) - if (time1-time0 > 15.d0) then +! if (time1-time0 > 15.d0) then 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) call wall_time(time0) - if ( iproc == 1 .or. i_generator < 100 .or. done) then - call omp_set_lock(global_selection_buffer_lock) - call push_pt2_results_async_recv(zmq_socket_push,b%mini,sending) - call push_pt2_results_async_send(zmq_socket_push, (/i_generator/), (/pt2_data/), global_selection_buffer, (/task_id/), 1,sending) - global_selection_buffer%cur = 0 - call omp_unset_lock(global_selection_buffer_lock) - endif +! endif + + call push_pt2_results_async_recv(zmq_socket_push,b%mini,sending) + if ( iproc == 1 .or. i_generator < 100 .or. done) then + call omp_set_lock(global_selection_buffer_lock) + call push_pt2_results_async_send(zmq_socket_push, (/i_generator/), (/pt2_data/), global_selection_buffer, (/task_id/), 1,sending) + global_selection_buffer%cur = 0 + call omp_unset_lock(global_selection_buffer_lock) else - call push_pt2_results_async_recv(zmq_socket_push,b%mini,sending) call push_pt2_results_async_send(zmq_socket_push, (/i_generator/), (/pt2_data/), b, (/task_id/), 1,sending) endif diff --git a/src/cipsi/selection_buffer.irp.f b/src/cipsi/selection_buffer.irp.f index 79899139..1f743e0e 100644 --- a/src/cipsi/selection_buffer.irp.f +++ b/src/cipsi/selection_buffer.irp.f @@ -92,38 +92,51 @@ subroutine merge_selection_buffers(b1, b2) allocate(val(sze), detmp(N_int, 2, sze)) i1=1 i2=1 - do i=1,nmwen - if ( (i1 > b1%cur).and.(i2 > b2%cur) ) then - exit - else if (i1 > b1%cur) then - val(i) = b2%val(i2) - detmp(1:N_int,1,i) = b2%det(1:N_int,1,i2) - detmp(1:N_int,2,i) = b2%det(1:N_int,2,i2) - i2=i2+1 - else if (i2 > b2%cur) then - val(i) = b1%val(i1) - detmp(1:N_int,1,i) = b1%det(1:N_int,1,i1) - detmp(1:N_int,2,i) = b1%det(1:N_int,2,i1) - i1=i1+1 - else - if (b1%val(i1) <= b2%val(i2)) then - val(i) = b1%val(i1) - detmp(1:N_int,1,i) = b1%det(1:N_int,1,i1) - detmp(1:N_int,2,i) = b1%det(1:N_int,2,i1) - i1=i1+1 + + select case (N_int) +BEGIN_TEMPLATE + case $case + do i=1,nmwen + if ( (i1 > b1%cur).and.(i2 > b2%cur) ) then + exit + else if (i1 > b1%cur) then + val(i) = b2%val(i2) + detmp(1:$N_int,1,i) = b2%det(1:$N_int,1,i2) + detmp(1:$N_int,2,i) = b2%det(1:$N_int,2,i2) + i2=i2+1 + else if (i2 > b2%cur) then + val(i) = b1%val(i1) + detmp(1:$N_int,1,i) = b1%det(1:$N_int,1,i1) + detmp(1:$N_int,2,i) = b1%det(1:$N_int,2,i1) + i1=i1+1 else - val(i) = b2%val(i2) - detmp(1:N_int,1,i) = b2%det(1:N_int,1,i2) - detmp(1:N_int,2,i) = b2%det(1:N_int,2,i2) - i2=i2+1 + if (b1%val(i1) <= b2%val(i2)) then + val(i) = b1%val(i1) + detmp(1:$N_int,1,i) = b1%det(1:$N_int,1,i1) + detmp(1:$N_int,2,i) = b1%det(1:$N_int,2,i1) + i1=i1+1 + else + val(i) = b2%val(i2) + detmp(1:$N_int,1,i) = b2%det(1:$N_int,1,i2) + detmp(1:$N_int,2,i) = b2%det(1:$N_int,2,i2) + i2=i2+1 + endif endif - endif - enddo + enddo + do i=nmwen+1,b2%N + val(i) = 0.d0 +! detmp(1:$N_int,1,i) = 0_bit_kind +! detmp(1:$N_int,2,i) = 0_bit_kind + enddo +SUBST [ case, N_int ] +(1); 1;; +(2); 2;; +(3); 3;; +(4); 4;; +default; N_int;; +END_TEMPLATE + end select deallocate(b2%det, b2%val) - do i=nmwen+1,b2%N - val(i) = 0.d0 - detmp(1:N_int,1:2,i) = 0_bit_kind - enddo b2%det => detmp b2%val => val b2%mini = min(b2%mini,b2%val(b2%N)) From 6b0162e22991717d21cab150b952b9fdb006d824 Mon Sep 17 00:00:00 2001 From: ydamour Date: Tue, 15 Mar 2022 16:25:21 +0100 Subject: [PATCH 13/17] dipole moments x,y,z --- src/determinants/dipole_moments.irp.f | 4 +++- src/tools/print_dipole.irp.f | 2 +- 2 files changed, 4 insertions(+), 2 deletions(-) diff --git a/src/determinants/dipole_moments.irp.f b/src/determinants/dipole_moments.irp.f index 8a5f1a2d..c3675e7d 100644 --- a/src/determinants/dipole_moments.irp.f +++ b/src/determinants/dipole_moments.irp.f @@ -55,11 +55,13 @@ END_PROVIDER - subroutine print_z_dipole_moment_only + subroutine print_dipole_moments implicit none print*, '' print*, '' print*, '****************************************' + print*, 'x_dipole_moment = ',x_dipole_moment + print*, 'y_dipole_moment = ',y_dipole_moment print*, 'z_dipole_moment = ',z_dipole_moment print*, '****************************************' end diff --git a/src/tools/print_dipole.irp.f b/src/tools/print_dipole.irp.f index 8351308e..1d095af9 100644 --- a/src/tools/print_dipole.irp.f +++ b/src/tools/print_dipole.irp.f @@ -1,5 +1,5 @@ program print_dipole implicit none - call print_z_dipole_moment_only + call print_dipole_moments end From cf96b74b52443e63d20c8d414315299ce26591ad Mon Sep 17 00:00:00 2001 From: Anthony Scemama Date: Wed, 16 Mar 2022 10:26:18 +0100 Subject: [PATCH 14/17] Remove debug --- external/qp2-dependencies | 2 +- src/ao_basis/spherical_to_cartesian.irp.f | 2 +- src/determinants/dipole_moments.irp.f | 22 +++++++++++----------- 3 files changed, 13 insertions(+), 13 deletions(-) diff --git a/external/qp2-dependencies b/external/qp2-dependencies index bc856147..90ee61f5 160000 --- a/external/qp2-dependencies +++ b/external/qp2-dependencies @@ -1 +1 @@ -Subproject commit bc856147f6e626a6616b20344e5b8e3f30f44a92 +Subproject commit 90ee61f5041c7c94a0c605625a264860292813a0 diff --git a/src/ao_basis/spherical_to_cartesian.irp.f b/src/ao_basis/spherical_to_cartesian.irp.f index 33a3bc89..336161f8 100644 --- a/src/ao_basis/spherical_to_cartesian.irp.f +++ b/src/ao_basis/spherical_to_cartesian.irp.f @@ -1,7 +1,7 @@ ! Spherical to cartesian transformation matrix obtained with ! Horton (http://theochem.github.com/horton/, 2015) -! First index is the index of the carteisan AO, obtained by ao_power_index +! First index is the index of the cartesian AO, obtained by ao_power_index ! Second index is the index of the spherical AO BEGIN_PROVIDER [ double precision, cart_to_sphe_0, (1,1) ] diff --git a/src/determinants/dipole_moments.irp.f b/src/determinants/dipole_moments.irp.f index 8a5f1a2d..733dd535 100644 --- a/src/determinants/dipole_moments.irp.f +++ b/src/determinants/dipole_moments.irp.f @@ -9,7 +9,7 @@ double precision :: weight, r(3) double precision :: cpu0,cpu1,nuclei_part_z,nuclei_part_y,nuclei_part_x - call cpu_time(cpu0) +! call cpu_time(cpu0) z_dipole_moment = 0.d0 y_dipole_moment = 0.d0 x_dipole_moment = 0.d0 @@ -26,10 +26,10 @@ enddo enddo - print*,'electron part for z_dipole = ',z_dipole_moment - print*,'electron part for y_dipole = ',y_dipole_moment - print*,'electron part for x_dipole = ',x_dipole_moment - +! print*,'electron part for z_dipole = ',z_dipole_moment +! print*,'electron part for y_dipole = ',y_dipole_moment +! print*,'electron part for x_dipole = ',x_dipole_moment +! nuclei_part_z = 0.d0 nuclei_part_y = 0.d0 nuclei_part_x = 0.d0 @@ -38,18 +38,18 @@ nuclei_part_y += nucl_charge(i) * nucl_coord(i,2) nuclei_part_x += nucl_charge(i) * nucl_coord(i,1) enddo - print*,'nuclei part for z_dipole = ',nuclei_part_z - print*,'nuclei part for y_dipole = ',nuclei_part_y - print*,'nuclei part for x_dipole = ',nuclei_part_x - +! print*,'nuclei part for z_dipole = ',nuclei_part_z +! print*,'nuclei part for y_dipole = ',nuclei_part_y +! print*,'nuclei part for x_dipole = ',nuclei_part_x +! do istate = 1, N_states z_dipole_moment(istate) += nuclei_part_z y_dipole_moment(istate) += nuclei_part_y x_dipole_moment(istate) += nuclei_part_x enddo - call cpu_time(cpu1) - print*,'Time to provide the dipole moment :',cpu1-cpu0 +! call cpu_time(cpu1) +! print*,'Time to provide the dipole moment :',cpu1-cpu0 END_PROVIDER From 7da10a66cfb6c6bba3913f669e3bd79b79fdcf3e Mon Sep 17 00:00:00 2001 From: ydamour Date: Wed, 16 Mar 2022 11:42:26 +0100 Subject: [PATCH 15/17] output print dipole for n_states >= 1 + read_wf true --- src/determinants/dipole_moments.irp.f | 15 ++++++++++++--- src/tools/print_dipole.irp.f | 2 ++ 2 files changed, 14 insertions(+), 3 deletions(-) diff --git a/src/determinants/dipole_moments.irp.f b/src/determinants/dipole_moments.irp.f index 932c1c51..06fca0cd 100644 --- a/src/determinants/dipole_moments.irp.f +++ b/src/determinants/dipole_moments.irp.f @@ -57,11 +57,20 @@ END_PROVIDER subroutine print_dipole_moments implicit none + integer :: i print*, '' print*, '' print*, '****************************************' - print*, 'x_dipole_moment = ',x_dipole_moment - print*, 'y_dipole_moment = ',y_dipole_moment - print*, 'z_dipole_moment = ',z_dipole_moment + write(*,'(A10)',advance='no') ' State : ' + do i = 1,N_states + write(*,'(i16)',advance='no') i + end do + write(*,*) '' + write(*,'(A17,100(1pE16.8))') 'x_dipole_moment = ',x_dipole_moment + write(*,'(A17,100(1pE16.8))') 'y_dipole_moment = ',y_dipole_moment + write(*,'(A17,100(1pE16.8))') 'z_dipole_moment = ',z_dipole_moment + !print*, 'x_dipole_moment = ',x_dipole_moment + !print*, 'y_dipole_moment = ',y_dipole_moment + !print*, 'z_dipole_moment = ',z_dipole_moment print*, '****************************************' end diff --git a/src/tools/print_dipole.irp.f b/src/tools/print_dipole.irp.f index 1d095af9..8db9aa09 100644 --- a/src/tools/print_dipole.irp.f +++ b/src/tools/print_dipole.irp.f @@ -1,5 +1,7 @@ program print_dipole implicit none + read_wf = .True. + TOUCH read_wf call print_dipole_moments end From e63fa3fdd5f9c1770fccad32d93bdd045b7ccf6f Mon Sep 17 00:00:00 2001 From: Anthony Scemama Date: Fri, 25 Mar 2022 09:30:43 +0100 Subject: [PATCH 16/17] Fixed tooth_width=0.0 --- src/cipsi/pt2_stoch_routines.irp.f | 3 +-- src/cipsi/selection.irp.f | 2 +- src/dressing/alpha_factory.irp.f | 2 +- 3 files changed, 3 insertions(+), 4 deletions(-) diff --git a/src/cipsi/pt2_stoch_routines.irp.f b/src/cipsi/pt2_stoch_routines.irp.f index 5019c957..db0b1527 100644 --- a/src/cipsi/pt2_stoch_routines.irp.f +++ b/src/cipsi/pt2_stoch_routines.irp.f @@ -842,9 +842,8 @@ END_PROVIDER do t=1, pt2_N_teeth tooth_width = tilde_cW(pt2_n_0(t+1)) - tilde_cW(pt2_n_0(t)) if (tooth_width == 0.d0) then - tooth_width = sum(tilde_w(pt2_n_0(t):pt2_n_0(t+1))) + tooth_width = max(1.d-15,sum(tilde_w(pt2_n_0(t):pt2_n_0(t+1)))) endif - ASSERT(tooth_width > 0.d0) do i=pt2_n_0(t)+1, pt2_n_0(t+1) pt2_w(i) = tilde_w(i) * pt2_W_T / tooth_width end do diff --git a/src/cipsi/selection.irp.f b/src/cipsi/selection.irp.f index acb91fb5..d4d44c2d 100644 --- a/src/cipsi/selection.irp.f +++ b/src/cipsi/selection.irp.f @@ -1550,7 +1550,7 @@ subroutine bitstring_to_list_in_selection( string, list, n_elements, Nint) use bitmasks implicit none BEGIN_DOC - ! Gives the inidices(+1) of the bits set to 1 in the bit string + ! Gives the indices(+1) of the bits set to 1 in the bit string END_DOC integer, intent(in) :: Nint integer(bit_kind), intent(in) :: string(Nint) diff --git a/src/dressing/alpha_factory.irp.f b/src/dressing/alpha_factory.irp.f index 5eeeb1a6..c7adffe3 100644 --- a/src/dressing/alpha_factory.irp.f +++ b/src/dressing/alpha_factory.irp.f @@ -1179,7 +1179,7 @@ subroutine bitstring_to_list_in_selection( string, list, n_elements, Nint) use bitmasks implicit none BEGIN_DOC - ! Gives the inidices(+1) of the bits set to 1 in the bit string + ! Gives the indices(+1) of the bits set to 1 in the bit string END_DOC integer, intent(in) :: Nint integer(bit_kind), intent(in) :: string(Nint) From 1a890a78df6d90d73a633cccc1c0b9cc35e373be Mon Sep 17 00:00:00 2001 From: Anthony Scemama Date: Fri, 25 Mar 2022 09:32:56 +0100 Subject: [PATCH 17/17] dsqrt --- src/cipsi/pt2_stoch_routines.irp.f | 6 +++--- 1 file changed, 3 insertions(+), 3 deletions(-) diff --git a/src/cipsi/pt2_stoch_routines.irp.f b/src/cipsi/pt2_stoch_routines.irp.f index db0b1527..c7cee1ac 100644 --- a/src/cipsi/pt2_stoch_routines.irp.f +++ b/src/cipsi/pt2_stoch_routines.irp.f @@ -523,15 +523,15 @@ subroutine pt2_collector(zmq_socket_pull, E, relative_error, pt2_data, pt2_data_ ! 1/(N-1.5) : see Brugger, The American Statistician (23) 4 p. 32 (1969) if(c > 2) then eqt = dabs((pt2_data_S2(t) % pt2(pt2_stoch_istate) / c) - (pt2_data_S(t) % pt2(pt2_stoch_istate)/c)**2) ! dabs for numerical stability - eqt = sqrt(eqt / (dble(c) - 1.5d0)) + eqt = dsqrt(eqt / (dble(c) - 1.5d0)) pt2_data_err % pt2(pt2_stoch_istate) = eqt eqt = dabs((pt2_data_S2(t) % variance(pt2_stoch_istate) / c) - (pt2_data_S(t) % variance(pt2_stoch_istate)/c)**2) ! dabs for numerical stability - eqt = sqrt(eqt / (dble(c) - 1.5d0)) + eqt = dsqrt(eqt / (dble(c) - 1.5d0)) pt2_data_err % variance(pt2_stoch_istate) = eqt eqta(:) = dabs((pt2_data_S2(t) % overlap(:,pt2_stoch_istate) / c) - (pt2_data_S(t) % overlap(:,pt2_stoch_istate)/c)**2) ! dabs for numerical stability - eqta(:) = sqrt(eqta(:) / (dble(c) - 1.5d0)) + eqta(:) = dsqrt(eqta(:) / (dble(c) - 1.5d0)) pt2_data_err % overlap(:,pt2_stoch_istate) = eqta(:)