mix guess for UHF

input/methods

@@ -1,5 +1,5 @@
 # RHF UHF KS  MOM
-  T   T   F   F
+  F   T   F   F
 # MP2* MP3 MP2-F12
   F   F   F
 # CCD CCSD CCSD(T) 

input/options

@@ -1,5 +1,5 @@
-# HF: maxSCF thresh   DIIS n_diis guess_type ortho_type
+# HF: maxSCF thresh   DIIS n_diis guess_type ortho_type mix_guess
-       128    0.0000001  T    10      1          1
+       128    0.0000001  T    10      1          1      T
 # MP: 
       
 # CC: maxSCF thresh   DIIS n_diis

mol/h2.xyz

@@ -1,4 +1,4 @@
 2
 
 H 0.0 0.0 0.0
-H 0.0 0.0 1.4
+H 0.0 0.0 3.0

src/HF/UHF.f90

@@ -1,4 +1,4 @@
-subroutine UHF(maxSCF,thresh,max_diis,guess_type,nNuc,ZNuc,rNuc,ENuc,nBas,nO,S,T,V,Hc,ERI,dipole_int,X,EUHF,e,c,P)
+subroutine UHF(maxSCF,thresh,max_diis,guess_type,mix,nNuc,ZNuc,rNuc,ENuc,nBas,nO,S,T,V,Hc,ERI,dipole_int,X,EUHF,e,c,P)
 
 ! Perform unrestricted Hartree-Fock calculation
 
@@ -10,6 +10,7 @@ subroutine UHF(maxSCF,thresh,max_diis,guess_type,nNuc,ZNuc,rNuc,ENuc,nBas,nO,S,T
   integer,intent(in)            :: maxSCF
   integer,intent(in)            :: max_diis
   integer,intent(in)            :: guess_type
+  logical,intent(in)            :: mix 
   double precision,intent(in)   :: thresh
   integer,intent(in)            :: nBas
 
@@ -138,6 +139,10 @@ subroutine UHF(maxSCF,thresh,max_diis,guess_type,nNuc,ZNuc,rNuc,ENuc,nBas,nO,S,T
       c(:,:,ispin) = matmul(X(:,:),cp(:,:,ispin))
     end do
 
+!   Mix guess for UHF solution in singlet states
+
+    if(mix .and. nSCF == 1) call mix_guess(nBas,nO,c)
+
 !   Compute density matrix 
 
     do ispin=1,nspin

src/HF/mix_guess.f90

@@ -0,0 +1,48 @@
+subroutine mix_guess(nBas,nO,c)
+
+!  Guess mixing for UHF calculations on singlet states
+
+  implicit none
+  include 'parameters.h'
+
+! Input variables
+
+  integer,intent(in)            :: nBas
+  integer,intent(in)            :: nO(nspin)
+
+! Local variables
+
+  double precision,parameter   :: th = 0.1d0
+
+  double precision,allocatable :: HOMO(:), LUMO(:)
+  double precision,allocatable :: HOMOa(:),LUMOa(:)
+  double precision,allocatable :: HOMOb(:),LUMOb(:)
+
+! Output variables
+
+  double precision,intent(inout):: c(nBas,nBas,nspin)
+
+! Memory allocation
+
+  allocate(HOMO(nBas), LUMO(nBas),  &
+           HOMOa(nBas),LUMOa(nBas), &
+           HOMOb(nBas),LUMOb(nBas))
+
+! Perform HOMO and LUMO rotation for guess mixing
+
+  HOMO(:) = c(:,nO(1),1)
+  LUMO(:) = c(:,nO(1)+1,1)
+
+  HOMOa(:) =  cos(th)*HOMO(:) + sin(th)*LUMO(:)
+  HOMOb(:) =  cos(th)*HOMO(:) - sin(th)*LUMO(:)
+
+  LUMOa(:) = -sin(th)*HOMO(:) + cos(th)*LUMO(:)
+  LUMOb(:) =  sin(th)*HOMO(:) + cos(th)*LUMO(:)
+
+  c(:,nO(1),1)   = HOMOa(:)
+  c(:,nO(1)+1,1) = LUMOa(:)
+
+  c(:,nO(2),2)   = HOMOb(:)
+  c(:,nO(2)+1,2) = LUMOb(:)
+
+end subroutine mix_guess

src/QuAcK/QuAcK.f90

@@ -106,7 +106,7 @@ program QuAcK
 
   integer                       :: maxSCF_HF,n_diis_HF
   double precision              :: thresh_HF
-  logical                       :: DIIS_HF,guess_type,ortho_type
+  logical                       :: DIIS_HF,guess_type,ortho_type,mix
 
   integer                       :: maxSCF_CC,n_diis_CC
   double precision              :: thresh_CC
@@ -167,7 +167,7 @@ program QuAcK
 
 ! Read options for methods
 
-  call read_options(maxSCF_HF,thresh_HF,DIIS_HF,n_diis_HF,guess_type,ortho_type, &
+  call read_options(maxSCF_HF,thresh_HF,DIIS_HF,n_diis_HF,guess_type,ortho_type,mix, &
                     maxSCF_CC,thresh_CC,DIIS_CC,n_diis_CC,                           &
                     TDA,singlet,triplet,spin_conserved,spin_flip,                    &
                     maxSCF_GF,thresh_GF,DIIS_GF,n_diis_GF,linGF,eta_GF,renormGF,     &
@@ -298,7 +298,7 @@ program QuAcK
     unrestricted = .true.
 
     call cpu_time(start_HF)
-    call UHF(maxSCF_HF,thresh_HF,n_diis_HF,guess_type,nNuc,ZNuc,rNuc,ENuc, & 
+    call UHF(maxSCF_HF,thresh_HF,n_diis_HF,guess_type,mix,nNuc,ZNuc,rNuc,ENuc, & 
              nBas,nO,S,T,V,Hc,ERI_AO,dipole_int_AO,X,EUHF,eHF,cHF,PHF)
     call cpu_time(end_HF)
 

src/QuAcK/read_options.f90

@@ -1,4 +1,4 @@
-subroutine read_options(maxSCF_HF,thresh_HF,DIIS_HF,n_diis_HF,guess_type,ortho_type, &
+subroutine read_options(maxSCF_HF,thresh_HF,DIIS_HF,n_diis_HF,guess_type,ortho_type,mix, &
                         maxSCF_CC,thresh_CC,DIIS_CC,n_diis_CC,                           &
                         TDA,singlet,triplet,spin_conserved,spin_flip,                    &
                         maxSCF_GF,thresh_GF,DIIS_GF,n_diis_GF,linGF,eta_GF,renormGF,     &
@@ -20,6 +20,7 @@ subroutine read_options(maxSCF_HF,thresh_HF,DIIS_HF,n_diis_HF,guess_type,ortho_t
   integer,intent(out)           :: n_diis_HF
   integer,intent(out)           :: guess_type
   integer,intent(out)           :: ortho_type
+  logical,intent(out)           :: mix
 
   integer,intent(out)           :: maxSCF_CC
   double precision,intent(out)  :: thresh_CC
@@ -85,11 +86,13 @@ subroutine read_options(maxSCF_HF,thresh_HF,DIIS_HF,n_diis_HF,guess_type,ortho_t
   n_diis_HF    = 5
   guess_type   = 1
   ortho_type   = 1
+  mix          = .false.
 
   read(1,*) 
-  read(1,*) maxSCF_HF,thresh_HF,answer1,n_diis_HF,guess_type,ortho_type
+  read(1,*) maxSCF_HF,thresh_HF,answer1,n_diis_HF,guess_type,ortho_type,answer2
 
   if(answer1 == 'T') DIIS_HF = .true.
+  if(answer2 == 'T') mix     = .true.
 
   if(.not.DIIS_HF) n_diis_HF = 1