From fc4a0623262bbba559477c4d6633d4c49bfe6b0d Mon Sep 17 00:00:00 2001 From: Pierre-Francois Loos Date: Tue, 27 Oct 2020 17:23:48 +0100 Subject: [PATCH] mix guess for UHF --- input/methods | 2 +- input/options | 4 ++-- mol/h2.xyz | 2 +- src/HF/UHF.f90 | 7 +++++- src/HF/mix_guess.f90 | 48 ++++++++++++++++++++++++++++++++++++++ src/QuAcK/QuAcK.f90 | 20 ++++++++-------- src/QuAcK/read_options.f90 | 23 ++++++++++-------- 7 files changed, 81 insertions(+), 25 deletions(-) create mode 100644 src/HF/mix_guess.f90 diff --git a/input/methods b/input/methods index 617d1fb..86c7c20 100644 --- a/input/methods +++ b/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) diff --git a/input/options b/input/options index 8e22114..4c47176 100644 --- a/input/options +++ b/input/options @@ -1,5 +1,5 @@ -# HF: maxSCF thresh DIIS n_diis guess_type ortho_type - 128 0.0000001 T 10 1 1 +# HF: maxSCF thresh DIIS n_diis guess_type ortho_type mix_guess + 128 0.0000001 T 10 1 1 T # MP: # CC: maxSCF thresh DIIS n_diis diff --git a/mol/h2.xyz b/mol/h2.xyz index a302682..dcaf1a5 100644 --- a/mol/h2.xyz +++ b/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 diff --git a/src/HF/UHF.f90 b/src/HF/UHF.f90 index 558bafd..4123dd5 100644 --- a/src/HF/UHF.f90 +++ b/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 diff --git a/src/HF/mix_guess.f90 b/src/HF/mix_guess.f90 new file mode 100644 index 0000000..3cca61a --- /dev/null +++ b/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 diff --git a/src/QuAcK/QuAcK.f90 b/src/QuAcK/QuAcK.f90 index 924f641..4ce5392 100644 --- a/src/QuAcK/QuAcK.f90 +++ b/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,14 +167,14 @@ program QuAcK ! Read options for methods - call read_options(maxSCF_HF,thresh_HF,DIIS_HF,n_diis_HF,guess_type,ortho_type, & - 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, & - maxSCF_GW,thresh_GW,DIIS_GW,n_diis_GW,linGW,eta_GW, & - COHSEX,SOSEX,TDA_W,G0W,GW0, & - doACFDT,exchange_kernel,doXBS, & - BSE,dBSE,dTDA,evDyn, & + 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, & + maxSCF_GW,thresh_GW,DIIS_GW,n_diis_GW,linGW,eta_GW, & + COHSEX,SOSEX,TDA_W,G0W,GW0, & + doACFDT,exchange_kernel,doXBS, & + BSE,dBSE,dTDA,evDyn, & nMC,nEq,nWalk,dt,nPrint,iSeed,doDrift) ! Weird stuff @@ -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) diff --git a/src/QuAcK/read_options.f90 b/src/QuAcK/read_options.f90 index 1abed8a..1492197 100644 --- a/src/QuAcK/read_options.f90 +++ b/src/QuAcK/read_options.f90 @@ -1,11 +1,11 @@ -subroutine read_options(maxSCF_HF,thresh_HF,DIIS_HF,n_diis_HF,guess_type,ortho_type, & - 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, & - maxSCF_GW,thresh_GW,DIIS_GW,n_diis_GW,linGW,eta_GW, & - COHSEX,SOSEX,TDA_W,G0W,GW0, & - doACFDT,exchange_kernel,doXBS, & - BSE,dBSE,dTDA,evDyn, & +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, & + maxSCF_GW,thresh_GW,DIIS_GW,n_diis_GW,linGW,eta_GW, & + COHSEX,SOSEX,TDA_W,G0W,GW0, & + doACFDT,exchange_kernel,doXBS, & + BSE,dBSE,dTDA,evDyn, & nMC,nEq,nWalk,dt,nPrint,iSeed,doDrift) ! Read desired methods @@ -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(answer1 == 'T') DIIS_HF = .true. + if(answer2 == 'T') mix = .true. if(.not.DIIS_HF) n_diis_HF = 1