# Convergence by energy and/or density?

The SCF is an iterative procedure to produce the ground state energy and wavefunction for Hartree-Fock, MCSCF, and DFT. Due to the iteration, one needs to establish convergence criteria. Depending on the program used, the default criteria could be just the change in energy (ORCA), the RMS/max change in the density (Gaussian, Q-Chem), or a combination of both (Psi4).

I would think you would always want to use both the energy and density in the convergence criteria. Since you need to compute both anyways, it should negligibly increase the cost to check both and would give more information about convergence than either energy/density alone. What would be the advantage of not including the energy or density in the SCF convergence criteria?

• If only aiming to get the SCF energy, then going for only convergence in the energy is fine, but if doing any post-SCF calculation like coupled cluster or CI, then one must make sure that the largest difference in the density is very small (10^-8 in my applications). This is absolutely crucial, because in my applications, the energy converges several iterations before the largest density difference, meaning that checking only the energy, is not enough and can be a red herring. I'm surprised that ORCA only checks the energy by default. But they are not focused on accuracy. May 10 '20 at 20:34
• It may also just be a programming accident. For example the devs just decided to not include that criteria... May 10 '20 at 20:46
• @NikeDattani that could be a mistake on my part. I'm mainly familiar with Gaussian and a little bit Psi4. I went looking for how other programs did it, ORCA's website made it seem like it was just energy, but that might just be for SCF calculations and I just noticed they say other convergence criteria are affected by tightness options, but they don't specify what those are.
– Tyberius
May 10 '20 at 21:02
• I'd also be interested if anyone can show an example where the density differences converge to some reasonable amount, like 10^-8 or even 10^-6, and the energy is still not converged in the micro-Hartree digit. May 10 '20 at 21:12
• @Nike Dattani Just an info about ORCA: They not only check total energy convergence but also 1-electron energy convergence, which necessarily changes if the density is still not converged due to the linear dependency. May 11 '20 at 6:44

Approaching the fixed-point problem an SCF solves from the direct minimisation (DM) point of view (for the relation of DM to SCF, see my previous post, the key quantities to check for convergence are actually the total energy and it's derivative, i.e. the occupied-virtual block of the Fock / Kohn-Sham matrix. The latter is after all the gradient and that's what you want to drive to zero in a fixed point. In an SCF context you could theoretically do the same and to the best of my knowledge some codes like PySCF actually use these two criteria for convergence.

Of course there are other choices, which can be done and using the density works the same, since a zero density change implies a zero change in the Fock matrix and vice versa. From that point of view the density is kind of a substitute for checking convergence in the gradient, but of course not the same threshold can be used for both checks and tolerances need to be adjusted between them.

Since clearly the gradient being zero is the whole aim of an SCF, I find it dubious to only check for convergence in the energy as it could accidentally stagnate without being at a fixed point. In turn if you check for the density / gradient change being zero you should be safe. However, as @Tyberius already pointed out, computing the energy is pretty much for free if you compute density and Fock matrix and therefore performing a check on the energy does not cost and shouldn't do harm either.

Talking about convergence speeds: In absolute numbers the energy always converges first, since it depends quadratically on the density (so error 1e-3 in the density usually translates to 1e-6 in the energy). Therefore convergence criteria on energy should be the square of the convergence tolerance in the density / gradient. But notice that many codes are a bit sloppy here and do not print the norm of the density / gradient change in their outputs, but the norm squared, such that both numbers are again on the same scale (e.g. ABINIT does this). So one has to be a little careful when looking at the numbers presented from codes.

Let me point out one last subtlety from a mathematical point of view. When checking the convergence in density / gradient one is faced with the problem that these quantities are vectors or matrices, but of course we want to check convergence in a number. So we need to take the norm of the difference. But which norm should we choose? Most codes to the best of my knowledge use the Frobenius or l2 norms (so just square the elements, add them up and take the square root), but this is just one choice. There are plenty of other norms, see the wikipedia articles on the Lp norms and the matrix norms to get an idea. Which norm is the best one to choose depends a little on which property you are after in your calculation: Total energy, forces / gradients wrt. nuclear positions, partial charges etc, just because these ask different questions about the mathematical properties of the wavefunction.

Now, the take-away from this is not that one should use a different norm depending on what property one is after. This is just impractical and also in calculations with their finite basis set sizes all norms are equivalent up to a constant. But this constant might not be small and usually depends on the size of the basis (it typically grows with larger bases). So depending on what property you are after in the end, the number of digits you can trust in the computed answer differs. For example, even if the density is converged to 6 digits in the Frobenius norm, the forces might only be correct to 5 and the partial charges to 4.

Note that Q-Chem and Psi4 are not actually measuring the change in the density, but rather the orbital gradient. There is a subtle difference: no change in the density does not mean the energy has been minimized, but if the orbital gradient is zero you're mathematically certain to be at an extremal point, which can be either a local minimum or a saddle-point solution.