4.5.7 Relaxed Constraint Algorithm (RCA)

The relaxed constraint algorithm (RCA) is an ingenious and simple means of minimizing the SCF energy that is particularly effective in cases where the initial guess is poor. The latter is true, for example, when employing a user-specified basis (when the “core” or GWH guess must be employed) or when near-degeneracy effects imply that the initial guess will likely occupy the wrong orbitals relative to the desired converged solution.

Briefly, RCA begins with the SCF problem as a constrained minimization of the energy as a function of the density matrix, $E(𝐏)$.^{136, 139} The constraint is that the density matrix be idempotent, $𝐏\cdot 𝐏=𝐏$, which basically forces the occupation numbers to be either zero or one. The fundamental realization of RCA is that this constraint can be relaxed to allow sub-idempotent density matrices, $𝐏\cdot 𝐏\le 𝐏$. This condition forces the occupation numbers to be between zero and one. Physically, we expect that any state with fractional occupations can lower its energy by moving electrons from higher energy orbitals to lower ones. Thus, if we solve for the minimum of $E(𝐏)$ subject to the relaxed sub-idempotent constraint, we expect that the ultimate solution will nonetheless be idempotent.

In fact, for Hartree-Fock this can be rigorously proven. For density functional theory, it is possible that the minimum will have fractional occupation numbers but these occupations have a physical interpretation in terms of ensemble DFT. The reason the relaxed constraint is easier to deal with is that it is easy to prove that a linear combination of sub-idempotent matrices is also sub-idempotent as long as the linear coefficients are between zero and one. By exploiting this property, convergence can be accelerated in a way that guarantees the energy will go down at every step.

The implementation of RCA in Q-Chem closely follows the “Energy DIIS” implementation of the RCA algorithm.⁵⁰⁹ Here, the current density matrix is written as a linear combination of the previous density matrices:

$$𝐏(x)=\sum _i{x}_i{𝐏}_i$$

(4.36)

To a very good approximation (exact for Hartree-Fock) the energy for $𝐏(x)$ can be written as a quadratic function of x:

$$𝐄(x)=\sum _i{E}_i{x}_i+\frac{1}{2}\sum _i{x}_i({𝐏}_i-{𝐏}_j)\cdot ({𝐅}_i-{𝐅}_j)x_j$$

(4.37)

At each iteration, $x$ is chosen to minimize $𝐄(x)$ subject to the constraint that all of the $x_i$ are between zero and one. The Fock matrix for $𝐏(x)$ is further written as a linear combination of the previous Fock matrices,

$$𝐅(x)=\sum _i{x}_i{𝐅}_i+\delta {𝐅}_{xc}(x)$$

(4.38)

where $\delta {𝐅}_{xc}(x)$ denotes a (usually quite small) change in the exchange-correlation part that is computed once $x$ has been determined. We note that this extrapolation is very similar to that used by DIIS. However, this procedure is guaranteed to reduce the energy $𝐄(x)$ at every iteration, unlike DIIS.

In practice, the RCA approach is ideally suited to difficult convergence situations because it is immune to the erratic orbital swapping that can occur in DIIS. On the other hand, RCA appears to perform relatively poorly near convergence, requiring a relatively large number of steps to improve the precision of a good approximate solution. It is thus advantageous in many cases to run RCA for the initial steps and then switch to DIIS either after some specified number of iterations or after some target convergence threshold has been reached. Finally, note that by its nature RCA considers the energy as a function of the density matrix. As a result, it cannot be applied to restricted open shell calculations which are explicitly orbital-based. Note: RCA interacts poorly with INCDFT, so INCDFT is disabled by default when an RCA or RCA_DIIS calculation is requested. To enable INCDFT with such a calculation, set INCDFT = 2 in the $rem section. RCA may also have poor interactions with incremental Fock builds; if RCA fails to converge, setting INCFOCK = FALSE may improve convergence in some cases.

Job-control variables for RCA are listed below, along with an example of its use.

$molecule
0 2
   H    1.004123   -0.180454    0.000000
   O   -0.246002    0.596152    0.000000
   O   -1.312366   -0.230256    0.000000
$end

$rem
   UNRESTRICTED         true
   METHOD               hf
   BASIS                cc-pVDZ
   SCF_GUESS            gwh
   SCF_ALGORITHM        RCA_DIIS
   THRESH               9
$end