Chebyshev polynomials are defined as follows
and like all orthogonal polynomials they obey a three terms recursion relation :
Any function f(x) with x in [-1,1] can be approximated by a truncated series of Chebyshev polynomials :
or alternatively
where the expansion coefficients are given by
The KPM uses (15) and (17) whereas the the FOE relies on (14) and (16).
Let us consider the the density of states of our system, defined as :
where 
are the eigenvalues of 
. The first step
is to define a rescaled Hamiltonian 
whose eigenvalues 
lie in the interval 
:
where
Equation (19) is the operator equivalent of the change of variables
which we will use thru out this paper.
This requires the knowledge of the find the lowest and highest
eigenvalues of ( 
and 
respectively). This can be done very efficiently (actually in an O(N)
way) using the Lanczos method.
One can define a scaled DOS which can be represented by a Chebyshev polynomial expansion
where the 
are called the moments and are defined by
(17) (in which f(x) had been substituted with D(x)). In
practice we only consider a finite number 
of moments in the
series (23). Such an abrupt truncation of (23) would
result in unwanted Gibbs oscillations, in analogy to Fourier series.
The important point made by Silver and co-workers[3] is to consider instead smooth truncations of the form :
where the 
are called Gibbs damping factors. Those factors
can be chosen in such a way as to impose certain properties to
, like for instance smoothness, positiveness, etc. ,
depending on the application. Many choices of Gibbs damping factors
are possible, for a more detail discussion of possible choices we direct
the reader to Ref. [3]. In the context of this work we
have restricted ourselves to the so-called Jackson damping factors
which guarantee that is non-negative [3].
To understand the nature of the approximation introduced
by the KPM one examine the relationship between the approximate DOS in
(24) and the exact DOS. This can more easily be done after a
change of variable to work in 
space where 
and the
mth Chebyshev polynomial is simply 
. In space
the approximated DOS is related to the exact DOS by a convolution :
where 
is the kernel polynomial, which is an
approximation to a Dirac function (normalized on 
)
which has a width proportional to 
. In Figure 2 we have
plotted 
with and without Gibbs damping factors.
Figure 2: The kernel polynomial for
. The dotted line represents the raw undamped
polynomial, which exhibits strong Gibbs oscillations. The solid line
is the kernel polynomial using Jackson damping factors which
strongly suppresses the oscillations and enforces positivity.
The approximate DOS can therefore be seen as a smearing of the exact DOS with a know resolution function.
It is clear that (25) is very general and can be used to
approximate any function f(x), by just replacing D and 
by
f and 
, respectively.
Assuming that we have calculated some moments (how this is done
will be addressed in the next section) we now show how they can be used
to determine the Fermi energy 
(or its rescaled equivalent
) and to calculate the electronic energy.
The number of electrons in the system can easily be found by integrating the DOS :
where 
is the -space equivalent of the
Heaviside function and 
. Equation (27)
can be seen as the definition of the Fermi energy .
The electronic energy of the system is defined as the energy integral over the occupied states of the DOS,
If we substitute the smeared DOS ( 
as defined in
(25) for 
in (27) we get an
approximation for the number of electrons (as a function of 
) :
with 
and 
.
In practice this last expression will be used to determine the
which gives the correct number of electrons :
. In this respect one should point out
the importance of choosing Gibbs damping factors which
guarantee that is non-negative which in turn
defines a unique solution for .
Similarly if we replace in (28) with
we get the smeared-DOS (SD) approximation to the electronic energy :
An alternative expression to calculate the electronic energy is found by taking the kernel polynomial approximation to the Fermi-function in (28) (see equ. 25) :
This yields the smeared-Fermi (SF) approximation to the electronic energy :
It can be shown that the substitution of 
for
in equation (27) yields the same expression
for 
as in (29)[4].
Figure 3: Approximate DOS as a function of energy for a)
64 Si atoms calculated with a 100 moments and b) 64 atoms of Cu. The
sharp vertical lines show the positions of the exact eigenvalues, with
a height proportional to their degeneracy.
To illustrate the method presented in the previous section, we plotted
in Figure 3a) and b) the approximate DOS at the
-point as a function of energy for a 64 Si atoms system
calculated with 100 moments and for 64 Cu atoms with 200 moments
respectively. For comparison we also show the exact eigenvalues. One
can clearly see that the approximated DOS is a smeared version of the
exact DOS.
Figure 4: Comparison of the energy from the KPM
vs the number of moments, calculated using the SD (circles) and the
SF (diamonds) energy expressions to the exact energy (solid line).
a) Electronic energy for 64 atoms of Si, b) Electronic energy for 64
atoms of Cu.
The KPM offers two expressions to compute the electronic energy given
the moments and some Gibbs damping factors . We
designate these as the smeared-DOS (SD) and the smeared-Fermi (SF)
expressions respectively (see section C.3). In Figure
4 we compare the convergence of the two expressions
with respect to the number of moments. Two observations can be
made. First we see that the SF expression converges much faster than
the SD expression. The reasons for this have been discussed in Ref.\
[3]. In light of this we will always use the SF expressions
in our calculations. Second we see that we seem to have 
, and it can be shown that this is in general
the case.
In Figure 5 we show the error in the total energy as a function of the number of moments for different systems, calculated using the KPM.
Figure 5: Convergence of the error in the electronic
energy as a function of the number of moments for different
systems. The error is defined as the difference between the exact
electronic energy (obtained from the eigenvalues of the Hamiltonian)
and the KPM approximation. All the energies where calculated using the
SF expression. We see that the convergence for metals is generally
slower.
The calculation for Si were done by using the Goodwin, Skinner and Pettifor TB model for Si, and the calculations for Al using the Harrison parameterization, both these models are orthogonal. The results for Cu and Mo were obtained using the non-orthogonal tight-binding model developed by Mehl and Papaconstantopoulos [5]. We see that for most systems the electronic energy is converged to within 10meV/atom (which is the generally accepted accuracy of the TB method) for less than 200 moments. We see that the convergence is slower in metallic systems. Although the convergence of the total energy is slow with respect to the number of moments, it should be stressed that the energy differences converge faster. To illustrate this point we have calculated the vacancy formation energy in bcc Mo. In Figure 6 we show how the total energy and the formation energy converge with respect to the number of moments.
Figure 6: Vacancy in bcc Mo. Convergence of energies
as a function of the number of moments. a) Total energy for the perfect
crystal, for the crystal with a defect and vacancy formation energy
(difference of the two others). The dotted lines show the exact results.
b) Energy differences for the same quantities as in a). We see that
the error in the formation energy is much smaller.
From a computational point of view, the most time consuming part of
the KPM is the calculation of the moments . Using equation
(17) the moments are defined by
where we used the fact that the DOS may be written as a trace :
If we introduce the following notation :
where the last equation is a consequence of (13), then one can express the moments as follow :
We see that the computationally intensive part of determining the
moments is the calculation of the 
using (36),
which involves one matrix-matrix-multiplication (MMM) for each value of m.
One very useful feature of the KPM is that to calculate moments
one requires 
MMMs.
