# Quantifying relativistic effects in DFT calculations

I am wondering what is the best approach to quantify the contribution of relativistic corrections to a result of an electronic structure calculation. As far as I understand, from a user perspective one can simply "turn on" some relativistic effects in most electronic structure codes, and the result should contain an approximation to relativistic effects. A small caveat is that a specially contracted/parametrized basis set should be used for the relativistic calculation (Gaussian-type orbitals).

What I am wondering is whether I can just take the difference between a result obtained with relativistic correction and one without, and assume that what I obtain corresponds to the magnitude of the relativistic correction, even if the basis set used in the calculation are different? Or, if I am interested in quantifying such a difference, should I use the same basis set for the relativistic and the non-relativistic calculation, and, if so, should it be the relativistic one? I could also imagine that the right approach should be to converge both the relativistic and the ''standard'' calculation with respect to their own basis sets, but that might not be computationally feasible for some larger calculations, so I am wondering about the approach that maximizes the error cancellation due to the finite basis.

As an example, calculating the dipole polarizability of a single Argon atom with PBE0 functional with aug-cc-pVDZ basis set (as implemented in NWChem) I get $$9.98$$ a.u. for the isotropic term, whereas if I add a DK2 correction with the aug-cc-pVDZ-DK basis set I obtain $$10.0$$ a.u. Can I say that the difference of 0.02 a.u. in the results is due to the DK2 correction, or is a large portion of the difference due to the difference in the aug-cc-pVDZ and the aug-cc-pVDZ-DK basis sets?

• Out of curiosity, how many atoms are in the systems of your interest? Surely for a single argon atom (which was the example you gave) you can get much better accuracy with a wavefunction-based method rather than DFT, so there must be a larger system you have in mind for which DFT is not an unwise method to choose. Jan 11 at 10:58
• The example of the single argon atom was just something I can run on my PC quickly to put some illustrative numbers for my question. My 'real' problem would be predicting the VCD spectrum of a small-ish oligopeptide. Jan 12 at 15:14
• The reason I asked that question is because it makes a big difference. I thought about mentioning uncontracted basis sets too, but you likely can't use them for an oligopeptide due to linear dependencies. Linear dependencies in aug-cc-pVZ-unc are not a problem for a single argon atom as in your example, because there won't be any linear dependencies in the single atom calculation, but there will be substantial overlap in a many-atom molecule. Since you're using DFT, I was certain that you would be doing calculations on a system that is too large to use an uncontracted aVDZ basis set properly. Jan 13 at 5:00

The straightforward way for quantifying relativistic effects with NWChem (and with any other software): with the same DECONTRACTED basis set do the nonrelativistic calculation of your property, then do the corresponding relativistic calculation.

Concerning relativistic calculations of properties you should be aware if there is the so-called picture change involved. For polarizabilities (valence property) this effect is small, for core properties this effect is more profound.

In reality, the contribution to an energy due to relativity, would be the difference between the "exact" energy with a model that captures all relativistic effects, versus a model that captures everything except any relativistic effects.

For an energy to be "exact", when using Gaussian-type orbitals as you described, you need the energy to be converged to its CBS (complete basis set) limit and to its full correlation limit (this means using "the exact density functional" or "full configuration interaction"), in extremely large basis sets, and extrapolating to the CBS limit.

However, no quantum mechanical model will be able to capture "all" relativistic effects because we don't yet have a quantum mechanical model that perfectly captures the relativistic effects associated with gravity (general relativistic effects). Ignoring gravity (its effect on the interaction energy between two atoms is anticipated to be something like $$10^{-40}$$ in atomic units anyway), we still don't have a simple model that captures all remaining relativistic effects (special relativity effects), because the Dirac equation only works for one electron, the Dirac-Coulomb-Breit equation works for many electrons but only approximates the true relativistic behavior, and QED (quantum electrodynamics) doesn't provide us with an easy way to calculate molecular energies.

The best that you can likely do today, is to calculate an approximation of the relativistic effects on an energy, by using something like the X2C (exact 2-component) Hamiltonian, which for quantum chemistry purposes is quite a good approximation to the effects of special relativity. You mentioned the DK2 (Douglas-Kroll) Hamiltonian, which is part of a series of more and more accurate Hamiltonians (DK2,DK3,DK4,etc.), but the X2C Hamiltonian is equivalent to a an "infinite-order" DK Hamiltonian, and is just as easy to use as any of the DK Hamiltonians, so the DK method is now obsolete and there is no reason to use it now that the X2C Hamiltonian is implemented in a lot of quantum chemistry packages and is superior to the DK Hamiltonian in every conceivable way.

Since it's unlikely for you to be able to get energies at the CBS limit and full correlation limit for a system that's large enough for you to need to use DFT instead of a wavefunction-based method, the best estimate you can get for the contribution of relativity to an energy, is likely by doing an energy calculation with the X2C model and with the non-relativistic model, using the same basis set in both cases and the same density functional in both cases. Although X2C basis sets such as aug-cc-pVDZ-X2C exist, you don't need to use them for this type of comparison. If you do want to use it, you can still compare the energy that you get with the X2C model and an X2C basis set, to the energy that you get with a non-relativistic model and an equivalent "standard" basis set, since this "effect of relativity on the energy" will only be approximate anyway.

• I finally caved in and logged into SE on my tablet, to answer this question. For almost two years I had not signed-in on this device, partially to help me reduce the number of hours I spend on this site :) This year I won't be able to spend much (or any) time on a computer during weekdays though, and I wanted to answer this question! Jan 8 at 13:33
• It is still a bit unclear to me whether a comparison of DFT+X2C/normal basis vs. DFT/normal basis OR a comparison of DFT+X2C/X2C basis vs. DFT/normal basis is more faithful. You mention both possibility, with the caveat that the result is approximate anyways, but is there a way to know which approach is less approximate, if it even makes sense? Jan 8 at 14:51
• @Szgoger I would need more space than a comment to answer that, and I don't get much keyboard time on weekdays now, so asking that in a new post would not only make it more likely to get an answer from one of the other ~7500 users, but would also allow the answerer to have more space to give a detailed answer. Jan 9 at 12:31