It is the intend of this section to introduce the subject to a point that the reader can understand the origin of the algorithms. We will be overly simplifying, schematic and sometimes strictly incorrect in the presentation to avoid cluttering details.
The induced magnetic 
is related to the external,
applied magnetic field 
by
where 
is the magnetic susceptibility tensor, which depends on
the position of the induced field r, and the position 
of the external field. Its spatial variation means that different
nuclei will feel a different induced magnetic field, and hence exhibit
a characterstic ``chemical shift'' when probed with NMR
techniques. Thus, NMR is widely used to identify substances, but is
also useful for structural analysis.
Our goal is to compute the tensor 
, but we do so in Fourier
space, since we will be dealing with infinite periodic systems such as
crystals. If they are not periodic (such as molecules and amorphous
systems), we can always make them periodic by a replication at
sufficiently large distances. Therefore, we work mostly in
Fourier-space. For ease of exposition, we will stay in 
-space
here. By definition, 
is computed as the second derivative of
the (quantum-mechanical) energy E as a function of the magnetic
field:
The energy E depends on the magnetic field via the equations of quantum mechanics, where we compute the second derivative via second order perturbation theory:
Here, the 
are the valence band wave
functions, i.e. those quantum-mechanical states occupied by
electrons, and 
are the conduction band
wave functions, which are unoccupied. The effect of the magnetic field
is contained in the ``perturbation'' operator 
. The wave
functions are solutions to a Schrödinger-type eigenvalue equation
for the hermitian Hamiltonian operator H:
which has many solutions, ascendingly ordered by their (energy)
eigenvalues 
. When expanding the wave functions
in a basis set of size N, (4) turns
into a large matrix equation for which the 
are
eigenvectors. We will henceforth call the wave functions ``vectors'',
because they are represented by a vector of expansion coefficients
with respect to a (Fourier) basis set. The size N of the matrix can
easily reach 50000 or more. But do we really need all its solutions
, or is it enough to just compute the
much fewer valence wave functions 
?
Typically, m is at most a few hundred, and thus much smaller than
N. A few lowest eigenvalues and eigenvectors can be efficiently
computed by means of an iterative algorithm. Computing all the
eigenvectors is very demanding both computationally and in terms of
memory (just storing all the eigenvectors would require 40 GB of
memory). The problem is solved by realizing that Eq. (3)
can be rewritten as:
with the definition:
Now the problem is to compute 
in
Eq. (6), which still requires the knowledge of all
. But there is an alternate route to computing the
's. In fact, they are solutions to the equation:
Equation (7) is a large system of linear equations of the form
, and has to be solved for each of the
. Since there are efficient iterative
schemes for solving linear equations, this is about as expensive as
computing the valence wave functions 
. However,
we now avoid the prohibitive cost of computing the rest of the
spectrum 
. This trick enables us to handle
much larger problems.
We summarize the section by presenting the two key algorithms which we use to compute the NMR chemical shift.
and eigenvalues
of the large 
matrix H in
Eq. (4). Typically m is a few hundred, but N is
several tens of thousands.
to the m large 
linear systems of
Eq.(7).
