From 1a12e7f30882684cdc8dddcf96df05449ee54613 Mon Sep 17 00:00:00 2001 From: Anthony Scemama Date: Mon, 9 Jan 2023 10:51:12 +0100 Subject: [PATCH] Fixed AO normalization problem --- src/ao_basis/EZFIO.cfg | 6 +++--- src/ao_basis/aos.irp.f | 15 ++++++++------- src/ao_one_e_ints/pot_ao_ints.irp.f | 4 ++-- src/ao_one_e_ints/pseudopot.f90 | 28 ++++++++++++++-------------- 4 files changed, 27 insertions(+), 26 deletions(-) diff --git a/src/ao_basis/EZFIO.cfg b/src/ao_basis/EZFIO.cfg index dd61b1be..3ac16446 100644 --- a/src/ao_basis/EZFIO.cfg +++ b/src/ao_basis/EZFIO.cfg @@ -57,13 +57,13 @@ default: false [ao_normalized] type: logical -doc: Use normalized basis functions +doc: Normalize the atomic orbitals interface: ezfio, provider -default: true +default: false [primitives_normalized] type: logical -doc: Use normalized primitive functions +doc: Normalize the primitive basis functions interface: ezfio, provider default: true diff --git a/src/ao_basis/aos.irp.f b/src/ao_basis/aos.irp.f index 3a9e9fb7..dafea9c4 100644 --- a/src/ao_basis/aos.irp.f +++ b/src/ao_basis/aos.irp.f @@ -63,15 +63,14 @@ END_PROVIDER ! Coefficients including the |AO| normalization END_DOC - do i=1,ao_num - l = ao_shell(i) - ao_coef_normalized(i,:) = shell_coef(l,:) * shell_normalization_factor(l) - end do double precision :: norm,overlap_x,overlap_y,overlap_z,C_A(3), c - integer :: l, powA(3), nz + integer :: l, powA(3) + integer, parameter :: nz=100 integer :: i,j,k - nz=100 + + ao_coef_normalized(:,:) = ao_coef(:,:) + C_A = 0.d0 do i=1,ao_num @@ -80,7 +79,7 @@ END_PROVIDER powA(2) = ao_power(i,2) powA(3) = ao_power(i,3) - ! Normalization of the primitives + ! GAMESS-type normalization of the primitives if (primitives_normalized) then do j=1,ao_prim_num(i) call overlap_gaussian_xyz(C_A,C_A,ao_expo(i,j),ao_expo(i,j), & @@ -91,6 +90,7 @@ END_PROVIDER ! Normalization of the contracted basis functions if (ao_normalized) then norm = 0.d0 + l = ao_shell(i) do j=1,ao_prim_num(i) do k=1,ao_prim_num(i) call overlap_gaussian_xyz(C_A,C_A,ao_expo(i,j),ao_expo(i,k),powA,powA,overlap_x,overlap_y,overlap_z,c,nz) @@ -98,6 +98,7 @@ END_PROVIDER enddo enddo ao_coef_normalization_factor(i) = 1.d0/dsqrt(norm) + ao_coef_normalized(i,:) *= ao_coef_normalization_factor(i) else ao_coef_normalization_factor(i) = 1.d0 endif 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 dc19f6c7..928053ad 100644 --- a/src/ao_one_e_ints/pot_ao_ints.irp.f +++ b/src/ao_one_e_ints/pot_ao_ints.irp.f @@ -18,6 +18,8 @@ BEGIN_PROVIDER [ double precision, ao_integrals_n_e, (ao_num,ao_num)] double precision :: A_center(3),B_center(3),C_center(3) double precision :: overlap_x,overlap_y,overlap_z,overlap,dx,NAI_pol_mult + ao_integrals_n_e = 0.d0 + if (read_ao_integrals_n_e) then call ezfio_get_ao_one_e_ints_ao_integrals_n_e(ao_integrals_n_e) @@ -36,8 +38,6 @@ BEGIN_PROVIDER [ double precision, ao_integrals_n_e, (ao_num,ao_num)] else - ao_integrals_n_e = 0.d0 - !$OMP PARALLEL & !$OMP DEFAULT (NONE) & !$OMP PRIVATE (i,j,k,l,m,alpha,beta,A_center,B_center,C_center,power_A,power_B,& diff --git a/src/ao_one_e_ints/pseudopot.f90 b/src/ao_one_e_ints/pseudopot.f90 index 7321dff7..e02dea3b 100644 --- a/src/ao_one_e_ints/pseudopot.f90 +++ b/src/ao_one_e_ints/pseudopot.f90 @@ -1950,26 +1950,26 @@ xq(17)=-3.34785456738322 xq(18)=-3.94476404011563 xq(19)=-4.60368244955074 xq(20)=-5.38748089001123 -wq(1)= 2.229393645534151E-013 -wq(2)= 4.399340992273176E-010 -wq(3)= 1.086069370769280E-007 -wq(4)= 7.802556478532063E-006 -wq(5)= 2.283386360163528E-004 -wq(6)= 3.243773342237853E-003 -wq(7)= 2.481052088746362E-002 +wq(1)= 2.229393645534151D-013 +wq(2)= 4.399340992273176D-010 +wq(3)= 1.086069370769280D-007 +wq(4)= 7.802556478532063D-006 +wq(5)= 2.283386360163528D-004 +wq(6)= 3.243773342237853D-003 +wq(7)= 2.481052088746362D-002 wq(8)= 0.109017206020022 wq(9)= 0.286675505362834 wq(10)= 0.462243669600610 wq(11)= 0.462243669600610 wq(12)= 0.286675505362834 wq(13)= 0.109017206020022 -wq(14)= 2.481052088746362E-002 -wq(15)= 3.243773342237853E-003 -wq(16)= 2.283386360163528E-004 -wq(17)= 7.802556478532063E-006 -wq(18)= 1.086069370769280E-007 -wq(19)= 4.399340992273176E-010 -wq(20)= 2.229393645534151E-013 +wq(14)= 2.481052088746362D-002 +wq(15)= 3.243773342237853D-003 +wq(16)= 2.283386360163528D-004 +wq(17)= 7.802556478532063D-006 +wq(18)= 1.086069370769280D-007 +wq(19)= 4.399340992273176D-010 +wq(20)= 2.229393645534151D-013 npts=20 ! call gauher(xq,wq,npts)