Quantum Chemistry of Solids / 15-Electron Correlations in Molecules and Crystals

Gνpσs (t) =

µ0λr

5.3 Atomic Orbital Laplace-transformed MP2 Theory

for Periodic Systems

177

The smooth decay of the energy contributions in the contraction can, of course, be used to safely truncate the lattice summation. For example, the central cell energy contribution e02r provides a good estimate for the contribution from distant cells. In most cases, the magnitude e02r is reasonably well approximated by e002 /|o − r|6. A system with a unit-cell translation length of 2.5 ˚A can, thus, be expected to have no significant contribution to the correlation energy contraction beyond approximately the fifth or sixth neighboring cell. However, although relatively few cells need to be considered in the energy contraction (rmax=4–6 cells for 1D systems such as polyacetylene, typically), the spatial framework of duplicated cells needs to be large enough so that all significant charge distributions, including the ones centered near the edge, are properly expanded.

5.3.2 Laplace MP2 for Periodic Systems: Bandgap

To correct the HF bandgap at the MP2 level the quasiparticle energy formalism was applied [185, 188]. The quasiparticle (MP2) energies for a given state g is given by (closed-shell case):

MP 2 =

g

+ U (g) + V (g) =

g

+ ∆ MP 2

(5.83)

g

−

−

g

U (g) =

−

(ia|gb) [2(ia|gb) − (ib|ga)]

(5.84)

i,ab

a + b

i

g

−

−

V (g) =

(ia|jb) [2(ia|jb) − (ig|ja)]

(5.85)

ij,a

a + g

i

j

where i and j are occupied, and a and b are unoccupied HF orbitals of the system. Consequently, i, j and a, b are occupied and unoccupied HF orbital energies, respectively. In a finite-basis calculation, the number of virtual orbitals is finite and determined by the basis set. The Laplace-transform of the energy denominators is applied to Eqs.(5.84) and (5.85) and using an AO-basis formulation one obtains

M P 2	=	∞			νpσs	νpσs	(5.86)
∆ g			dt		Gµ0λr	(t)vµ0λr

0µ0λrνpσs

where

νpσs	= 2(µ0λr|νpσs) − (µ0σs|νpλr)
vµ0λr
(µ0λr|νpσs) = dr1dr2|r1 − r2|−1µ0(r1)λr(r1)|νp(r2)σs(r2)		(5.87)

are two-electron integrals of AOs, for orbitals centered on cells with real-space coordinates 0, r, p, and s, respectively. All real-space integrations are performed over all space. ∆ Mg P 2 is the MP2 correction to the HF eigenvalue as defined in (5.83). The other required quantity in (5.86) is

Xµt 0γt Yνtpδu (γtδu|κvεw)

γtδuκvεw
× Wλtrκv Yσtsεw + Xλt rκv Zσt sεw		(5.88)

178 5 Electron Correlations in Molecules and Crystals

The elements of the X and Y matrices are given in (5.80). The elements of the W and Z matrices are defined by

Wµtpγs = C(k†)µg exp[− g(k†)t]C(k†)γg exp[k†(p − s)]	(5.89)
Zµt pγs = −C(k†)µg exp[ g (k†)t]C(k†)γg exp[k†(p − s)]	(5.90)

where for evaluating the correction to the bandgap, W and Z are constructed from the highest-occupied crystal orbitals (HOCO) and lowest-unoccupied crystal orbitals (LUCO) at point k at which the bandgap is minimum. Note that adding a constant value (e.g., the Fermi energy) to all occupied orbital energies, and subtracting the same value from all virtual orbital energies does not a ect the denominators of (5.84) and (5.85) nor the MP2 corrections, but they guarantee that all the integrals containing terms such as exp[− i(k)t] and exp[− a(k)t] are finite. Depending on whether g is the HOCO or the LUCO, U (g(k)) + V (g(k))(g(k)) corresponds to the MP2 correction to the ionization potential or the electron a nity. The subtraction of these two corrections, in the spirit of Koopman’s theorem, yields the correction to the bandgap. Indirect bandgaps are obtained in a similar fashion by considering two di erent k points.

Finally, (5.86) is approximated by numerical Gauss–Legendre quadrature after a logarithmic transform of the t variable

MP 2	≈ ωn		νpσs	νpσs	(5.91)
∆ g			Gµ0λr	(tn)vµ0λr

nµ0λrνpσs

where tn and ωn are the quadrature points and weights, respectively. Usually 3–5 points are enough to evaluate the gap corrections, although more are necessary for the energy [188].

The Laplace MP2 method is applied in [185] for the calculation of the correlation corrections to the bandgap in 1D and 2D systems. Let us consider the results for a 2D system – hexagonal boron nitride in a one-layer model. Like graphite, BN can exist in the form of nanotubes. Because these nanotubes can be viewed as hexagonal BN sheets rolled onto themselves, hexagonal boron nitride has been the object of renewed interest. The main absorption of a thin film of hexagonal BN is found at 6.2 eV and a sharp fall at about 5.8 eV was attributed to the direct bandgap on the basis of semiempirical calculations [190]. The quasiparticle band energy structure at both the Hartree–Fock and the correlated level was calculated in [191] using a method very similar to MP2, the second-order many-body Green-function approach. Hexagonal BN was shown to have a rather complex band structure. Points P and Q are two high-symmetry points within a Brillouin zone with hexagonal symmetry – the hexagon corner and center of edge, respectively. The HF bandgap is minimum at k-point P and was evaluated to be 12.46 eV and the quasiparticle gap to be 2.88 eV. At point Q, the HF gap is 16.45 eV and correlation was said to decrease it to 6.48 eV.

PBC MP2 calculations [185] have been unable to corroborate these results. In these calculations, the HF band structure was shown to have a minimum bandgap of 13.73 eV at a point that is neither P nor Q, but rather in between, at a point of higher symmetry still. After examination of the phases of the Bloch orbitals at these points, one can say that point P is best characterized as corresponding to a trans alternation,

5.3 Atomic Orbital Laplace-transformed MP2 Theory

for Periodic Systems

179

whereas Q corresponds to a cis alternation. Laplace MP2 calculation, thus, seems to indicate first, that the Hartree–Fock method favors a form that is a mixture of these two distinct alternations, and second, that for a meaningful calculation, great care should be taken in the convergence of the HF density by using an appropriate k-mesh and spatial framework. The correlation correction to the bandgap is very important, as suggested by the MP2 studies. Using 5 Laplace quadrature points and 49 cells to carry out the final contraction, the MP2 gap correction is –6.25 eV. The direct MP2 bandgap can be thus estimated to be at most 7.5 eV, approximately 1.5 eV above the experimental absorption peak. This constitutes a rather good agreement given that much larger basis sets would be needed to properly describe the conduction band that is thought to have a very large extent into the vacuum region.

The Laplace MP2 algorithm of a correlated band-structure calculation is discussed in [188] and applied to trans-polyacetilene.

The post-HF methods for periodic systems were discussed in Sections 5.2 and 5.3. The incremental technique and the Laplace MP2 approach can be mentioned as the most successful approaches to the correlation in solids [5, 192]. The former uses standard molecular codes with clusters of various shapes and sizes to calculate the correlation energy, the latter adopts a formulation of correlation directly in AO basis.

The incremental scheme can be realized using the di erent post-Hartree–Fock methods (MP2, coupled-cluster CCSD(T)). This scheme is applicable both for nonconducting crystals and metals. The use of finite clusters allows modern sophisticated molecular codes to be applied (like MOLPRO code [116]). The incremental scheme can be used not only for calculation of the energy but also for the energy gradients necessary in the atomic-structure optimization. Unfortunately, the application of the incremental scheme requires some “hand” work when choosing the cluster size and shape for any new system. To the best of our knowledge the code WANNIER (for crystals) [158] allows consideration only of the crystals with a small number of atoms in the unit cell. If clusters are used instead of crystals – the finite nature of the system under consideration a ects the results. It is not so easy to set a reasonable cluster to approximate a crystal, especially when there are many atoms in the unit cell. We do not think that these disadvantages are intrinsic features of the incremental scheme, but are mainly connected with the necessity of some standard code generation.

The essential advantage of the Laplace MP2 method is its incorporation in publicly available Gaussian03 [107] computer code with PBC introduction for HF calculations. Due to the Laplace-transform the dense grids of k-points can be a orded, the virtual space is not truncated so the method gives the true MP2 energy. The Laplace MP2 band structure can be calculated, including the correlation correction to the bandgap. As no results of its application to 3D systems were published it is di cult to judge how reliable this method is for three-dimensional crystals.

In the next section we consider a local MP2 scheme for periodic systems that is now incorporated in more or less “standard” CRYSCOR correlation computer code for nonconducting crystals [109].

180 5 Electron Correlations in Molecules and Crystals

5.4 Local MP2 Electron-correlation Method for Nonconducting Crystals

5.4.1 Local MP2 Equations for Periodic Systems

The local MP2 electron-correlation method for nonconducting crystals [109] is an extension to crystalline solids of the local correlation MP2 method for molecules (see Sect. 5.1.5), starting from a local representation of the occupied and virtual HF subspaces. The localized HF crystalline orbitals of the occupied states are provided in the LCAO approximation by the CRYSTAL program [23] and based on a Boys localization criterion. The localization technique was considered in Sect. 3.3.3. The label im of the occupied localized Wannier functions (LWF) Wim = Wi(r − Rm) includes the type of LWF and translation vector Rm, indicating the primitive unit cell, in which the LWF is centered (m = 0 for the reference cell). The index i runs from 1 to Nb, the number of filled electron bands used for the localization procedure; the correlation calculation is restricted usually to valence bands LWFs. The latter are expressed as a linear combination of the Gaussian-type atomic orbitals (AOs) χµ(r − Rn) = χµn numbered by index µ = 1, . . . , M (M is the number of AOs in the reference cell) and the cell n translation vector Rn:

Wim(r) = lµn,imχµn(r)	(5.92)

µn

AOs constitute the nonorthogonal basis set used for the solution of the HF crystalline equations. To exploit the translation invariance property of LWFs the matrix Liµn = lµn,i0 has been introduced.

As in the molecular case (Sect. 5.1.5) for the virtual space the projected AOs

ˆ		occ	|l l| onto the occupied
(PAOs) are used. To generate PAOs the projector P =		l
HF space is expressed in terms of the density matrix	in an AO representation, P ,

which is calculated very accurately by the CRYSTAL program via Brillouin-zone integration over all canonical occupied crystalline orbitals |l :

ˆ	|ν Pµν µ|		(5.93)
P =
	µν
	ˆ	ˆ	ˆ
The projector onto the virtual space is hence obtained: Q = 1		− P . By acting with Q

on any one of the local functions χµn ≡ χβ of the original AO set, a Q-projected AO is obtained, or simply PAO:

ˆ	(5.94)
Qχβ (r) = Qνµχµn = χ˜β (r)

µn

In what follows, PAOs are identified with the Latin letter corresponding to the Greek letter labeling the AOs from which they are generated, χ˜β (r) ↔ |b ≡ |bn .

The set of all PAOs so generated (the standard PAO set) constitutes a nonorthogonal, linearly dependent, incomplete set of local functions, strictly orthogonal to all WFs. The problem of the adequacy of the standard set to describe excitations is a delicate issue, because the basis sets usually adopted to solve the HF problem for

5.4 Local MP2 Electron-correlation Method for Nonconducting Crystals

181

crystals are calibrated so as to describe accurately only the occupied subspace. Indications from molecular experience could not be easily transferable to the periodic case, especially in the case of dense crystals.

An essential feature of O(N ) local techniques is the introduction of some metric expressing the “distance“ between local functions, and the approximate treatment or total neglect of terms entering the exact equations, according to the distance between the functions on which they depend. The wealth of experience gained in molecular studies can be exploited to this end.

A preliminary step is the definition of WF domains. According to the standard convention, the domain of the ith WF, Di ≡ {χα}i is the set of all AOs associated with the smallest and most compact set of atoms such that Di “spans” Wi, it is clearly su cient to define this set for WFs in the reference cell (im = i0). The criterion adopted to select the spanning set is based on the one proposed by Boughton and Pulay [193]. Briefly, one includes in the set, step by step, all AOs belonging to a group of symmetry-equivalent atoms (that is, atoms that are carried into each other by a symmetry operation that reproduces the LWF into itself). At each step, the norm N BP of the di erence between Wi and the best possible linear combination of AOs in the set is determined. The process is stopped when N BP is below a certain threshold T BP , for instance, 0.02; T BP will be called the Boughton-Pulay (BP) parameter. The sequence according to which atoms are introduced in the domain is as follows. First,

a “kernel” is selected including atoms A that contribute a Mulliken population QA
		i
larger than a certain threshold, e.g., 0.2 (QiA = |	α	Pi,βn Pi,αn (νn|µm)|; β ≡ νm);
β A

next, groups of symmetry-related atoms are introduced in order of increasing distance from the kernel (if the distance is the same, atoms with larger QAi are included first).

After defining the domains, it is expedient to introduce, for any two WFs, i and j, a measure dij of the distance between them: the standard definition was adopted, namely, the minimum distance between any two atoms belonging to the respective kernels. A classification of LWF pairs by distance is then possible, based on the idea (strongly supported by the results of molecular calculations) that excitations from two local occupied orbitals become progressively less important with increasing distance d between the two centers. Strong (d = 0), weak (0 < d ≤ d1), distant (d1 < d < d2) and very distant (d2 < d) pairs can be so termed, and the respective excitations are treated in a progressively less accurate way; in molecular calculations a standard choice is d1 = 8 Bohr, d2 = 15 Bohr (see Sect. 5.1.5). The total neglect of contributions from very distant pairs is crucial to the realization of O(N ) computational schemes. In the following, two WFs will be said to be “close-by” each other, if they do not form a very distant pair. For each WF pair, a pair-excitation domain or simply, “pair domain”, can be defined, which is the union of the respective domains (Di,j Di Dj ). Given the biunivocal correspondence between AOs and PAOs (see (5.94), χ˜β will be said to “belong to” Dij (briefly, |b Dij ), if χβ does. Excitations from i, j are usually restricted to PAOs belonging to Dij .

The use of distance criteria in the definition of the shape of the LWF domains and for controlling the level of the approximation is certainly not the only possibility, even if it is the simplest one; for instance, some “topological” index might be used instead, as is already possible in the MOLPRO program [116].

182 5 Electron Correlations in Molecules and Crystals

Definitions and conventions introduced are not specific of any given level of local correlation theory and can be used both in CC and MP theory. In [109] the MP2, the simplest local correlation method, was considered.

The first-order perturbative correction to the HF wavefunction (5.37) is written as a combination of contravariant doubly excited configurations

		|Ψ (1) =		|Ψ (ij)	(5.95)
			ij
where				Tabij |Φabij
		|Ψ (ij) =			(5.96)
			ab
The exact equations that determine the excitation amplitudes T ij					are
				ab
			FacTcdij Sdb + SacTcdij Fdb
	0 = Rabij	= Kabij + cd
			FilTcdlj + Tcdil Flj # Sdb
		−Sac " l			(5.97)
Here, Kij	is the two-electron “exchange” integral,
ab
		Kabij = (ia|jb)			(5.98)

Fil and Fdb are elements of the WF or PAO representation of the Fock operator, respectively, and Sab is the overlap between PAOs. Equation (5.97) must be solved iteratively in the unknown amplitudes until all residuals Rabij are zero. For any given set of amplitudes, the Hylleraas functional can be calculated:

(
E2 = E2ij =	(Kabij + Rabij )(2Tabij − Tbaij )	(5.99)
ij)	ij ab

It is easily shown that E2 ≥ E(2) (the MP2 estimate of the correlation energy), the equality sign holding true only at convergence, when all residuals are zero. In the above equations, the locality assumption is introduced according to the vicinity criteria, by restricting the sums over (ij) to close-by WF pairs, and those over (ab) to PAOs Dij .

For the case of a crystalline application, the fact that all two-index quantities (Fock matrix in the PAO and in the WF representation, overlap matrix in the PAO representation) and all four-index quantities (amplitudes, residuals, two-electron integrals) are translationally invariant was used. For instance,

T im,j(n+m)	= T	0,jn	(5.100)
an +m,bm +n		an ,bm

Furthermore, the total correction energy can be written as a sum of correction energies per cell: E2 = LE2cell, where L is the total number of cells in the crystal with periodic boundary conditions, and the expression for E2cell is the same as in (5.99), but with the sum restricted to “zero-cell” WF pairs, that is, with the first WF in the reference zero cell:

5.4 Local MP2 Electron-correlation Method for Nonconducting Crystals		183
Ecell =
	Ei0,jn
2	2

ijn

So, finally, it is confined to a local problem, in the vicinity of the zero cell. The computation becomes O(0) with respect to the size L of the macroscopic crystal and (asymptotically) O(N ) with respect to the size N of its unit cell.

In order to reformulate (5.97) for a given zero-cell WF pair (i0, jn) in the periodic case, the following notation is introduced,

Rab = Kab + Aab + Bab

(5.101)

where a, b ≡ |am , |bn and the apex i0, jn is implicit in all four terms. This allows an “internal” contribution, Aab (second and third terms in (5.97), which depends on the amplitudes referring to the same i0, jn pair, and an “external” Bab contribution (last two terms in (5.97)), which depends on amplitudes involving other WF pairs, to be distinguished.

The most computationally demanding step in the implementation of the local MP2 method is the evaluation of the two-electron integrals K12 = (i1a1|i2a2). The straightforward solution is to express them as a linear combination of analytical KA integrals over AOs, with coe cients depending on the four local functions. Other techniques are possible, however. In [157] the technique of fitted LWF for local electron-correlation methods was suggested. A very promising one is the so-called density-fitting approximation [146]. Exchange integral K is obtained in terms of two-center repulsion integrals JF1F2 between two AO-like functions F1(r) and F2(r), after expressing each product distribution as a linear combination of a suitable set of “fitting functions” Fn(r).

In the problem of the four-index transformation from {KA} to {K} advantage can be taken, as is customary in molecular calculations, of the fact that the same KA is needed for many K integrals. When the involved LWF are su ciently distant the approximate multipolar evaluation of K integrals is used. For an e cient use of this technique, a number of parameters must be accurately calibrated: the distance between WFs beyond which the multipolar technique substitutes the exact one; the value of the parameter defining the truncation of the AO expansion of WFs and PAOs, according to the value of the coe cient or of the distance of the AO from the LWF center. The last two aspects are interrelated: almost negligible tails at very large distance may give enormous high-order multipoles, leading to unphysical values of the interaction.

The possibility to eliminate the WFs’ tails may have a significant impact on the e ciency of LMP2 methods. For these methods the most expensive step is the evaluation and transformation of four-index two-electron integrals, where the computational cost is governed by the spatial extent of the AO support of the individual WFs. The details of the evaluation of these integrals can be found in [109].

To improve the e ciency of the local correlation methods in nonconducting crystals the technique for constructing compact (fitted) WFs in AO basis sets was suggested [157]. Nevertheless, the results are general and can also be applied to other cases where the LCAO expansion of WFs is a critical issue.

184 5 Electron Correlations in Molecules and Crystals

5.4.2 Fitted Wannier Functions for Periodic Local Correlation Methods

Although WFs are spatially localized, the localization region can be relatively large. As an illustration one can consider the case of diamond [157]. In the LCAO approximation, the WFs are represented by the AO expansion coe cients. The coe cients of the four translationally inequivalent WFs (which correspond to the four upper valence bands of the HF reference) larger (in absolute value) than a threshold of 10−4 comprise the atomic functions centered in as much as 435(!) unit cells. Such a wide extent of WFs makes accurate calculations, in particular local MP2 calculations, for solids, computationally di cult.

The slow decay of WFs is often related to their mutual orthogonality, which gives rise to long-range oscillating tails, containing no physically relevant information. For quite sometime now, the possibility to avoid these tails by sacrificing orthogonality has been discussed. Many important results concerning the localization properties of WFs and nonorthogonal localized orbitals (NOLO) (following most of the authors, we do not use the term “Wannier functions” for nonorthogonal localized functions) were obtained for the one-dimensional (1D) case, [194–197]. It was shown that in the 1D case the spread of the NOLOs can be lower than that of the WFs and the term “ultralocalized” was introduced for NOLO [195]. The LCAO expansion of the WF Wi0(r) is defined as

Wi(r) =	CµR,iχµ(r − R) ≡ Cµiχµ(r)	(5.102)
µR	µ

where χµ(r − R) represents the atomic function centered in the Rth unit cell, CµR,i is the corresponding AO-WF coe cient and index µ combines the indices µ and R. Equation (5.102) defines WFs for the reference unit cell. In [157] it was suggested to maximize the variational functional that defines the L¨owdin population of the given localized function on a set of atoms surrounding the center of the localization:

no	(no−Lw)		2	(no)	ˆ(LW )	(no)		(5.103)
Ii	= i Si(loc) Cµi			= Wi	|PSi(loc)	|Wi

where Cµi(no−Lw) are the expansion coe cients of the NOLO Wi(no)(r) over the L¨owdinorthogonolized atomic basis set functions φ(µLW )(r)

	Wi(no)(r) =	µ	Cµi(no−Lw)φµ(LW )(r)				(5.104)
	φµ(LW )(r) =	µ	Sµ(AOµ ) −		1/2		(5.105)
					φµ (r)
			ˆ(LW )
(AO)					is an operator projecting onto the
Here Sµ µ	is the AO overlap matrix, and P (loc)
			Si
space of the L¨owdin-orthogonalized AOs from the given set of atoms:
	ˆ(LW )		(LW )		(LW )	|
	PSi(loc) =		|φµ	φµ			(5.106)
		S(loc)
	µ i

5.4 Local MP2 Electron-correlation Method for Nonconducting Crystals

185

The set of atomic orbitals Si(loc) determines the region (or the set of atoms), the corresponding localized functions Wi(no)(r) should concentrate within.

Maximization of the functional (5.103) with a constraint of normalization of the functions Wi(no)(r) can be carried out by diagonalizing the matrix of the operator

ˆ(LW )

PSi(loc) (5.106). The eigenvector of this matrix, corresponding to the highest eigen-

value, holds the coe cients of the expansion (5.104), providing the most localized NOLO in the sense of the functional (5.103). The basis for this matrix could be formed by the Bloch functions of the occupied bands, WFs or even atomic orbitals, projected onto the occupied space.

During generation of the nonorthogonal functions one should not lose the rank of the space spanned by these functions (i.e. the linear independency of the functions should be kept). For the crystals of su ciently high symmetry the rank can be preserved automatically, if the basis used for the variational procedure is symmetrized [198]. The localized functions are basis functions of irreducible representations (irreps) of the site groups related to the centroids of these functions [13]. Thus, if the basis for the functional (5.103) is symmetrized according to these site group and irrep, the centering site of the localized function cannot change [13,198]. So, when the localized functions are centered on an isolated symmetry point, the symmetrization of the basis keeps the functions linearly independent. Moreover, even if the centroid of the localized function is not an isolated symmetry point (lies on the symmetry axis or in the symmetry plane), but the symmetry of the crystal is still high enough, the rank still may be kept by just taking symmetry into account. Yet, in this case the latter is not guaranteed. If preservation of rank cannot be assured by applying symmetry, one may have to apply some additional constraints – e.g., as suggested in [199], the constraint of the nonsingularity for the localized functions overlap matrix, or in [200], the constraint of freezing the centroids of the localized functions. However, the latter constraint might collide with the requirement of higher localization, if the centers of the functions are not fixed by the symmetry.

The variational procedure should be performed in the space of the first vectors of the corresponding irreps. In other words, only symmetry-unique localized functions are to be obtained variationally. The remaining functions can be generated by applying the appropriate symmetry operators [13, 198]. If there is more than one independent set of localized functions centered on the same site and transforming according to the same irrep, one can obtain all of them in the variational procedure either by taking the eigenfunctions of the diagonalized matrix corresponding to several of the highest eigenvalues instead of just one, or by performing independent procedures for their generation. In the latter case, the variational procedure should be run in the WF-basis, and only one out of several WFs with the considered symmetry, which are centered on the same site, should be included in this basis. The symmetry properties of the WFs and localized functions can be determined in high-symmetry cases without the actual construction of WFs, by analyzing the symmetry of the Bloch functions corresponding to the valence bands. The symmetry information on the NOLOs (the site groups of the WF centroids and the irreps of the site groups) can be extracted using the induced representation theory. Let the symmetry of the WFs be defined unumbiguosly by the symmetry of Bloch states forming the occupied crystalline orbitals. For example, as was seen in Chap. 3, the four upper valence-band

186 5 Electron Correlations in Molecules and Crystals

Bloch states in diamond-type crystals are connected with four WFs centered at the center of the atom–atom bond with the site-symmetry group D3d and transforming over the one-dimensional a1g -irrep of this group. The coset representatives in the decomposition of space group Oh7 in cosets over subgroup D3d transform this WF into another three equivalent WFs directed along the interatomic bonds. In MgO crystal four upper valence bands, nondegenerated and triple degenerated, correspond to four WFs centered at the oxygen site (site symmetry group Oh coincides with the point group of a crystal) and transforming over irreps a1g and t1u of Oh group. This is the case when NOLOs di erent by symmetry are centered at the same atom.

In low-symmetry cases one needs to know also the centering points of the initial orthogonal symmetry-adapted WFs. Let instead of the diamond structure one consider ZnS (sphalerite) structure. The centering of four WFs corresponding to four upper valence bands, is defined by symmetry as C3v , i.e. WFs centers are at the atom–atom bond. The actual centering can not be found only from symmetry considerations. In this case, the additional information is necessary to fix the WFs centers. One possibility – to use the nonsymmetry-adapted WFs generated by CRYSTAL code [23] and a posteriori to apply the technique of Casassa et al. [201] to nonsymmetrized WFs. This technique is considered in the next section. If the WFs centers are not fixed the method under consideration can give highly localized, but linearly dependent functions, which of course are of no use.

In some special low-symmetry cases, when there are several independent WFs of the same symmetry (like s- and pz -WFs of oxygen in rutile structure TiO2), the method also requires the orthogonal WFs themselves to use them as a basis for construction NOLOs. In other cases one can also use Bloch functions as such a basis.

Let the symmetry-adapted localized orthogonal WF be constructed and used furthermore in the post-Hartree–Fock calculations. In numerical local MP2 calculations the sum (5.102) has to be truncated according to some threshold. It is possible to make a truncation based on the magnitude of the calculated WFs coe cients in this sum. In e ect, all contributions due to WF coe cients below some threshold could be disregarded. Unfortunately, only relatively loose thresholds can presently be a orded in practical calculations [109] and the numbers of translation vectors R taken into account may become a critical parameter significantly influencing the accuracy of the calculation. For example, the four-index transformation of the two-electron integrals in the CRYSCOR computer code [117] uses tolerances for the LCAO coe cients of the WFs of 10−2–10−3. Such loose thresholds can hardly guarantee su cient accuracy for reliable results, yet tighter thresholds render the calculations too expensive.

In order to treat the truncation of the WFs’ tails in an improved way the fitting of WFs was suggested [157]. The localization criterion is introduced

	Ii(f it) =		Wi(r) − Wi(r)		2	(5.107)
					dr
		!
where the function W	(r) is the truncated (fitted) localized function					(5.108)
!i	Wi(r) =		Cµiφµ(r)
	!

µ Si(f it)

The sum (5.108) di ers from the infinite sum (5.102) because the summation is made over Si(f it) containing an incomplete set of atomic functions centered only on atoms,