4.6 Large Molecules and Linear Scaling Methods

4.6.9 occ-RI-K Exchange Algorithm

The occupied orbital RI-K (occ-RI-K) algorithm620 is a new scheme for building the exchange matrix (K) partially in the MO basis using the RI approximation. occ-RI-K typically matches current alternatives in terms of both the accuracy (energetics identical to standard RI-K) and convergence (essentially unchanged relative to conventional methods). On the other hand, this algorithm exhibits significant speedups over conventional integral evaluation (14x) and standard RI-K (3.3x) for a test system, a graphene fragment (C68H22) using cc-pVQZ basis set (4400 basis functions), whereas the speedup increases with the size of the AO basis set. Thus occ-RI-K helps to make larger basis set hybrid DFT calculations more feasible, which is quite desirable for achieving improved accuracy in DFT calculations with modern functionals.

The idea of the occ-RI-K formalism comes from a simple observation that the exchange energy EK and its gradient can be evaluated from the diagonal elements of the exchange matrix in the occupied-occupied block Kii, and occupied-virtual block Kia, respectively, rather than the full matrix in the AO representation, Kμν. Mathematically,

EK = -μνPμνKμν (4.50)
= -μνcμiKμνcνi
= -iKii


EKΔai=2Kai (4.51)

where Δ is a skew-symmetric matrix used to parameterize the unitary transformation U, which represents the variations of the MO coefficients as follows:

U=e(Δ-ΔT). (4.52)

From Eq. 4.50 and 4.51 it is evident that the exchange energy and gradient need just Kiν rather than Kμν.

In regular RI-K one has to compute two quartic terms,1022 whereas there are three quartic terms for the occ-RI-K algorithm. The speedup of the latter with respect to former can be explained from the following ratio of operations; refer to Ref. 620 for details.

# of RI-K quartic operations# of occ-RI-K quartic operations=oNX2+oN2Xo2X2+o2NX+o2NX=N(X+N)o(X+2N) (4.53)

With a conservative approximation of X2N, the speedup is 34(N/o). The occ-RI-K algorithm also involves some cubic steps which should be negligible in the very large molecule limit. Tests in the Ref. 620 suggest that occ-RI-K for small systems with large basis will gain less speed than a large system with small basis, because the cubic terms will be more dominant for the former than the latter case.

In the course of SCF iteration, the occ-RI-K method does not require us to construct the exact Fock matrix explicitly. Rather, kiν contributes to the Fock matrix in the mixed MO and AO representations (Fiν) and yields orbital gradient and DIIS error vectors for converging SCF. On the other hand, since occ-RI-K does not provide exactly the same unoccupied eigenvalues, the diagonalization updates can differ from the conventional SCF procedure. In Ref. 620, occ-RI-K was found to require, on average, the same number of SCF iterations to converge and to yield accurate energies.

The occ-RI-K can be invoked by either setting OCC_RI_K to be true, or (since Q-Chem 5.2) specifying auxiliary basis set for K using AUX_BASIS_K.

       Controls the use of the occ-RI-K approximation for constructing the exchange matrix
       False Do not use occ-RI-K.
       True Use occ-RI-K.
       Larger the system, better the performance occ-RI-K for exchange energy gradient evaluation

A very attractive feature of occ-RI-K framework is that one can compute the exchange energy gradient with respect to nuclear coordinates with the same leading quartic-scaling operations as the energy calculation.

The occ-RI-K formulation yields the following formula for the gradient of exchange energy in global Coulomb-metric RI:

EKx = (ij|ij)x (4.54)
= μνPijcμicνjCijP(μν|P)x-RSijCijRCijS(R|S)x.

The superscript x represents the derivative with respect to a nuclear coordinate. Note that the derivatives of the MO coefficients cμi are not included here, because they are already included in the total energy derivative calculation by Q-Chem via the derivative of the overlap matrix.

In Eq. 4.54, the construction of the density fitting coefficients (CμνP) has the worst scaling of 𝒪(M4) because it involves MO to AO back transformations:

CμνP=ijcμicνjCijP (4.55)

where the operation cost is o2NX+o[NB2]X.

       Turn on the nuclear gradient calculations
       FALSE Do not invoke occ-RI-K based gradient
       TRUE Use occ-RI-K based gradient
       Use "RI_J false"

Example 4.19  Q-Chem input for a energy and gradient calculations with occ-RI-K method.

0   1
  C    0.0000000    0.3057430    5.7138876
  C    0.0000000   -0.5442831    4.4256275
  C    0.0000000    0.3675825    3.1787857
  C    0.0000000   -0.4853210    1.8908707
  C    0.0000000    0.4264823    0.6439435
  C    0.0000000   -0.4264823   -0.6439435
  C    0.0000000    0.4853210   -1.8908707
  C    0.0000000   -0.3675825   -3.1787857
  C    0.0000000    0.5442831   -4.4256275
  C    0.0000000   -0.3057430   -5.7138876
  H    0.8973039    0.9418895    5.7433715
  H   -0.8973039    0.9418895    5.7433715
  H   -0.8965960   -1.1827891    4.4155828
  H    0.8965960   -1.1827891    4.4155828
  H    0.8966512    1.0057416    3.1940146
  H   -0.8966512    1.0057416    3.1940146
  H   -0.8966512   -1.1234851    1.8761493
  H    0.8966512   -1.1234851    1.8761493
  H    0.8966507    1.0646443    0.6587590
  H   -0.8966507    1.0646443    0.6587590
  H   -0.8966507   -1.0646443   -0.6587590
  H    0.8966507   -1.0646443   -0.6587590
  H    0.8966512    1.1234851   -1.8761493
  H   -0.8966512    1.1234851   -1.8761493
  H   -0.8966512   -1.0057416   -3.1940146
  H    0.8966512   -1.0057416   -3.1940146
  H    0.8965960    1.1827891   -4.4155828
  H   -0.8965960    1.1827891   -4.4155828
  H    0.8973039   -0.9418895   -5.7433715
  H    0.0000000    0.3580832   -6.5913113
  H   -0.8973039   -0.9418895   -5.7433715
  H    0.0000000   -0.3580832    6.5913113

   JOBTYPE         force
   EXCHANGE        HF
   BASIS           cc-pVTZ
   AUX_BASIS       cc-pVTZ-JK
   OCC_RI_K        1
   RI_K_GRAD       1
   INCFOCK         0
   PURECART        1111