From 9c1f81180076f5721d53123013437e2374e4c390 Mon Sep 17 00:00:00 2001 From: Anthony Scemama Date: Thu, 20 Mar 2025 12:14:31 +0100 Subject: [PATCH] Fixed representations in spherical AOs --- src/ao_one_e_ints/ao_ortho_canonical.irp.f | 5 ++-- src/ao_one_e_ints/ao_overlap.irp.f | 4 ++-- src/ao_one_e_ints/kin_ao_ints.irp.f | 4 ++-- src/ao_one_e_ints/pot_ao_ints.irp.f | 4 ++-- src/mo_basis/mos.irp.f | 27 ++++++++++++++++++---- src/scf_utils/roothaan_hall_scf.irp.f | 4 ---- 6 files changed, 31 insertions(+), 17 deletions(-) diff --git a/src/ao_one_e_ints/ao_ortho_canonical.irp.f b/src/ao_one_e_ints/ao_ortho_canonical.irp.f index bdfa2096..523e49f7 100644 --- a/src/ao_one_e_ints/ao_ortho_canonical.irp.f +++ b/src/ao_one_e_ints/ao_ortho_canonical.irp.f @@ -99,8 +99,9 @@ BEGIN_PROVIDER [ double precision, ao_cart_to_sphe_inv, (ao_sphe_num,ao_num) ] R(:,:) = ao_cart_to_sphe_coef(:,:) call dgemm('N','T', m, m, k, 1.d0, R, k, R, k, 0.d0, S, m) - call get_pseudo_inverse(S, k, k, m, Sinv, k, 1.d-20) - call dgemm('T','T', m, m, k, 1.d0, R, k, Sinv, k, 0.d0, Rinv, m) + call get_pseudo_inverse(S, k, k, m, Sinv, k, 1.d-12) + call dgemm('T','N', m, m, k, 1.d0, R, k, Sinv, k, 0.d0, Rinv, m) + integer :: i do i=1,ao_num diff --git a/src/ao_one_e_ints/ao_overlap.irp.f b/src/ao_one_e_ints/ao_overlap.irp.f index cc676c4b..17eb7cbc 100644 --- a/src/ao_one_e_ints/ao_overlap.irp.f +++ b/src/ao_one_e_ints/ao_overlap.irp.f @@ -317,12 +317,12 @@ BEGIN_PROVIDER [ double precision, ao_sphe_overlap, (ao_sphe_num,ao_sphe_num) ] double precision, allocatable :: tmp(:,:) allocate (tmp(ao_sphe_num,ao_num)) - call dgemm('T','N',ao_sphe_num,ao_num,ao_num, 1.d0, & + call dgemm('N','N',ao_sphe_num,ao_num,ao_num, 1.d0, & ao_cart_to_sphe_inv,size(ao_cart_to_sphe_inv,1), & ao_overlap,size(ao_overlap,1), 0.d0, & tmp, size(tmp,1)) - call dgemm('N','N',ao_sphe_num,ao_sphe_num,ao_num, 1.d0, & + call dgemm('N','T',ao_sphe_num,ao_sphe_num,ao_num, 1.d0, & tmp, size(tmp,1), & ao_cart_to_sphe_inv,size(ao_cart_to_sphe_inv,1), 0.d0, & ao_sphe_overlap,size(ao_sphe_overlap,1)) diff --git a/src/ao_one_e_ints/kin_ao_ints.irp.f b/src/ao_one_e_ints/kin_ao_ints.irp.f index f450721d..1e3e684d 100644 --- a/src/ao_one_e_ints/kin_ao_ints.irp.f +++ b/src/ao_one_e_ints/kin_ao_ints.irp.f @@ -199,12 +199,12 @@ BEGIN_PROVIDER [ double precision, ao_sphe_kinetic_integrals, (ao_sphe_num,ao_sp double precision, allocatable :: tmp(:,:) allocate (tmp(ao_sphe_num,ao_num)) - call dgemm('T','N',ao_sphe_num,ao_num,ao_num, 1.d0, & + call dgemm('N','N',ao_sphe_num,ao_num,ao_num, 1.d0, & ao_cart_to_sphe_inv,size(ao_cart_to_sphe_inv,1), & ao_kinetic_integrals,size(ao_kinetic_integrals,1), 0.d0, & tmp, size(tmp,1)) - call dgemm('N','N',ao_sphe_num,ao_sphe_num,ao_num, 1.d0, & + call dgemm('N','T',ao_sphe_num,ao_sphe_num,ao_num, 1.d0, & tmp, size(tmp,1), & ao_cart_to_sphe_inv,size(ao_cart_to_sphe_inv,1), 0.d0, & ao_sphe_kinetic_integrals,size(ao_sphe_kinetic_integrals,1)) diff --git a/src/ao_one_e_ints/pot_ao_ints.irp.f b/src/ao_one_e_ints/pot_ao_ints.irp.f index c1544a5d..b4b4aa02 100644 --- a/src/ao_one_e_ints/pot_ao_ints.irp.f +++ b/src/ao_one_e_ints/pot_ao_ints.irp.f @@ -618,12 +618,12 @@ BEGIN_PROVIDER [ double precision, ao_sphe_integrals_n_e, (ao_sphe_num,ao_sphe_n double precision, allocatable :: tmp(:,:) allocate (tmp(ao_sphe_num,ao_num)) - call dgemm('T','N',ao_sphe_num,ao_num,ao_num, 1.d0, & + call dgemm('N','N',ao_sphe_num,ao_num,ao_num, 1.d0, & ao_cart_to_sphe_inv,size(ao_cart_to_sphe_inv,1), & ao_integrals_n_e,size(ao_integrals_n_e,1), 0.d0, & tmp, size(tmp,1)) - call dgemm('N','N',ao_sphe_num,ao_sphe_num,ao_num, 1.d0, & + call dgemm('N','T',ao_sphe_num,ao_sphe_num,ao_num, 1.d0, & tmp, size(tmp,1), & ao_cart_to_sphe_inv,size(ao_cart_to_sphe_inv,1), 0.d0, & ao_sphe_integrals_n_e,size(ao_sphe_integrals_n_e,1)) diff --git a/src/mo_basis/mos.irp.f b/src/mo_basis/mos.irp.f index 1eecca6c..4da8ac9b 100644 --- a/src/mo_basis/mos.irp.f +++ b/src/mo_basis/mos.irp.f @@ -352,14 +352,31 @@ BEGIN_PROVIDER [ double precision, mo_sphe_coef, (ao_sphe_num, mo_num) ] BEGIN_DOC ! MO coefficients in the basis of spherical harmonics AOs. END_DOC - double precision, allocatable :: tmp(:,:) - allocate (tmp(ao_sphe_num,ao_num)) + double precision, allocatable, dimension (:,:) :: C0, S, tmp + allocate(C0(ao_sphe_num,mo_num), S(mo_num,mo_num), tmp(mo_num,ao_sphe_num)) - call dgemm('T','N',ao_sphe_num,ao_num,ao_num, 1.d0, & + call dgemm('T','N',ao_sphe_num,mo_num,ao_num, 1.d0, & ao_cart_to_sphe_coef,ao_num, & mo_coef,size(mo_coef,1), 0.d0, & - mo_sphe_coef, size(mo_sphe_coef,1)) + C0, size(C0,1)) - deallocate (tmp) + ! C0^T S S^0 + call dgemm('T','N',mo_num,ao_sphe_num,ao_sphe_num, 1.d0, & + C0,size(C0,1), & + ao_sphe_overlap,size(ao_sphe_overlap,1), 0.d0, & + tmp, size(tmp,1)) + + call dgemm('N','N',mo_num,mo_num,ao_sphe_num, 1.d0, & + tmp, size(tmp,1), & + C0,size(C0,1), 0.d0, & + S, size(S,1)) + + integer :: m + m = ao_sphe_num + mo_sphe_coef = C0 + call ortho_lowdin(S,size(S,1),mo_num,mo_sphe_coef,ao_sphe_num,m,1.d-10) + + + deallocate (tmp, S, C0) END_PROVIDER diff --git a/src/scf_utils/roothaan_hall_scf.irp.f b/src/scf_utils/roothaan_hall_scf.irp.f index 66bedb07..7a73121b 100644 --- a/src/scf_utils/roothaan_hall_scf.irp.f +++ b/src/scf_utils/roothaan_hall_scf.irp.f @@ -239,10 +239,6 @@ END_DOC endif i = j enddo - TOUCH mo_coef - call mo_as_eigvectors_of_mo_matrix(Fock_matrix_mo,size(Fock_matrix_mo,1), & - size(Fock_matrix_mo,2),mo_label,1,.true.) - call restore_symmetry(ao_num, mo_num, mo_coef, size(mo_coef,1), 1.d-10) if(do_mom)then call reorder_mo_max_overlap