\documentclass[prl,showpacs,superscriptaddress,twocolumn,preprintnumbers,amssymb,amsfonts]{revtex4}
%\documentclass[prl,showpacs,superscriptaddress,galley,preprintnumbers,amssymb,amsfonts]{revtex4}

\usepackage{amsmath}
\usepackage{bm}
\usepackage{epsfig}
\usepackage{graphics}
%\usepackage{subfigure}

\def\beq{\begin{equation}}
\def\eeq{\end{equation}}
\def\beqn{\begin{eqnarray}}
\def\eeqn{\end{eqnarray}}

\def\E{{\bf E}}
\def\B{{\bf B}}
\def\P{{\bf P}}
\def\M{{\bf M}}
\def\K{{\bf K}}
\def\k{{\bf k}}
\def\r{{\bf r}}
\def\del{\partial}
\def\F{{\cal F}}
\def\A{{\cal A}}
\def\ket#1{\vert #1 \rangle}
\def\ev#1{\langle #1 \rangle}
\def\ip#1#2{\langle #1 \vert #2 \rangle}
\def\me#1#2#3{\langle #1 \vert #2 \vert #3 \rangle}
\def\inh{^{\rm (in)}}
\def\intbz{\int_{\rm BZ}\frac{d^3k}{(2\pi)^3}}

\begin{document}
\title{Magnetoelectric polarizability and axion electrodynamics in crystalline insulators}
\author{Andrew M. Essin}
\affiliation{Department of Physics, University of California,
Berkeley, CA 94720} 
\author{Joel~E.~Moore}
\affiliation{Department of Physics, University of California,
Berkeley, CA 94720} \affiliation{Materials Sciences Division,
Lawrence Berkeley National Laboratory, Berkeley, CA 94720}
\author{David Vanderbilt}
\affiliation{Department of Physics and Astronomy, Rutgers University, Piscataway, NJ 08854}
\date{October 16, 2008}

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% \marginparwidth 2.7in
% \marginparsep 0.5in
% \def\dvm#1{\marginpar{\small DV: #1}}
% \def\jmm#1{\marginpar{\small JM: #1}}
% \def\aem#1{\marginpar{\small AE: #1}}
%uncomment next lines to have no commentaries
% \def\dvm#1{}
% \def\jjm#1{}
% \def\aem#1{}
%
\def\scr{\scriptsize}
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%


\begin{abstract}
Spin-orbit coupling can lead in two- and three-dimensional solids to
time-reversal-invariant insulating phases that are ``topological''
in the same sense as the integer quantum Hall effect and similarly
have protected edge or surface states.  The three-dimensional
topological insulator is known to have unusual magnetoelectric
properties referred to as ``axion electrodynamics'': it supports an
electromagnetic coupling $\Delta{\cal L}_{EM} = (\theta e^2 / 2 \pi h)
{\bf E} \cdot {\bf B}$ with $\theta=\pi$, giving a half-integer
surface Hall conductivity $\sigma_{xy}=(n + 1/2) e^2 / h$.  An approach to $\theta$
in any three-dimensional crystal is developed based on the Berry-phase theory of
polarization: $\theta e^2/ 2 \pi h$ is the bulk 
orbital magnetoelectric polarizability (the polarization induced by
an applied magnetic field).  We compute the orbital magnetoelectric polarizability for a simple model and show that it predicts the fractional part of surface $\sigma_{xy}$, computed
using a slab geometry.  Although $\theta$ is not
quantized once time-reversal and inversion symmetries are broken,
it remains a bulk quantity for the same reasons as ordinary polarization.
\end{abstract}

\pacs{73.43.-f, 85.75.-d, 73.20.At, 03.65.Vf, 75.80.+q}

\maketitle

%\dvm{\scr Abstract:\\
%1. The order in which things are mentionedin the
%abstract doesn't quite seem to match the order in which things are
%presented in the paper.  (Maybe it doesn't matter?...)
%{\bf Note: In my 10/15 revision I at least removed the word ``first''
%in one place.}\\
%2. I'm not sure what is intended by the word ``universal'' in the last
%sentence.  Remove it?\\
%3. The abstract strikes me as being a bit too long.}

The behavior of electric and magnetic fields can be strongly
modified inside an insulating solid.  Changes
to the dielectric constant or magnetic permeability could be
regarded as just numerical renormalizations of vacuum properties,
but the solid environment can also yield entirely new behavior.
For example, it has been known for two decades that
three-dimensional (3D) insulators can generate an
electromagnetic coupling of the form ($c=1$)
%
\beq
\Delta {\cal L}_{EM} = \frac{\theta e^2}{2\pi h} {\bf E} \cdot {\bf B} =
\frac{\theta e^2}{16\pi h} \epsilon^{\alpha \beta \gamma \delta} 
	F_{\alpha \beta} F_{\gamma \delta}.
\label{axioncoupling}
\eeq
%
This coupling was labeled ``axion
electrodynamics''~\cite{wilczekaxion} since a Lagrangian density of
this form can describe the interaction of a dynamical ``axion field''
$\theta(\bm{x},t)$ with the electromagnetic field. 
It plays no role in
vacuum electrodynamics since the parameter $\theta$ 
couples to a total derivative, $\epsilon^{\alpha \beta \gamma \delta}
 F_{\alpha \beta} F_{\gamma \delta} = 2 \epsilon^{\alpha \beta \gamma \delta} 
 \partial_\alpha (A_\beta F_{\gamma \delta})$, and so does not modify
the classical equations of motion if $\theta$ is constant.
%
%

This theory is invariant under $\theta \rightarrow \theta+2 \pi$,
and the time-reversal operator $T$ maps $\theta\rightarrow-\theta$ (since ${\bf E} \rightarrow {\bf E}$ and ${\bf B} \rightarrow -{\bf B}$), so $T$ invariance is
consistent with $\theta=\pi$ as well as $\theta=0$.~\cite{wilczekaxion}
%
% \dvm{Close parenthesis missing.  I didn't fix it because I wasn't
% sure if you wanted the citation inside it our outside it.}
%
Thus $\theta$ is a ``topological invariant'' of $T$-invariant
insulators, that is, a property stable under continuous changes
as long as the material remains insulating and $T$-symmetric.
The prototypical invariant of this type is the TKNN integer $C$
that gives the quantum Hall conductance, $\sigma_{xy} = C e^2/h$.~\cite{tknn}

%\aem{\scr My suggestions for reworking: From symmetry considerations
%one might expect $\theta$ to be non-zero only when both inversion and
%time-reversal (TR) symmetries are broken in the material, but surprisingly
%TR symmetry need not be broken.
%Since the theory is invariant under $\theta \rightarrow \theta+2 \pi$
%(a nonobvious result),
%and TR maps $\bm{B}\rightarrow-\bm{B}$, TR invariance is
%consistent with $\theta=\pi$ as well as $\theta=0$.}

% Until recently it was unclear how to determine whether a given
% $T$-invariant insulator realized the $\theta=0$ or $\theta=\pi$
% case.  $T$-invariant band insulators turn out to have ``topological
% invariants'' (properties stable under continuous changes as
% long as the material remains insulating) beyond those familiar
% in the general case, such as the TKNN integer~\cite{tknn}
% that distinguishes integer quantum Hall states on a lattice.
% These topological invariants result from the Berry connection
% of the Bloch wavefunctions of a band structure, which behaves in
% many respects like a vector potential in electromagnetism.

%\aem{\scr  Thus $\theta$ is a ``topological invariant'' of TR-invariant
%insulators, that is, a property stable under continuous changes
%as long as the material remains insulating and TR-symmetric.
%The prototypical invariant is the TKNN integer~\cite{tknn}
%that distinguishes integer quantum Hall states on a lattice.}

% \dvm{I was also going to cut some of the redundant wording in the
% next paragraph, although I was going to keep some transition wording
% like:}
%
$T$-invariant band insulators turn out to have topological invariants
beyond those of the TKNN type.  In 2D there is a $\mathbb{Z}_2$
invariant~\cite{km2} that distinguishes ``ordinary'' from
``$\mathbb{Z}_2$-odd'' insulators, with ``quantum spin Hall'' states~\cite{kane&mele-2005,zhangscience1} providing examples of the latter.
In 3D there is a similar invariant~\cite{moore&balents-2006,rroy3d,fu&kane&mele-2007}
that can be computed either from the 2D invariant on certain
planes~\cite{moore&balents-2006} or, with additional symmetry,
from a certain index at the eight time-reversal-invariant
momenta~\cite{fu&kane&mele-2007}.  The resulting ``strong
topological insulator'' in 3D has Dirac surface modes
that have recently been observed in photoemission on
Bi$_{0.9}$Sb$_{0.1}$~\cite{hsieh}.  The $\theta=\pi$ value in Eq.~(\ref{axioncoupling})
%
% \dvm{We should be consistent about how we refer to equation numbers.
% I thought PR preferred ``Eq.~(\ref{axioncoupling})'' and not just
% ``(\ref{axioncoupling})''.  But whatever we decide, be consistent
% throughout the paper...}
%
is related
to these surface Dirac modes~\cite{pankratov,fradkinaxion,wilczekaxion,fu&kane2-2006,qilong}: the fact
that each gapless 2D Dirac mode, under a weak applied field, becomes a
gapped Hall insulator with $|\sigma_{xy}| = e^2/2 h$ is the
microscopic origin of the $\theta$ coupling.

%\aem{\scr Recent work on two- and three- dimensional TR-invariant band 
%insulators has revealed for each dimensionality a $\mathbb{Z}_2$
%invariant ~\cite{km2,moore&balents-2006,rroy3d,fu&kane&mele-2007}
 %that distinguishes ``ordinary'' from
%``$\mathbb{Z}_2$-odd'' insulators, with spin-Hall states providing
%an example of the latter in two dimensions.
%The 3D invariant can be computed either from the 2D invariant on certain
%planes~\cite{moore&balents-2006} or, with additional symmetry,
%from a certain index at the eight TR-invariant
%momenta~\cite{fu&kane&mele-2007}.  The resulting ``strong
%topological insulator'' in 3D has Dirac surface modes
%that have recently been observed in photoemission on
%Bi$_{0.9}$Sb$_{0.1}$~\cite{hsieh}.}

%\dvm{New paragraph.}
%
Insulators that instead have broken time-reversal $T$ and also inversion $P$---so-called ``magnetoelectric'' (or ``multiferroic'') materials---have been the topic of intense experimental and theoretical investigations
in recent years \cite{spaldin-s05,fiebig-jpdp05}.
%
% \dvm{Refs added.  If I forget to send you the bibtex file, remind me...}
%
In general the first-order magnetoelectric coupling
$\alpha$ is a tensor giving a contribution of the form
$\alpha_{ij}E_iB_j$, as indeed expected for
magnetoelectric couplings arising from lattice displacements
and Zeeman interactions.  Here we discuss only the frozen-lattice
coupling via the orbital magnetization, which we denote
as the orbital magnetoelectic polarizability (OMP).
Remarkably the OMP is indeed only a {\it scalar}
as in Eq.~(1) and generates ``axion electrodynamics.''

% \dvm{We have a notation problem.  Here I used $m,n$ for band indices,
% while later $m,n$ are used for surface branch-choice integers, and
% even later $n$ is used for layer number.  I don't like $\alpha,\beta$
% for band indices.  How about $\mu,\nu$?}
%
In the present letter, we provide a simple derivation of this result via an expression for the OMP as an integral of the Chern-Simons 3-form over the Brillouin zone.  In particular, defining the Berry connection
$\A^{\mu\nu}_j=i\me{u_\mu}{\del_j}{u_\nu}$ where $\ket{u_\nu}$ is the
cell-periodic Bloch function of occupied band $\nu$ and
$\del_j=\partial/\partial k_j$, we obtain
%
\beq
\theta=\frac{1}{2\pi}\int_{\rm BZ} d^3k \; \epsilon_{ijk}\,
{\rm Tr}[\A_i\del_j\A_k-i\frac{2}{3}\A_i\A_j\A_k] 
\label{theta}
\eeq
%
where the trace is over occupied bands.  This is essentially Eq.~(78) of recent work by Qi, Hughes, and Zhang~\cite{qilong} that, among other results, importantly made the connection between 3D $\mathbb{Z}_2$ topological insulators
and axion electrodynamics with $\theta=\pi$.  Our derivation shows that this integral describes the orbital magnetoelectric polarizability (rather than a polarization~\cite{qilong}) and follows directly from an extension \cite{xiao} of the Berry-phase theory of polarization \cite{ksv} to the case of slow spatial variations of the Hamiltonian.  The OMP angle $\theta$ is a bulk
property in exactly the same sense as electric
%
% \dvm{Perhaps better: ``...in exactly the same sense as electric...'' ?}
%
polarization~\cite{ortizmartin,ksv}.  
% 
% %Rather than working from 2D invariants as in previous definitions, the 3D topological insulator was redefined in terms of the second Chern class invariant in 4D.
% %A similar mathematical object has also recently been obtained by semiclassical arguments in a study of inhomogeneity-induced polarization~\cite{xiao}.
% 
% %\dvm{My strategy for PRL-style articles is to try to frame the
% %context with 3-4 intro paragraphs (above), then give the ``power
% %punch'' paragraph saying briefly what is done in the present paper
% %and why it is the best thing since the invention of the wheel,
% %and only then start on details.  This is that paragraph.
% %I am very open to improvements, modifications,
% %etc.  I decided to move Eq.~(10) to become Eq.~(1) in this paragraph.
% %It is against my usual principles to put an equation in this
% %paragraph, but rules are made to be broken...}
% 
% %In the present Letter, we build on these developments by providing a transparent derivation of an expression for the OMP in terms of an integral of a Chern-Simons 3-form over the Brillouin zone. 
% %
% %\dvm{Next part modified relative to the previous draft in an
% %effort to merge in some of the previous introductory material.}
% %
% %This is essentially Eq.~(79) of Ref.~\cite{qilong}. However, it is derived here using 
% 
% When $T$ is present, the
% classification of normal ($\theta=0$) and strong topological
% ($\theta=\pi$) insulators is recovered.  Even when it is not,
%
% \dvm{Note: Here I tentatively cut ``When $T$ is present, the
% classification of normal ($\theta=0$) and strong topological
% ($\theta=\pi$) insulators is recovered.  Even when it is not,'' and
% joined up the paragraphs.}
%
We show that Eq.~(\ref{theta}) is amenable to direct calculation
using established band-structure methods and discuss its symmetry
properties.  Explicit numerical calculations on model
crystals are presented in order to validate the theory
and to illustrate how a non-zero $\theta$ corresponds to a
``fractionalization'' of the quantum Hall effect (QHE) at the
surfaces of magnetoelectric insulators~\cite{wilczekaxion,fu&kane2-2006,qilong}.

From Eq.~(\ref{axioncoupling}) it follows that
%
\beq
\theta \frac{e^2}{2\pi h} =
\frac{\partial M}{\partial E} =
\frac{\partial}{\partial E}\,\frac{\partial}{\partial B}\,\ev{H} =
\frac{\partial P}{\partial B}
\label{eq:theta}
\eeq
%
evaluated at $E=B=0$, where $\langle H \rangle$ is the ground-state
energy and we have used the commutativity of partial derivatives.
Thus the OMP can be approached from several points of view.
(i) It can be regarded as describing the {\it electric} polarization
arising from the application of a small {\it magnetic} field.
(ii) It can be regarded as describing the orbital {\it magnetization}
arising from the application of a small {\it electric} field.
(iii) It also gives the (dissipationless) {\it surface Hall conductivity}
at the surface of the crystal, provided that the surface
is insulating.  
% \aem{For symmetry, shouldn't this list also include
% surface charge?  JEM: True, like (iii) that follows from (i) and (ii), but personally I can live without it.}
  Note that (iii) follows from (ii): for a surface
with unit normal $\hat{\bf n}$ and electric field along $\E$, the
resulting surface current ${\bf K}=\M\times\hat{\bf n}$ is proportional
to $\E\times\hat{\bf n}$.  Since $\theta$ is only well-defined modulo
$2\pi$, one sees that the surface conductivity of an insulating surface is determined,
modulo the quantum $e^2/h$, by the bulk bandstructure alone.
%
%\dvm{\scr Clarify that the quantum corresponds to the possibility that
%there is a surface band with nonzero Chern number?  Maybe in an
%explanatory footnote; I prefer not to interrupt the flow of the
%text too much...}
%
There is an elegant analogy here to the case of electric polarization,
where the surface charge of an insulating surface is determined,
modulo the quantum $e/S$, by the bulk bandstructure alone
($S$ is the surface cell area).

%\dvm{\scr Orphan text: ``When an electrical field is applied, the
%surface Hall conductivity gives currents that produce an orbital
%magnetic moment parallel to the applied field.  When all surfaces
%have conductivity $\theta e^2 / (2 \pi h)$ relative to the outward
%normal, so that there are no discontinuities around an outer
%surface, (Hence $\theta$ seems more logically described as a
%magnetoelectric polarizability than a polarization~\cite{qilong}.)''\\
%Also `` ... integrating out the electrons in a one-loop gradient
%expansion~\cite{qilong} to generate the effective action for the EM
%fields (\ref{axioncoupling}).''\\
%Also ``We focus here on the second as it can be adapted to
%established DFT-based methods to compute polarization in
%interacting systems and clarifies the relationship between
%magnetoelectric polarizability and ordinary polarization.\\
%Drop these?}

The above discussion suggests several possible approaches
to the derivation of a bulk formula for the OMP $\theta$.
One would be to follow (ii) and compute the orbital
magnetization~\cite{thonhauser,niugroup} in an applied electrical
field.  We choose to focus here on (i) instead, working via
$dP/dB$.  Our derivation starts from the semiclassical wavepacket
analysis of Xiao {\it et al.}~\cite{xiao}, who analyzed how a slow
spatial variation of the local Hamiltonian in an insulator may
generate an additional polarization $\P\inh$ associated with the
inhomogeneity, beyond the spatially integrated
periodic-crystal Berry-phase polarization. For the case of an
orthorhombic 3D crystal with $M$ occupied bands in which the slow
spatial variation occurs along the $y$ direction in a supercell
of length $l_y$, they obtain
%
\beq
\ev{\Delta P_x\inh}
=  \frac{e}{4} \int_0^1 d\lambda \intbz \int_0^{l_y}\frac{dy}{l_y}
   \;\epsilon_{ijkl}\, {\rm Tr}\,[\F_{ij}\F_{kl}]
\label{eq:DPin}
\eeq
%
% \aem{The only ref I know that does this stuff without differential
% forms is Weinberg's QFT, vol. 2, but his notation has its own
% hurdles.  If we can reference to discussions in terms of forms,
% Nakahara would be standard, and the recent Fukui et al preprint has
% a nice compact discussion.}
% \dvm{I leave it to your judgment.}
%
for the change in the supercell-averaged inhomogeneously induced
polarization that occurs as a global parameter $\lambda$ evolves
adiabatically from 0 to 1.  Here indices $ijkl$ run over
$(k_x,k_y,y,\lambda)$,
%
%\beq
${\cal F}_{ij} = \partial_i {\cal A}_j - \partial_j
{\cal A}_i -i [{\cal A}_i, {\cal A}_j]$
%\eeq
%
is the Berry field-strength tensor ($\A_\lambda=i\me{u}{\del_\lambda}{u}$),
and the trace and commutator refer to band indices.

Because $\F$ is gauge-covariant,
the integrand in Eq.~(\ref{eq:DPin}) is explicitly gauge-invariant.
In fact it is just the non-Abelian second Chern class \cite{nakahara},
so that Eq.~(\ref{eq:DPin}) is path-invariant modulo a quantum
$e/a_zl_y$, where $a_z$ is the lattice constant in the $z$ direction.
Moreover, the $\lambda$ integral can be performed to obtain an
expression in terms of the non-Abelian Chern-Simons 3-form \cite{nakahara}.  Thus,
%
\beq
\ev{P_x\inh}
=  e \!\intbz \int_0^{l_y}\!\frac{dy}{l_y}\,
\epsilon_{ijk}{\rm Tr}[\A_i\del_j\A_k-\frac{2 i}{3}\A_i\A_j\A_k]
\label{eq:Pin}
\eeq
%
where $ijk$ now run only over $(k_x,k_y,y)$.  The integrand is not
gauge-invariant, and while it is not obvious that the integral is
gauge-invariant, the Chern-Simons theory implies that it is,
modulo the quantum $e/a_z l_y$.

We now apply this result to study the polarization induced by a
magnetic field, which may be regarded as arising from a spatial
inhomogeneity in the electromagnetic vector potential ${\bf A}(\r)$.
In particular, we compute $\ev{P_x\inh}$ for the case of
${\bf A}=By\hat{\bf z}$ with $B=h/ea_zl_y$, i.e., a $B$-field along
$\hat{x}$ with one flux quantum threading the supercell.  This
has the effect of taking $k_z\rightarrow k_z+eB y/\hbar$, and this is
the only way $y$ will enter into the Hamiltonian,
%\aem{slight modifications to wording}
so
that $\ket{\del_yu}=(Be/\hbar)\ket{\del_{k_z} u}$ and
%
\beq \label{PofB}
\ev{P_x\inh}=\frac{Be^2}{\hbar}\intbz \; \epsilon_{ijk}\,
{\rm Tr}[\A_i\del_j\A_k-i\frac{2}{3}\A_i\A_j\A_k]
\eeq
%
where $ijk$ now run over $(k_x,k_y,k_z)$.
%
% \dvm{Minor change to next sentence.}
%
Using Eq.~(\ref{eq:theta})
we easily arrive at Eq.~(\ref{theta}) as one of our
main results.

% \dvm{I think it might be good to have a paragraph about symmetries.
% Possibly it might go better elsewhere, but I put it here for
% the time being.\\ But there is an aspect I did not mention.
% Suppose one has a crystal with TR symmetry and $\theta=\pi$ and
% one places it in an electric field.  Then the system still has
% TR symmetry (the E-field doesn't break that) so you would think
% that the resulting magnetization and surface current must remain
% zero, in violation of $\theta=\pi$.  The escape hatch is that the
% crystal with $\theta=\pi$ has metallic surfaces unless TR is broken
% at the surfaces.  I didn't get into this here.  Should I have?
% Can we pack an explanation into a short explanatory footnote and
% attach it?  Or should we perhaps expand this paragraph and put it
% nearer the end of the paper? JEM: I agree with your point, but
% I'm not sure there is anything else I would cut to put this in.
% I modified the symmetry paragraph somewhat.  Although P and
% T act in similar ways on the EM fields, their action on the
% electronic system is different, mostly because T is antiunitary
% and P is unitary.  Probably theta=pi is associated with pairs of
% bands linked by inversion, and an inversion-trivial band can't
% contribute to theta.}

Before proceeding, we briefly discuss the symmetry properties of
the OMP.  Clearly the combination $\E\cdot\B$ in Eq.~(\ref{axioncoupling})
is odd under $T$ and also odd under inversion $P$ (although it is
even under the combination $PT$~\footnote{Both our tight-binding
model and the model in~\cite{qilong} have $PT$ symmetry, leading to twofold degeneracies even when $P$ and $T$ are
separately broken.  We have verified that the $\theta$ term is
robust to breaking of $PT$.}). It is also
odd under any improper rotation, such as a simple mirror reflection.
This implies that $\theta=-\theta$ if the crystal has {\it any}
%
% \dvm{I put the thing about $PT$ in parentheses so that this ``any''
% would not look like it could refer to it.}
%
of the
above symmetries.  This would force an ordinary coupling to vanish,
but since $\theta$ is only well-defined modulo $2\pi$, it actually
only forces $\theta=0$ or $\pi$.  Thus, one can obtain an insulator with quantized $\theta=\pi$ not only for $T$-invariant systems (regardless of whether they obey inversion
symmetry), but also for inversion- and mirror-symmetric crystals
regardless of $T$ symmetry~\cite{fu&kane2-2006}.  When none of these
symmetries is present, one generically has a non-zero (and non-$\pi$)
value of $\theta$, but restricted to the surprisingly simple scalar
form of Eq.~(\ref{axioncoupling}).

%\newpage
%\dvm{Temporary newpage}

% \dvm{OK, I stopped here.  My suggestion would actually be to
% introduce the TB model first, then use it to carry out several
% tests:\\
% 1. Calculation of $\theta$ via surface Hall conductivity.\\
% 2. Calculation of $\theta$ via 4D integral of 2nd Chern form.\\
% 3. Calculation of $\theta$ direction from Eq.~(2) using a smooth
%    gauge.\\
% Perhaps also show tests demonstrating that breaking P but not T,
% or T but not P, still leaves $\theta=\pi$ exactly, and resolving
% the question of whether a single band can have a nonzero $\theta$.\\
% I think we have a bit of a space problem, so I suspect that much
% of what was previously in the draft will have to be reduced
% significantly, especially if we want to include s.}

In the remainder of this Letter, we validate and illustrate the
above theory via numerical calculations on a tight-binding
Hamiltonian that generates non-zero values of $\theta$.  We start with
the model of Fu, Kane, and Mele~\cite{fu&kane&mele-2007} for
a 3D topological insulator on the diamond lattice,
%
% \dvm{Sentence streamlined; OK?}
%
and add a staggered Zeeman term to break $T$:
%
\begin{multline} \label{modelham}
 H = \sum_{\langle i j \rangle} t_{ij} c_i^\dagger c_j
  + i \frac{4 \lambda_{SO}}{a^2} \sum_{\langle \langle i j \rangle \rangle}
		 c_i^\dagger \bm{\sigma} \cdot (\bm{d}_{ij}^1 \times \bm{d}_{ij}^2) c_j \\
		+ \bm{h}\cdot \left( \sum_{i\in A} c_i^\dagger \bm{\sigma} c_i 
		- \sum_{i\in B} c_i^\dagger \bm{\sigma} c_i \right).
\end{multline}
%
In the first term, the nearest-neighbor hopping amplitude depends
on the bond direction; we take $t_{[111]} = 3t + m \cos \beta$,
$t_{\left[\bar{1}11\right]} = t_{\left[1\bar{1}1\right]} =
t_{\left[11\bar{1}\right]} = t$.
%
% \dvm{Next sentence rewritten for more clarity (I found the notation
% quite confusing on a first reading, although admittedly the same
% confusion was there in the FKM paper).  I don't like calling
% $\bm{s}$ an ``operator'' because the $c^\dagger_i$ and $c_j$ are
% really the operators.  They are really just Pauli spin matrices,
% right?  Any reason not to use $\bm{sigma}$ insteadof $\bm{s}$?}
%
% \aem{My notation was trying to be consistent with FKM; any changes you
% want to make are fine with me.  I tend to think of Pauli spin matrices
% as specifying a spinor basis and a real-space basis, both of which seem
% unnecessary, but if it's clearer then that's what we should do.}
% \dvm{Maybe you're right.  I guess in the original notation `s' was
% an OPERATOR in the spin sector, while $c_i^\dagger$ was the 
% operator in the spatial sector only.  So you can revert it to the
% way it was, and call `s' an operator, if you prefer.}
%
The second term describes spin-dependent hopping between pairs
of second neighbors $\langle\langle ij\rangle\rangle$, where
$\bm{d}_{ij}^1$ and $\bm{d}_{ij}^2$ are the connecting first-neighbor
legs and $\bm{\sigma}$ are the Pauli spin matrices.
%
The linear size of the face-centered cubic (fcc) conventional
cell is $a$.  The last term is the new
ingredient and describes a staggered Zeeman field with opposite
signs on the two fcc sublattices $A$ and $B$.

At half filling, with $0<m<2t$ and $\lambda_{SO}$ sufficiently
large, the original model has a
direct band gap of $2m$ at $\beta = \pi$.  We take $|\bm{h}| = m
\sin \beta$ and choose $\bm{h}$ in the [111] direction; varying the
single parameter $\beta$ keeps the gap constant and interpolates
between the ordinary ($\beta = 0$) and the topological ($\beta =
\pi$) insulator.
%
% \aem{this could go if space is short.}
%
% The Hamiltonian linearized
% about the band minimum has the form of a Dirac Hamiltonian with a
% ``complex mass'' $m (\cos \beta + i \sin \beta \gamma^5)$, which generates a coupling of the form
% (\ref{axioncoupling})~\cite{wilczekaxion}.
%Alternatively, $\beta$
%can be seen a fixed value of momentum in a fourth spatial
%dimension.~\cite{qilong}

We have calculated $\theta$ using four different methods and found that
they agree (Figure~\ref{polplot}).  Using 
established methods~\cite{MarzariVanderbilt97}, we can obtain $\theta$
from Eq.~(\ref{theta}), as shown for $\beta = \pi/4$ and $\beta = \pi/2$ (filled squares).
This computation requires finding a smooth gauge for $\A$, which is
equivalent to finding localized Wannier functions.
%width=0.45\textwidth
\begin{figure}[!ht]
 \includegraphics[width=3.2in]{thetaplot2.eps} 
 \caption{ The magnetoelectric polarizability $\theta$
 (in units of $e^2/2\pi h$).  The curve is obtained from the second Chern
 integral, Eq.~(\ref{secondchernapp}).  The filled squares are computed by 
 the Chern-Simons form, Eq.~(\ref{theta}).  The open squares are the slopes of 
 $P$ \emph{vs}. $B$, Eq.~(\ref{polarization}).  The remaining points are obtained
 by layer-resolved $\sigma_{xy}$ calculations using Eq.~(\ref{conductance}).
		\label{polplot} } 
\end{figure}

Since a tight-binding model is not computationally demanding, we have also
calculated the polarization resulting from one flux quantum
in a large supercell, using~\cite{ksv}
%
\begin{equation} \label{polarization}
 P_i = e \int_{BZ} \! \frac{d^3 k}{(2\pi)^3} \mathrm{Tr}\,\mathcal{A}_i \,.
\end{equation}
Plotting the values for a number of supercell sizes against the corresponding $B$
allows us to approximate $dP/dB$, shown as open squares in Figure~\ref{polplot}.
Further, we have
calculated the orbital magnetization response in a slab geometry, as we
describe shortly.  We have also exploited the simple
nature of the tight-binding model to compute $\theta(\beta)$ from the 
formula~\cite{xiao,qilong} (see also Eq.~(\ref{eq:DPin}))
%
% \dvm{Do you mean ``(see also Eq.~(\ref{eq:DPin}))'' ?}
%
\beq
\theta(\beta) = \frac{1}{16\pi} \int_0^\beta\,d\beta^\prime \int\,d^3{\bf k}\,\epsilon_{ijkl} \mathrm{Tr}[{\cal F}_{ij}({\bf k},\beta^\prime) {\cal F}_{kl} ({\bf k},\beta^\prime)],
\label{secondchernapp}
\eeq
%
%\aem{If we want to eliminate F, we could rewrite this.}
plotted as the curve in Figure~\ref{polplot}.
%This works because the Hamiltonian Eq.~(\ref{modelham}) has  $\theta=0$ for $\beta=0$.
%
% \dvm{Say somewhere that this is plotted as the solid curve in the fig?}
% %
The sign of $\theta$ is fixed here 
by the choice of a particular path
in Hamiltonian space from the initial $\theta=0$ value at $\beta=0$.

We now discuss the surface Hall conductivity,
possibly the most accessible experimental realization of
$\theta$. 
% \aem{I have just made a bold statement; should we tone it down?} 
% \dvm{I like it.}
Its fractional
part in units of $e^2/h$ is determined purely by bulk properties.
Consider a material with coupling $\theta$ in a slab geometry
that is finite in the ${\bf \hat z}$
direction and surrounded by $\theta=0$ vacuum.
%
% \dvm{\tiny Small changes: (1) removed 'only' in front of 'finite in z'
% since if it is infinite in $xy$ there will not be any boundaries
% between surfaces as mentioned at the end of the paragraph.
% (2) Too early to refer to Fig.~1.  From the way you cited Fig.~1
% previously, I expected it to be a sketch of the slab geometry.}
%
This leads~\cite{wilczekaxion} to surface Hall
conductivity $(\theta+2 \pi r) e^2 / (2 \pi h)$ at the top surface
and $(-\theta - 2 \pi s) e^2/(2 \pi h)$ at the bottom surface,
for some surface-dependent
integers $r$ and $s$, consistent with the total Hall
conductivity of the slab being an integer multiple of $e^2/h$.
%
% \aem{Is this any better?}
% \dvm{Fine.}
%
In a finite geometry, if the total Hall conductance $(r-s) e^2/h$ is nonzero,
there must be conducting 1D edge channels at the boundaries between
surfaces.

The 3D topological insulator has $\theta=\pi$ so that each surface
should have half-integer Hall conductance.  The spatial
contributions to the Hall conductance in the slab
geometry can be resolved as follows.  The unit cell is a supercell containing
some number $N$ of original unit cells in the ${\bf \hat z}$
direction, with translational invariance remaining in the ${\bf
\hat x}$ and ${\bf \hat y}$ directions.  The TKNN integer is~\cite{tknn,ass}
%
\beqn
C &=& \frac{1}{2\pi} \int\!d^2 k\,\sum_\nu\,{\cal F}^{\nu \nu}_{xy}
= \frac{i}{2\pi} \int\!d^2 k\,\sum_\nu \epsilon_{ij}
\ip{\partial_i u_\nu}{\partial_j u_\nu} \cr
&=&
\frac{i}{2\pi} \int\!d^2 k\,{\rm Tr} \left[\mathcal{P} \epsilon_{ij}
\partial_i \mathcal{P} \partial_j \mathcal{P} \right].
\eeqn
%
Here $\nu$ runs over occupied bands, $|u_\nu \rangle$ are the
Bloch wavefunctions, $i$ and $j$ take the values $k_x$ and $k_y$,
and $\mathcal{P}=\sum_\nu |u_\nu\rangle \langle u_\nu|$ is
the projection operator onto the occupied subspace, which is
manifestly gauge-invariant.
%The Berry field strength tensor is
%
%\beq
%{\cal F}_{ij} = \partial_i {\cal A}_j - \partial_j {\cal A}_i -i [{\cal A}_i, {\cal A}_j],\quad {\cal A}^{\alpha %\beta}_i = i \langle u_\alpha | \partial_i | u_\beta \rangle.
%\eeq
%
%Script letters distinguish the Berry connection and field strength
%${\cal A}$ and ${\cal F}$, which are in general nonabelian, from
%electromagnetic $A$ and $F$. 
To know how different ${\bf
\hat z}$ layers contribute to the TKNN number, we define
a real-space projection $\tilde{\mathcal{P}}_n$ onto layer
$n$ within the slab supercell, and compute
%
\beq
C(n) = \frac{i}{2\pi} \int\,d^2 k\,{\rm Tr} \left[\mathcal{P} \epsilon_{ij}
(\partial_i \mathcal{P}) \tilde{\mathcal{P}}_n (\partial_j \mathcal{P})
\right]. \label{conductance}
\eeq



The results, presented in Fig.~\ref{currplot}, confirm that the surface
layers have half-integer Hall conductance when $\beta = \pi$ in
(\ref{modelham}) and that the sign on each surface is
switched by local $T$-breaking perturbations.
This suggests a possible experiment on the
3D topological insulator: applying
a weak uniform magnetic field in almost any direction to a crystal
or film leads to an integer quantum Hall state with the sign
depending on the field direction, $\sigma_{xy} = \pm e^2/h$.  The
state with $r-s=0$ does not occur because the magnetic field acts
oppositely at the two surfaces.

% \aem{I've hacked away here pretty extensively.  Let me know if
% it hangs together any better than before.}
% \dvm{Good.}
%
The ambiguity in the surface conductance (the integer $r$ above) indicates
that there may be one or more integer quantum-Hall layers at each surface.
It is the same as the ambiguity of specifying the bulk property $\theta$, 
discussed in connection with Eq.~(\ref{eq:Pin}).  This ambiguity also
appears in the supercell polarization in small field:
%  These surface layers result from bulk physics since for each $\beta$, the fractional part of the surface Hall conductance in the calculation above is equal to the bulk magnetoelecric polarizability.  The bulk polarization is~\cite{ksv}
% %
% \begin{equation} \label{polarization}
%  P_i = e \int_{BZ} \! \frac{d^3 k}{(2\pi)^3} \mathrm{Tr}\,\mathcal{A}_i \,.
% \end{equation}
%
Eq.~(\ref{polarization}) only gives the polarization
(and hence differences in polarization) modulo the ``quantum of
polarization'' \cite{ksv} which, for the case of the
flux-threaded supercell of Eq.~(\ref{PofB}),
is $\Delta P_x = e /a_z l_y$.
%  This ambiguity is connected
%with that of the OMP as follows.
 Since the magnetic field
is $B_x = h/e a_z l_y$,
%
\begin{equation}
 \Delta \frac{P_x}{B_x} = \frac{e^2}{h} = 2\pi \frac{e^2}{2\pi h}.
\end{equation}
%
Therefore, $\theta$ has an associated ambiguity of $2\pi$,
in agreement with the predictions of axion electrodynamics.

\begin{figure}
  \includegraphics[width=3.2in]{currentpi.eps}
  \caption{(Color online)
   The layer-resolved Hall conductivity (in units of $e^2/h$) at
   $\beta = \pi$ in a slab of twenty layers, with $m=t/2$ and
   $\lambda_{SO} = t/4$, terminated in $(\bar{1}11)$ planes.
	 Since there are gapless surface states at
   $\beta = \pi$, a surface term was added to Eq.~(\ref{modelham}) that breaks $T$ locally and
   gaps the surface states.  Switching the sign of this
   perturbation at one surface (the two curves) switches the current direction
   at that surface, causing a transition between integer quantum Hall
   states.
		}  \label{currplot}
\end{figure}


% This confirms the equivalence of the second Chern definition to
% previous definitions of the 3D topological insulator.

%Showing this equivalence for one model is essentially sufficient to
%show equivalence in general, since both the $\mathbb{Z}_2$ integral
%in $T$-invariant planes and second Chern integral are homotopy
%invariants.

Finally, using Eq.~(\ref{theta}) we have calculated $\theta$ for a
Hamiltonian that breaks $PT$ (as well as $P$ and $T$) by adding a weak,
uniform (\emph{i.e.}, not staggered) Zeeman coupling.  For some values
of $\beta$ this lifts all degeneracies, allowing separation of the 
intraband and interband contributions to $\theta$.  A single band can have 
nonzero $\theta$ as long as there are more than two bands in total~\cite{hopfinsulator}.
%\aem{ref. Hopf paper for comparison?}
Because
interband contributions are nonzero in general, $\theta$ is 
a property of the whole occupied spectrum and not a sum of individual band
contributions, unlike the polarization.
%\aem{In this sense, it's different from
%polarization, which only depends on the diagonal elements of A and so
%is additive.}
%
% \dvm{I find the next sentence to be very confusing.  What, exactly, do
% you mean by a band structure that separatesinto independent two-band
% blocks?  You mean a band structure generated by a tight-binding model
% that is separable in some sense?  But what is the correct general
% statement that you want to make, outside the restrictions of a
% tight-binding framework?  Either make it more precise, or drop it?}
%
% Since $\theta = 0$ for any band structure that separates into
% independent two-band blocks with one filled band in each block,
%The magnetoelectric polarizability can be viewed as arising from twisting of occupied bands around each other in the Bloch Hilbert space.

%To understand this, note
%that any two-band insulator with one occupied band has $\theta = 2
%\pi n$ and hence zero magnetoelectric polarizability by our
%definition~\cite{hopfinsulator}.   The same applies to any band
%structure that can be decomposed into independent two-band blocks
%with each block having one occupied band.  While symmetry
%alone would suggest that any T-breaking insulator without
%fine-tuning will have a nonzero $\theta$, occupied bands must twist
%around each other in the Bloch Hilbert space to generate a nonzero
%magnetoelectric polarizability.  In more precise language, the
%magnetoelectric polarizability depends strongly on the nonabelian
%nature of the $U(N)$ vector bundle with $N$ occupied bands, and is
%zero for the case of independent Abelian $U(1)$ bundles.

The weak-field-driven integer quantum Hall
transition in a film of 3D topological insulator is one
experimental signature that $\theta = \pi$.
For general $\theta$, optical conductivity measurement of a single
surface's Hall conductivity may be feasible.
A more
exotic signature is that a magnetic monopole acquires
electrical charge in a $\theta$-vacuum~\cite{witteneffect}.  Recently ``spin ice'' materials such as the
pyrochlore magnet Dy$_2$Ti$_2$O$_7$ have been argued to
support monopole excitations~\cite{castelnovo}; if the background
material has nonzero $\theta$, such monopoles
carry fractional electrical charge $\theta e / 2 \pi$.

%A question for future work is whether other three-dimensional topological invariants~\cite{schnyder} can also show magnetoelectric or related effects.

% \dvm{Added my grant too.}

The authors acknowledge useful discussions with A.~Selem and I.~Souza.
The work was supported by the Western Institute of Nanoelectronics
(AME), NSF DMR-0804413 (JEM), and NSF DMR-0549198 (DV).


%\bibliography{../bigbib}
\begin{thebibliography}{26}
\expandafter\ifx\csname natexlab\endcsname\relax\def\natexlab#1{#1}\fi
\expandafter\ifx\csname bibnamefont\endcsname\relax
  \def\bibnamefont#1{#1}\fi
\expandafter\ifx\csname bibfnamefont\endcsname\relax
  \def\bibfnamefont#1{#1}\fi
\expandafter\ifx\csname citenamefont\endcsname\relax
  \def\citenamefont#1{#1}\fi
\expandafter\ifx\csname url\endcsname\relax
  \def\url#1{\texttt{#1}}\fi
\expandafter\ifx\csname urlprefix\endcsname\relax\def\urlprefix{URL }\fi
\providecommand{\bibinfo}[2]{#2}
\providecommand{\eprint}[2][]{\url{#2}}

\bibitem[{\citenamefont{Wilczek}(1987)}]{wilczekaxion}
\bibinfo{author}{\bibfnamefont{F.}~\bibnamefont{Wilczek}},
  \bibinfo{journal}{Phys. Rev. Lett.} \textbf{\bibinfo{volume}{58}},
  \bibinfo{pages}{1799} (\bibinfo{year}{1987}).

\bibitem[{\citenamefont{Thouless et~al.}(1982)\citenamefont{Thouless, Kohmoto,
  Nightingale, and den Nijs}}]{tknn}
\bibinfo{author}{\bibfnamefont{D.~J.} \bibnamefont{Thouless}},
  \bibinfo{author}{\bibfnamefont{M.}~\bibnamefont{Kohmoto}},
  \bibinfo{author}{\bibfnamefont{M.~P.} \bibnamefont{Nightingale}},
  \bibnamefont{and} \bibinfo{author}{\bibfnamefont{M.}~\bibnamefont{den Nijs}},
  \bibinfo{journal}{Phys. Rev. Lett.} \textbf{\bibinfo{volume}{49}},
  \bibinfo{pages}{405} (\bibinfo{year}{1982}).

\bibitem[{\citenamefont{Kane and Mele}(2005{\natexlab{a}})}]{km2}
\bibinfo{author}{\bibfnamefont{C.~L.} \bibnamefont{Kane}} \bibnamefont{and}
  \bibinfo{author}{\bibfnamefont{E.~J.} \bibnamefont{Mele}},
  \bibinfo{journal}{Phys. Rev. Lett.} \textbf{\bibinfo{volume}{95}},
  \bibinfo{pages}{226801} (\bibinfo{year}{2005}{\natexlab{a}}).

\bibitem[{\citenamefont{Kane and Mele}(2005{\natexlab{b}})}]{kane&mele-2005}
\bibinfo{author}{\bibfnamefont{C.~L.} \bibnamefont{Kane}} \bibnamefont{and}
  \bibinfo{author}{\bibfnamefont{E.~J.} \bibnamefont{Mele}},
  \bibinfo{journal}{Phys. Rev. Lett.} \textbf{\bibinfo{volume}{95}},
  \bibinfo{pages}{146802} (\bibinfo{year}{2005}{\natexlab{b}}).

\bibitem[{\citenamefont{{Bernevig} et~al.}(2006)\citenamefont{{Bernevig},
  {Hughes}, and {Zhang}}}]{zhangscience1}
\bibinfo{author}{\bibfnamefont{B.~A.} \bibnamefont{{Bernevig}}},
  \bibinfo{author}{\bibfnamefont{T.~L.} \bibnamefont{{Hughes}}},
  \bibnamefont{and} \bibinfo{author}{\bibfnamefont{S.-C.}
  \bibnamefont{{Zhang}}}, \bibinfo{journal}{Science}
  \textbf{\bibinfo{volume}{314}}, \bibinfo{pages}{1757} (\bibinfo{year}{2006}).

\bibitem[{\citenamefont{Moore and Balents}(2007)}]{moore&balents-2006}
\bibinfo{author}{\bibfnamefont{J.~E.} \bibnamefont{Moore}} \bibnamefont{and}
  \bibinfo{author}{\bibfnamefont{L.}~\bibnamefont{Balents}},
  \bibinfo{journal}{Phys. Rev. B} \textbf{\bibinfo{volume}{75}},
  \bibinfo{pages}{121306} (\bibinfo{year}{2007}).

\bibitem[{\citenamefont{Roy}()}]{rroy3d}
\bibinfo{author}{\bibfnamefont{R.}~\bibnamefont{Roy}},
  \bibinfo{howpublished}{cond-mat/0607531}.

\bibitem[{\citenamefont{Fu et~al.}(2007)\citenamefont{Fu, Kane, and
  Mele}}]{fu&kane&mele-2007}
\bibinfo{author}{\bibfnamefont{L.}~\bibnamefont{Fu}},
  \bibinfo{author}{\bibfnamefont{C.~L.} \bibnamefont{Kane}}, \bibnamefont{and}
  \bibinfo{author}{\bibfnamefont{E.~J.} \bibnamefont{Mele}},
  \bibinfo{journal}{Phys. Rev. Lett.} \textbf{\bibinfo{volume}{98}},
  \bibinfo{pages}{106803} (\bibinfo{year}{2007}).

\bibitem[{\citenamefont{Hsieh et~al.}(2008)\citenamefont{Hsieh, Qian, Wray,
  Xia, Hor, Cava, and Hasan}}]{hsieh}
\bibinfo{author}{\bibfnamefont{D.}~\bibnamefont{Hsieh}},
  \bibinfo{author}{\bibfnamefont{D.}~\bibnamefont{Qian}},
  \bibinfo{author}{\bibfnamefont{L.}~\bibnamefont{Wray}},
  \bibinfo{author}{\bibfnamefont{Y.}~\bibnamefont{Xia}},
  \bibinfo{author}{\bibfnamefont{Y.~S.} \bibnamefont{Hor}},
  \bibinfo{author}{\bibfnamefont{R.~J.} \bibnamefont{Cava}}, \bibnamefont{and}
  \bibinfo{author}{\bibfnamefont{M.~Z.} \bibnamefont{Hasan}},
  \bibinfo{journal}{Nature} \textbf{\bibinfo{volume}{452}},
  \bibinfo{pages}{970} (\bibinfo{year}{2008}).

\bibitem[{\citenamefont{Volkov and Pankratov}(1985)}]{pankratov}
\bibinfo{author}{\bibfnamefont{B.~A.} \bibnamefont{Volkov}} \bibnamefont{and}
  \bibinfo{author}{\bibfnamefont{O.~A.} \bibnamefont{Pankratov}},
  \bibinfo{journal}{JETP Lett.} \textbf{\bibinfo{volume}{42}},
  \bibinfo{pages}{178} (\bibinfo{year}{1985}).

\bibitem[{\citenamefont{Fradkin et~al.}(1986)\citenamefont{Fradkin, Dagotto,
  and Boyanovsky}}]{fradkinaxion}
\bibinfo{author}{\bibfnamefont{E.}~\bibnamefont{Fradkin}},
  \bibinfo{author}{\bibfnamefont{E.}~\bibnamefont{Dagotto}}, \bibnamefont{and}
  \bibinfo{author}{\bibfnamefont{D.}~\bibnamefont{Boyanovsky}},
  \bibinfo{journal}{Phys. Rev. Lett.} \textbf{\bibinfo{volume}{57}},
  \bibinfo{pages}{2967} (\bibinfo{year}{1986}).

\bibitem[{\citenamefont{Fu and Kane}(2007)}]{fu&kane2-2006}
\bibinfo{author}{\bibfnamefont{L.}~\bibnamefont{Fu}} \bibnamefont{and}
  \bibinfo{author}{\bibfnamefont{C.~L.} \bibnamefont{Kane}},
  \bibinfo{journal}{Phys. Rev. B} \textbf{\bibinfo{volume}{76}},
  \bibinfo{pages}{045302} (\bibinfo{year}{2007}).

\bibitem[{\citenamefont{Qi et~al.}()\citenamefont{Qi, Hughes, and
  Zhang}}]{qilong}
\bibinfo{author}{\bibfnamefont{X.}~\bibnamefont{Qi}},
  \bibinfo{author}{\bibfnamefont{T.~L.} \bibnamefont{Hughes}},
  \bibnamefont{and} \bibinfo{author}{\bibfnamefont{S.-C.} \bibnamefont{Zhang}},
  \bibinfo{howpublished}{arxiv:0802.3537}.

\bibitem[{\citenamefont{Spaldin and Fiebig}(2005)}]{spaldin-s05}
\bibinfo{author}{\bibfnamefont{N.~A.} \bibnamefont{Spaldin}} \bibnamefont{and}
  \bibinfo{author}{\bibfnamefont{M.}~\bibnamefont{Fiebig}},
  \bibinfo{journal}{Science} \textbf{\bibinfo{volume}{309}},
  \bibinfo{pages}{391} (\bibinfo{year}{2005}).

\bibitem[{\citenamefont{Fiebig}(2005)}]{fiebig-jpdp05}
\bibinfo{author}{\bibfnamefont{M.}~\bibnamefont{Fiebig}}, \bibinfo{journal}{J.
  Phys. D-Applied Phys.} \textbf{\bibinfo{volume}{38}}, \bibinfo{pages}{R123}
  (\bibinfo{year}{2005}).

\bibitem[{\citenamefont{Xiao et~al.}()\citenamefont{Xiao, Shi, Clougherty, and
  Niu}}]{xiao}
\bibinfo{author}{\bibfnamefont{D.}~\bibnamefont{Xiao}},
  \bibinfo{author}{\bibfnamefont{J.}~\bibnamefont{Shi}},
  \bibinfo{author}{\bibfnamefont{D.~P.} \bibnamefont{Clougherty}},
  \bibnamefont{and} \bibinfo{author}{\bibfnamefont{Q.}~\bibnamefont{Niu}},
  \bibinfo{howpublished}{arxiv:0711.1855}.

\bibitem[{\citenamefont{King-Smith and Vanderbilt}(1993)}]{ksv}
\bibinfo{author}{\bibfnamefont{R.~D.} \bibnamefont{King-Smith}}
  \bibnamefont{and}
  \bibinfo{author}{\bibfnamefont{D.}~\bibnamefont{Vanderbilt}},
  \bibinfo{journal}{Phys. Rev. B} \textbf{\bibinfo{volume}{47}},
  \bibinfo{pages}{1651} (\bibinfo{year}{1993}).

\bibitem[{\citenamefont{Ortiz and Martin}(1994)}]{ortizmartin}
\bibinfo{author}{\bibfnamefont{G.}~\bibnamefont{Ortiz}} \bibnamefont{and}
  \bibinfo{author}{\bibfnamefont{R.~M.} \bibnamefont{Martin}},
  \bibinfo{journal}{Phys. Rev. B} \textbf{\bibinfo{volume}{49}},
  \bibinfo{pages}{14202} (\bibinfo{year}{1994}).

\bibitem[{\citenamefont{Thonhauser et~al.}(2005)\citenamefont{Thonhauser,
  Ceresoli, Vanderbilt, and Resta}}]{thonhauser}
\bibinfo{author}{\bibfnamefont{T.}~\bibnamefont{Thonhauser}},
  \bibinfo{author}{\bibfnamefont{D.}~\bibnamefont{Ceresoli}},
  \bibinfo{author}{\bibfnamefont{D.}~\bibnamefont{Vanderbilt}},
  \bibnamefont{and} \bibinfo{author}{\bibfnamefont{R.}~\bibnamefont{Resta}},
  \bibinfo{journal}{Phys. Rev. Lett.} \textbf{\bibinfo{volume}{95}},
  \bibinfo{pages}{137205} (\bibinfo{year}{2005}).

\bibitem[{\citenamefont{Xiao et~al.}(2005)\citenamefont{Xiao, Shi, and
  Niu}}]{niugroup}
\bibinfo{author}{\bibfnamefont{D.}~\bibnamefont{Xiao}},
  \bibinfo{author}{\bibfnamefont{J.}~\bibnamefont{Shi}}, \bibnamefont{and}
  \bibinfo{author}{\bibfnamefont{Q.}~\bibnamefont{Niu}},
  \bibinfo{journal}{Phys. Rev. Lett.} \textbf{\bibinfo{volume}{95}},
  \bibinfo{pages}{137204} (\bibinfo{year}{2005}).

\bibitem[{\citenamefont{Nakahara}(1998)}]{nakahara}
\bibinfo{author}{\bibfnamefont{M.}~\bibnamefont{Nakahara}},
  \emph{\bibinfo{title}{{Geometry, Topology and Physics}}}
  (\bibinfo{publisher}{Institute of Physics}, \bibinfo{year}{1998}).

\bibitem[{\citenamefont{Marzari and Vanderbilt}(1997)}]{MarzariVanderbilt97}
\bibinfo{author}{\bibfnamefont{N.}~\bibnamefont{Marzari}} \bibnamefont{and}
  \bibinfo{author}{\bibfnamefont{D.}~\bibnamefont{Vanderbilt}},
  \bibinfo{journal}{Phys. Rev. B} \textbf{\bibinfo{volume}{56}},
  \bibinfo{pages}{12847} (\bibinfo{year}{1997}).

\bibitem[{\citenamefont{Avron et~al.}(1983)\citenamefont{Avron, Seiler, and
  Simon}}]{ass}
\bibinfo{author}{\bibfnamefont{J.~E.} \bibnamefont{Avron}},
  \bibinfo{author}{\bibfnamefont{R.}~\bibnamefont{Seiler}}, \bibnamefont{and}
  \bibinfo{author}{\bibfnamefont{B.}~\bibnamefont{Simon}},
  \bibinfo{journal}{Phys. Rev. Lett.} \textbf{\bibinfo{volume}{51}},
  \bibinfo{pages}{51} (\bibinfo{year}{1983}).

\bibitem[{\citenamefont{Moore et~al.}()\citenamefont{Moore, Ran, and
  Wen}}]{hopfinsulator}
\bibinfo{author}{\bibfnamefont{J.~E.} \bibnamefont{Moore}},
  \bibinfo{author}{\bibfnamefont{Y.}~\bibnamefont{Ran}}, \bibnamefont{and}
  \bibinfo{author}{\bibfnamefont{X.-G.} \bibnamefont{Wen}},
  \bibinfo{howpublished}{arXiv:0804.4527, to appear in Phys. Rev. Lett.}

\bibitem[{\citenamefont{Witten}(1979)}]{witteneffect}
\bibinfo{author}{\bibfnamefont{E.}~\bibnamefont{Witten}},
  \bibinfo{journal}{Phys. Lett. B} \textbf{\bibinfo{volume}{86}},
  \bibinfo{pages}{283} (\bibinfo{year}{1979}).

\bibitem[{\citenamefont{Castelnovo et~al.}(2008)\citenamefont{Castelnovo,
  Moessner, and Sondhi}}]{castelnovo}
\bibinfo{author}{\bibfnamefont{C.}~\bibnamefont{Castelnovo}},
  \bibinfo{author}{\bibfnamefont{R.}~\bibnamefont{Moessner}}, \bibnamefont{and}
  \bibinfo{author}{\bibfnamefont{S.~L.} \bibnamefont{Sondhi}},
  \bibinfo{journal}{Nature} \textbf{\bibinfo{volume}{451}}, \bibinfo{pages}{42}
  (\bibinfo{year}{2008}).

\end{thebibliography}

\end{document}
\begin{table}
\begin{tabular*}{3 in}{@{\extracolsep{\fill}}  c  c  c }
&Polarization&Magnetoelectric\cr
&&polarizability\cr
\hline\hline\\
$d_{\rm min}$&1&3\\
%2 \log (2) / 3 \approx
Observable & ${\bf P} = {\partial \langle H \rangle / \partial E}$ & $M_{ij} = \partial \langle H \rangle / \partial E_i \partial B_j$\cr
&&$= \delta_{ij} \theta e^2/(2 \pi h)$\cr
Quantum& $\Delta {\bf P} = e {\bf R}/ \Omega$ &$\Delta M = e^2/h$ \cr
Surface&$q=({\bf P}_1 - {\bf P}_2) \cdot {\bf \hat n}$&$\sigma_{xy} = (M_1 - M_2)$\cr
EM coupling & ${\bf P} \cdot {\bf E}$ & $M {\bf E} \cdot {\bf B} $\cr
CS form&${\cal A}_i$&
$\epsilon_{ijk} ({\cal A}_i {\cal F}_{jk} + i {\cal A}_i {\cal A}_j {\cal A}_k/3)$\cr
Chern form& $\epsilon_{ij} \partial_i {\cal A}_j $&$ \epsilon_{ijkl} {\cal F}_{ij} {\cal F}_{kl}$\cr

\end{tabular*}
\caption{Comparison of Berry-phase theories of polarization and magnetoelectric polarizability.}
\end{table}


