A Least-Squares Commutator in the Iterative Subspace Method for

Oct 6, 2016 - subspace (LCIIS) approach is explored for accelerating self- consistent field ... inversion of the iterative subspace (DIIS) methods in ...
0 downloads 0 Views 2MB Size
Subscriber access provided by CORNELL UNIVERSITY LIBRARY

Article

A least-squares commutator in the iterative subspace (LCIIS) method for accelerating self-consistent field convergence Haichen Li, and David J. Yaron J. Chem. Theory Comput., Just Accepted Manuscript • DOI: 10.1021/acs.jctc.6b00666 • Publication Date (Web): 06 Oct 2016 Downloaded from http://pubs.acs.org on October 18, 2016

Just Accepted “Just Accepted” manuscripts have been peer-reviewed and accepted for publication. They are posted online prior to technical editing, formatting for publication and author proofing. The American Chemical Society provides “Just Accepted” as a free service to the research community to expedite the dissemination of scientific material as soon as possible after acceptance. “Just Accepted” manuscripts appear in full in PDF format accompanied by an HTML abstract. “Just Accepted” manuscripts have been fully peer reviewed, but should not be considered the official version of record. They are accessible to all readers and citable by the Digital Object Identifier (DOI®). “Just Accepted” is an optional service offered to authors. Therefore, the “Just Accepted” Web site may not include all articles that will be published in the journal. After a manuscript is technically edited and formatted, it will be removed from the “Just Accepted” Web site and published as an ASAP article. Note that technical editing may introduce minor changes to the manuscript text and/or graphics which could affect content, and all legal disclaimers and ethical guidelines that apply to the journal pertain. ACS cannot be held responsible for errors or consequences arising from the use of information contained in these “Just Accepted” manuscripts.

Journal of Chemical Theory and Computation is published by the American Chemical Society. 1155 Sixteenth Street N.W., Washington, DC 20036 Published by American Chemical Society. Copyright © American Chemical Society. However, no copyright claim is made to original U.S. Government works, or works produced by employees of any Commonwealth realm Crown government in the course of their duties.

Page 1 of 52

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Journal of Chemical Theory and Computation

A least-squares commutator in the iterative subspace (LCIIS) method for accelerating self-consistent field convergence Haichen Li and David J. Yaron∗ Department of Chemistry, Carnegie Mellon University, Pittsburgh, PA 15213 E-mail: [email protected]

Abstract A least-squares commutator in the iterative subspace (LCIIS) approach is explored for accelerating self-consistent field (SCF) calculations. LCIIS is similar to direct inversion of the iterative subspace (DIIS) methods in that the next iterate of the density matrix is obtained as a linear combination of past iterates. However, whereas DIIS methods find the linear combination by minimizing a sum of error vectors, LCIIS minimizes the Frobenius norm of the commutator between the density matrix and the Fock matrix. This minimization leads to a quartic problem that can be solved iteratively through a constrained Newton’s method. The relationship between LCIIS and DIIS is discussed. Numerical experiments suggest that LCIIS leads to faster convergence than other SCF convergence accelerating methods in a statistically significant sense, and in a number of cases LCIIS leads to stable SCF solutions that are not found by other methods. The computational cost involved in solving the quartic minimization problem is small compared to the typical cost of SCF iterations and the approach is easily integrated into existing codes. LCIIS can therefore serve as a powerful addition to SCF convergence accelerating method in computational quantum chemistry packages.

1

ACS Paragon Plus Environment

Journal of Chemical Theory and Computation

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

1

Introduction

Self-consistent field (SCF) calculations are pervasive in modern quantum chemistry studies. Both Hartree-Fock (HF) theory and Kohn-Sham density functional theory (KS-DFT) require SCF convergence, and correlated methods typically start with an SCF solution. KS-DFT is widely acknowledged as the preferred quantum chemical tool in a variety of fields. 1,2 As a result, SCF solutions often consume a considerable portion of the compute time in quantum chemical studies. A na¨ıve implementation of SCF, sometimes referred to as a Roothaan-Hall process, 3 is based on repeated diagonalization of the Fock matrix. (In this work, we use Fock matrix to refer to both the Fock matrix in HF theory and the Kohn-Sham matrix in KS-DFT theory.) This iterative process does not have convergence guarantees and often leads to oscillatory or diverging behavior. 4 This work advances the long-standing goal of accelerating the Roothaan-Hall process by reducing the number of iterations without noticeably increasing the per-iteration computational cost. A number of techniques have been explored to improve the robustness and convergence rate of the Roothaan-Hall process. These include the direct inversion of the iterative subspace (DIIS) 5 family of methods, 6–11 the linear-expansion shooting techniques (LIST) 12,13 which also fall into the framework of DIIS methods, 14 the optimal-damping algorithm, 15 the levelshifting method, 16–18 and the fractional occupation method. 19 Within the DIIS family, early work identified the commutator-DIIS (CDIIS, sometimes abbreviated as just DIIS) method as particularly effective 6 and CDIIS remains among the most commonly used approaches. The energy-DIIS (EDIIS) 9 and the augmented Roothaan-Hall DIIS (ADIIS) 10 methods lead to more stable convergence and, in difficult cases, may converge to a lower final energy. This enhanced stability arises primarily in the early SCF iterations, during which EDIIS or ADIIS can lower the total energy more rapidly than CDIIS. But as the density matrix approaches convergence, CDIIS becomes significantly faster. 9,10,20 This behavior motivates methods such as EDIIS+CDIIS or ADIIS+CDIIS, that switch from EDIIS/ADIIS to CDIIS as the solution approaches convergence. DIIS methods generate a density matrix for the next iteration by 2

ACS Paragon Plus Environment

Page 2 of 52

Page 3 of 52

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Journal of Chemical Theory and Computation

taking a linear combination of the density matrices resulting from previous iterations. The coefficients for this linear combination are obtained by minimizing an error vector, where the vector has the property of going to zero at SCF convergence. 6 The initial letter C in CDIIS refers to the use of the commutator between the Fock matrix and the density matrix as the error vector. The LIST family of methods can also be viewed as a DIIS procedure with an alternative error vector. 14 This family includes LIST-direct (LISTd), LIST-indirect (LISTi), and LIST-better (LISTb), with LISTb leading to the smoothest and fastest SCF convergence. 13 The above methods are sometimes referred to as first-order SCF convergence accelerating methods because each is a modification to the Roothaan-Hall process. In the Roothaan-Hall process, diagonalizing the Fock matrix is equivalent to direct minimization of the first-order Taylor expansion for the total HF/KS-DFT energy, with the expansion being in terms of the density matrix. 3 Higher-order SCF convergence techniques have also been developed. Bacskay’s quadratically convergent SCF procedure (QC-SCF) 21,22 is essentially an adaptation of Newton’s method to the HF/KS-DFT problem. QC-SCF has a much higher chance of reaching SCF convergence; however, due to the implicit evaluation of the Hessian matrix, the computational cost is substantially greater than first-order SCF methods. As a result, it is typically invoked only when first-order SCF methods fail to converge. Methods that go beyond first order but maintain low computational cost have also been developed. A quasi-Newton method, denoted as second-order SCF (SOSCF), has been found to be competitive to DIIS for spin-restricted closed shell and spin-restricted open shell calculations. 23,24 However, for spin-unrestricted calculations on transition metal complexes, SOSCF requires expert knowledge of SCF convergence as the startup condition has to be adjusted carefully. 24 Trust region methods 3,25,26 enjoy strong convergence properties. However, they require full energy evaluation at every iteration with high numerical accuracy and are thus much more computationally expensive than first-order methods under certain scenarios. 27 Moreover, due to fundamental differences from the Roothaan-Hall process, the

3

ACS Paragon Plus Environment

Journal of Chemical Theory and Computation

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

corresponding SCF code must be substantially rewritten. As a result, many existing SCF convergence aiding techniques such as dynamic damping and level-shifting would need to be redesigned to couple to SOSCF or trust region methods. Advantages of DIIS methods include easy generalization to spin-unrestricted HF/KSDFT, easy addition as a module to existing Roothaan-Hall implementations and, in most cases, easy integration with other techniques. Moreover, the flexibility of the DIIS formalism extends its use to other methods such as coupled cluster theory, 28 multireference wavefunctions, 29 and geometry optimization. 30–32 Generally speaking, DIIS methods are still considered as the most well-rounded SCF convergence methods, and they are the default convergence method in many modern quantum chemistry packages. For example, Gaussian 09 33 uses EDIIS+CDIIS by default, while Q-Chem, 34 NWChem 35 and Psi4 36 use CDIIS as their default. GAMESS 37,38 and ORCA 39 use SOSCF as the default SCF convergence method, with a recommendation to switch to DIIS methods or QC-SCF on problematic SCF convergence cases, as may arise with transition metal complexes. In this work, we develop and test a new first-order SCF convergence accelerating method which falls partly under the framework of DIIS methods. We follow the DIIS approach of writing a prediction of the density matrix as a linear combination of density matrices obtained from past iterates. We also follow the CDIIS approach of determining the coefficients for the linear combination from the requirement that the Fock matrix should commute with the density matrix at convergence. However, instead of minimizing a sum of error vectors, we find the linear combination of density matrices that, when used to construct the Fock matrix, minimizes this commutator. We denote this scheme as least-squares commutator in the iterative subspace (LCIIS). Similar to CDIIS, LCIIS requires only that the coefficients for the linear combination of density matrices sum to 1. This is in contrast to those interpolation schemes, such as EDIIS and ADIIS, which also require each linear expansion coefficient to be greater than 0. In Section 2, we discuss the LCIIS formalism and contrast it with CDIIS. Section 3

4

ACS Paragon Plus Environment

Page 4 of 52

Page 5 of 52

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Journal of Chemical Theory and Computation

examines the additional computational cost introduced by LCIIS. Section 4 discusses the computational details. In Section 5 and Section 6, we present numerical experiments to explore the performance of LCIIS in comparison with other DIIS methods. To explore stability, we consider cases that are particularly difficult to converge and have been used as test cases in previous studies. 9,10,12–14,20,40 We also evaluate the convergence properties of various first-order SCF methods by doing statistical comparisons on a data set of 1,035 transition metal complexes.

2

Description of LCIIS

In this section, we state the objective function to be minimized in LCIIS and discuss the approach used to search for that minimum. LCIIS is similar to DIIS methods in that, at ˜ n is given as a linear the nth SCF iteration, a prediction for the converged density matrix D P i ˜n = combination of M retained density matrices {D n−M +1 , . . . , D n } by D i ci D . The coefficients for this linear combination are also used to construct a prediction for the converged P Fock matrix F˜ n from M retained Fock matrices {F n−M +1 , . . . , F n } by F˜ n = i ci F i , where F i is the Fock matrix built from density matrix D i as part of the previous SCF iterations. In HF theory, the Fock matrix is linear in the density matrix and so F˜ n is equivalent to the ˜ n . In KS-DFT, the assumption of linearity is also Fock matrix that would be built from D a reasonable approximation. 20,41 The linear expansion coefficients cn,LCIIS are obtained by ˜ n ] to vanish, which is equivalent to the first-order Karushrequiring the commutator [F˜ n , D Kuhn-Tucker (KKT) conditions for the HF/KS-DFT problem. This commutator can be written as a linear combination of commutators between retained Fock matrices and density matrices. ˜ n] = [F˜ n , D

X

ci cj [F i , D j ].

ij

5

ACS Paragon Plus Environment

(1)

Journal of Chemical Theory and Computation

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 6 of 52

In LCIIS, the linear coefficients cn,LCIIS for SCF iteration n are obtained by minimizing the ˜ n ], subject to the trace constraint 3 that the linear coefficients sum Frobenius norm of [F˜ n , D to 1. A prediction of the converged Fock matrix F˜ n,LCIIS is then constructed from the linear coefficient vector and retained Fock matrices according to Equation (2) below. Note that in ˜ n is not idempotent and thus is sometimes referred to as a pseudodensity matrix. 9 general D F˜ n,LCIIS is then diagonalized to give the next idempotent density matrix D n+1 as in a regular Roothaan-Hall process. F˜ n,LCIIS =

X

cn,LCIIS F i, i

(2)

i

cn,LCIIS = arg min f LCIIS (c) s. t. c

X

ci = 1,

(3)

i

2



˜ n ]T · [F˜ n , D ˜ n] . ˜ n ] f LCIIS (c) = [F˜ n , D

= tr [F˜ n , D F

(4)

Introducing the four dimensional tensor,  T Tijkl = tr [F i , D j ] · [F k , D l ] ,

(5)

allows the LCIIS target function f LCIIS (c) to be written as a homogeneous quartic form, f LCIIS (c) =

X

ci cj ck cl Tijkl .

(6)

ijkl

Note that the size of the tensor depends on the number of retained Fock matrices M which is practically bounded by a small constant integer. The gradient g and Hessian H of f LCIIS (c) may be evaluated as

gi (c) =

X  ∂ LCIIS f (c) = 2 cj ck cl Tijkl + Tjikl , ∂ci jkl

Hij (c) =

X  ∂2 f LCIIS (c) = 2 ck cl Tijkl + Tikjl + Tiklj + Tjikl + Tkijl + Tkilj . ∂ci ∂cj kl

6

ACS Paragon Plus Environment

(7) (8)

Page 7 of 52

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Journal of Chemical Theory and Computation

We find the LCIIS linear combination coefficients cn,LCIIS iteratively using Equation (10) below which is obtained by applying Newton’s method to the Lagrangian defined in Equation (9) , L(c, λ) = f LCIIS (c) − λ(1 −

X

ci ),

(9)

i



(t+1)





−1 

(t)

c  H(c ) 1   = λ(t+1) 1T 0

(t)

(t)



g(c ) + λ   , 0

(10)

where superscript (t) indicates the tth step in the iterative process. The initial coefficients c(0) are taken as the CDIIS coefficients and the initial Lagrangian multiplier λ(0) is set to 0. Although there is no guarantee that the resulting solution is a global minimum, our numerical results in Section 5 and Section 6 show satisfactory performance with the minimum that is located. The numerical stability of LCIIS is comparable to that of CDIIS, as measured by the condition number of the corresponding Hessians (see Section S2 of the Supporting Information). The relation between LCIIS and CDIIS is made apparent by comparison with Pulay’s CDIIS scheme, 6

c

n,CDIIS

2

X X

i i ci [F , D ] s. t. ci = 1, = arg min c

F

i

(11)

i

In LCIIS, the quantity being minimized is that of Equation (3) which may be written as,

X

2

˜ n ˜ n 2

i j ci cj [F , D ]

[F , D ] = F

ij

F

(12)

LCIIS becomes equivalent to CDIIS when the density matrix D j in Equation (12) is replaced

7

ACS Paragon Plus Environment

Journal of Chemical Theory and Computation

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 8 of 52

with D i ,

X

2

X

2

i i i i c c [F , D ] = c [F , D ]



, i j i ij

F

i

F

(13)

where we have used the condition that the coefficients ci sum to 1. CDIIS is based on minimization of an error vector and may be regarded as a projected quasi-Newton’s method. 3,42 LCIIS is based on an alternative approach in which the Frobenius norm of the commutator between the predicted converged density matrix and the Fock matrix is minimized. The fact that optimizing the LCIIS target function leads to a quartic problem guarantees that LCIIS is not a form of DIIS with a specific choice of error vector, in contrast to the LIST methods. 14 LCIIS and CDIIS become equivalent when the changes in density matrix between iterations are small, such that D j ≈ D i for any j. This provides an alternative explanation for the good performance of CDIIS as the SCF iterations approach convergence. While both LCIIS and CDIIS are based on the condition that the Fock and density matrices commute at convergence, explicit minimization of this commutator in LCIIS likely generates pseudodensity matrices that best commute with the corresponding Fock matrix at each iteration. It is therefore plausible that LCIIS will lead to enhanced performance. In the results presented below, the density and Fock matrices are expressed in an orthogonal basis. Such a basis may be readily obtained from a non-orthogonal atomic orbital basis through L¨owdin orthogonalization. 43 LCIIS may also be implemented in a non-orthogonal basis by using an appropriate form for the commutator (see Section S5 of the Supporting Information).

8

ACS Paragon Plus Environment

Page 9 of 52

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Journal of Chemical Theory and Computation

3

Computational cost of LCIIS

3.1

Time complexity analysis

The computational cost of LCIIS is dominated by the construction of the tensor T defined in Equation (5). Suppose there are N basis functions in the system, and we retain at most M density/Fock matrix pairs. With proper programming, SCF iterations beyond iteration M require evaluation of 2M − 1 new commutators and M 4 − (M − 1)4 new tensor elements. In the scenario of classical SCF implementations where the density matrix, the Fock matrix, and the two-electron integrals are treated as dense arrays, the time complexity of evaluating the commutators is O(M N 3 ), and the time complexity of evaluating the tensor elements from these commutators is O(M 3 N 2 ). In practice, M is a small constant integer with a value of about 20 such that N is much greater than M 2 for chemical systems of modern interest. The computational cost of LCIIS is therefore dominated by the 2M − 1 commutator evaluations, with time complexity of O(M N 3 ). In comparison, the CDIIS algorithm requires evaluation of only one additional commutator for each SCF iteration, such that the overall time complexity of CDIIS’s cost-determining step is O(N 3 ). On the other hand, in this scenario, the time complexity of the construction of the Fock matrix is O(N 4 ) when the exact Hartree-Fock exchange is required, 44 or O(N 3 ) but with a large prefactor when a pure, density-fitted KSDFT calculation is performed; 45 either of which is much more computationally expensive than the cost of an LCIIS step. In the scenario of SCF implementations for large systems, the computational cost may be greatly reduced by taking advantage of the sparsity of both the Fock/density matrices and the two-electron integrals. 34 It has been shown that the time complexity of the construction of the Fock matrix can be routinely reduced to O(N 2 ) ∼ O(N 3 ), 34,46 and that the Fock matrix and density matrix are sparse unless the chemical system is metallic. 47 (Asymptotic scaling of O(N ) may also be achieved on very large systems, but with substantial prefactors. 48 ) Note that because the number of significant elements in these sparse matrices grows

9

ACS Paragon Plus Environment

Journal of Chemical Theory and Computation

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

linearly with system size asymptotically, the time complexity of matrix manipulations in this scenario is O(N ). 49 As a result, the time complexity of evaluating the commutator becomes O(M N ), and the time complexity of evaluating the tensor elements of T becomes O(M 3 N ). Construction of the Fock matrix, which scales as O(N 2 ) ∼ O(N 3 ), will dominate the cost for large N . For intermediate N , reducing the value of M provides a means to keep the cost of evaluating T far below that of constructing the Fock matrix. We also note that taking advantage of the sparse character of the Hamiltonian may require use of a non-orthogonal atomic orbital basis. Section S5 of the Supporting Information shows that LCIIS behaves similarly with and without transformation to an orthogonal basis. Given the time complexity analysis in the above two scenarios, we conclude that the computational cost of LCIIS is negligible compared with the construction of the Fock matrix, provided the corresponding matrix manipulation code is implemented optimally.

3.2

Relative CPU time results

To support the above time complexity analysis, Table 1 and Table 2 provide relative CPU times of LCIIS accelerated SCF compared against Gaussian 09 SCF with default settings, for cases that vary both the density functional and the size of the basis set. The timing results are obtained using our implementation of LCIIS in a developmental version of Gaussian 09, based on Gaussian 09 D.01. From the above analysis, we expect the per iteration CPU time of LCIIS accelerated SCF to be not significantly higher than the per iteration CPU time of Gaussian 09 default SCF. CPU times were obtained by performing the same calculation 10 times and then taking the minimum CPU time, as reported at the end of the Gaussian 09 output file. For calculations with less than 400 basis functions, Gaussian 09 precomputes and stores two-electron integrals in the main memory (which greatly reduces the computational cost). In such cases, calculations were run on a single core. Above 400 functions, direct SCF is used where two-electron integrals are computed on-the-fly. In such cases, calculations were run using 10 cores. Because parallel efficiency is almost 100% for direct SCF, the resulting 10

ACS Paragon Plus Environment

Page 10 of 52

Page 11 of 52

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Journal of Chemical Theory and Computation

estimates of CPU time are reliable. Table 1 shows cases where SCF convergence can be readily achieved. Because the number of retained density/Fock matrix pairs is always less than the number of SCF iterations, the change in CPU time is less than 1.5%, which is comparable to the accuracy of the timing estimates. These results demonstrate that in such cases LCIIS introduces a negligible additional computational cost. Table 1: Relative CPU time results on CH3 CHO, which represents cases where SCF convergence is readily achieved and LCIIS does not lead to a reduced number of iterations. NBF is Number of Basis Functions. All results are based on minimum CPU time over 10 runs. (†Using 10 cores.) DFT b3lyp

m06

Basis set aug-cc-pvdz aug-cc-pvtz †aug-cc-pvqz aug-cc-pvdz aug-cc-pvtz †aug-cc-pvqz

NBF 105 230 424 105 230 424

Number of iterations LCIIS G09 default 11 11 10 10 10 10 11 11 12 12 12 12

Change in CPU time 1.4% 0.7% 0.2% 0.0% 0.7% 0.3%

Table 2: Relative CPU time results on Ru4 CO which represents cases where LCIIS leads to a saving in number of iterations. †See Table 1 for details. DFT b3lyp

m06

Basis set def2svp def2tzvp †def2qzvp def2svp def2tzvp †def2qzvp

NBF 152 222 402 152 222 402

Number of iterations LCIIS G09 default Change 43 54 -20.4% 44 51 -17.0% 51 53 -3.8% 43 49 -12.2% 45 47 -4.3% 44 69 -36.2%

Change in CPU time -15.0% -7.7% -4.5% -9.0% -0.9% -34.9%

Results shown in Table 2 are for a case where LCIIS converges in fewer iterations than Gaussian 09 default, and both methods require a considerable number of iterations to reach convergence. For smaller bases, where two-electron integrals are precomputed and stored, the results find that LCIIS increases per iteration CPU time by ≈ 5%. However, the saving in number of iterations outweighs this small additional computational cost. For the larger basis 11

ACS Paragon Plus Environment

Journal of Chemical Theory and Computation

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

set, where direct SCF is used, the change in number of iterations and the change in CPU time become very close. Therefore, for large basis sets with direct SCF, the computational cost of LCIIS is totally negligible. In summary, these timing results demonstrate that the additional computations of LCIIS, relative to those of CDIIS, lead to a negligible increase in CPU time per SCF iteration. The per-iteration computational cost is instead dominated by evaluation of the Fock matrix and other computations required by all first-order SCF methods. Comparisons of number of iterations across first-order methods therefore accurately reflects the SCF computation time.

4

Computational methods

In Section 5, we compare the performance of six first-order SCF methods: EDIIS+CDIIS, CDIIS, CDIIS with damping and level-shifting (CDIIS-DL), ADIIS+CDIIS, LISTb, and LCIIS. We invoked the first three methods directly from Gaussian 09. EDIIS+CDIIS is the default method used in Gaussian 09’s SCF routine, with dynamic damping performed only at the first iteration, level-shifting turned off, and retention of at most 20 Fock matrices. When CDIIS is requested, Gaussian 09 automatically turns on both dynamic damping and level-shifting, therefore we consider the default implementation of CDIIS in Gaussian 09 as CDIIS-DL. On the other hand, (plain) CDIIS in this work refers to requesting CDIIS in Gaussian 09 with dynamic damping only at the first iteration and with level-shifting turned off. We have implemented the LCIIS method in a developmental version of Gaussian 09 based on Gaussian 09 D.01. Our LCIIS implementation also performs dynamic damping only at the first iteration, turns off level-shifting, and retains at most 20 Fock matrices following the standard setting of DIIS methods in Gaussian 09. We have implemented ADIIS+CDIIS and LISTb in MATLAB R2015a (with Gaussian 09 being used to generate the initial guess molecular orbitals and to build Fock matrices from density matrices) using the settings reported in the original publications. Specifically, our LISTb implementation retains at

12

ACS Paragon Plus Environment

Page 12 of 52

Page 13 of 52

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Journal of Chemical Theory and Computation

most 5 Fock matrices as recommended in reference 13. Our ADIIS+CDIIS implementation retains at most 6 Fock matrices and switches from ADIIS to CDIIS when the absolute value of the energy change in successive steps is less than 10−2 Hartree, as recommended in reference 10. In all calculations, we preserve Gaussian 09 default settings as much as possible. SCF convergence criteria were set to Gaussian 09 default values: 10−6 and 10−8 for maximum and root-mean-squares change in density matrix elements, respectively. The SCF procedure was terminated, and convergence considered as failed, if the iteration number went above 128. Similar results are obtained with the maximum number of iterations increased to 256 (see Section S3 of the Supporting Information). The benchmark systems in Section 5 include all computationally feasible cases previously studied in numerous past works, 9,10,12,13,20 corresponding to 6 molecules which exhibit various levels of difficulty regarding SCF convergence. The geometries of CH3 CHO, [Cd(Im)]2+ , and Ru4 (CO) were taken from reference 10; the geometry of UF4 was taken from reference 20; CrC has a bond length of 2.00 ˚ A; SiH3 H is tetrahedral with three 1.48 ˚ A bonds and one elongated 4.00 ˚ A bond. A commonly-used strategy to verify that a wavefunction is a local minimum with regard to second-order variations in orbital coefficients is the wavefunction stability test. 50 In all cases, we successfully located low energy solutions that passed the corresponding RHF → RHF or UHF → UHF wavefunction stability tests 14 in Gaussian 09. To systematically compare the performance of LCIIS with that of EDIIS+CDIIS and CDIIS, we performed B3LYP/lanl2dz single point calculations on the lowest two spin states of 1,035 transition metal complexes with EDIIS+CDIIS, CDIIS, CDIIS-DL, and LCIIS. We refer to these complexes as the MITM-COD data set, for distinct Molecules and Ions with Transition Metals from the Crystallography Open Database 51 (COD). Each molecule in the MITM-COD data set contains exactly one first or second row (Sc-Zn, Y-Cd) transition metal atom, and has less than 30 atoms in its chemical formula. To avoid biasing the data, we used Mercury 52,53 to extract all molecules/ions satisfying these two conditions from all crystal structures currently available in COD. Occasionally, we added hydrogen atoms manually as

13

ACS Paragon Plus Environment

Journal of Chemical Theory and Computation

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 14 of 52

some crystallographic information files (CIFs) do not record hydrogen atom coordinates.

5

Benchmarks on well studied SCF convergence cases

We first report benchmark results on 8 individual SCF convergence cases which have been extensively studied in numerous past works. 9,10,12–14,20,40 We pay particular attention to the comparison of LCIIS with EDIIS+CDIIS, because past work comparing methods across a variety of cases has identified EDIIS+CDIIS as the best overall first-order SCF method. 14,20 The acetaldehyde molecule in Figure 1(a) represents the majority of closed-shell organic molecules where SCF convergence is readily achieved at either HF or KS-DFT level. All 6 methods lead to the same final SCF energy in almost the same number of iterations. The results validate the conclusion that the performance of ADIIS+CDIIS should be identical to the performance of EDIIS+CDIIS at the HF level. 20 SCF convergence on SiH3 H is moderately difficult, with all methods converging to the same RHF → RHF stable solution in under 35 iterations. The convergence rates do, however, differ between methods with LCIIS leading to the fastest convergence. The results are also in agreement with past work that found EDIIS+CDIIS to perform better than LISTb and CDIIS. 20 Interestingly, it seems that here damping and level-shifting together have a significant negative impact on the performance of CDIIS. In calculations on [Cd(Im)]2+ , all 6 methods converge to the same energy and in each case, the final solution is unstable with regards to the RHF → RHF stability test. This hints at a commonality between all 6 of these first-order SCF methods. The QC-SCF approach finds a solution that passes the RHF → RHF stability test and has an energy of −5667.009041 Hartree, which is lower than the energy of the unstable solution by 0.35 kCal/mol. LCIIS and EDIIS+CDIIS performed similarly in this case, converging in 16 and 15 iterations respectively. Past work 10,13 has found solutions that differ from those reported here, which may be attributed to small differences in the implementation as detailed in Section S1 of the

14

ACS Paragon Plus Environment

15

10

10 log 10 |E i+1 - E i |

15

5 RHF/6-31g* Harris guess E* = -152.914403

0 -5 -10

5 RSVWN5/6-31g* Harris guess E* = -290.457688

0 -5 -10

0

5

10

15

0

5

10

15 20 Iteration number i

Iteration number i

15

15

10

10

5 RB3LYP/3-21g Harris guess E* = -5667.009041

0 -5

30

35

5 RB3LYP/6-31g Harris guess E* = -1082.114808

0 -5 -10

-10 0

5

10 15 Iteration number i

20

25

0

30

10

20

30 40 50 Iteration number i

(c) [Cd(Im)]2+

60

70

80

(d) CrC

15

15

10 #

5

$ $

10

! " ! "

$ $ $ $

log 10 |E i+1 - E i |

log 10 |E i+1 - E i |

25

(b) SiH3 H

log 10 |E i+1 - E i |

log 10 |E i+1 - E i |

(a) CH3 CHO

RB3LYP/lanl2dz Harris guess E* = -488.704767

0 -5 -10

"

5

! !

#

UB3LYP/lanl2dz Triplet guess E* = -488.740761

0 -5 -10

0

20

40

60 80 Iteration number i

100

120

0

(e) Ru4 (CO) singlet restricted 15

15

10

10

5 RB3LYP/lanl2dz Harris guess E* = -451.219602

0 -5

20

40

60 80 Iteration number i

100

120

(f) Ru4 (CO) singlet broken-symmetry

log 10 |E i+1 - E i |

log 10 |E i+1 - E i |

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Journal of Chemical Theory and Computation

log 10 |E i+1 - E i |

Page 15 of 52

-10

5 UB3LYP/lanl2dz Triplet guess E* = -451.259876

0 -5 -10

0

5

10

15 20 25 Iteration number i

30

35

40

(g) UF4 singlet restricted

0

20

40 60 Iteration number i

80

100

120

(h) UF4 singlet broken-symmetry

Figure 1: Convergence curves on well studied SCF convergence cases. E* is the energy where the corresponding wavefunction has passed the proper stability test: RHF → RHF for (a), (b), (c), (d), (e), and (g); UHF → UHF for (f) and (h). All energies are in Hartrees. 15

ACS Paragon Plus Environment

Journal of Chemical Theory and Computation

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Supporting Information. CrC is a case where LCIIS shows its potential strength in locating low energy, stable SCF solutions. CDIIS-DL, ADIIS+CDIIS, and LCIIS were able to find the RHF → RHF stable solution with energy −1082.114808 Hartree. This is 0.50 kCal/mol lower than the solution found by EDIIS+CDIIS, CDIIS, and LISTb. (The energy from our EDIIS+CDIIS calculation agrees exactly with that from previous work. 20 ) Among the three methods which converged to the stable solution, LCIIS outperforms CDIIS-DL by about 20 iterations, and outperforms ADIIS+CDIIS by about 50 iterations. This is an example of a case where LCIIS leads to a lower-energy solution that is not found by EDIIS+CDIIS. At B3LYP/lanl2dz level of theory, Ru4 (CO) and UF4 are challenging SCF convergence cases. When starting from commonly used initial guess schemes such as Harris guess or superposition of atomic densities (SAD) guess, the global minimum solutions for these systems had not been previously achieved through first-order SCF methods. 14 Following the recommendations in reference 20, we performed both restricted calculations starting from Harris guess and broken-symmetry unrestricted calculations starting from converged triplet calculations. We first discuss the four methods that always converge on these systems: EDIIS+CDIIS, CDIIS, CDIIS-DL, and LCIIS. These four methods lead to the same converged solution, as reflected in the energy and stability, in all but one situation. The exception is restricted calculation of Ru4 (CO), where LCIIS is the only method to find a stable solution. The energy of this stable solution of −488.704767 Hartree is 8.39 kCal/mol lower than the unstable solution of −488.691402 Hartree and equal to the energy obtained previously using a triplet state as the initial guess. 20 The ADIIS+CDIIS and LISTb methods converge more slowly than the other methods such that in Figures 1(e) and 1(f) the convergence criteria used here have not yet been satisfied by 128 iterations. In the restricted calculation of Ru4 (CO), the energy for ADIIS+CDIIS at iteration 128 is equivalent to that of EDIIS+CDIIS, CDIIS, and CDIIS-DL, and the energy for LISTb at iteration 128 is in agreement with that of the stable solution. This

16

ACS Paragon Plus Environment

Page 16 of 52

Page 17 of 52

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Journal of Chemical Theory and Computation

supports the motivation of LIST methods, which is to approximate the corrected HohnbergKohn-Sham energy in the iterative subspace, 12,13,54 and so obtain a reasonable estimate of the converged KS-DFT energy. However, the failure to converge in 128 iterations points to the ill-conditioned nature of the LISTb equations. 14,20 The numerical results shown here for ADIIS and LISTb differ somewhat from the published values 10,12,13 due to differences in settings regarding the basis set and integral grid, as detailed in Section S1 of the Supporting Information. In summary, convergence was attained for all 8 cases studied using LCIIS, EDIIS+CDIIS, CDIIS, and CDIIS-DL. Among these 4 methods, EDIIS+CDIIS and CDIIS each found 4 stable solutions, and CDIIS-DL found 5 stable solutions. ADIIS+CDIIS and LISTb each converged in 6 cases and found 3 or 2 stable solutions. LCIIS, on the other hand, converged in all 8 cases and found 6 stable solutions. These results suggest that LCIIS may serve as a useful alternative to either EDIIS+CDIIS or CDIIS, the two most commonly used default first-order SCF methods in modern quantum chemistry packages. In terms of convergence rate, LCIIS and EDIIS+CDIIS are similar in that LCIIS slightly outperforms EDIIS+CDIIS for SiH3 H and UF4 singlet restricted, while EDIIS+CDIIS slightly outperforms LCIIS for [Cd(Im)]2+ . Meanwhile, in 2 of these 8 cases, LCIIS leads to lower energy solutions than EDIIS+CDIIS.

6 6.1

Benchmarks on the MITM-COD data set Convergence of various first-order SCF methods

The previous section examined 8 cases in detail because those cases have been a primary focus of past work on SCF convergence acceleration. Here, we compare methods on the MITM-COD data set described in Section 4. We have chosen transition metal complexes because SCF convergence problems often arise in this class of compounds, possibly due to low HOMO-LUMO gaps. 40,55–57 To systematically evaluate the chance of convergence 17

ACS Paragon Plus Environment

Journal of Chemical Theory and Computation

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

CDIIS-DL

CDIIS-DL

Unconverged in 135 cases Confidence interval: [113, 158]

500

Page 18 of 52

100

0

0 20

40

60

80

100

120

20

40

60

100

120

100

120

100

120

60 80 100 Number of iterations

120

CDIIS Unconverged in 125 cases Confidence interval: [105, 148]

500

80

CDIIS 100

0

0 20

40

60

80

100

120

20

40

EDIIS+CDIIS 500

60

80

EDIIS+CDIIS

Unconverged in 115 cases Confidence interval: [96, 137]

100

0

0 20

40

60

80

100

120

20

40

Unconverged in 83 cases Confidence interval: [66, 102]

500

60

80

LCIIS

LCIIS 100

0

0 20

40

60 80 Number of iterations

100

120

20

(a) On all 2,070 cases

40

(b) On the selected 558 “hard-to-converge” cases

Figure 2: Histograms of number of iterations needed to reach SCF convergence on the MITM-COD data set, across four first-order SCF methods. The dot in the middle of the horizontal bar indicates the mean of iteration numbers, the length of the bar indicates the standard deviation, and the dot with a vertical bar indicates the median. The number of unconverged cases are listed along with 95% confidence intervals constructed using nonparametric bootstrapping. and convergence rate of different first-order SCF methods, we report statistical benchmark results on this data set with four first-order SCF methods implemented in Gaussian 09, namely CDIIS-DL, CDIIS, EDIIS+CDIIS, and LCIIS. There are 1,035 complexes in the MITM-COD data set, and we performed B3LYP/lanl2dz single point calculations on the two lowest spin states of each complex. Some statistics of these 2,070 calculations are shown in Figure 2(a), including median and mean number of iterations required to reach convergence along with the number of cases that failed to converge. We also selected a “hard-to-converge” data set composed of 558 cases where any of the four methods required more than 60 iterations to converge. The corresponding statistics on this selected data set is shown in Figure 2(b). The results indicate that, of these four first-order SCF methods, LCIIS has the fastest average rate of convergence and so is the most computationally efficient. Moreover, LCIIS has a higher chance of converging within 128 iterations. Figure 2(a) shows the number of unconverged cases along with 95% confidence intervals. The number of cases for which LCIIS fails to converge is lower than, and outside the confidence intervals of, the other three DIIS methods, and vice versa. Therefore, for the chemical systems represented 18

ACS Paragon Plus Environment

Page 19 of 52

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Journal of Chemical Theory and Computation

by the MITM-COD data set, LCIIS leads to the highest chance of SCF convergence in a statistically significant sense. When the maximum number of iterations allowed to reach convergence is increased to 256, the number of cases that converge increases for all methods. Nevertheless, the faster convergence rate of LCIIS still gives it an advantage over the other methods regarding the likelihood of reaching convergence (see Section S3 of the Supporting Information). Table 3: Statistical comparison of convergence rate of first-order SCF methods on the entire MITM-COD data set. †Conf. int. refers to confidence interval constructed using the same method as in Figure 2. Number of cases where both methods lead to the same solution LCIIS CDIIS-DL LCIIS CDIIS LCIIS EDIIS+CDIIS CDIIS CDIIS-DL EDIIS+CDIIS CDIIS

1711 1795 1774 1713 1760

Statistics of iteration number on cases with same solutions Median Mean †Conf. int. of mean 22.0 26.1 (25.6, 26.8) 26.0 34.8 (33.7, 35.8) 23.0 26.7 (26.1, 27.3) 25.0 29.6 (28.8, 30.4) 23.0 26.5 (25.9, 27.1) 24.0 28.2 (27.5, 28.9) 24.0 28.4 (27.7, 29.2) 26.0 34.4 (33.4, 35.4) 24.0 27.8 (27.1, 28.4) 25.0 28.9 (28.2, 29.7)

To better evaluate the convergence rate of the first-order methods, Tables 3 and 4 list oneon-one comparisons of pairs of first-order SCF methods. For these direct comparisons, we include only those cases where both methods converge to the same solution. We consider two methods to lead to the same solution if the absolute difference of the final energies is less than 10−6 Hartree. The tables report the median and mean of the number of iterations needed to reach convergence. Bootstrapping confidence intervals of the mean are also given. The nonoverlapping bootstrapping confidence intervals of LCIIS with those of other methods indicate that LCIIS’s convergence rate is faster than any of the other methods in a statistically significant sense. It is also interesting to note that, because the mean number of iterations for EDIIS+CDIIS falls outside the confidence interval of CDIIS (and vice versa), EDIIS+CDIIS 19

ACS Paragon Plus Environment

Journal of Chemical Theory and Computation

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Table 4: Statistical comparison of first-order SCF methods on the selected 558 “hard-toconverge” cases. †See Table 3 for details. Number of cases where both methods lead to the same solution LCIIS CDIIS-DL LCIIS CDIIS LCIIS EDIIS+CDIIS CDIIS CDIIS-DL EDIIS+CDIIS CDIIS

268 323 318 261 307

Statistics of iteration number on cases with same solutions Median Mean †Conf. int. of mean 42.0 45.5 (43.8, 47.4) 73.0 75.1 (72.8, 77.5) 41.0 44.4 (43.0, 46.0) 49.0 54.3 (52.0, 56.7) 41.0 44.1 (42.6, 45.8) 43.0 49.4 (47.4, 51.6) 47.0 52.7 (50.5, 55.3) 71.0 73.3 (71.1, 75.7) 43.0 47.8 (46.1, 49.8) 47.0 52.0 (50.0, 54.3)

convergence may also be considered faster than CDIIS in a statistically significant sense. Table 5: Number of lowest-energy solutions found by different first-order SCF methods for the 2,070 cases in MITM-COD data set.

CDIIS-DL CDIIS EDIIS+CDIIS LCIIS

Number 1804 1820 1827 1833

Confidence interval [1773, 1832] [1790, 1848] [1797, 1855] [1803, 1860]

Although LCIIS tends to converge faster, none of the four first-order methods considered here has a statistically significant advantage regarding the likelihood of finding the lowestenergy solution. Table 5 lists the number of cases for which each first-order method finds the lowest-energy solution. For this analysis, we use the lowest-energy solution obtained from any of these four first-order methods as an estimate of the lowest-energy SCF solution. LCIIS finds the greatest number of lowest-energy solutions. However, because the confidence interval of LCIIS includes the point estimates of the other methods and vice versa, the differences between the methods is not statistically significant. The advantages of LCIIS therefore relate only to the speed of finding a converged solution.

20

ACS Paragon Plus Environment

Page 20 of 52

Page 21 of 52

6.2

Performance of LCIIS and CDIIS during geometry optimizations

To further examine LCIIS’s use in practical quantum chemistry calculations, we also evaluate LCIIS’s performance in the common use case of geometry optimization. In this situation, the Harris guess is typically used only for the initial starting geometry. As the geometry is updated, during the search for a minimum-energy structure, the converged density matrix from the previous geometry is often used as the starting density for the SCF iterations. At the final stages of geometry optimization, where the change in geometry is small, this starting density is very close to the converged density of the current geometry. Here, we compare LCIIS to CDIIS for this situation, by examining the number of SCF iterations required in the final step of a geometry optimization. We also note that more advanced initial guess schemes have been developed for geometry optimizations, 58 but these are not explored here. CDIIS

600 400 200

15 cases in the right most bin Confidence interval: [9, 24]

*

0 5

10

15

20

25

30

35

LCIIS

600

Mean = 12.51 Median = 11 Standard deviation = 5.12

400 200

LCIIS vs CDIIS iteration number Reference line with slope 1

60

Mean = 12.63 Median = 11 Standard deviation = 5.75

SCF iteration number with LCIIS

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Journal of Chemical Theory and Computation

50

40

30

20

10

5 cases in the right most bin Confidence interval: [2, 11]

*

0 5

10

15 20 25 Number of iterations

30

0

35

0

10

20

30

40

50

60

SCF iteration number with CDIIS

(a) Histograms of SCF iteration numbers

(b) Scatter plot of number of iterations required for SCF convergence in LCIIS versus CDIIS

Figure 3: Statistical results of CDIIS and LCIIS’s performance for the final step of geometry optimization, where the density matrix from the previous geometry provides a good guess. Figure 3(a) shows results from 1,705 molecules in the MITM-COD data set whose geometries can be successfully optimized with Gaussian 09’s default setting. In a geometry optimization, the converged density matrix from the penultimate geometry can be consid21

ACS Paragon Plus Environment

Journal of Chemical Theory and Computation

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

ered as an extremely good initial guess for the last geometry. Nevertheless, the histograms in Figure 3(a) show that the average number of required SCF iterations is still greater than 10. Although LCIIS shows a slight advantage in terms of the average iteration number, this result is not statistically significant. LCIIS does show a statistically significant advantage on cases that require a larger number of SCF iterations. Only 5 cases required more than 35 SCF iterations with LCIIS (bootstrapping confidence interval [2, 11]), while with CDIIS there are 15 such cases (bootstrapping confidence interval [9, 24]). The tendency of LCIIS to converge faster when larger number of iterations are required is further illustrated in Figure 3(b). LCIIS and CDIIS perform similarly when less than 20 iterations are required for convergence. When more than 20 iterations are required, LCIIS tends to converge more quickly than CDIIS.

6.3

LCIIS as a potential complement to existing first-order SCF methods

Another criterion along which to judge the usefulness of a new first-order SCF method is whether there exist cases where the new method leads to a converged stable solution while existing methods do not. If such cases exist, the new method provides a useful alternative to try before resorting to expensive higher-order techniques. Beyond the 8 well-studied cases in Section 5, we report another 8 cases from the MITM-COD data set where LCIIS leads to stable SCF solutions while the other 3 methods fail to do so (Figure 4). These calculations were performed at B3LYP/lanl2dz level starting from the default Harris guess, with other settings at Gaussian 09 defaults. Figures 4(b), 4(e), and 4(g) are cases where LCIIS successfully achieved convergence while other methods failed to converge in 128 iterations. The remainder of Figure 4 are cases where the other methods converged to unstable solutions. These latter cases are of interest because, in cases where SCF convergence is problematic, different first-order SCF methods sometimes lead to different SCF solutions. Comparison across multiple methods is therefore a reasonable approach to test the reliabil22

ACS Paragon Plus Environment

Page 22 of 52

Page 23 of 52

10 log 10 |E i+1 - E i |

log 10 |E i+1 - E i |

10 5 COD entry ID 7221858 0 -5

!"#

5 COD entry ID 2006447 0 -5 -10

-10 0

10

20

30 40 50 Iteration number i

60

70

0

80

20

(a) Co(C7 H3 NO5 )(H2 O)3 quadruplet

40

60 80 Iteration number i

100

120

(b) Co(C5 H5 )(C3 S5 ) triplet

10

10

5

5

log 10 |E i+1 - E i |

log 10 |E i+1 - E i |

! COD entry ID 1540998 0 -5

" COD entry ID 7027308

0 -5 -10

-10 0

10

20

30 40 50 Iteration number i

60

70

0

80

20

(c) [Cr(H2 O)6 ]3+ doublet

60 80 Iteration number i

100

120

10

5

log 10 |E i+1 - E i |

log 10 |E i+1 - E i |

40

(d) Fe(C5 H5 )(C6 H8 P) triplet

10

COD entry ID 2009698 0 -5 -10

!"# !

5

COD entry ID 8102691 0 -5 -10

0

20

40

60 80 Iteration number i

100

120

0

20

(e) FeCl3 (C6 H12 S3 ) doublet

40

60 80 Iteration number i

100

120

(f) PdI2 (C12 H8 N2 ) triplet

10

10

!"#

5

log 10 |E i+1 - E i |

log 10 |E i+1 - E i |

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Journal of Chemical Theory and Computation

COD entry ID 4318245 0 -5 -10

5 COD entry ID 1535632 0 -5 -10

0

20

40

60 80 Iteration number i

100

120

0

(g) RuCl3 (NO)(C3 H7 NO)(C2 H6 OS) triplet

10

20 30 Iteration number i

40

50

60

(h) TiBr4 triplet

Figure 4: Instances in the MITM-COD data set where only LCIIS converges to a stable solution. All calculations were performed at B3LYP/lanl2dz level starting from Harris guess with other settings set to Gaussian 09 default. All energies are in Hartrees. 23

ACS Paragon Plus Environment

Journal of Chemical Theory and Computation

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

ity of the solutions. The faster average convergence rate of LCIIS, along with the potential to find low-energy solutions not found by other methods, suggest it is a useful addition to first-order SCF methods. The above relates to the extent to which LCIIS complements the other first-order methods. Section S4 of the Supporting Information considers the extent to which other methods can complement LCIIS. The results indicate that any one of these methods can find the lowest-energy solution for at most 57% of those cases where LCIIS fails to do so. The availability of a variety of first-order methods can reduce, but not eliminate, the need for costly higher-order SCF solution methods.

7

Conclusions

A new first-order method is introduced for accelerating SCF convergence. This LCIIS approach has some important fundamental distinctions from existing DIIS methods. LCIIS is similar to DIIS methods, in that the density matrix used for each iteration in the RoothaanHall process is written as a linear combination of past density matrices. It also has overlap with CDIIS in that both are based on the condition that the density matrix commutes with the Fock matrix, which is equivalent to the first-order KKT conditions of the HF/KS-DFT problem. However, most DIIS methods are essentially quadratic minimization schemes which differ from one another regarding the choice of the error vector and/or constraints placed on the linear coefficients. LCIIS instead involves a quartic formalism and so cannot be recast as DIIS with an alternate error vector. Rather than minimizing a sum of error vectors based on the commutator between the Fock and density matrices, LCIIS directly minimizes the commutator in the iterative subspace. Therefore, if the current iterative subspace consists of only idempotent density matrices and incorporates an SCF solution, LCIIS may be able to directly find that solution. Practically, at each SCF iteration, LCIIS searches for the desired solution within a limited iterative subspace consisting of both idempotent and non-

24

ACS Paragon Plus Environment

Page 24 of 52

Page 25 of 52

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Journal of Chemical Theory and Computation

idempotent matrices. Nevertheless, the empirical results indicate this significantly pushes forward the progress of SCF iterations. The results presented here suggest that LCIIS has advantages over other DIIS techniques. On eight difficult SCF convergence cases that have been the subject of extensive past work, LCIIS performs well and finds stable, lowest-energy solutions more often than the other first-order SCF methods. To obtain a more reliable measure of performance, LCIIS is also compared to three commonly used DIIS methods on 1,035 transition metal complexes. Statistics across this data set show that none of these four first-order techniques is more likely than the others to find the lowest-energy solution. However, LCIIS is found to be more computationally efficient, in that it converges more rapidly in a statistically significant sense. A related finding is that LCIIS has a higher chance of reaching SCF convergence than the other techniques. Finally, in some cases, LCIIS can find stable solutions not found by other first-order SCF methods and so is a useful complement to existing methods. We therefore consider LCIIS as a powerful alternative to EDIIS+CDIIS or CDIIS, the two most commonly used first-order SCF methods in modern quantum chemistry packages.

Acknowledgement This work was supported by the National Science Foundation under Grant CHE-1027985. We thank Prof. Gustavo E. Scuseria and Dr. Alejandro J. Garza for insightful discussion regarding singlet, broken-symmetry calculations on Ru4 (CO) and UF4 .

Supporting Information Supporting Information includes the following: detailed computational settings of [Cd(Im)]2+ and spin-restricted calculations of Ru4 (CO) and UF4 , investigation of the numerical stability of solving the LCIIS constrained quartic problem using Newton’s method, statistical benchmarks of first-order SCF methods with a different setting on the maximum allowed number of 25

ACS Paragon Plus Environment

Journal of Chemical Theory and Computation

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

SCF iterations, other first-order SCF methods as complements to LCIIS, and performance of LCIIS with the commutator defined in a non-orthogonal basis. This information is available free of charge via the Internet at http://pubs.acs.org.

References (1) Scalmani, G.; Frisch, M. J. J. Chem. Theory Comput. 2012, 8, 2193–2196. (2) Bulik, I. W.; Scalmani, G.; Frisch, M. J.; Scuseria, G. E. Phys. Rev. B 2013, 87, 035117. (3) Høst, S.; Olsen, J.; Jans´ık, B.; Thøgersen, L.; Jørgensen, P.; Helgaker, T. J. Chem. Phys. 2008, 129, 124106. (4) Natiello, M. A.; Scuseria, G. E. Int. J. Quantum Chem. 1984, 26, 1039–1049. (5) Pulay, P. Chem. Phys. Lett. 1980, 73, 393–398. (6) Pulay, P. J. Comput. Chem. 1982, 3, 556–560. (7) Sellers, H. Int. J. Quantum Chem. 1993, 45, 31–41. (8) Hutter, J.; L¨ uthi, H. P.; Parrinello, M. Comput. Mater. Sci. 1994, 2, 244–248. (9) Kudin, K. N.; Scuseria, G. E.; Cances, E. J. Chem. Phys. 2002, 116, 8255–8261. (10) Hu, X.; Yang, W. T. J. Chem. Phys. 2010, 132, 054109. (11) Banerjee, A. S.; Suryanarayana, P.; Pask, J. E. Chem. Phys. Lett. 2016, (12) Wang, Y. A.; Yam, C. Y.; Chen, Y. K.; Chen, G. J. Chem. Phys. 2011, 134, 241103. (13) Chen, Y. K.; Wang, Y. A. J. Chem. Theory Comput. 2011, 7, 3045–3048. (14) Garza, A. J.; Scuseria, G. E. J. Chem. Phys. 2015, 142, 164104. 26

ACS Paragon Plus Environment

Page 26 of 52

Page 27 of 52

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Journal of Chemical Theory and Computation

(15) Cances, E.; Le Bris, C. Int. J. Quantum Chem. 2000, 79, 82–90. (16) Saunders, V.; Hillier, I. Int. J. Quantum Chem. 1973, 7, 699–705. (17) Mitin, A. J. Comput. Chem. 1988, 9, 107–110. (18) D¨om¨ot¨or, G.; B´an, M. Comput. Chem. 1989, 13, 53–57. (19) Rabuck, A. D.; Scuseria, G. E. J. Chem. Phys. 1999, 110, 695–700. (20) Garza, A. J.; Scuseria, G. E. J. Chem. Phys. 2012, 137, 054110. (21) Bacskay, G. B. Chem. Phys. 1981, 61, 385–404. (22) Bacskay, G. B. Chem. Phys. 1982, 65, 383–396. (23) Chaban, G.; Schmidt, M. W.; Gordon, M. S. Theor. Chem. Acc. 1997, 97, 88–95. (24) Neese, F. Chem. Phys. Lett. 2000, 325, 93–98. (25) Thøgersen, L.; Olsen, J.; Yeager, D.; Jørgensen, P.; Salek, P.; Helgaker, T. J. Chem. Phys. 2004, 121, 16–27. (26) Høst, S.; Jans´ık, B.; Olsen, J.; Jørgensen, P.; Reine, S.; Helgaker, T. Phys. Chem. Chem. Phys. 2008, 10, 5344–5348. (27) ADF2016. Scientific Computing & Modelling (SCM), Theoretical Chemistry, Vrije Universiteit, Amsterdam, The Netherlands, http://www.scm.com/. (28) Scuseria, G. E.; Lee, T. J.; Schaefer, H. F. Chem. Phys. Lett. 1986, 130, 236–239. (29) Hamilton, T. P.; Pulay, P. J. Chem. Phys. 1986, 84, 5728–5734. (30) Cs´asz´ar, P.; Pulay, P. J. Mol. Struct. 1984, 114, 31–34. ¨ Schlegel, H. B. Phys. Chem. Chem. Phys. 2002, 4, 11–15. (31) Farkas, O.;

27

ACS Paragon Plus Environment

Journal of Chemical Theory and Computation

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

(32) Li, X.; Frisch, M. J. J. Chem. Theory Comput. 2006, 2, 835–839. (33) Frisch, M. J.; Trucks, G. W.; Schlegel, H. B.; Scuseria, G. E.; Robb, M. A.; Cheeseman, J. R.; Scalmani, G.; Barone, V.; Mennucci, B.; Petersson, G. A.; Nakatsuji, H.; Caricato, M.; Li, X.; Hratchian, H. P.; Izmaylov, A. F.; Bloino, J.; Zheng, G.; Sonnenberg, J. L.; Hada, M.; Ehara, M.; Toyota, K.; Fukuda, R.; Hasegawa, J.; Ishida, M.; Nakajima, T.; Honda, Y.; Kitao, O.; Nakai, H.; Vreven, T.; Montgomery, J. A., Jr.; Peralta, J. E.; Ogliaro, F.; Bearpark, M.; Heyd, J. J.; Brothers, E.; Kudin, K. N.; Staroverov, V. N.; Kobayashi, R.; Normand, J.; Raghavachari, K.; Rendell, A.; Burant, J. C.; Iyengar, S. S.; Tomasi, J.; Cossi, M.; Rega, N.; Millam, J. M.; Klene, M.; Knox, J. E.; Cross, J. B.; Bakken, V.; Adamo, C.; Jaramillo, J.; Gomperts, R.; Stratmann, R. E.; Yazyev, O.; Austin, A. J.; Cammi, R.; Pomelli, C.; Ochterski, J. W.; Martin, R. L.; Morokuma, K.; Zakrzewski, V. G.; Voth, G. A.; Salvador, P.; Dannenberg, J. J.; Dapprich, S.; Daniels, A. D.; Farkas, O.; Foresman, J. B.; Ortiz, J. V.; Cioslowski, J.; Fox, D. J. Gaussian 09 ; Revision D.01; Gaussian Inc.: Wallingford CT 2009. (34) Shao, Y.; Molnar, L. F.; Jung, Y.; Kussmann, J.; Ochsenfeld, C.; Brown, S. T.; Gilbert, A. T.; Slipchenko, L. V.; Levchenko, S. V.; O’Neill, D. P.; DiStasio Jr, R. A.; Lochan, R. C.; Wang, T.; Beran, G. J.; Besley, N. A.; Herbert, J. M.; Yeh Lin, C.; Van Voorhis, T.; Hung Chien, S.; Sodt, A.; Steele, R. P.; Rassolov, V. A.; Maslen, P. E.; Korambath, P. P.; Adamson, R. D.; Austin, B.; Baker, J.; Byrd, E. F. C.; Dachsel, H.; Doerksen, R. J.; Dreuw, A.; Dunietz, B. D.; Dutoi, A. D.; Furlani, T. R.; Gwaltney, S. R.; Heyden, A.; Hirata, S.; Hsu, C.-P.; Kedziora, G.; Khalliulin, R. Z.; Klunzinger, P.; Lee, A. M.; Lee, M. S.; Liang, W.; Lotan, I.; Nair, N.; Peters, B.; Proynov, E. I.; Pieniazek, P. A.; Min Rhee, Y.; Ritchie, J.; Rosta, E.; David Sherrill, C.; Simmonett, A. C.; Subotnik, J. E.; Lee Woodcock III, H.; Zhang, W.; Bell, A. T.; Chakraborty, A. K.; Chipman, D. M.; Keil, F. J.; Warshel, A.; Hehre, W. J.; Schae-

28

ACS Paragon Plus Environment

Page 28 of 52

Page 29 of 52

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Journal of Chemical Theory and Computation

fer III, H. F.; Kong, J.; Krylov, A. I.; Gill, P. M. W.; Head-Gordon, M. Phys. Chem. Chem. Phys. 2006, 8, 3172–3191. (35) Valiev, M.; Bylaska, E. J.; Govind, N.; Kowalski, K.; Straatsma, T. P.; Van Dam, H. J.; Wang, D.; Nieplocha, J.; Apra, E.; Windus, T. L.; de Jong, W. A. Comput. Phys. Commun. 2010, 181, 1477–1489. (36) Turney, J. M.; Simmonett, A. C.; Parrish, R. M.; Hohenstein, E. G.; Evangelista, F. A.; Fermann, J. T.; Mintz, B. J.; Burns, L. A.; Wilke, J. J.; Abrams, M. L.; Russ, N. J.; Leininger, M. L.; Janssen, C. L.; Seidl, E. T.; Allen, W. D.; Schaefer, H. F.; King, R. A.; Valeev, E. F.; Sherrill, C. D.; Crawford, T. D. WIREs Comput. Mol. Sci. 2012, 2, 556– 565. (37) Schmidt, M. W.; Baldridge, K. K.; Boatz, J. A.; Elbert, S. T.; Gordon, M. S.; Jensen, J. H.; Koseki, S.; Matsunaga, N.; Nguyen, K. A.; Su, S.; Windus, T. L.; Dupuis, M.; Montgomery Jr, J. A. J. Comput. Chem. 1993, 14, 1347–1363. (38) Gordon, M. S.; Michael, W. S. 2005, (39) Neese, F. WIREs Comput. Mol. Sci. 2012, 2, 73–78. (40) Sulzer, D.; Iuchi, S.; Yasuda, K. Chem. Phys. Lett. 2015, 635, 201–204. (41) Lipparini, F.; Scalmani, G.; Mennucci, B.; Frisch, M. J. J. Chem. Theory Comput. 2011, 7, 610–617. (42) Rohwedder, T.; Schneider, R. J. Math. Chem. 2011, 49, 1889–1914. (43) L¨owdin, P.-O. J. Chem. Phys. 1950, 18, 365–375. (44) Weigend, F. Phys. Chem. Chem. Phys. 2002, 4, 4285–4291. (45) Weigend, F. Phys. Chem. Chem. Phys. 2006, 8, 1057–1065. (46) Strout, D. L.; Scuseria, G. E. J. Chem. Phys. 1995, 102, 8448–8452. 29

ACS Paragon Plus Environment

Journal of Chemical Theory and Computation

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

(47) Echenique, P.; Alonso, J. L. Mol. Phys. 2007, 105, 3057–3098. (48) Hern´andez, E.; Gillan, M.; Goringe, C. Phys. Rev. B 1996, 53, 7147. (49) Saravanan, C.; Shao, Y.; Baer, R.; Ross, P. N.; Head-Gordon, M. J. Comput. Chem. 2003, 24, 618–622. (50) Sharada, S. M.; St¨ uck, D.; Sundstrom, E. J.; Bell, A. T.; Head-Gordon, M. Mol. Phys. 2015, 113, 1802–1808. (51) Graˇzulis, S.; Daˇskeviˇc, A.; Merkys, A.; Chateigner, D.; Lutterotti, L.; Quir´os, M.; Serebryanaya, N. R.; Moeck, P.; Downs, R. T.; Le Bail, A. Nucleic Acids Res. 2012, 40, D420–D427. (52) Macrae, C. F.; Bruno, I. J.; Chisholm, J. A.; Edgington, P. R.; McCabe, P.; Pidcock, E.; Rodriguez-Monge, L.; Taylor, R.; Streek, J. v.; Wood, P. A. J. Appl. Crystallogr. 2008, 41, 466–470. (53) Macrae, C. F.; Edgington, P. R.; McCabe, P.; Pidcock, E.; Shields, G. P.; Taylor, R.; Towler, M.; Streek, J. v. d. J. Appl. Crystallogr. 2006, 39, 453–457. (54) Zhang, Y. A.; Wang, Y. A. J. Chem. Phys. 2009, 130, 144116. (55) Dong, H.; Lin, B.; Gilmore, K.; Hou, T.; Lee, S.-T.; Li, Y. J. Power Sources 2015, 299, 371–379. (56) Chaban, V. V.; Fileti, E. E. Phys. Chem. Chem. Phys. 2015, 17, 15739–15745. (57) Colherinhas, G.; Fileti, E. E.; Chaban, V. V. Phys. Chem. Chem. Phys. 2015, 17, 17413–17420. (58) Muhlbach, A. H.; Vaucher, A. C.; Reiher, M. J. Chem. Theory Comput. 2016, 12, 1228–1235.

30

ACS Paragon Plus Environment

Page 30 of 52

Page 31 of 52

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Journal of Chemical Theory and Computation

CDIIS 100 0

Error vector 20

60

100

140

LCIIS

Norm minimization

100 0

20

60 100 Number of iterations

140

Figure 5: TOC graphic.

31

ACS Paragon Plus Environment

Journal of Chemical Theory and Computation

E(EDIIS+CDIIS) E(CDIIS) E(CDIIS-DL) E(ADIIS+CDIIS) E(LISTb) E(LCIIS)

15

log 10 |E i+1 - E i |

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22

10 5

= = = = = =

-152.914403 -152.914403 -152.914403 -152.914403 -152.914403 -152.914403

Page 32 of 52 (stable) (stable) (stable) (stable) (stable) (stable)

RHF/6-31g* Harris guess E* = -152.914403

0 -5 -10 0

ACS Paragon Plus Environment

5

Iteration number i

10

15

Page 33 of 52

Journal of Chemical Theory and Computation

E(EDIIS+CDIIS) E(CDIIS) E(CDIIS-DL) E(ADIIS+CDIIS) E(LISTb) E(LCIIS)

15

log 10 |E i+1 - E i |

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22

10 5

= = = = = =

-290.457688 -290.457688 -290.457688 -290.457688 -290.457688 -290.457688

(stable) (stable) (stable) (stable) (stable) (stable)

RSVWN5/6-31g* Harris guess E* = -290.457688

0 -5 -10 0

5

ACS Paragon Plus Environment

10

15 20 Iteration number i

25

30

35

Journal of Chemical Theory and Computation

E(EDIIS+CDIIS) E(CDIIS) E(CDIIS-DL) E(ADIIS+CDIIS) E(LISTb) E(LCIIS)

15

log 10 |E i+1 - E i |

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22

10 5

= = = = = =

-5667.008480 -5667.008480 -5667.008480 -5667.008480 -5667.008480 -5667.008480

Page 34 of 52 (unstable) (unstable) (unstable) (unstable) (unstable) (unstable)

RB3LYP/3-21g Harris guess E* = -5667.009041

0 -5 -10 0

5

ACS Paragon Plus Environment

10 15 Iteration number i

20

25

30

Page 35 of 52

Journal of Chemical Theory and Computation

E(EDIIS+CDIIS) E(CDIIS) E(CDIIS-DL) E(ADIIS+CDIIS) E(LISTb) E(LCIIS)

15

log 10 |E i+1 - E i |

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22

10 5

= = = = = =

-1082.114016 -1082.114016 -1082.114808 -1082.114808 -1082.114016 -1082.114808

(unstable) (unstable) (stable) (stable) (unstable) (stable)

RB3LYP/6-31g Harris guess E* = -1082.114808

0 -5 -10 0

10

20

ACS Paragon Plus Environment

30 40 50 Iteration number i

60

70

80

Journal of Chemical Theory and Computation

E(EDIIS+CDIIS) E(CDIIS) E(CDIIS-DL) E(ADIIS+CDIIS) E(LISTb) E(LCIIS)

15

log 10 |E i+1 - E i |

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22

10 5

= = = = = =

-488.691402 -488.691402 -488.691402 -488.691402 -488.704767 -488.704767

Page 36 of 52 (unstable) (unstable) (unstable) (unconverged) (unconverged) (stable)

RB3LYP/lanl2dz Harris guess E* = -488.704767

0 -5 -10 0

20

ACS Paragon Plus Environment

40

60 80 Iteration number i

100

120

Page 37 of 52

Journal of Chemical Theory and Computation

E(EDIIS+CDIIS) E(CDIIS) E(CDIIS-DL) E(ADIIS+CDIIS) E(LISTb) E(LCIIS)

15

log 10 |E i+1 - E i |

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22

10 5 0 -5

= = = = = =

-488.740761 (stable) -488.740761 (stable) -488.740761 (stable) -488.740760 (unconverged) -488.740756 (unconverged) -488.740761 (stable) UB3LYP/lanl2dz Triplet guess E* = -488.740761

-10 0

20

ACS Paragon Plus Environment

40

60 80 Iteration number i

100

120

Journal of Chemical Theory and Computation

E(EDIIS+CDIIS) E(CDIIS) E(CDIIS-DL) E(ADIIS+CDIIS) E(LISTb) E(LCIIS)

15

log 10 |E i+1 - E i |

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22

10 5

= = = = = =

-451.178080 -451.178080 -451.178080 -451.178080 -451.178080 -451.178080

Page 38 of 52 (unstable) (unstable) (unstable) (unstable) (unstable) (unstable)

RB3LYP/lanl2dz Harris guess E* = -451.219602

0 -5 -10 0

5

10

ACS Paragon Plus Environment

15 20 25 Iteration number i

30

35

40

Page 39 of 52

Journal of Chemical Theory and Computation

E(EDIIS+CDIIS) E(CDIIS) E(CDIIS-DL) E(ADIIS+CDIIS) E(LISTb) E(LCIIS)

15

log 10 |E i+1 - E i |

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22

10 5

= = = = = =

-451.259876 -451.259876 -451.259876 -451.251014 -451.251014 -451.259876

(stable) (stable) (stable) (unstable) (unstable) (stable)

UB3LYP/lanl2dz Triplet guess E* = -451.259876

0 -5 -10 0

20

ACS Paragon Plus Environment

40 60 Iteration number i

80

100

120

Journal of Chemical Theory and Computation

Page 40 of 52

CDIIS-DL 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33

Unconverged in 135 cases Confidence interval: [113, 158]

500 0 20

40

60

80

100

120

CDIIS Unconverged in 125 cases Confidence interval: [105, 148]

500 0 20

40

60

80

100

120

EDIIS+CDIIS 500

Unconverged in 115 cases Confidence interval: [96, 137]

0 20

40

60

80

100

120

LCIIS Unconverged in 83 cases Confidence interval: [66, 102]

500 0 20

40

60 Plus Environment 80 ACS Paragon Number of iterations

100

120

Page 41 of 52

Journal of Chemical Theory and Computation

CDIIS-DL 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33

100 0 20

40

60

80

100

120

100

120

100

120

60Paragon Plus80 ACS Environment 100

120

CDIIS 100 0 20

40

60

80

EDIIS+CDIIS 100 0 20

40

60

80

LCIIS 100 0 20

40

Number of iterations

Journal of Chemical Theory and Computation

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33

CDIIS

600

Page 42 of 52

Mean = 12.63 Median = 11 Standard deviation = 5.75

400 200 0 5

10

15

20

15 cases in the right most bin Confidence interval: [9, 24] * 25 30 35

LCIIS

600

Mean = 12.51 Median = 11 Standard deviation = 5.12

400 200 0 5

10

5 cases in the right most bin Confidence interval: [2, 11] * ACS Paragon Plus Environment 15 20 25 30 35 Number of iterations

Page 43 Journal of 52 ofLCIIS Chemical Theory andnumber Computation vs CDIIS iteration 60 Reference line with slope 1

SCF iteration number with LCIIS

1 2 50 3 4 5 40 6 7 8 9 30 10 11 1220 13 14 15 10 16 17 18 19 0 20 0 21

ACS Paragon Plus Environment 10 20 30 40 50

SCF iteration number with CDIIS

60

Journal of Chemical Theory and Computation

E(EDIIS+CDIIS) E(CDIIS) E(CDIIS-DL) E(LCIIS)

10

log 10 |E i+1 - E i |

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22

5

= = = =

-1073.592830 -1073.592830 -1073.592830 -1073.595671

Page 44 of 52 (unstable) (unstable) (unstable) (stable)

COD entry ID 7221858 0 -5 -10 0

10

20

ACS Paragon Plus Environment

30 40 50 Iteration number i

60

70

80

Page 45 of 52

Journal of Chemical Theory and Computation

E(EDIIS+CDIIS) E(CDIIS) E(CDIIS-DL) E(LCIIS)

10

log 10 |E i+1 - E i |

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22

5

= = = =

-503.192940 -503.192645 -503.192914 -503.193717

(unconverged) (unconverged) (unconverged) (stable)

COD entry ID 2006447 0 -5 -10 0

20

ACS Paragon Plus Environment

40

60 80 Iteration number i

100

120

Journal of Chemical Theory and Computation

E(EDIIS+CDIIS) E(CDIIS) E(CDIIS-DL) E(LCIIS)

10

log 10 |E i+1 - E i |

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22

5

= = = =

-543.835793 -543.835793 -543.835793 -543.836417

Page 46 of 52 (unstable) (unstable) (unstable) (stable)

COD entry ID 1540998 0 -5 -10 0

10

20

ACS Paragon Plus Environment

30 40 50 Iteration number i

60

70

80

Page 47 of 52

Journal of Chemical Theory and Computation

E(EDIIS+CDIIS) E(CDIIS) E(CDIIS-DL) E(LCIIS)

10

log 10 |E i+1 - E i |

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22

5

= = = =

-556.473254 -556.473472 -556.472044 -556.474241

(unconverged) (unstable) (unstable) (stable)

COD entry ID 7027308 0 -5 -10 0

20

ACS Paragon Plus Environment

40

60 80 Iteration number i

100

120

Journal of Chemical Theory and Computation

E(EDIIS+CDIIS) E(CDIIS) E(CDIIS-DL) E(LCIIS)

10

log 10 |E i+1 - E i |

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22

5

= = = =

-434.411648 -434.411647 -434.411648 -434.411648

Page 48 of 52 (unconverged) (unconverged) (unconverged) (stable)

COD entry ID 2009698 0 -5 -10 0

20

ACS Paragon Plus Environment

40

60 80 Iteration number i

100

120

Page 49 of 52

Journal of Chemical Theory and Computation

E(EDIIS+CDIIS) E(CDIIS) E(CDIIS-DL) E(LCIIS)

10

log 10 |E i+1 - E i |

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22

5

= = = =

-720.947798 -720.947798 -720.947798 -720.953685

(unstable) (unstable) (unconverged) (stable)

COD entry ID 8102691 0 -5 -10 0

20

ACS Paragon Plus Environment

40

60 80 Iteration number i

100

120

Journal of Chemical Theory and Computation

E(EDIIS+CDIIS) E(CDIIS) E(CDIIS-DL) E(LCIIS)

10

log 10 |E i+1 - E i |

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22

5

= = = =

-681.965188 -681.965923 -681.965932 -681.977610

Page 50 of 52 (unconverged) (unconverged) (unconverged) (stable)

COD entry ID 4318245 0 -5 -10 0

20

ACS Paragon Plus Environment

40

60 80 Iteration number i

100

120

Page 51 of 52

Journal of Chemical Theory and Computation

E(EDIIS+CDIIS) E(CDIIS) E(CDIIS-DL) E(LCIIS)

10

log 10 |E i+1 - E i |

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22

5

= = = =

-110.870147 -110.870147 -110.870147 -110.871883

(unstable) (unstable) (unstable) (stable)

COD entry ID 1535632 0 -5 -10 0

10

ACS Paragon Plus Environment

20 30 Iteration number i

40

50

60

Journal of Chemical Theory and Computation

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

CDIIS 100 0

Error vector 20

60

100

140

LCIIS

Norm minimization

100 0

20

60 100 Number of iterations

140

ACS Paragon Plus Environment

Page 52 of 52