LCPQ/quack
4
1
Fork 0
mirror of https://github.com/pfloos/quack synced 2025-01-03 10:05:59 +01:00
Code Releases Activity

cleanup qsGW

Browse Source
This commit is contained in:
Pierre-Francois Loos 2023-11-27 13:30:53 +01:00
parent 39b8b8d6f1
commit cdb9220872
4 changed files with 41 additions and 46 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

49
src/GW/qsRGW.f90
View File

@ -158,12 +158,9 @@ subroutine qsRGW(dotest,maxSCF,thresh,max_diis,doACFDT,exchange_kernel,doXBS,dop
nSCF = nSCF + 1 nSCF = nSCF + 1
! Buid Hartree matrix ! Build Hartree-exchange matrix
call Hartree_matrix_AO_basis(nBas,P,ERI_AO,J) call Hartree_matrix_AO_basis(nBas,P,ERI_AO,J)
! Compute exchange part of the self-energy
call exchange_matrix_AO_basis(nBas,P,ERI_AO,K) call exchange_matrix_AO_basis(nBas,P,ERI_AO,K)
! AO to MO transformation of two-electron integrals ! AO to MO transformation of two-electron integrals
@ -206,29 +203,6 @@ subroutine qsRGW(dotest,maxSCF,thresh,max_diis,doACFDT,exchange_kernel,doXBS,dop
if(nSCF > 1) Conv = maxval(abs(err)) if(nSCF > 1) Conv = maxval(abs(err))
! DIIS extrapolation
if(max_diis > 1) then
n_diis = min(n_diis+1,max_diis)
call DIIS_extrapolation(rcond,nBasSq,nBasSq,n_diis,err_diis,F_diis,err,F)
end if
! Diagonalize Hamiltonian in AO basis
Fp = matmul(transpose(X),matmul(F,X))
cp(:,:) = Fp(:,:)
call diagonalize_matrix(nBas,cp,eGW)
c = matmul(X,cp)
call AOtoMO(nBas,c,SigCp,SigC)
! Compute new density matrix in the AO basis
P(:,:) = 2d0*matmul(c(:,1:nO),transpose(c(:,1:nO)))
!------------------------------------------------------------------------
! Compute total energy
!------------------------------------------------------------------------
! Kinetic energy ! Kinetic energy
ET = trace_matrix(nBas,matmul(P,T)) ET = trace_matrix(nBas,matmul(P,T))
@ -249,6 +223,27 @@ subroutine qsRGW(dotest,maxSCF,thresh,max_diis,doACFDT,exchange_kernel,doXBS,dop
EqsGW = ET + EV + EJ + EK EqsGW = ET + EV + EJ + EK
! DIIS extrapolation
if(max_diis > 1) then
n_diis = min(n_diis+1,max_diis)
call DIIS_extrapolation(rcond,nBasSq,nBasSq,n_diis,err_diis,F_diis,err,F)
end if
! Diagonalize Hamiltonian in AO basis
Fp = matmul(transpose(X),matmul(F,X))
cp(:,:) = Fp(:,:)
call diagonalize_matrix(nBas,cp,eGW)
c = matmul(X,cp)
call AOtoMO(nBas,c,SigCp,SigC)
! Density matrix
P(:,:) = 2d0*matmul(c(:,1:nO),transpose(c(:,1:nO)))
! Print results ! Print results
call dipole_moment(nBas,P,nNuc,ZNuc,rNuc,dipole_int_AO,dipole) call dipole_moment(nBas,P,nNuc,ZNuc,rNuc,dipole_int_AO,dipole)

6
src/GW/ufG0W0.f90
View File

@ -209,7 +209,7 @@ subroutine ufG0W0(dotest,TDA_W,nBas,nC,nO,nV,nR,nS,ENuc,ERHF,ERI,eHF)
timing = end_timing - start_timing timing = end_timing - start_timing
write(*,*) write(*,*)
write(*,'(A65,1X,F9.3,A8)') 'Total wall time for construction of supermatrix = ',timing,' seconds' write(*,'(A65,1X,F9.3,A8)') 'Total CPU time for construction of supermatrix = ',timing,' seconds'
write(*,*) write(*,*)
else else
@ -324,7 +324,7 @@ subroutine ufG0W0(dotest,TDA_W,nBas,nC,nO,nV,nR,nS,ENuc,ERHF,ERI,eHF)
timing = end_timing - start_timing timing = end_timing - start_timing
write(*,*) write(*,*)
write(*,'(A65,1X,F9.3,A8)') 'Total wall time for construction of supermatrix = ',timing,' seconds' write(*,'(A65,1X,F9.3,A8)') 'Total CPU time for construction of supermatrix = ',timing,' seconds'
write(*,*) write(*,*)
end if end if
@ -341,7 +341,7 @@ subroutine ufG0W0(dotest,TDA_W,nBas,nC,nO,nV,nR,nS,ENuc,ERHF,ERI,eHF)
timing = end_timing - start_timing timing = end_timing - start_timing
write(*,*) write(*,*)
write(*,'(A65,1X,F9.3,A8)') 'Total wall time for diagonalization of supermatrix = ',timing,' seconds' write(*,'(A65,1X,F9.3,A8)') 'Total CPU time for diagonalization of supermatrix = ',timing,' seconds'
write(*,*) write(*,*)
!-----------------! !-----------------!

6
src/GW/ufGW.f90
View File

@ -213,7 +213,7 @@ subroutine ufGW(dotest,TDA_W,nBas,nC,nO,nV,nR,nS,ENuc,ERHF,ERI,eHF)
timing = end_timing - start_timing timing = end_timing - start_timing
write(*,*) write(*,*)
write(*,'(A65,1X,F9.3,A8)') 'Total wall time for construction of supermatrix = ',timing,' seconds' write(*,'(A65,1X,F9.3,A8)') 'Total CPU time for construction of supermatrix = ',timing,' seconds'
write(*,*) write(*,*)
else else
@ -338,7 +338,7 @@ subroutine ufGW(dotest,TDA_W,nBas,nC,nO,nV,nR,nS,ENuc,ERHF,ERI,eHF)
timing = end_timing - start_timing timing = end_timing - start_timing
write(*,*) write(*,*)
write(*,'(A65,1X,F9.3,A8)') 'Total wall time for construction of supermatrix = ',timing,' seconds' write(*,'(A65,1X,F9.3,A8)') 'Total CPU time for construction of supermatrix = ',timing,' seconds'
write(*,*) write(*,*)
end if end if
@ -355,7 +355,7 @@ subroutine ufGW(dotest,TDA_W,nBas,nC,nO,nV,nR,nS,ENuc,ERHF,ERI,eHF)
timing = end_timing - start_timing timing = end_timing - start_timing
write(*,*) write(*,*)
write(*,'(A65,1X,F9.3,A8)') 'Total wall time for diagonalization of supermatrix = ',timing,' seconds' write(*,'(A65,1X,F9.3,A8)') 'Total CPU time for diagonalization of supermatrix = ',timing,' seconds'
write(*,*) write(*,*)
!-----------------! !-----------------!

26
src/HF/RHF.f90
View File

@ -104,43 +104,43 @@ subroutine RHF(dotest,maxSCF,thresh,max_diis,guess_type,level_shift,nNuc,ZNuc,rN
do while(Conv > thresh .and. nSCF < maxSCF) do while(Conv > thresh .and. nSCF < maxSCF)
! Increment ! Increment
nSCF = nSCF + 1 nSCF = nSCF + 1
! Build Fock matrix ! Build Fock matrix
call Hartree_matrix_AO_basis(nBas,P,ERI,J) call Hartree_matrix_AO_basis(nBas,P,ERI,J)
call exchange_matrix_AO_basis(nBas,P,ERI,K) call exchange_matrix_AO_basis(nBas,P,ERI,K)
F(:,:) = Hc(:,:) + J(:,:) + 0.5d0*K(:,:) F(:,:) = Hc(:,:) + J(:,:) + 0.5d0*K(:,:)
! Check convergence ! Check convergence
err = matmul(F,matmul(P,S)) - matmul(matmul(S,P),F) err = matmul(F,matmul(P,S)) - matmul(matmul(S,P),F)
if(nSCF > 1) Conv = maxval(abs(err)) if(nSCF > 1) Conv = maxval(abs(err))
! Kinetic energy ! Kinetic energy
ET = trace_matrix(nBas,matmul(P,T)) ET = trace_matrix(nBas,matmul(P,T))
! Potential energy ! Potential energy
EV = trace_matrix(nBas,matmul(P,V)) EV = trace_matrix(nBas,matmul(P,V))
! Hartree energy ! Hartree energy
EJ = 0.5d0*trace_matrix(nBas,matmul(P,J)) EJ = 0.5d0*trace_matrix(nBas,matmul(P,J))
! Exchange energy ! Exchange energy
EK = 0.25d0*trace_matrix(nBas,matmul(P,K)) EK = 0.25d0*trace_matrix(nBas,matmul(P,K))
! Total energy ! Total energy
ERHF = ET + EV + EJ + EK ERHF = ET + EV + EJ + EK
! DIIS extrapolation ! DIIS extrapolation
if(max_diis > 1) then if(max_diis > 1) then
@ -149,22 +149,22 @@ subroutine RHF(dotest,maxSCF,thresh,max_diis,guess_type,level_shift,nNuc,ZNuc,rN
end if end if
! Level-shifting ! Level shift
if(level_shift > 0d0 .and. Conv > thresh) call level_shifting(level_shift,nBas,nO,S,c,F) if(level_shift > 0d0 .and. Conv > thresh) call level_shifting(level_shift,nBas,nO,S,c,F)
! Diagonalize Fock matrix ! Diagonalize Fock matrix
Fp = matmul(transpose(X),matmul(F,X)) Fp = matmul(transpose(X),matmul(F,X))
cp(:,:) = Fp(:,:) cp(:,:) = Fp(:,:)
call diagonalize_matrix(nBas,cp,eHF) call diagonalize_matrix(nBas,cp,eHF)
c = matmul(X,cp) c = matmul(X,cp)
! Density matrix ! Density matrix
P(:,:) = 2d0*matmul(c(:,1:nO),transpose(c(:,1:nO))) P(:,:) = 2d0*matmul(c(:,1:nO),transpose(c(:,1:nO)))
! Dump results ! Dump results
write(*,'(1X,A1,1X,I3,1X,A1,1X,F16.10,1X,A1,1X,F16.10,1X,A1,1X,F16.10,1X,A1,1X,E10.2,1X,A1,1X)') & write(*,'(1X,A1,1X,I3,1X,A1,1X,F16.10,1X,A1,1X,F16.10,1X,A1,1X,F16.10,1X,A1,1X,E10.2,1X,A1,1X)') &
'|',nSCF,'|',ERHF + ENuc,'|',EJ,'|',EK,'|',Conv,'|' '|',nSCF,'|',ERHF + ENuc,'|',EJ,'|',EK,'|',Conv,'|'