update DIIS in HF and qsGW

input/options

@@ -9,7 +9,7 @@
 # GF: maxSCF thresh  DIIS n_diis lin eta renorm reg
       256    0.00001 T    5      T   0.0 3      F
 # GW: maxSCF thresh  DIIS n_diis lin eta COHSEX SOSEX TDA_W G0W GW0 reg
-      256    0.0000001 F    5      T   0.0 F      F     F     F   F   F
+      25     0.00001  T    5      T   0.0 F      F     F     F   F   F
 # GT: maxSCF thresh  DIIS n_diis lin eta TDA_T reg
       10    0.00001  T    5      T   0.0 F     F  
 # ACFDT: AC Kx  XBS

mol/h2.xyz

@@ -1,4 +1,4 @@
 2
 
 H 0. 0. 0.
-H 0. 0. 0.740848
+H 0. 0. 1.2

src/GW/SRG_qsGW.f90

@@ -142,6 +142,8 @@ subroutine SRG_qsGW(maxSCF,thresh,max_diis,doACFDT,exchange_kernel,doXBS,BSE,BSE
   error_diis(:,:) = 0d0
   rcond           = 0d0
 
+   print*,max_diis
+
 !------------------------------------------------------------------------
 ! Main loop
 !------------------------------------------------------------------------
@@ -190,11 +192,11 @@ subroutine SRG_qsGW(maxSCF,thresh,max_diis,doACFDT,exchange_kernel,doXBS,BSE,BSE
 
     ! DIIS extrapolation 
 
-    n_diis = min(n_diis+1,max_diis)
+    if(max_diis > 1) then
-    if(abs(rcond) > 1d-7) then 
+
+      n_diis = min(n_diis+1,max_diis)
       call DIIS_extrapolation(rcond,nBasSq,nBasSq,n_diis,error_diis,F_diis,error,F)
-    else
+
-      n_diis = 0
     end if
 
     ! Diagonalize Hamiltonian in AO basis
@@ -211,7 +213,7 @@ subroutine SRG_qsGW(maxSCF,thresh,max_diis,doACFDT,exchange_kernel,doXBS,BSE,BSE
 
     ! Save quasiparticles energy for next cycle
 
-    Conv = maxval(abs(eGW - eOld))
+    Conv = maxval(abs(error))
     eOld(:) = eGW(:)
 
     !------------------------------------------------------------------------

src/GW/qsGW.f90

@@ -239,11 +239,11 @@ subroutine qsGW(maxSCF,thresh,max_diis,doACFDT,exchange_kernel,doXBS,COHSEX,BSE,
 
     ! DIIS extrapolation 
 
-    n_diis = min(n_diis+1,max_diis)
+    if(max_diis > 1) then
-    if(abs(rcond) > 1d-7) then 
+
+      n_diis = min(n_diis+1,max_diis)
       call DIIS_extrapolation(rcond,nBasSq,nBasSq,n_diis,error_diis,F_diis,error,F)
-    else
+
-      n_diis = 0
     end if
 
     ! Diagonalize Hamiltonian in AO basis
@@ -260,7 +260,7 @@ subroutine qsGW(maxSCF,thresh,max_diis,doACFDT,exchange_kernel,doXBS,COHSEX,BSE,
 
     ! Save quasiparticles energy for next cycle
 
-    Conv = maxval(abs(eGW - eOld))
+    Conv = maxval(abs(error))
     eOld(:) = eGW(:)
 
     !------------------------------------------------------------------------

src/HF/RHF.f90

@@ -122,15 +122,11 @@ subroutine RHF(maxSCF,thresh,max_diis,guess_type,level_shift,nNuc,ZNuc,rNuc,ENuc
 
 !   DIIS extrapolation
 
-    n_diis = min(n_diis+1,max_diis)
+    if(max_diis > 1) then
-    if(abs(rcond) > 1d-15) then
 
+      n_diis = min(n_diis+1,max_diis)
       call DIIS_extrapolation(rcond,nBasSq,nBasSq,n_diis,error_diis,F_diis,error,F)
 
-    else
-
-      n_diis = 0
-
     end if
 
 !   Level-shifting

src/HF/UHF.f90

@@ -170,19 +170,14 @@ subroutine UHF(maxSCF,thresh,max_diis,guess_type,mix,level_shift,nNuc,ZNuc,rNuc,
     
 !   DIIS extrapolation
 
-    n_diis = min(n_diis+1,max_diis)
+    if(max_diis > 1) then
-
-    if(minval(rcond(:)) > 1d-15) then
 
+      n_diis = min(n_diis+1,max_diis)
       do ispin=1,nspin
         if(nO(ispin) > 1) call DIIS_extrapolation(rcond(ispin),nBasSq,nBasSq,n_diis,err_diis(:,1:n_diis,ispin), &
                                                 F_diis(:,1:n_diis,ispin),err(:,:,ispin),F(:,:,ispin))
       end do
 
-    else
-
-      n_diis = 0
-
     end if
 
 !   Level-shifting

src/utils/DIIS_extrapolation.f90

@@ -45,10 +45,17 @@ subroutine DIIS_extrapolation(rcond,n_err,n_e,n_diis,error,e,error_in,e_inout)
 
 ! Solve linear system
 
-  call linear_solve(n_diis+1,A,b,w,rcond)
+  call easy_linear_solve(n_diis+1,A,b,w)
+! call linear_solve(n_diis+1,A,b,w,rcond)
+
+  call matout(n_diis+1,1,w)
 
 ! Extrapolate
 
-    e_inout(:) = matmul(w(1:n_diis),transpose(e(:,1:n_diis)))
+! print*,e_inout
+
+  e_inout(:) = matmul(w(1:n_diis),transpose(e(:,1:n_diis)))
+
+! print*,e_inout
 
 end subroutine DIIS_extrapolation

src/utils/wrap_lapack.f90

@@ -153,7 +153,7 @@ subroutine linear_solve(N,A,b,x,rcond)
   implicit none
 
   integer,intent(in)             :: N
-  double precision,intent(out)    :: A(N,N),b(N),rcond
+  double precision,intent(out)   :: A(N,N),b(N),rcond
   double precision,intent(out)   :: x(N)
 
   integer                        :: info,lwork