diff --git a/input/basis b/input/basis index eb58543..05207b0 100644 --- a/input/basis +++ b/input/basis @@ -1,29 +1,24 @@ -1 6 -S 8 1.00 - 17880.0000000 0.0007380 - 2683.0000000 0.0056770 - 611.5000000 0.0288830 - 173.5000000 0.1085400 - 56.6400000 0.2909070 - 20.4200000 0.4483240 - 7.8100000 0.2580260 - 1.6530000 0.0150630 -S 8 1.00 - 17880.0000000 -0.0001720 - 2683.0000000 -0.0013570 - 611.5000000 -0.0067370 - 173.5000000 -0.0276630 - 56.6400000 -0.0762080 - 20.4200000 -0.1752270 - 7.8100000 -0.1070380 - 1.6530000 0.5670500 +1 10 +S 4 1.00 + 528.5000000 0.0009400 + 79.3100000 0.0072140 + 18.0500000 0.0359750 + 5.0850000 0.1277820 S 1 1.00 - 0.4869000 1.0000000 -P 3 1.00 - 28.3900000 0.0460870 - 6.2700000 0.2401810 - 1.6950000 0.5087440 + 1.6090000 1.0000000 +S 1 1.00 + 0.5363000 1.0000000 +S 1 1.00 + 0.1833000 1.0000000 P 1 1.00 - 0.4317000 1.0000000 + 5.9940000 1.0000000 +P 1 1.00 + 1.7450000 1.0000000 +P 1 1.00 + 0.5600000 1.0000000 D 1 1.00 - 2.2020000 1.0000000 + 4.2990000 1.0000000 +D 1 1.00 + 1.2230000 1.0000000 +F 1 1.00 + 2.6800000 1.0000000 diff --git a/input/molecule b/input/molecule index edeba31..c78e87e 100644 --- a/input/molecule +++ b/input/molecule @@ -1,4 +1,4 @@ # nAt nEla nElb nCore nRyd - 1 5 5 0 0 + 1 1 1 0 0 # Znuc x y z - Ne 0.0 0.0 0.0 + He 0.0 0.0 0.0 diff --git a/input/weight b/input/weight index eb58543..05207b0 100644 --- a/input/weight +++ b/input/weight @@ -1,29 +1,24 @@ -1 6 -S 8 1.00 - 17880.0000000 0.0007380 - 2683.0000000 0.0056770 - 611.5000000 0.0288830 - 173.5000000 0.1085400 - 56.6400000 0.2909070 - 20.4200000 0.4483240 - 7.8100000 0.2580260 - 1.6530000 0.0150630 -S 8 1.00 - 17880.0000000 -0.0001720 - 2683.0000000 -0.0013570 - 611.5000000 -0.0067370 - 173.5000000 -0.0276630 - 56.6400000 -0.0762080 - 20.4200000 -0.1752270 - 7.8100000 -0.1070380 - 1.6530000 0.5670500 +1 10 +S 4 1.00 + 528.5000000 0.0009400 + 79.3100000 0.0072140 + 18.0500000 0.0359750 + 5.0850000 0.1277820 S 1 1.00 - 0.4869000 1.0000000 -P 3 1.00 - 28.3900000 0.0460870 - 6.2700000 0.2401810 - 1.6950000 0.5087440 + 1.6090000 1.0000000 +S 1 1.00 + 0.5363000 1.0000000 +S 1 1.00 + 0.1833000 1.0000000 P 1 1.00 - 0.4317000 1.0000000 + 5.9940000 1.0000000 +P 1 1.00 + 1.7450000 1.0000000 +P 1 1.00 + 0.5600000 1.0000000 D 1 1.00 - 2.2020000 1.0000000 + 4.2990000 1.0000000 +D 1 1.00 + 1.2230000 1.0000000 +F 1 1.00 + 2.6800000 1.0000000 diff --git a/src/QuAcK/QuAcK.f90 b/src/QuAcK/QuAcK.f90 index 2db13b5..b222703 100644 --- a/src/QuAcK/QuAcK.f90 +++ b/src/QuAcK/QuAcK.f90 @@ -623,9 +623,10 @@ program QuAcK call AO_values_grid(nBas,nShell,CenterShell,TotAngMomShell,KShell,DShell,ExpShell,nGrid,root,AO,dAO) call density(nGrid,nBas,PHF(:,:,1),AO(:,:),rho(:)) call MO_values_grid(nBas,nGrid,cHF(:,:,1),AO,dAO,MO,dMO) - call f_grid(nBas,nO(1),nGrid,MO,ERI_MO_basis,f) + call f_grid(nBas,nO(1),nGrid,weight,MO,ERI_MO_basis,f) call mu_grid(nGrid,rho,f,mu) call ec_srlda(nGrid,weight,rho,mu) + call fc_srlda(nEl(1),nBas,nGrid,weight,MO,rho,mu) !------------------------------------------------------------------------ ! End of QuAcK diff --git a/src/QuAcK/ec_srlda.f90 b/src/QuAcK/ec_srlda.f90 index 85b653a..7bcac27 100644 --- a/src/QuAcK/ec_srlda.f90 +++ b/src/QuAcK/ec_srlda.f90 @@ -18,7 +18,7 @@ subroutine ec_srlda(nGrid,weight,rho,mu) double precision :: r double precision :: rs double precision :: ecsr - double precision :: ec + double precision :: ec,vcup,vcdw double precision,parameter :: thres = 1d-15 ! Initialization @@ -33,9 +33,10 @@ subroutine ec_srlda(nGrid,weight,rho,mu) rs = (4d0*pi*r/3d0)**(-1d0/3d0) - call srlda(rs,mu(iG),ecsr) +! call srlda(rs,mu(iG),ecsr) + call lsdsr(rs,0d0,mu(iG),ecsr,vcup,vcdw) - ec = ec + weight(iG)*ecsr*rho(iG) + ec = ec + weight(iG)*ecsr*r end if diff --git a/src/QuAcK/f_grid.f90 b/src/QuAcK/f_grid.f90 index 63bd3f8..edf476f 100644 --- a/src/QuAcK/f_grid.f90 +++ b/src/QuAcK/f_grid.f90 @@ -1,4 +1,4 @@ -subroutine f_grid(nBas,nO,nGrid,MO,ERI,f) +subroutine f_grid(nBas,nO,nGrid,weight,MO,ERI,f) ! Compute f @@ -10,6 +10,7 @@ subroutine f_grid(nBas,nO,nGrid,MO,ERI,f) integer,intent(in) :: nBas integer,intent(in) :: nO integer,intent(in) :: nGrid + double precision,intent(in) :: weight(nGrid) double precision,intent(in) :: MO(nBas,nGrid) double precision,intent(in) :: ERI(nBas,nBas,nBas,nBas) @@ -18,6 +19,7 @@ subroutine f_grid(nBas,nO,nGrid,MO,ERI,f) integer :: p,q integer :: i,j integer :: iG + double precision :: toto ! Output variables @@ -27,6 +29,26 @@ subroutine f_grid(nBas,nO,nGrid,MO,ERI,f) f(:) = 0d0 + do p=1,nBas + do i=1,nO + do j=1,nO + do iG=1,ngrid + + f(iG) = f(iG) + MO(i,iG)*MO(p,iG)*ERI(i,j,p,j) + + end do + end do + end do + end do + + toto=0d0 + do iG=1,nGrid + toto = toto + weight(iG)*f(iG) + end do + print*,'toto=',toto + + f(:) = 0d0 + do p=1,nBas do q=1,nBas do i=1,nO diff --git a/src/QuAcK/mu_grid.f90 b/src/QuAcK/mu_grid.f90 index 23a6ff2..b86425e 100644 --- a/src/QuAcK/mu_grid.f90 +++ b/src/QuAcK/mu_grid.f90 @@ -27,7 +27,7 @@ subroutine mu_grid(nGrid,rho,f,mu) do iG=1,ngrid - n2 = rho(iG)**2 + n2 = 0.25d0*rho(iG)**2 if(abs(n2) > thres) then