% \section{WANNIER FUNCTIONS IN OTHER CONTEXTS}
% \label{sec:other}

As described in Sec.~\ref{sec:bloch-wannier}, Wannier
functions provide an alternative, localized, description of a
manifold of states spanned by the eigenstates (energy bands) of a
single-particle Hamiltonian that describes electrons in a periodic
potential.
%%
The equivalence of the Wannier representation and the
eigenstate representation may be expressed in terms of the
band projection operator $\hat{P}$, see \equ{mbandproj}.
This operator satisfies the idempotency condition
$P^2=P$, which embodies simultaneously the
requirements of orthogonality and Pauli exclusion.

From their conception, and until relatively recently, Wannier
functions have been used almost exclusively in this context, namely to
represent a manifold of single-particle orbitals for
electrons. Furthermore, as discussed in Sec.~\ref{sec:theory}, we need
not restrict ourselves to an isolated group of states, such as the
occupied manifold: the disentanglement procedure enables a subspace of
a larger manifold, e.g., of occupied and unoccupied states, to be
selected which may then be ``wannierized''. This has, for example,
opened up areas of application in which Wannier functions are used as
tight-binding basis functions for electronic structure and transport
calculations, as described in Sec.~\ref{sec:winterp} and
Sec.~\ref{sec:basis}.

From a general mathematical point of view, howewer, the set of
orthogonal eigenfunctions of any self-adjoint (Hermitian) operator may
be ``rotated'' by unitary transformation to another orthogonal 
basis that spans the same space.
%%
As we have seen, the unitary transformation is arbitrary and may be
chosen to render the new basis set maximally-localized, which has
computational advantages when it comes to representing physical
quantities that are short-ranged. When the operator in question
has translational symmetry, the maximally-localized functions thus
obtained are reminiscent of the Wannier functions familiar from
electronic structure theory. Often, such a basis is also preferable to
using another localized basis because information regarding the
symmetries of the self-adjoint operator from which the basis is
derived is encoded within it.


These ideas have led to the appropriation of the MLWF formalism
described in Sec.~\ref{sec:theory} for contexts other than the
description of electrons: the single-particle electronic Hamiltonian
is replaced by another suitable periodic self-adjoint operator, and
the Bloch eigenstates are replaced by the eigenfunctions of that
operator, which may then be transformed to give a MLWF-like
representation that may be used as an optimal and compact basis for
the desired calculation, for example, to analyze the eigenmodes of the
static dielectric matrix in ice and liquid water~\cite{lu-prl08}.


Below we review the three most prominent of these applications, namely
to the study of phonons, photonic crystals, and cold atom lattices.


\subsection{Phonons}
\label{sec:phonons}

Lattice vibrations in periodic crystals are usually described in
terms of normal modes, which constitute a delocalized orthonormal
basis for the space of vibrations of the lattice such that an arbitrary
displacement of the atoms in the crystal may be expressed in terms of
a linear combination of its normal modes.
%%%
By analogy with the electronic case, \textcite{kohn-prb73} first
showed (for isolated phonon branches in one-dimensional) that a
similar approach could be used for constructing a localized
orthonormal basis for lattice vibrations that span the same space as
the delocalized normal modes. The approach was subsequently
generalized to isolated manifolds in three-dimensions by
\textcite{tindemansvaneijndhoven-jpcsp75}. The localized modes are now
generally referred to as \emph{lattice Wannier functions}
(LWFs)~\cite{rabe-prb95, iniguez-prb00}.
%%%

Following the notation of Sec.~\ref{sec:el-ph}, we denote by 
$\qq$ the phonon wavevector, and by $e_\qq$ the matrix whose 
columns are the eigenvectors of the dynamical matrix. As in
case of electronic Wannier functions, the phases of these eigenvectors
are undetermined. A unitary transformation of the form
%%%
\beq
\left[ \tilde{e}_{\qq} \right]_{\mu\nu} 
 = \left[ M_{\qq} e_{\qq}\right]_{\mu\nu},
\label{eq:lwf-gauge}
\eeq
%%%
performed within a subspace of dispersion branches that is
invariant with respect to the space group of the crystal, results in
an equivalent representation of generalized extended modes
$\left[\tilde{e}_{\qq}\right]_{\mu\nu}$ that are also orthonormal. 
%
LWFs may then be defined by
%%%
\beq
\left[w_{\RR}\right]_{\mu\nu} = \frac{1}{N_p} \sum_{\qq}
e^{-i\qq\cdot\RR} \left[\tilde{e}_{\qq}\right]_{\mu\nu}, 
\label{eq:lwf-transform}
\eeq
%%%
with the associated inverse transform
%%%
\beq
\left[\tilde{e}_{\qq}\right]_{\mu\nu} = \sum_{\RR} e^{i\qq\cdot\RR}
\left[w_{\RR}\right]_{\mu\nu}.
\eeq
%%%
By construction, the LWFs are periodic according to $w_{\RR+\bf{t}} =
w_{\RR}$, where $\bf{t}$ is a translation vector of the Born-von
K\'{a}rm\'{a}n supercell.

%%%
The freedom inherent in \equ{lwf-gauge} allows very localized LWFs to
be constructed, by suitable choice of the transformation matrix $M_{\qq}$. 
As noted by \textcite{kohn-prb73}, the proof of exponential
localization of LWFs follows the same reasoning as for electronic
Wannier functions (Sec.~\ref{sec:wfprop}).


%%%
The formal existence of LWFs was first invoked in order to justify
the construction of approximate so-called \emph{local modes} of
vibration which were used in effective Hamiltonians for the study of
systems exhibiting strong coupling between electronic states and
lattice instabilities, such as perovskite
ferroelectrics~\cite{thomas-prl68, pytte-pr69}. 

\textcite{zhong-prl94} used first-principles methods in order to
calculate the eigenvector associated with a soft mode at $\qq=\0$ in
BaTiO$_3$. A localized displacement pattern, or local mode, of the
atoms in the cell was then parametrized, taking account of the
symmetries associated with the soft mode, and the parameters were
fitted to reproduce the calculated soft mode eigenvector at
$\qq=\0$. The degree of localization of the local mode was determined
by setting to zero all displacement parameters beyond the second shell
of atoms surrounding the central atom. Although this spatial
truncation results in the local modes being non-orthogonal, it does
not hamper the accuracy of practical calculations. As the local modes
are constructed using information only from the eigenvector at
$\qq=\0$, they do not correspond to a particular phonon branch in the
Brillouin zone.
%%%
\textcite{rabe-prb95} generalized the approach to allow fitting to
more than just $\qq=\0$, but rather to a small set of, usually
high-symmetry, $\qq$-points. The phase-indeterminacy of the
eigenvectors is exploited in order to achieve optimally rapid decay of
the local modes. 
%%%
Another approach, introduced by \textcite{iniguez-prb00}, constructs
local modes via a projection method that preserves the correct
symmetry. The procedure is initiated from simple atomic displacements
as trial functions. The quality of the local modes thus obtained may
be improved by systematically densifying the $\qq$-point mesh that is
used in \equ{lwf-transform}. Although there is no formal criterion of
maximal-localization in the approach, it also results in
non-orthogonal local modes that decay exponentially.

These ideas for generating local modes from first-principles
calculations have been particularly successful for the study of
structural phase transitions in ferroelectrics such as
BaTiO$_3$~\cite{zhong-prl94, zhong-prb95},
PbTiO$_3$~\cite{waghmare-prb97}, KbNiO3~\cite{krakauer-jpm99},
Pb$_3$GeTe$_4$~\cite{cockayne-prb97} and perovskite
superlattices~\cite{lee-jap08}.
%%% 

The use of maximal-localization as an exclusive criterion for
determining LWFs was first introduced by \textcite{giustino-prl06}. In
this work, a real-space periodic position operator for non-interacting
phonons was defined, by analogy with the periodic position operator for
non-interacting electrons (\equ{rcenter}). 
%%%

%%% 
The problem of minimizing the total spread of a set of WFs in
real-space is equivalent to the problem of simultaneously 
diagonalizing the three non-commuting matrices corresponding to the
three components of the position operator represented in the WF
basis, and \textcite{giustino-prl06} use the method outlined by
\textcite{gygi-cpc03} to achieve this. It is worth noting that
\textcite{giustino-prl06} furthermore define a generalized spread
functional that, with a single parameter, allows a trade-off between
localization in energy (the eigenstate or Bloch limit) and
localization in space (the Wannier limit), resulting in so-called
\emph{mixed Wannier-Bloch functions} which may be obtained for
the electrons as well as the phonons.

%%%
Finally, as first pointed out by
\textcite{kohn-prb73}, and subsequently by \textcite{giustino-prb07},
maximally-localized lattice Wannier functions correspond to
displacements of individual atoms. This may be
seen by considering a vibrational eigenmode, $\hat{\bf e}^{\nu}_{\qq
s} \equiv {\bf e}^{\nu}_{\qq s} e^{i\qq\cdot\RR}$, and noting that it
may be expressed as 
%%%
\beq
\label{eq:eigmode}
\hat{\bf e}^{\nu}_{\qq s} = 
    \sum_{s'\RR'}\,e^{i\qq\cdot\RR'}
\delta_{\RR\RR'}\delta_{s s'} 
    {\bf e}^{\nu}_{\qq s'}.
\eeq
%%%
\equ{eigmode} stands in direct correspondence to the electronic
analogue given by the inverse of \equ{wannier}, from which we conclude
that the LWFs do indeed correspond to individual atomic displacements
$\delta_{\RR\RR'}\delta_{ss'}$ and, furthermore, that the required
unitary transformation is the matrix of eigenvectors
$\left[e_{\qq}\right]_{\mu\nu}$.
%
As discussed in Sec.~\ref{sec:el-ph}, \textcite{giustino-prb07}
exploit this property for the efficient interpolation of dynamical
matrices and calculation of electron-phonon couplings. 



\subsection{Photonic crystals}
\label{sec:phot-crys}

Photonic crystals are periodic arrangements of dielectric materials
that are designed and fabricated in order to control the flow of
light~\cite{john-prl87, yablonovitch-prl87}. They are very much to light 
what semiconductors are to  
electrons and, like semiconductors that exhibit an electronic band
gap in which an electron may not propagate, photonic
crystals can be engineered to exhibit photonic band gaps: ranges of
frequencies in which light is forbidden to propagate in the
crystal. In the electronic case, a band gap results from scattering
from the periodic potential due to the ions in the crystal; in the
photonic case, it arises from scattering from the periodic dielectric 
interfaces of the crystal. Again by analogy with electronic
materials, localized defect states can arise in the gap by the
deliberate introduction of defects into a perfect photonic crystal
structure. The ability to control the nature of these states promises
to lead to entirely light-based integrated circuits, which would have
a number of advantages over their electronic counterparts, including
greater speeds of prpagation, greater bandwidth, and smaller energy
losses~\cite{joannopoulos-n97}.

In SI units, Maxwell's equations in source-free regions of space are
\begin{align}
&\nabla\cdot{\bf E} = 0,& &\nabla\cdot\BB = 0,& \\
&\nabla\times{\bf E} = -\frac{\partial\BB}{\partial t},&
&\nabla\times\HH = \frac{\partial\DD}{\partial t},& 
\end{align}
where the constitutive relations between the fields are
\beq
\DD = \eps_{\rm r}\eps_0{\bf E}, \quad
\BB = \mu_{\rm r}\mu_0\HH.
\eeq
Considering non-magnetic materials ($\mu_{\rm r}=1$) with a position
dependent dielectric constant $\eps_{\rm r}(\rr)$, and fields that
vary with a sinusoidal dependence $e^{-i\omega t}$, it is
straightforward to derive electromagnetic wave equations in terms of
either the electric field ${\bf E}$ or the magnetic field $\HH$, 
\bea
\nabla\times\left(\nabla\times{\bf E}(\rr)\right) &=&
      \frac{\omega^2}{c^2}\eps_{\rm r}(\rr){\bf E}(\rr) ,
      \label{eq:wave-eqn-E}\\
\nabla\times\left(\eps^{-1}_{\rm r}(\rr) \nabla\times\HH(\rr)\right) &=&
      \frac{\omega^2}{c^2}\HH(\rr), 
      \label{eq:wave-eqn-H}
\eea
where $c=(\mu_0\eps_0)^{-1/2}$ is the speed of light.
%%

For a perfect periodic dielectric structure, $\eps_{\rm
r}(\rr)=\eps_{\rm r}(\rr + \RR)$, where $\RR$ is a lattice
vector. Application of Bloch's theorem leads to solutions that are
indexed by wavevector $\kk$, which may be chosen to lie in the first
Brillouin zone, and a band index $n$. For example 
$$
\HH_{\bnk}(\rr)=e^{i\kk\cdot\rr}{\bf u}_{\bnk}(\rr) ,
$$
where ${\bf u}_{\bnk}(\rr)={\bf u}_{\bnk}(\rr+\RR)$ is the periodic
part of the magnetic field Bloch function. The electromagnetic wave
equations can be solved, and hence the Bloch functions obtained, by a
number of methods including finite-difference time
domain~\cite{yee-itoap66, taflove-book00}, transfer
matrix~\cite{pendry-prl92, pendry-jpm96}, empirical tight-binding
methods~\cite{lidorikis-prl98, yariv-ol99}, and Galerkin techniques in
which the field is expanded in a set of orthogonal basis
functions~\cite{mogilevtsev-jlt99}. Within the latter
class, use of a plane-wave basis set is particularly
common~\cite{ho-prl90, johnson-oe01}.


The operators $\nabla\times\nabla$ and $\nabla\times\eps^{-1}_{\rm
r}(\rr)\nabla$ are self-adjoint and, therefore, the fields satisfy 
orthogonality relations given by\footnote{The notation 
$\bf{A}\cdot\bf{B} = \sum_{i=1}^{3}A_i B_i$, and denotes the 
scalar product of the vectors $\bf{A}$ and $\bf{B}$, with 
Cartesian components $\{A_{i}\}$ and $\{B_{i}\}$, respectively.}
\bea
\int d\rr\, \HH^{\ast}_{n\kk}(\rr) \cdot \HH_{n'\kk'}(\rr) &=&
         \delta_{nn'}\delta(\kk-\kk'), \\
\int d\rr\, \epsilon_{\rm r}(\rr)\: {\bf E}^{\ast}_{n\kk}(\rr)\cdot  
         {\bf E}_{n'\kk'}(\rr) &=&
         \delta_{nn'}\delta(\kk-\kk').
\eea

\textcite{leung-josabp93} first suggested that transforming to a basis
of Wannier functions localized in real space would be advantageous for
computational efficiency, especially when dealing with defects in
photonic crystals which, using conventional methods, require very
large supercells for convergence. Although of great formal importance
for justifying the existence of a suitable localized basis, and hence
the tight-binding approach, the non-uniqueness of the transformation
between Bloch states and Wannier functions caused difficulties. As a
result, early work was limited to the case of single, isolated
bands~\cite{leung-josabp93, konotop-josabp97} or composite bands in
which the matrix elements $U_{mn}^{(\kk)}$ were treated as parameters
to fit the tight-binding band structure to the plane-wave result.

The formalism for obtaining maximally-localized Wannier functions,
however, removed this obstacle and several applications of MLWFs to
calculating the properties of photonic crystals have been reported
since, in both
two-dimensional~\cite{garcia-martin-n03, whittaker-prb03, jiao-ijqe06}
and three-dimensional~\cite{takeda-prb06} photonic crystal structures,
as well as for the case of entangled bands~\cite{hermann-josabp08}
(see \onlinecite{busch-jpcm03} for an early review).

Typically one chooses to solve for either the electric field ${\bf E}_{\bnk}$
or the magnetic field $\HH_{\bnk}$. Once the Bloch states for the
periodic crystal are obtained, a basis of magnetic or electric field
Wannier functions may be constructed using the usual definition, e.g.,
for the magnetic field 
%
\beq
\WW_{n\RR}^{(\HH)}(\rr) = 
         \ibz d\kk\,e^{-i\kk\cdot\RR}
\sum_{m} U_{mn}^{(\kk)}
         \HH_{m\kk}(\rr), 
\eeq
satisfying orthogonality relations
%
\beq
\int d\rr\, \WW_{n\RR}^{(\HH)\ast}\cdot\WW_{n'\RR'}^{(\HH)} 
   = \delta_{nn'}\delta_{\RR\RR'} ,
\eeq
where the unitary transformation $U_{mn}^{(\kk)}$ is chosen in the
same way described in Sec.~\ref{sec:theory} such that the sum of the
quadratic spreads of the Wannier functions is minimsed, i.e., such
that the Wannier functions are maximally localized. 

Concentrating on the magnetic field, it may be expanded in the basis of
Wannier functions with some expansion coefficients $c_{n\RR}$,
\beq
\HH(\rr) = \sum_{n\RR} c_{n\RR} \WW_{n\RR}^{(\HH)}(\rr),
\label{eq:H-expansion}
\eeq
which on substitution into \equ{wave-eqn-H} gives the
tight-binding representation of the wave-equation for the magnetic
field in the Wannier function basis. 

The utility of the approach becomes evident when considering the
presence of a defect in the dielectric lattice such that $\eps_{\rm
r}(\rr)\rightarrow \eps_{\rm r}(\rr) + \delta\eps(\rr)$. The magnetic
field wave equations become
\beq
\nabla\times\left( \left[\eps^{-1}_{\rm r}(\rr)+\Delta^{-1}(\rr)\right]
      \nabla\times\HH(\rr)\right) = \frac{\omega^2}{c^2}\HH(\rr),
\eeq
where
\beq
\Delta(\rr) = \frac{-\delta\eps(\rr)}
        {\eps_{\rm r}(\rr)[\eps_{\rm r}(\rr)+\delta\eps(\rr)]}.
\eeq
Using the Wannier functions from the defect-free
calculation as a basis in which to expand $\HH(\rr)$, as in
\equ{H-expansion}, the wave equations may be written in matrix form,
\beq
\sum_{n'\RR'}\left(A_{nn'}^{\RR\RR'} + B_{nn'}^{\RR\RR'} \right)
   c_{n'\RR'} = \frac{\omega^2}{c^2}c_{n\RR},
\eeq
where 
%
\beq
A_{nn'}^{\RR\RR'} = \ibz d\kk\, e^{i\kk\cdot(\RR-\RR')}\sum_{m}
   U_{nm}^{(\kk)\dagger} \left(\frac{\omega_{m\kk}}{c}\right)^2
   U_{mn'}^{(\kk)},
\eeq
and
\beq
B_{nn'}^{\RR\RR'}=\int d\rr\,\Delta(\rr)
                  \left[\nabla\times\WW_{n\RR}(\rr)\right]^{\ast}
             \cdot\left[\nabla\times\WW_{n'\RR'}(\rr)\right].
\eeq
Due to the localization and compactness of the basis, these matrix
equations may be solved efficiently to find frequencies of localized
cavity modes, dispersion relations for waveguides, and the
transmission and reflection properties of complex waveguide
structures. Fig.~\ref{fig:photonic}, for example, shows the photonic
band structure for a three-dimensional photonic crystal structure 
with a two-dimensional defect. 
%
\begin{figure}[h!]
\begin{center}
\includegraphics[width=1.6in]{figs/takeda-prb06-fig17}
\includegraphics[width=1.6in]{figs/takeda-prb06-fig19b}\\
\vspace{3mm}
\includegraphics[angle=90,width=2.5in]{figs/takeda-prb06-fig18}
\caption{(Color online) Photonic band structure (bottom) of the 3D Si woodpile structure
intercalated with a 2D layer consisting of a square lattice of
square rods (top left). Solid lines indicate the photonic band structure
calculated by the plane-wave expansion (PWE) method, and black points
indicate that reproduced by the MLWFs. Shaded regions indicate
the photonic band structure of the woodpile projected onto the 2D
$\kk_\parallel$ space. The square rods in the 2D layer are chosen
to structurally match the woodpile. The thickness of the layer is $0.8a$, 
where $a$ is the lattice parameter of the woodpile structure. Top right: absolute
value of the 17th MLWF of the magnetic field in the $yz$ plane.
Adapted from \onlinecite{takeda-prb06}.}
\label{fig:photonic}
\end{center}
\end{figure}
%



\subsection{Cold atoms in optical lattices}
\label{sec:cold-atoms}

A good 70 years after Albert Einstein predicted that a system of
non-interacting bosons would undergo a phase transition to a new state
of matter in which there is macroscopic occupation of the lowest
energy quantum state, major achievements in laser cooling and
evaporative cooling of atoms enabled the first experimental
realizations of Bose-Einstein condensation~\cite{anderson-s95,
davis-prl95, bradley-prl95} and the award of Nobel prizes in 1997 and
2001. Since then, the study of cold atoms trapped in optical lattices
has flourished. For a reviews see \textcite{morsch-rmp06,
bloch-rmp08}.

Ultracold atoms trapped in optical lattices provide a versatile
alternative to electrons in crystal lattices for the study of quantum
phenomena. Indeed, they have a number of advantages over
the solid state in this respect, such as the absence of lattice
defects, the absence of a counterpart to the electron-phonon
interaction, and the ability to control precisely both the nature of
the inter-atomic interactions and the depth and periodicity of the
optical lattice potential.

The second quantized Hamiltonian for a system of $N$
weakly interacting bosons of zero spin and mass $m$ in a (periodic)
external potential $V_{0}(\rr)=V_{0}(\rr + \RR)$ is
given by~\cite{yukalov-lp09}
%
\bea
\hat{\mathcal{H}} &=& \intr \hat{\Psi}^{\dagger}(\rr) \left[
-\frac{\hbar^2}{2m}\nabla^2 + V_{0}(\rr) \right]
\hat{\Psi}(\rr) \nonumber \\ 
& & + \: \frac{g}{2} \intr
\hat{\Psi}^{\dagger}(\rr) \hat{\Psi}^{\dagger}(\rr)
\hat{\Psi}(\rr)\hat{\Psi}(\rr),
\label{eq:cold-atom-ham}
\eea
where $g=4\pi a_s \hbar^2 / m$, it is assumed that the atoms
interact via a short-range pseudopotential with $a_s$ as the
$s$-wave scattering length, and $\hat{\Psi}(\rr)$ and
$\hat{\Psi}^{\dagger}(\rr)$ are bosonic field operators obeying 
canonical commutation relations~\cite{fetter-book-ch1}.
%

In a Bose-Einstein condensate, wherein the condensate particle
densities are typically of the order of $10^{14}$\:cm$^{-3}$ or more,
the mean-field limit of this Hamiltonian is usually taken, which leads
to the Gross-Pitaevskii equation, also known as the non-linear
Schr\"{o}dinger equation (NLSE),
\beq
i\hbar \frac{\partial}{\partial t} \varphi (\rr, t) = 
\left[ -\frac{\hbar^2}{2m} + V_{0}(\rr) + 
g \left|\varphi(\rr,t)\right|^2 \right] \varphi(\rr,t).
\label{eq:gp}
\eeq
$\varphi(\rr,t)$ is the condensate wavefunction, the squared
norm of which gives the condensate density. The Gross-Pitaevskii
equation has been used with remarkable success in the study of
BEC~\cite{leggett-rmp01}.

As shown, for example, by \textcite{alfimov-pre02}, a basis of Wannier
functions, localized to each site $\alpha$ of the optical lattice, may
be used to expand the condensate wavefunction,
\beq
\varphi(x,t) = \sum_{n,\alpha} c_{n\alpha}(t)w_{n\alpha}(x).
\label{eq:condensate-exp}
\eeq
The Wannier functions are related to the Bloch eigenstates
$\psi_{nk}(x)$ of the eigenvalue equation
\beq
\left[ -\frac{d^2}{dx^2} + V_{0}(x)\right] \psi_{nk}(x) =
\epsilon_{nk}\psi_{nk}(x)
\eeq
by the usual Wannier transformation 
\beq
w_{n\alpha}(x) = \frac{L}{2\pi} \int_{\rm BZ} dk\,
e^{i\phi_{n}(k)} e^{-ik\alpha L}\psi_{nk}(x)
.
\eeq
Substituting \equ{condensate-exp} into \equ{gp}, leads to a tight-binding
formulation known as the discrete non-linear Schr\"{o}dinger equation
(DNLSE),
\beq
i\frac{d}{d t}c_{n\alpha} = \sum_{\beta} c_{n\beta}
\tilde{\epsilon}_{n,\alpha-\beta} 
 + \: \sum_{\beta,\gamma,\eta} \sum_{i, j, k}
c^{\ast}_{i\beta} c_{j\gamma} c_{k\eta} U_{nijk}^{\alpha\beta\gamma\eta},
\label{eq:dnlse}
\eeq
where
%
\beq
\tilde{\epsilon}_{n\alpha} = \frac{L}{2\pi} 
\int_{\rm BZ} dk\, e^{-ik\alpha L}
\epsilon_{nk},
\eeq
and the interaction matrix is given by
\beq
U_{nijk}^{\alpha\beta\gamma\eta} = g \int dx\,
w_{n\alpha}(x) w_{i\beta}(x)
w_{j\gamma}(x) w_{k\eta}(x).
\eeq

Truncating the first term on the right-hand side of \equ{dnlse} to
nearest-neighbors only, and the second term to on-site
($\alpha=\beta=\gamma=\eta$) terms within a single band ($n=i=j=k$)
results in the usual tight-binding description~\cite{chiofalo-epjd00,
trombettoni-prl01}, 
\beq
i\frac{d}{d t}c_{n\alpha} = c_{n\alpha} \tilde{\epsilon}_{n0} 
+ \tilde{\epsilon}_{n1}(c_{n,\alpha-1} + c_{n,\alpha+1})
+ U_{nnnn}^{\alpha\alpha\alpha\alpha}
\left|c_{n\alpha}\right|^2 c_{n\alpha}.
\label{eq:dnlse-simple}
\eeq

As pointed out by \textcite{alfimov-pre02}, using a WF basis enables
the range and type of interactions encapsulated in the DNLSE to be
systematically controlled and improved. On the most part, however, WFs
have been used in the context of the NLSE in order to carry out
formal derivations and to justify the use of empirical or
semi-empirical tight-binding models.

An interesting analogy with electrons in atomic lattices manifests
itself when the filling of sites in the optical lattice is low and
hence particle correlations need to be accounted for more
rigorously. This is done via the Bose-Hubbard model, developed by
\textcite{fisher-prb89} in the context of He-4, and first applied to
cold atoms in optical lattices by \textcite{jaksch-prl98}. The
Bose-Hubbard Hamiltonian is derived from \equ{cold-atom-ham} by
expanding the boson field operator in terms of WFs of a single band,
localized at the lattice sites,
\beq
\hat{\Psi} = \sum_\alpha \hat{b}_\alpha w_{\alpha}(\rr),
\eeq
where the bosonic particle creation and annihilation operators,
$\hat{b}_\alpha$ and $\hat{b}^{\dagger}_{\alpha}$, respectively, satisfy
the usual commutation rules~\cite{fetter-book-ch1}.
%%
This, on approximation to nearest-neighbour coupling,
and on-site only interactions, results in the standard Bose-Hubbard
Hamiltonian~\cite{jaksch-prl98}
\beq
\hat{H}_{\rm BH} = -J \sum_{<\alpha,\beta>} \hat{b}^{\dagger}_{\alpha}
\hat{b}_{\beta} + 
\frac{U}{2} \sum_\alpha \hat{n}_\alpha (\hat{n}_\alpha - 1),
\eeq
where 
$\hat{n}_\alpha=\hat{b}^{\dagger}_{\alpha}\hat{b}_{\alpha}$ is the
number operator for lattice site $\alpha$, and the nearest-neighbor
hopping and on-site repulsion parameters are given by
\beq
J = -\intr w_{0}(\rr) 
\left[ -\frac{\hbar^2}{2m}\nabla^2 + V_{0}(\rr) \right] 
w_{1}(\rr),
\eeq
and
%
\beq
U = g \int d\rr\, \left|w(\rr)\right|^4,
\eeq
which may be calculated explicitly using WFs constructed from Bloch
eigenstates~\cite{vaucher-njp07, shotter-pra08}. The Bose-Hubbard
model is the bosonic analogue to the Hubbard model for fermions. 
As in the latter case, the behaviour of the model
depends on the competition between hopping ($J$) and on-site
($U$) energies which determines whether the system is in a superfluid or
a Mott insulator phase.

Finally, we note that in work that is closely related to, and combines
elements from both ideas developed in photonic crystals and cold
atoms, WFs have also been used to represent polaritons in coupled
cavity arrays, a class of systems that serves as another experimental
realization of the Bose-Hubbard model~\cite{hartmann-np06,
hartmann-lpr08}.



