8.4: Perturbative Approaches
( \newcommand{\kernel}{\mathrm{null}\,}\)
The situation becomes much more difficult if we need to account for direct interactions between the particles. Let us assume that the interaction may be reduced to that between their pairs (as it is the case at their Coulomb interaction and most other interactions 30 ), so that it may be described by the following "pair-interaction" Hamiltonian ˆUint =12N∑k,k′=1k≠k′ˆuint (rk,rk′), with the front factor of 1/2 compensating the double-counting of each particle pair. The translation of this Pairoperator to the second-quantization form may be done absolutely similarly to the derivation of Eq. (87), Hamiltonian and gives a similar (though naturally more involved) result ˆUint =12∑j,j′,l,l′ujj′′l′ˆa†jˆa†j′ˆal′ˆal, where the two-particle matrix elements are defined similarly to Eq. (82): uij′′l′≡⟨βjβj′|ˆuint|βlβl′⟩. The only new feature of Eq. (114) is a specific order of the indices of the creation operators. Note the mnemonic rule of writing this expression, similar to that for Eq. (87): each term corresponds to moving a pair of particles from states l and l′ to states j ’ and j (in this order!) factored with the corresponding two-particle matrix element (115).
However, with the account of such term, the resulting Heisenberg equations of the time evolution of the creation/annihilation operators are nonlinear, so that solving them and calculating observables from the results is usually impossible, at least analytically. The only case when some general results may be obtained is the weak interaction limit. In this case, the unperturbed Hamiltonian contains only single-particle terms such as (79), and we can always (at least conceptually :-) find such a basis of orthonormal single-particle states βj in which that Hamiltonian is diagonal in the Dirac representation:
ˆH(0)=∑jε(0)jˆa†jˆaj. Now we can use Eq. (6.14), in this basis, to calculate the interaction energy as a first-order perturbation: E(1)int=⟨N1,N2,…|ˆUint|N1,N2,…⟩=12⟨N1,N2,…|∑j,j′,l,l′uj′′l′ˆa†jˆa†j′ˆal′ˆal|N1,N2,…⟩≡12∑j,j′,l,l′uijll′⟨N1,N2,…|ˆa†jˆa†j′ˆal′ˆal|N1,N2,…⟩. Since, according to Eq. (63), the Dirac states with different occupancies are orthogonal, the last long bracket is different from zero only for three particular subsets of its indices:
(i) j≠j′,l=j, and l′=j′. In this case, the four-operator product in Eq. (117) is equal to ˆa†jˆa†j′ˆaj′ˆaj, and applying the commutation rules twice, we can bring it to the so-called normal ordering, with each creation operator standing to the right of the corresponding annihilation operator, thus forming the particle number operator (72) : ˆa†jˆa†j′ˆaj′ˆaj=±ˆa†jˆa†j′ˆajˆaj′=±ˆa†j(±ˆajˆa†j′)ˆaj′=ˆa†jˆajˆa†j′ˆaj′=ˆNjˆNj′, with a similar sign of the final result for bosons and fermions.
(ii) j≠j′,l=j′, and l′=j. In this case, the four-operator product is equal to ˆa†jˆa†j′ˆajˆaj′, and bringing it to the form ˆNjˆNj′, requires only one commutation: ˆa†jˆa†j′ˆajˆaj′=ˆa†j(±ˆajˆa†j′)ˆaj′=±ˆa†jˆajˆa†j′ˆaj′=±ˆNjˆNj′, with the upper sign for bosons and the lower sign for fermions.
(iii) All indices are equal to each other, giving ˆa†jˆa†j′ˆal′ˆal=ˆa†jˆa†jˆajˆaj. For fermions, such an operator (that "tries" to create or to kill two particles in a row, in the same state) immediately gives the null-vector. In the case of bosons, we may use Eq. (74) to commute the internal pair of operators, getting ˆa†jˆa†jˆajˆaj=ˆa†j(ˆajˆa†j−ˆI)ˆaj=ˆNj(ˆNj−ˆI). Note, however, that this expression formally covers the fermion case as well (always giving zero). As a result, Eq. (117) may be rewritten in the following universal form: E(1)int=12∑j,j′j≠j′NjNj′(ujjij′′′±ujjj′j)+12∑jNj(Nj−1)ujiji The corollaries of this important result are very different for bosons and fermions. In the former case, the last term usually dominates, because the matrix elements (115) are typically the largest when all basis functions coincide. Note that this term allows a very simple interpretation: the number of the diagonal matrix elements it sums up for each state (j) is just the number of interacting particle pairs residing in that state. In contrast, for fermions the last term is zero, and the interaction energy is proportional to the difference of the two terms inside the first parentheses. To spell them out, let us consider the case when there is no direct spin-orbit interaction. Then the vectors |β⟩j of the single-particle state basis may be represented as direct products |oj⟩⊗|mj⟩ of their orbital and spin-orientation parts. (Here, for the brevity of notation, I am using m instead of ms.) For spin- 1/2 particles, including electrons, mj may equal only either +1/2 or −1/2; in this case the spin part of the first matrix element, proportional to ujjjjj, equals ⟨m|⊗⟨m′‖m⟩⊗|m′⟩, where, as in the general Eq. (115), the position of a particular state vector in each direct product is encoding the particle’s number. Since the spins of different particles are defined in different Hilbert spaces, we may move their state vectors around to get ⟨m|⊗⟨m′‖m⟩⊗|m′⟩=(⟨m∣m⟩)1×(⟨m′∣m′⟩)2=1, for any pair of j and j ’. On the other hand, the second matrix element, ujj′j′j, is factored by ⟨m|⊗⟨m′‖m′⟩⊗|m⟩=(⟨m∣m′⟩)1×(⟨m′∣m⟩)2=δmm′. In this case, it is convenient to rewrite Eq. (121) in the coordinate representation, using singleparticle wavefunctions called spin-orbitals ψj(r)≡⟨r∣βj⟩=(⟨r∣o⟩⊗|m⟩)j. They differ from the spatial parts of the usual orbital wavefunctions of the type (4.233) only in that their index j should be understood as the set of the orbital-state and the spin-orientation indices. 31 Also, due to the Pauli-principle restriction of numbers Nj to either 0 or 1, Eq. (121) may be also rewritten without the explicit occupancy numbers, with the understanding that the summation is extended only over the pairs of occupied states. As a result, it becomes E(l)int=12∑j,j′j≠j′∫d3r∫d3r′[ψ∗j(r)ψ∗j′(r′)uint(r,r′)ψj(r)ψj′(r′)−ψ∗j(r)ψ∗j′(r′)uint(r,r′)ψj′(r)ψj(r′)] In particular, for a system of two electrons, we may limit the summation to just two states (j,j′= 1,2). As a result, we return to Eqs. (39)-(41), with the bottom (minus) sign in Eq. (39), corresponding to the triplet spin states. Hence, Eq. (126) may be considered as the generalization of the direct and exchange interaction balance picture to an arbitrary number of orbitals and an arbitrary total number N of electrons. Note, however, that this formula cannot correctly describe the energy of the singlet spin states, corresponding to the plus sign in Eq. (39), and also of the entangled triplet states. 32 The reason is that the description of entangled spin states, given in particular by Eqs. (18) and (20), requires linear superpositions of different Dirac states. (A proof of this fact is left for the reader’s exercise.)
Now comes a very important fact: the approximate result (126), added to the sum of unperturbed energies ε(0)j, equals the sum of exact eigenenergies of the so-called Hartree-Fock equation: 33 (−ℏ22m∇2+u(r))ψj(r)+∑j′≠j∫[ψ∗j′(r′)uint(r,r′)ψj(r)ψj′(r′)−ψ∗j′(r′)uint(r,r′)ψj′(r)ψj(r)]d3r′=εjψj(r), where u(r) is the external-field potential acting on each particle separately - see the second of Eqs. (79). An advantage of this equation in comparison with Eq. (126) is that it allows the (approximate) calculation of not only the energy spectrum of the system, but also the corresponding spin-orbitals, taking into account their electron-electron interaction. Of course, Eq. (127) is an integro-differential rather than just a differential equation. There are, however, efficient methods of numerical solution of such equations, typically based on iterative methods. One more important practical trick is the exclusion of the filled internal electron shells (see Sec. 3.7) from the explicit calculations, because the shell states are virtually unperturbed by the valence electron effects involved in typical atomic phenomena and chemical reactions. In this approach, the Coulomb field of the shells, described by fixed, pre-calculated, and tabulated pseudo-potentials, is added to that of the nuclei. This approach dramatically cuts the computing resources necessary for systems of relatively heavy atoms, enabling a pretty accurate simulation of electronic and chemical properties of rather complex molecules, with thousands of electrons. 34 As a result, the Hartree-Fock approximation has become the de-facto baseline of all socalled ab-initio ("first-principle") calculations in the very important field of quantum chemistry. 35
In departures from this baseline, there are two opposite trends. For larger accuracy (and typically smaller systems), several "post-Hartree-Fock methods", notably including the configuration interaction method, 36 that are more complex but may provide higher accuracy, have been developed.
There is also a strong opposite trend of extending such ab initio ("first-principle") methods to larger systems while sacrificing some of the results’ accuracy and reliability. The ultimate limit of this trend is applicable when the single-particle wavefunction overlaps are small and hence the exchange interaction is negligible. In this limit, the last term in the square brackets in Eq. (127) may be ignored, and the multiplier ψi(r) taken out of the integral, which is thus reduced to a differential equation formally just the Schrödinger equation for a single particle in the following self-consistent effective potential: uef(r)=u(r)+udir(r),udir(r)=∑j′≠j∫ψ∗j′(r′)uint(r,r′)ψj′(r′)d3r′ This is the so-called Hartree approximation - that gives reasonable results for some systems, 37 especially those with low electron density.
However, in dense electron systems (such as typical atoms, molecules, and condensed matter) the exchange interaction, described by the second term in the square brackets of Eqs. (126)-(127), may be as high as ∼30% of the direct interaction, and frequently cannot be ignored. The tendency of taking this interaction in the simplest possible form is currently dominated by the Density Functional Theory, 38 universally known by its acronym DFT. In this approach, the equation solved for each eigenfunction ψj(r) is a differential, Schrödinger-like Kohn-Sham equation [−ℏ22m∇2+u(r)+uKSdir(r)+uxc(r)]ψj(r)=εjψj(r), where uKSdir(r)=−eϕ(r),ϕ(r)=14πε0∫d3r′ρ(r′)|r−r′|,ρ(r)=−en(r) and n(r) is the total electron density in a particular point, calculated as n(r)≡∑jψ∗j(r)ψj(r). The most important feature of the Kohn-Sham Hamiltonian is the simplified description of the exchange and correlation effects by the effective exchange-correlation potential uxc(r). This potential is calculated in various approximations, most of them valid only in the limit when the number of electrons in the system is very high. The simplest of them (proposed by Kohn et al. in the 1960s) is the Local Density Approximation (LDA) in which the effective exchange potential at each point is a function only of the electron density n at the same point r, taken from the theory of a uniform gas of free electrons. 39 However, for many tasks of quantum chemistry, the accuracy given by the LDA is insufficient, because inside molecules the density n typically changes very fast. As a result, DFT has become widely accepted in that field only after the introduction, in the 1980 s, of more accurate, though more cumbersome models for uxc(r), notably the so-called Generalized Gradient Approximations (GGAs). Due to its relative simplicity, DFT enables the calculation, with the same computing resources and reasonable precision, some properties of much larger systems than the methods based on the Hartree-Fock theory. As the result, is has become a very popular tool of ab initio calculations. This popularity is enhanced by the availability of several advanced DFT software packages, some of them in the public domain.
Please note, however, that despite this undisputable success, this approach has its problems. From my personal point of view, the most offensive of them is the implicit assumption of the unphysical Coulomb interaction of an electron with itself (by dropping, on the way from Eq. (128) to Eq. (130), the condition j′≠j at the calculation of uKSdir ). As a result, for a reasonable description of some effects, the available DFT packages are either inapplicable at all or require substantial artificial tinkering. 40
Unfortunately, because of lack of time/space, for details I have to refer the interested reader to specialized literature. 41
30 A simple but important example from the condensed matter theory is the so-called Hubbard model, in which particle repulsion limits their number on each of localized sites to either 0 , or 1 , or 2 , with negligible interaction of the particles on different sites - though the next-neighbor sites are still connected by tunneling - as in Eq. (112).
31 The spin-orbitals (125) are also close to spinors (13), besides that the former definition takes into account that the spin s of a single particle is fixed, so that the spin-orbital may be indexed by the spin’s orientation m≡ms only. Also, if an orbital index is used, it should be clearly distinguished from j, i.e. the set of the orbital and spin indices. This is why I believe that the frequently met notation of spin-orbitals as ψj,s(r) may lead to confusion.
32 Indeed, due to the condition j′≠j, and Eq. (124), the calculated negative exchange interaction is limited to electron state pairs with the same spin direction - such as the factorable triplet states (\uparrow\uparrow and ↓↓ ) of a twoelectron system, in which the contribution of Eex, given by Eq. (41), to the total energy is also negative.
33 This equation was suggested in 1929 by Douglas Hartree for the direct interaction and extended to the exchange interaction by Vladimir Fock in 1930. To verify its compliance with Eq. (126), it is sufficient to multiply all terms of Eq. (127) by ψ∗(r), integrate them over all r-space (so that the right-hand side would give εj), and then sum these single-particle energies over all occupied states j.
34 For condensed-matter systems, this and other computational methods are applied to single elementary spatial cells, with a limited number of electrons in them, using cyclic boundary conditions.
35 See, e.g., A. Szabo and N. Ostlund, Modern Quantum Chemistry, Revised ed., Dover, 1996.
36 That method, in particular, allows the calculation of proper linear superpositions of the Dirac states (such as the entangled states for N=2, discussed above) which are missing in the generic Hartree-Fock approach - see, e.g., the just-cited monograph by Szabo and Ostlund.
37 An extreme example of the Hartree approximation is the Thomas-Fermi model of heavy atoms (with Z>>1 ), in which atomic electrons, at each distance r from the nucleus, are treated as an ideal, uniform Fermi gas, with a certain density n(r) corresponding to the local value uef(r), but a global value of their highest full single-particle energy, ε=0, to ensure the equilibrium. (The analysis of this model is left for the reader’s exercise.)
38 It had been developed by Walter Kohn and his associates (notably Pierre Hohenberg) in 1965-66, and eventually (in 1998) was marked with a Nobel Prize in Chemistry for W. Kohn.
39 Just for the reader’s reference: for a uniform, degenerate Fermi-gas of electrons (with the Fermi energy εF >> kBT), the most important, exchange part ux of uxc may be calculated analytically: ux=−(3/4π)e2kF/4πε0, where the Fermi momentum kF=(2meεF)1/2/ℏ is defined by the electron density: n=2(4π/3)k3F/(2π)3≡k3F/3π2.
40 As just a few examples, see N. Simonian et al., J. Appl. Phys. 113, 044504 (2013); M. Medvedev et al., Science 335,49 (2017); A. Hutama et al., J. Phys. Chem. C 121, 14888 (2017).
41 See, e.g., either the monograph by R. Parr and W. Yang, Density-Functional Theory of Atoms and Molecules, Oxford U. Press, 1994, or the later textbook J. A. Steckel and D. Sholl, Density Functional Theory: Practical Introduction, Wiley, 2009. For a popular review and references to more recent work in this still-developing field, see A. Zangwill, Phys. Today 68,34 (July 2015).