From 4342169da802673f16b746fa9481b2fa67b0a4f9 Mon Sep 17 00:00:00 2001 From: Emmanuel Giner Date: Thu, 22 Jul 2021 23:29:13 +0200 Subject: [PATCH 1/4] added no_ov_natorb --- etc/openmp.rc | 1 + src/determinants/density_matrix.irp.f | 43 +++++++++++++++++++++++---- src/tools/save_natorb_no_ov_rot.irp.f | 25 ++++++++++++++++ 3 files changed, 63 insertions(+), 6 deletions(-) create mode 100644 etc/openmp.rc create mode 100644 src/tools/save_natorb_no_ov_rot.irp.f diff --git a/etc/openmp.rc b/etc/openmp.rc new file mode 100644 index 00000000..7f71c9b8 --- /dev/null +++ b/etc/openmp.rc @@ -0,0 +1 @@ +export OMP_NESTED=True diff --git a/src/determinants/density_matrix.irp.f b/src/determinants/density_matrix.irp.f index 7c4a7fec..c72dbc9e 100644 --- a/src/determinants/density_matrix.irp.f +++ b/src/determinants/density_matrix.irp.f @@ -262,17 +262,48 @@ subroutine set_natural_mos iorb = list_virt(i) do j = 1, n_core_inact_act_orb jorb = list_core_inact_act(j) - if(one_e_dm_mo(iorb,jorb).ne. 0.d0)then - print*,'AHAHAH' - print*,iorb,jorb,one_e_dm_mo(iorb,jorb) - stop - endif enddo enddo call mo_as_svd_vectors_of_mo_matrix_eig(one_e_dm_mo,size(one_e_dm_mo,1),mo_num,mo_num,mo_occ,label) soft_touch mo_occ - end + +subroutine set_natorb_no_ov_rot + implicit none + BEGIN_DOC + ! Set natural orbitals, obtained by diagonalization of the one-body density matrix + ! in the |MO| basis + END_DOC + character*(64) :: label + double precision, allocatable :: tmp(:,:) + allocate(tmp(mo_num, mo_num)) + label = "Natural" + tmp = one_e_dm_mo + integer :: i,j,iorb,jorb + do i = 1, n_virt_orb + iorb = list_virt(i) + do j = 1, n_core_inact_act_orb + jorb = list_core_inact_act(j) + tmp(iorb, jorb) = 0.d0 + tmp(jorb, iorb) = 0.d0 + enddo + enddo + call mo_as_svd_vectors_of_mo_matrix_eig(tmp,size(tmp,1),mo_num,mo_num,mo_occ,label) + soft_touch mo_occ +end + +subroutine save_natural_mos_no_ov_rot + implicit none + BEGIN_DOC + ! Save natural orbitals, obtained by diagonalization of the one-body density matrix in + ! the |MO| basis + END_DOC + call set_natorb_no_ov_rot + call nullify_small_elements(ao_num,mo_num,mo_coef,size(mo_coef,1),1.d-10) + call orthonormalize_mos + call save_mos +end + subroutine save_natural_mos implicit none BEGIN_DOC diff --git a/src/tools/save_natorb_no_ov_rot.irp.f b/src/tools/save_natorb_no_ov_rot.irp.f new file mode 100644 index 00000000..e5b69fbf --- /dev/null +++ b/src/tools/save_natorb_no_ov_rot.irp.f @@ -0,0 +1,25 @@ +program save_natorb + implicit none + BEGIN_DOC +! Save natural |MOs| into the |EZFIO|. +! +! This program reads the wave function stored in the |EZFIO| directory, +! extracts the corresponding natural orbitals and setd them as the new +! |MOs|. +! +! If this is a multi-state calculation, the density matrix that produces +! the natural orbitals is obtained from an average of the density +! matrices of each state with the corresponding +! :option:`determinants state_average_weight` + END_DOC + read_wf = .True. + touch read_wf + call save_natural_mos_no_ov_rot + call save_ref_determinant + call ezfio_set_mo_two_e_ints_io_mo_two_e_integrals('None') + call ezfio_set_mo_one_e_ints_io_mo_one_e_integrals('None') + call ezfio_set_mo_one_e_ints_io_mo_integrals_kinetic('None') + call ezfio_set_mo_one_e_ints_io_mo_integrals_n_e('None') + call ezfio_set_mo_one_e_ints_io_mo_integrals_pseudo('None') +end + From cbb8d4aa6dce42fe8c41d9f9f85596b3c0b4e144 Mon Sep 17 00:00:00 2001 From: Emmanuel Giner Date: Mon, 13 Sep 2021 14:38:00 +0200 Subject: [PATCH 2/4] added the possibility to have a pure active space interaction --- src/basis_correction/print_routine.irp.f | 2 +- src/mu_of_r/EZFIO.cfg | 2 +- src/mu_of_r/basis_def.irp.f | 18 ++++++++++++++---- src/mu_of_r/mu_of_r_conditions.irp.f | 2 +- 4 files changed, 17 insertions(+), 7 deletions(-) diff --git a/src/basis_correction/print_routine.irp.f b/src/basis_correction/print_routine.irp.f index 05fbbf60..67c5c6c2 100644 --- a/src/basis_correction/print_routine.irp.f +++ b/src/basis_correction/print_routine.irp.f @@ -38,7 +38,7 @@ subroutine print_basis_correction write(*, '(A29,X,I3,X,A3,X,F16.10)') ' ECMD PBE-UEG , state ',istate,' = ',ecmd_pbe_ueg_mu_of_r(istate) enddo - else if(mu_of_r_potential.EQ."cas_ful")then + else if(mu_of_r_potential.EQ."cas_ful".or.mu_of_r_potential.EQ."cas_truncated".or.mu_of_r_potential.EQ."pure_act")then print*, '' print*,'Using a CAS-like two-body density to define mu(r)' print*,'This assumes that the CAS is a qualitative representation of the wave function ' diff --git a/src/mu_of_r/EZFIO.cfg b/src/mu_of_r/EZFIO.cfg index 5677b3ab..c774ec82 100644 --- a/src/mu_of_r/EZFIO.cfg +++ b/src/mu_of_r/EZFIO.cfg @@ -6,7 +6,7 @@ size: (becke_numerical_grid.n_points_final_grid,determinants.n_states) [mu_of_r_potential] type: character*(32) -doc: type of potential for the mu(r) interaction: can be [ hf| cas_ful | cas_truncated] +doc: type of potential for the mu(r) interaction: can be [ hf| cas_ful | cas_truncated | pure_act] interface: ezfio, provider, ocaml default: hf diff --git a/src/mu_of_r/basis_def.irp.f b/src/mu_of_r/basis_def.irp.f index 4da27cb0..fff9f581 100644 --- a/src/mu_of_r/basis_def.irp.f +++ b/src/mu_of_r/basis_def.irp.f @@ -76,7 +76,11 @@ BEGIN_PROVIDER [integer, n_basis_orb] ! ! It corresponds to all MOs except those defined as "deleted" END_DOC - n_basis_orb = n_all_but_del_orb + if(mu_of_r_potential == "pure_act")then + n_basis_orb = n_act_orb + else + n_basis_orb = n_all_but_del_orb + endif END_PROVIDER BEGIN_PROVIDER [integer, list_basis, (n_basis_orb)] @@ -89,9 +93,15 @@ BEGIN_PROVIDER [integer, list_basis, (n_basis_orb)] ! It corresponds to all MOs except those defined as "deleted" END_DOC integer :: i - do i = 1, n_all_but_del_orb - list_basis(i) = list_all_but_del_orb(i) - enddo + if(mu_of_r_potential == "pure_act")then + do i = 1, n_act_orb + list_basis(i) = list_act(i) + enddo + else + do i = 1, n_all_but_del_orb + list_basis(i) = list_all_but_del_orb(i) + enddo + endif END_PROVIDER BEGIN_PROVIDER [double precision, basis_mos_in_r_array, (n_basis_orb,n_points_final_grid)] diff --git a/src/mu_of_r/mu_of_r_conditions.irp.f b/src/mu_of_r/mu_of_r_conditions.irp.f index 148c65b3..5c41acdc 100644 --- a/src/mu_of_r/mu_of_r_conditions.irp.f +++ b/src/mu_of_r/mu_of_r_conditions.irp.f @@ -26,7 +26,7 @@ do ipoint = 1, n_points_final_grid if(mu_of_r_potential.EQ."hf")then mu_of_r_prov(ipoint,istate) = mu_of_r_hf(ipoint) - else if(mu_of_r_potential.EQ."cas_ful".or.mu_of_r_potential.EQ."cas_truncated")then + else if(mu_of_r_potential.EQ."cas_ful".or.mu_of_r_potential.EQ."cas_truncated".or.mu_of_r_potential.EQ."pure_act")then mu_of_r_prov(ipoint,istate) = mu_of_r_psi_cas(ipoint,istate) else print*,'you requested the following mu_of_r_potential' From eaa0764f556370be460a93f0974a9b80b7aad788 Mon Sep 17 00:00:00 2001 From: Emmanuel Giner Date: Mon, 20 Sep 2021 14:04:09 +0200 Subject: [PATCH 3/4] added one_e_dm_ao --- src/ao_two_e_ints/map_integrals.irp.f | 2 ++ src/determinants/density_matrix.irp.f | 8 ++++++++ 2 files changed, 10 insertions(+) diff --git a/src/ao_two_e_ints/map_integrals.irp.f b/src/ao_two_e_ints/map_integrals.irp.f index 55b2d5e2..c3b206e1 100644 --- a/src/ao_two_e_ints/map_integrals.irp.f +++ b/src/ao_two_e_ints/map_integrals.irp.f @@ -327,6 +327,8 @@ double precision function get_ao_two_e_integral(i,j,k,l,map) result(result) implicit none BEGIN_DOC ! Gets one AO bi-electronic integral from the AO map + ! + ! i,j,k,l in physicist notation END_DOC integer, intent(in) :: i,j,k,l integer(key_kind) :: idx diff --git a/src/determinants/density_matrix.irp.f b/src/determinants/density_matrix.irp.f index c72dbc9e..6a474662 100644 --- a/src/determinants/density_matrix.irp.f +++ b/src/determinants/density_matrix.irp.f @@ -415,6 +415,14 @@ END_PROVIDER END_PROVIDER +BEGIN_PROVIDER [ double precision, one_e_dm_ao, (ao_num, ao_num)] + implicit none + BEGIN_DOC + ! one_e_dm_ao = one_e_dm_ao_alpha + one_e_dm_ao_beta + END_DOC + one_e_dm_ao = one_e_dm_ao_alpha + one_e_dm_ao_beta +END_PROVIDER + subroutine get_occupation_from_dets(istate,occupation) implicit none From 91937d5346852e521796da471df93e320feaa342 Mon Sep 17 00:00:00 2001 From: Emmanuel Giner Date: Wed, 22 Sep 2021 17:39:06 +0200 Subject: [PATCH 4/4] removed zero in b/src/dft_utils_func/on_top_from_ueg.irp.f --- external/qp2-dependencies | 2 +- src/dft_utils_func/on_top_from_ueg.irp.f | 4 +++- 2 files changed, 4 insertions(+), 2 deletions(-) diff --git a/external/qp2-dependencies b/external/qp2-dependencies index ae9397c3..0e1ca913 160000 --- a/external/qp2-dependencies +++ b/external/qp2-dependencies @@ -1 +1 @@ -Subproject commit ae9397c3b4e689a487fdd4b1425af5f519d3ea82 +Subproject commit 0e1ca91313e4b6ba3ea042b6378c3ff483781fb1 diff --git a/src/dft_utils_func/on_top_from_ueg.irp.f b/src/dft_utils_func/on_top_from_ueg.irp.f index 717081a7..5b964a03 100644 --- a/src/dft_utils_func/on_top_from_ueg.irp.f +++ b/src/dft_utils_func/on_top_from_ueg.irp.f @@ -37,13 +37,15 @@ double precision function g0_UEG_mu_inf(rho_a,rho_b) rs = (3d0 / (4d0*pi*rho))**(1d0/3d0) ! JT: serious bug fixed 20/03/19 x = -d2*rs if(dabs(x).lt.50.d0)then - g0_UEG_mu_inf= 0.5d0 * (1d0- B*rs + C*rs**2 + D*rs**3 + E*rs**4)*dexp(x) +! g0_UEG_mu_inf= 0.5d0 * (1d0- B*rs + C*rs**2 + D*rs**3 + E*rs**4)*dexp(x) + g0_UEG_mu_inf= 0.5d0 * (1d0+ rs* (-B + rs*(C + rs*(D + rs*E))))*dexp(x) else g0_UEG_mu_inf= 0.d0 endif else g0_UEG_mu_inf= 0.d0 endif + g0_UEG_mu_inf = max(g0_UEG_mu_inf,1.d-14) end