From 31028f8979189d1cc2822d4dc60eb527ed639932 Mon Sep 17 00:00:00 2001 From: eginer Date: Tue, 16 Jul 2024 17:44:48 +0200 Subject: [PATCH] fixed some weird dependencies in TC, introduced an AO cholesky 2e function --- .../cipsi_tc_bi_ortho/stochastic_cipsi.irp.f | 32 +++++++++++++++++++ plugins/local/fci_tc_bi/fci_tc_bi_ortho.irp.f | 20 ++++++++++-- .../local/non_h_ints_mu/total_tc_int.irp.f | 18 +++++++---- plugins/local/slater_tc/slater_tc_opt.irp.f | 2 -- plugins/local/slater_tc/tc_hmat.irp.f | 4 ++- src/ao_two_e_ints/cholesky.irp.f | 15 ++++++++- 6 files changed, 79 insertions(+), 12 deletions(-) diff --git a/plugins/local/cipsi_tc_bi_ortho/stochastic_cipsi.irp.f b/plugins/local/cipsi_tc_bi_ortho/stochastic_cipsi.irp.f index e363830d..6b8f3b42 100644 --- a/plugins/local/cipsi_tc_bi_ortho/stochastic_cipsi.irp.f +++ b/plugins/local/cipsi_tc_bi_ortho/stochastic_cipsi.irp.f @@ -1,4 +1,36 @@ + +! --- + +subroutine run_pouet + + BEGIN_DOC + ! Selected Full Configuration Interaction with Stochastic selection and PT2. + END_DOC + + use selection_types + implicit none + integer :: i, j, k, ndet + integer :: to_select + logical :: has + type(pt2_type) :: pt2_data, pt2_data_err + double precision :: rss + double precision :: correlation_energy_ratio + double precision :: hf_energy_ref + double precision :: relative_error + double precision, allocatable :: zeros(:),E_tc(:), norm(:) + + logical, external :: qp_stop + double precision, external :: memory_of_double + + PROVIDE mo_l_coef mo_r_coef + PROVIDE H_apply_buffer_allocated distributed_davidson + + print*, ' Diagonal elements of the Fock matrix ' + do i = 1, mo_num + write(*,*) i, Fock_matrix_tc_mo_tot(i,i) + enddo +end ! --- subroutine run_stochastic_cipsi diff --git a/plugins/local/fci_tc_bi/fci_tc_bi_ortho.irp.f b/plugins/local/fci_tc_bi/fci_tc_bi_ortho.irp.f index 1c1c0411..f1de0fe3 100644 --- a/plugins/local/fci_tc_bi/fci_tc_bi_ortho.irp.f +++ b/plugins/local/fci_tc_bi/fci_tc_bi_ortho.irp.f @@ -65,7 +65,15 @@ subroutine run_cipsi_tc() if (.not. is_zmq_slave) then - PROVIDE psi_det psi_coef mo_bi_ortho_tc_two_e mo_bi_ortho_tc_one_e + if(.True.)then! DO NOT REMOVE THE IF(.TRUE.) !! + ! this has to be provided before mo_bi_ortho_tc_two_e to avoid twice the computation of ao_two_e_tc_tot + PROVIDE Fock_matrix_tc_mo_tot + ! because Fock_matrix_tc_mo_tot depends on ao_two_e_tc_tot + ! and that mo_bi_ortho_tc_two_e erase ao_two_e_tc_tot after being provided + endif + if(.True.)then ! DO NOT REMOVE THE IF(.TRUE.) !! + PROVIDE psi_det psi_coef mo_bi_ortho_tc_two_e mo_bi_ortho_tc_one_e + endif if((elec_alpha_num+elec_beta_num) .ge. 3) then if(three_body_h_tc) then @@ -90,8 +98,16 @@ subroutine run_cipsi_tc() call json_close else + if(.True.)then! DO NOT REMOVE THE IF(.TRUE.) !! + ! this has to be provided before mo_bi_ortho_tc_two_e to avoid twice the computation of ao_two_e_tc_tot + PROVIDE Fock_matrix_tc_mo_tot + ! because Fock_matrix_tc_mo_tot depends on ao_two_e_tc_tot + ! and that mo_bi_ortho_tc_two_e erase ao_two_e_tc_tot after being provided + endif - PROVIDE mo_bi_ortho_tc_one_e mo_bi_ortho_tc_two_e pt2_min_parallel_tasks + if(.True.)then! DO NOT REMOVE THE IF(.TRUE.) !! + PROVIDE mo_bi_ortho_tc_one_e mo_bi_ortho_tc_two_e pt2_min_parallel_tasks + endif if((elec_alpha_num+elec_beta_num) .ge. 3) then if(three_body_h_tc) then diff --git a/plugins/local/non_h_ints_mu/total_tc_int.irp.f b/plugins/local/non_h_ints_mu/total_tc_int.irp.f index 656f5f16..fb09168e 100644 --- a/plugins/local/non_h_ints_mu/total_tc_int.irp.f +++ b/plugins/local/non_h_ints_mu/total_tc_int.irp.f @@ -288,25 +288,31 @@ BEGIN_PROVIDER [double precision, ao_two_e_tc_tot, (ao_num, ao_num, ao_num, ao_n !$OMP END DO !$OMP END PARALLEL else - print*, ' ao_integrals_map will be used' - PROVIDE ao_integrals_map +! print*, ' ao_integrals_map will be used' +! PROVIDE ao_integrals_map + print*,'Cholesky vectors will be used ' + double precision :: get_ao_integ_chol,eri + eri = get_ao_integ_chol(1,1,1,1) ! FOR OPENMP !$OMP PARALLEL DEFAULT(NONE) & - !$OMP SHARED(ao_num, ao_two_e_tc_tot, ao_integrals_map) & - !$OMP PRIVATE(i, j, k, l) +!!! !$OMP SHARED(ao_num, ao_two_e_tc_tot, ao_integrals_map) & + !$OMP SHARED(ao_num, ao_two_e_tc_tot) & + !$OMP PRIVATE(i, j, k, l,eri) !$OMP DO COLLAPSE(3) do j = 1, ao_num do l = 1, ao_num do i = 1, ao_num do k = 1, ao_num ! < 1:i, 2:j | 1:k, 2:l > - ao_two_e_tc_tot(k,i,l,j) = ao_two_e_tc_tot(k,i,l,j) + get_ao_two_e_integral(i, j, k, l, ao_integrals_map) +! eri = get_ao_two_e_integral(i, j, k, l, ao_integrals_map) + eri = get_ao_integ_chol(i,k,j,l) + ao_two_e_tc_tot(k,i,l,j) = ao_two_e_tc_tot(k,i,l,j) + eri enddo enddo enddo enddo !$OMP END DO !$OMP END PARALLEL - FREE ao_integrals_map +! FREE ao_integrals_map endif if((tc_integ_type .eq. "numeric") .and. (.not. tc_save_mem)) then diff --git a/plugins/local/slater_tc/slater_tc_opt.irp.f b/plugins/local/slater_tc/slater_tc_opt.irp.f index 5651a299..3c4421f8 100644 --- a/plugins/local/slater_tc/slater_tc_opt.irp.f +++ b/plugins/local/slater_tc/slater_tc_opt.irp.f @@ -10,8 +10,6 @@ subroutine provide_all_three_ints_bi_ortho() implicit none double precision :: t1, t2 - PROVIDE ao_two_e_integrals_in_map - print *, ' start provide_all_three_ints_bi_ortho' call wall_time(t1) diff --git a/plugins/local/slater_tc/tc_hmat.irp.f b/plugins/local/slater_tc/tc_hmat.irp.f index cc780364..6323d129 100644 --- a/plugins/local/slater_tc/tc_hmat.irp.f +++ b/plugins/local/slater_tc/tc_hmat.irp.f @@ -30,7 +30,9 @@ BEGIN_PROVIDER [double precision, htilde_matrix_elmt_bi_ortho, (N_det,N_det)] print *, ' PROVIDING htilde_matrix_elmt_bi_ortho ...' call wall_time(t1) - call provide_all_three_ints_bi_ortho() + if(three_body_h_tc)then + call provide_all_three_ints_bi_ortho() + endif i = 1 j = 1 diff --git a/src/ao_two_e_ints/cholesky.irp.f b/src/ao_two_e_ints/cholesky.irp.f index a2d9d043..bfa6bd0a 100644 --- a/src/ao_two_e_ints/cholesky.irp.f +++ b/src/ao_two_e_ints/cholesky.irp.f @@ -1,3 +1,15 @@ +double precision function get_ao_integ_chol(i,j,k,l) + implicit none + BEGIN_DOC + ! CHOLESKY representation of the integral of the AO basis or (ij|kl) + ! i(r1) j(r1) 1/r12 k(r2) l(r2) + END_DOC + integer, intent(in) :: i,j,k,l + double precision, external :: ddot + get_ao_integ_chol = ddot(cholesky_ao_num, cholesky_ao_transp(1,i,j), 1, cholesky_ao_transp(1,k,l), 1) + +end + BEGIN_PROVIDER [ double precision, cholesky_ao_transp, (cholesky_ao_num, ao_num, ao_num) ] implicit none BEGIN_DOC @@ -162,7 +174,8 @@ END_PROVIDER np = int(np8,4) if (np <= 0) stop 'np<=0' - rank_max = min(np,20*elec_num*elec_num) +! rank_max = min(np,20*elec_num*elec_num) + rank_max = np call mmap(trim(ezfio_work_dir)//'cholesky_ao_tmp', (/ ndim8, rank_max /), 8, fd(1), .False., .True., c_pointer(1)) call c_f_pointer(c_pointer(1), L, (/ ndim8, rank_max /))