# DIIS method to accelerate SCF convergence for stretched geometries

I am implementing from scratch an Hartree-Fock calculation in the STO-3G basis set to perform Born-Oppenheimer molecular dynamics. I have a Restricted Hartree-Fock procedure that can reproduce very well the total energies of Ref. [1] and of Gaussian09 with pre-optimized molecular geometries.

Now, during a Born-Oppenhaimer MD simulation, I will have geometries that can be quite far from equilibrium. I saw that my SCF cycle does not converge for stretched geometries: the SCF cycle for the $CO$ molecule with a bond length of 2.132 bohr converges, while for a bond length of 3.132 bohr I reach the maximal number of allowed iterations (500).

In Table 10.4 of Ref. [2] it is clear that an implementation of the DIIS (Direct Inversion of the iterative Subspace) method speedup the SCF convergence. In addition to this speedup, it seems that the convergence for a stretched molecule can be easily reached with this method. For this reason I implemented the DIIS algorithm following Ref. [3].

I believe I implemented the DIIS method quite correctly, since I was able to obtain a speedup of SCF convergence for pre-optimized molecular geometries. In the following table you can find the number of simple SCF and SCF-DIIS steps I have, for different molecules:

Molecule    SCF     DIIS
H2O         16      9
CO          53      19
HeH+        8       7
CH4         11      8
FH          10      7
O2          41      18
N2          110     22


In every molecule I tested I get a speedup, which is quite impressive where a lot of normal SCF cycles were needed. This is a good result, but unfortunately the convergence of stretched molecules is not improved. By doubling the distances (in bohr) of the pre-optimized structure, the DIIS algorithm (as well as the normal SCF) fails to converge.

There is some step to add to the original DIIS algorithm of Ref. [3] that I am not considering? How I can make the SCF cycle convergent for stretched molecules, provided that these calculations converge in Gaussian09?

EDIT

Gaussian09 brakes down as my program does with the option SCF=NoDIIS. With the option SCF=DIIS the SCF converges even for very distorted molecules.

[1] A. Szabo and N. Ostlund, Modern Quantum Chemistry, Dover, 1996.

[2] T. Helgaker, P. Jørgensen and J. Olsen, Molecular Electronic-Structure Theory, Wiley, 2000.

[3] P. Pulay, Improved SCF Convergence Acceleration, Journal of Computational Chemistry, 1982.

• Are you sure the RHF is doing any good for you in this situation? This neither describes bond breaking nor intermolecular interactions very well. Also, as a rule of thumb if you hit SCF limits like 500-1000 cycle without convergence, the problem/solution is generally not along the line of "let it run for a little more". In your situation the wave-function itself is screwed; you can get a little bit more juice out of it with tricks, but the lemon is essentially rotten. – Greg Dec 12 '15 at 17:25
• @Greg No, I am not. I also implemented the UHF method, but I didn't find how to properly differentiate the initial guess for the Fock matrices for an even number of electrons (see this disscussion). However, Gaussian09 can converge also with RHF, this is why I am asking if I am missing something. It is true, however, that in Ref. [2] the method is not specified, therefore it could be UHF. – user18279 Dec 12 '15 at 17:32
• don't worry about the lack of polish; no one else is doing it. Well done to you and thanks for sharing. – user1945827 Dec 13 '15 at 14:48
• You can also easily implement "level shift" (shifting the virtual orbital energies by a constant). It may help in your case, especially if your initial guess is the converged solution of the previous MD step (and your step is not too large). – Anthony Scemama Dec 18 '15 at 1:47
• Hi have you got it to work? If not, and you are still interested I can take a look at it. – Erik Kjellgren Aug 23 '17 at 21:06