4. Density-Matrix Formulation

Density-functional theory together with the pseudopotential approximation
has established itself as the method of choice for performing large-scale
* ab initio* quantum-mechanical calculations. In particular, the fact
that the kinetic energy operator is diagonal in momentum-space, that
the Hartree and local pseudopotential contributions are straightforward to
calculate in
momentum-space, the development of fast Fourier transforms (FFTs)
to efficiently switch between momentum-space and real-space and their
natural relation to periodic boundary conditions (3.47)
has led the plane-wave basis set and momentum-space formalism to become
the most widely accepted method for performing such calculations
[79,80].

Efficient methods to solve the Kohn-Sham equations have been developed [81,82,83] which iteratively diagonalise the Hamiltonian. All of these plane-wave methods require a computational effort which scales as the cube of the system-size, and an amount of memory which scales as the square. Although these methods have made first-principles quantum-mechanical calculations available as a tool to a wide range of scientists in a variety of disciplines, this scaling ultimately limits the maximum system-size which can be treated now and in the near future, despite the rapid development of computer technology. A method which requires an effort and amount of memory which scales linearly with system-size would push the boundaries back much further, and so it is to the development of such a method that we now turn our attention.

There has been a great deal of success in developing linear-scaling tight-binding methods [84,85,86,87,88], but full density-functional methods have proved far more elusive. The first approach was the ``divide and conquer'' method [89,90,91] in which the large system is partitioned into overlapping subsystems. The original formulation divides the electronic density between these subsystems, solving for each in turn until self-consistency is reached. More recently, the density-matrix has been divided instead of the density[92] and these methods have been applied to large molecular systems [93,94]. The recursion method [95,96] has been used in a linear-scaling scheme to determine the electronic density by calculating diagonal elements of the Green's function [97]. Still other methods focus on the density of states, which could be used to obtain total energies by integrating up to the chemical potential [98,99,100,101].

As already mentioned, attempts to construct approximations for the density-functional itself have led to linear-scaling methods for metals [29,31,32]. Electronic structure methods suitable for metallic systems based upon the multiple scattering method have also been proposed [102,103,104].

The Fermi operator expansion method calculates the density-matrix at finite temperature in terms of a Chebychev expansion [105,106,107,108] or a rational representation [109,106,110,111]. The closely related kernel polynomial method [112,113] expands the zero temperature Fermi distribution by integrating an expansion of the delta function in which damping factors are used to suppress Gibbs oscillations.

One method related to the density-matrix schemes discussed next is based on a formulation of Kohn-Sham theory in terms of non-orthogonal localised orbitals [114]. A new energy functional of these localised orbitals which naturally leads to orthogonal orbitals at its minimum (which is the same as the ground-state minimum of the conventional functional) is introduced [115,116,117,118,119,120,121] and allows a linear-scaling method to be constructed, which has been used to study extended defects in silicon [122]. For a review of these methods and an explanation of their relationship to the density-matrix schemes, see [123].

The method most closely related to the approach introduced in this dissertation is the density-matrix minimisation method. The total energy is minimised with respect to the density-matrix and the purifying transformation (section 4.4.4) is used to impose the idempotency constraint [124,125,126,127,128,129].

In this chapter we introduce the concept of the density-matrix and show
how it can be applied to the system of non-interacting particles in
density-functional theory. We then show how a linear-scaling method
results naturally from the short-ranged nature of the density-matrix in
real-space, and discuss the constraints which must be imposed in order
to find the ground-state.

- 4.1 The density-matrix
- 4.2 Partial occupation of the Kohn-Sham orbitals
- 4.3 Density-matrix DFT
- 4.4 Constraints on the density-matrix
- 4.5 Requirements for linear-scaling methods
- 4.6 Non-orthogonal orbitals