10
1
mirror of https://github.com/pfloos/quack synced 2024-12-22 20:34:46 +01:00
This commit is contained in:
Pierre-Francois Loos 2019-04-18 23:22:23 +02:00
parent 31e3707ce6
commit deb5144161
8 changed files with 178 additions and 28 deletions

8
GoSph
View File

@ -9,9 +9,9 @@ if [ $# = 2 ]
then
cp examples/molecule.Sph_"$1" input/molecule
cp examples/basis.Sph.Ylm"$2" input/basis
cp int/Sph_ERI_"$2".dat int/ERI.dat
cp int/Sph_Kin_"$2".dat int/Kin.dat
cp int/Sph_Nuc_"$2".dat int/Nuc.dat
cp int/Sph_Ov_"$2".dat int/Ov.dat
cp ~/Integrals/QuAcK_Sph/Sph_ERI_"$2".dat int/ERI.dat
cp ~/Integrals/QuAcK_Sph/Sph_Kin_"$2".dat int/Kin.dat
cp ~/Integrals/QuAcK_Sph/Sph_Nuc_"$2".dat int/Nuc.dat
cp ~/Integrals/QuAcK_Sph/Sph_Ov_"$2".dat int/Ov.dat
./bin/QuAcK
fi

View File

@ -2,7 +2,7 @@
integer,parameter :: nspin = 2
integer,parameter :: nsp = 3
integer,parameter :: maxEns = 10
integer,parameter :: maxShell = 50
integer,parameter :: maxShell = 1024
integer,parameter :: n1eInt = 3
integer,parameter :: n2eInt = 4
integer,parameter :: n3eInt = 3

View File

@ -1,9 +1,129 @@
1 3
S 3 1.00
38.3600000 0.0238090
5.7700000 0.1548910
1.2400000 0.4699870
1 64
S 1 1.00
0.2976000 1.0000000
P 1 1.00
1.2750000 1.0000000
1.0000000 1.0000000
S 1 1.00
1.0000000 1.0000000
S 1 1.00
1.0000000 1.0000000
S 1 1.00
1.0000000 1.0000000
S 1 1.00
1.0000000 1.0000000
S 1 1.00
1.0000000 1.0000000
S 1 1.00
1.0000000 1.0000000
S 1 1.00
1.0000000 1.0000000
S 1 1.00
1.0000000 1.0000000
S 1 1.00
1.0000000 1.0000000
S 1 1.00
1.0000000 1.0000000
S 1 1.00
1.0000000 1.0000000
S 1 1.00
1.0000000 1.0000000
S 1 1.00
1.0000000 1.0000000
S 1 1.00
1.0000000 1.0000000
S 1 1.00
1.0000000 1.0000000
S 1 1.00
1.0000000 1.0000000
S 1 1.00
1.0000000 1.0000000
S 1 1.00
1.0000000 1.0000000
S 1 1.00
1.0000000 1.0000000
S 1 1.00
1.0000000 1.0000000
S 1 1.00
1.0000000 1.0000000
S 1 1.00
1.0000000 1.0000000
S 1 1.00
1.0000000 1.0000000
S 1 1.00
1.0000000 1.0000000
S 1 1.00
1.0000000 1.0000000
S 1 1.00
1.0000000 1.0000000
S 1 1.00
1.0000000 1.0000000
S 1 1.00
1.0000000 1.0000000
S 1 1.00
1.0000000 1.0000000
S 1 1.00
1.0000000 1.0000000
S 1 1.00
1.0000000 1.0000000
S 1 1.00
1.0000000 1.0000000
S 1 1.00
1.0000000 1.0000000
S 1 1.00
1.0000000 1.0000000
S 1 1.00
1.0000000 1.0000000
S 1 1.00
1.0000000 1.0000000
S 1 1.00
1.0000000 1.0000000
S 1 1.00
1.0000000 1.0000000
S 1 1.00
1.0000000 1.0000000
S 1 1.00
1.0000000 1.0000000
S 1 1.00
1.0000000 1.0000000
S 1 1.00
1.0000000 1.0000000
S 1 1.00
1.0000000 1.0000000
S 1 1.00
1.0000000 1.0000000
S 1 1.00
1.0000000 1.0000000
S 1 1.00
1.0000000 1.0000000
S 1 1.00
1.0000000 1.0000000
S 1 1.00
1.0000000 1.0000000
S 1 1.00
1.0000000 1.0000000
S 1 1.00
1.0000000 1.0000000
S 1 1.00
1.0000000 1.0000000
S 1 1.00
1.0000000 1.0000000
S 1 1.00
1.0000000 1.0000000
S 1 1.00
1.0000000 1.0000000
S 1 1.00
1.0000000 1.0000000
S 1 1.00
1.0000000 1.0000000
S 1 1.00
1.0000000 1.0000000
S 1 1.00
1.0000000 1.0000000
S 1 1.00
1.0000000 1.0000000
S 1 1.00
1.0000000 1.0000000
S 1 1.00
1.0000000 1.0000000
S 1 1.00
1.0000000 1.0000000
S 1 1.00
1.0000000 1.0000000

View File

@ -9,6 +9,6 @@
# GF2 GF3
F F
# G0W0 evGW qsGW
T T F
T F F
# MCMP2
F

View File

@ -1,4 +1,4 @@
# nAt nEla nElb nCore nRyd
1 1 1 0 0
1 25 25 0 0
# Znuc x y z
He 0.0 0.0 0.0
X 0.0 0.0 0.0

View File

@ -165,7 +165,7 @@ program QuAcK
! Read integrals
call read_integrals(nBas,S,T,V,Hc,ERI_AO_basis)
call read_integrals(nEl(:),nBas,S,T,V,Hc,ERI_AO_basis)
! Compute orthogonalization matrix

View File

@ -16,12 +16,25 @@ subroutine RHF(maxSCF,thresh,max_diis,guess_type,nBas,nO,S,T,V,Hc,ERI,X,ENuc,ERH
! Local variables
integer :: nSCF,nBasSq,n_diis
double precision :: ET,EV,EJ,EK,Conv,Gap
integer :: nSCF
integer :: nBasSq
integer :: n_diis
double precision :: ET,EV,EJ,EK
double precision :: Conv
double precision :: Gap
double precision :: rcond
double precision,external :: trace_matrix
double precision,allocatable :: error(:,:),error_diis(:,:),F_diis(:,:)
double precision,allocatable :: J(:,:),K(:,:),cp(:,:),F(:,:),Fp(:,:)
double precision,allocatable :: error(:,:)
double precision,allocatable :: error_diis(:,:)
double precision,allocatable :: F_diis(:,:)
double precision,allocatable :: J(:,:)
double precision,allocatable :: K(:,:)
double precision,allocatable :: cp(:,:)
double precision,allocatable :: F(:,:)
double precision,allocatable :: Fp(:,:)
double precision,allocatable :: ON(:)
integer :: i
! Output variables
@ -45,7 +58,7 @@ subroutine RHF(maxSCF,thresh,max_diis,guess_type,nBas,nO,S,T,V,Hc,ERI,X,ENuc,ERH
! Memory allocation
allocate(J(nBas,nBas),K(nBas,nBas),error(nBas,nBas), &
cp(nBas,nBas),Fp(nBas,nBas),F(nBas,nBas), &
cp(nBas,nBas),Fp(nBas,nBas),F(nBas,nBas),ON(nBas), &
error_diis(nBasSq,max_diis),F_diis(nBasSq,max_diis))
! Guess coefficients and eigenvalues
@ -65,6 +78,13 @@ subroutine RHF(maxSCF,thresh,max_diis,guess_type,nBas,nO,S,T,V,Hc,ERI,X,ENuc,ERH
P(:,:) = 2d0*matmul(c(:,1:nO),transpose(c(:,1:nO)))
! ON(:) = 0d0
! do i=1,nO
! ON(i) = 1d0
! ON(i) = dble(2*i-1)
! end do
! call density_matrix(nBas,ON,c,P)
! Initialization
@ -123,6 +143,8 @@ subroutine RHF(maxSCF,thresh,max_diis,guess_type,nBas,nO,S,T,V,Hc,ERI,X,ENuc,ERH
P(:,:) = 2d0*matmul(c(:,1:nO),transpose(c(:,1:nO)))
! call density_matrix(nBas,ON,c,P)
! Compute HF energy
ERHF = trace_matrix(nBas,matmul(P,Hc)) &

View File

@ -1,8 +1,9 @@
subroutine read_integrals(nBas,S,T,V,Hc,G)
subroutine read_integrals(nEl,nBas,S,T,V,Hc,G)
! Read one- and two-electron integrals from files
implicit none
include 'parameters.h'
! Input variables
@ -11,9 +12,10 @@ subroutine read_integrals(nBas,S,T,V,Hc,G)
! Local variables
logical :: debug
integer :: nEl(nspin)
integer :: mu,nu,la,si
double precision :: Ov,Kin,Nuc,ERI
double precision :: scale
double precision :: rs
! Output variables
@ -23,7 +25,13 @@ subroutine read_integrals(nBas,S,T,V,Hc,G)
debug = .false.
scale = 1d0
open(unit=1,file='input/sph')
read(1,*)
read(1,*) rs
rs = sqrt(dble(sum(nEl(:))))/2d0*rs
print*, 'Scaling integrals by ',rs
open(unit=8 ,file='int/Ov.dat')
open(unit=9 ,file='int/Kin.dat')
@ -44,7 +52,7 @@ subroutine read_integrals(nBas,S,T,V,Hc,G)
T = 0d0
do
read(9,*,end=9) mu,nu,Kin
T(mu,nu) = Kin/scale**2
T(mu,nu) = Kin/rs**2
enddo
9 close(unit=9)
@ -67,7 +75,7 @@ subroutine read_integrals(nBas,S,T,V,Hc,G)
do
read(11,*,end=11) mu,nu,la,si,ERI
ERI = ERI/scale
ERI = ERI/rs
! <12|34>
G(mu,nu,la,si) = ERI
! <32|14>