%% SELECT THE LAYOUT
%% The class supports further options.
%% See aipguide.pdf for details.
%% FRONTMATTER


\documentclass[
    ,final              ]{aipproc}
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
\usepackage{graphicx}
\usepackage{amsmath}
\usepackage{amsfonts}
\usepackage{amssymb}

\setcounter{MaxMatrixCols}{10}
%TCIDATA{OutputFilter=Latex.dll}
%TCIDATA{Version=5.50.0.2953}
%TCIDATA{<META NAME="SaveForMode" CONTENT="1">}
%TCIDATA{BibliographyScheme=Manual}
%TCIDATA{LastRevised=Thursday, January 08, 2009 15:40:31}
%TCIDATA{<META NAME="GraphicsSave" CONTENT="32">}
%TCIDATA{Language=American English}

\layoutstyle{8x11double}
\input{tcilatex}

\begin{document}

\title{Implementation of decoherence-protected logical qubits in medium-size
Josephson-junction arrays.}
\keywords{Topological protection, Josephson junction arrays, qubits}
\author{B. Dou\c{c}ot}
{address={Laboratoire de Physique Th\'{e}orique et Hautes \'Energies, CNRS
UMR 7589, Universit\'{e}s Paris 6 et 7, 4, place Jussieu, 75252 Paris Cedex
05 France} }
\author{L. B. Ioffe}
{address={Center for Materials Theory, Department of Physics and
Astronomy, Rutgers University 136 Frelinghuysen Rd., Piscataway NJ 08854 USA}
} 
\begin{abstract}
We outline the physical principles underlying the design of finite size
superconducting arrays that may serve as protected qubits. These systems are
built from fully frustrated rhombi consisting of four Josephson junctions;
their low energy spectrum is composed of localized $Z_{2}$ charge
excitations and $\pi$ phase vortices. Protection from both the phase and
flip errors requires a balance between charging and Josephson energies;
optimization of the array size results in the hierarchical construction of
these arrays. We have developed computational methods, based on a
combination of perturbation theory and numerical diagonalization, which
allow us to compute the global phase stiffness and the energy of localized
charge excitations in various geometries. We present the most promising
designs, which are currently under experimental investigation.
\end{abstract}

\classification{03.67.Lx,74.81.Fa}


\maketitle





%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%% MAINMATTER
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

\section{Introduction}

The goal of this paper is to identify and to compute the properties of the
small to medium sized Josephson junction arrays which spectrum contains two
low energy states that are protected from the external noise and can be used
as protected qubits for quantum computation. The elementary building block
of all arrays considered here is the four junction superconducting loop
threaded by the magnetic flux $\Phi =\Phi _{0}/2$ that we shall call rhombus
in the following, see Fig.~\ref{TwoRhombi}. \ The effective Josephson energy 
$V(\theta )$ of the \emph{ideal} rhombus is a $\pi $- periodic function of
the phase difference across it; most of our discussion below is applicable
to any element that has \emph{approximately} this property.  In the case of
rhombus, $V(\theta )$ has two minima at $\theta =\pm \pi /2$, separated by $%
\pi $ modulo $2\pi $. In the absence of transitions between these states
induced by quantum fluctuations the ground state of each element is doubly
degenerate.  The $\pi $- periodicity for $V(\theta )$ can be interpreted in
simple physical terms by noting that Fourier decomposition of $V(\theta
)=\sum_{n}A_{n}e^{in\theta }$ contains only even harmonics, that is $A_{n}=0$
for all odd values of $n$. The absence of all odd harmonics and, in
particular, the suppression of the fundamental oscillation ($n=\pm 1$) means
that single Cooper pairs cannot tunnel across the individual element and
thus are localized. For the element realized as a rhombus threaded by half a
flux quantum, the localization occurs due to the cancellation between the
probability amplitudes for a Cooper pair to follow the upper and the lower
arms of the rhombus while tunneling across. A single rhombus is \emph{not}
protected against local noise: flux deviations from the $\Phi _{0}/2$,
scatter of junction parameters lead to non-zero tunneling amplitude of a
single pair (odd harmonics $A_{n}\neq 0$) which lifts the degeneracy of the
two states. Fluctuations of this amplitude lead to the fast dephasing of
such elementary qubit.

This description of a single $\pi $- periodic element and its dephasing
suggests that a stronger protection against flux noise may be achieved by
composing $N$ such elements in a series which makes the single pair
tunneling exponentially small. Indeed, the global state of the system is
described by the phase difference $\phi $ between the ends of the chain
which is equal to the sum $\phi =\sum_{j=1}^{N}\theta _{j}$ \ of the phase
differences, $\theta _{j}$ across individual elements. The energy of the
system is also a $\pi $- periodic function of $\phi $ for the chain composed
of ideal elements. Let us first assume that the quantum fluctuations of $%
\phi $ are small, which can be achieved in practice by introducing a large
capacitance between the two ends of the chain. In this case, the system has
two global states corresponding to two minima of the total energy as a
function of the phase $\phi $. Because the absolute value of the phase is
not measurable these two states become ideally protected from the
decoherence in the limit $N\rightarrow \infty $ provided that the energy
barrier separating them remains finite. For a finite $L$ this protection is
only as good as a suppression of a single pair tunneling along the whole
chain. This suppression is possible thanks to a remarkable property of arrays
composed of $\pi $- periodic elements: in these arrays the parity of the
number of Cooper pairs on each node joining two or more elements becomes a
good quantum number in the absence of external noise or any static
imperfection. As will be explained later, these local conserved quantities
can be viewed as the eigenvalues of local $Z_{2}$ symmetries generators that
translate the phase at a given node by $\pi $ while leaving phases on other
nodes unchanged.

The protection provided by a single chain depends crucially on the ratio $E_{%
\mathrm{C}}/E_{\mathrm{J}}$ of the charging energy $E_{\mathrm{C}}$ to
Josephson energy $E_{\mathrm{J}}$. In the limit when the charging energy $E_{%
\mathrm{C}}$ is much smaller than the Josephson energy $E_{\mathrm{J}}$, the
quantum fluctuations of phase variables are small and therefore, the energy
difference between the even and odd parity states on a given node is small
as well. In this limit the noise and imperfections have large effect on
these almost degenerate states and can transfer the Cooper pair across. The
increase in $E_{\mathrm{C}}/E_{\mathrm{J}}$ increases quantum phase
fluctuations, and yields a finite energy gap $\Delta $ between the states
with different parities. In the presence of flux noise of typical magnitude $%
\delta \Phi $ per rhombus, a single Cooper pair tunneling amplitude across
each rhombus of order $\frac{\delta \Phi }{\Phi _{0}}E_{\mathrm{J}}$ is
generated. The dangerous process that transfers a single Cooper pair across
the whole system arises only at order $N$ in perturbation theory, therefore
its amplitude is proportional to $(\frac{\delta \Phi }{\Phi _{0}}E_{\mathrm{J%
}})^{N}/\Delta ^{N-1}$, which can be made arbitrarily small by choosing $N$
large enough as long as \mbox{$\frac{\delta\Phi}{\Phi_{0}}\frac{E_{\mathrm{J}}}{%
\Delta}<1$}. This basic protection mechanism has been recently demonstrated
experimentally by the observation of a strong reduction of the single Cooper
pair tunneling amplitude for a flux per rhombus $\Phi $ close to $\Phi _{0}/2
$ in a system with $N=4$, compared to one with $N=2$.~\cite{Gladchenko08}

So far, we discussed how to achieve protection against processes that lift
the degeneracy between the two classical ground states identified by the
value of the phase variable $\phi $. We now turn to the tunneling processes
between the classical ground states which are made possible by the finite
capacitance across the chain. Their effect is to lift the degeneracy in the
ground-state doublet and to select the basis in which the operator giving
parity of the number of Cooper pairs transferred across the whole chain is
diagonal. The amplitude of this process is sensitive to the electrostatic
potential across the whole chain which might fluctuate. This would lead to
the fluctuations of the energy splitting and thus dephasing. As we already
pointed out, the tunneling might be inhibited by large capacitance between
the ends of the chain. However, such capacitance would also imply that
corresponding charging energy is small which would make it difficult to
avoid thermal activation in the realistic setup. It would be very desirable
to reduce the tunneling between classical ground states exponentially fast
in the system size, as for the single Cooper pair tunneling along the whole
chain discussed above. Note that a process that translates $\phi $ by $\pi $
may be viewed as a single $\pi $ vortex crossing the chain. A natural way to
suppress such events is to couple $M$ identical $Z_{2}$ chains in parallel.
As we shall see, there is to some extent a conflict between the two
requirements of imposing simultaneously a good localization of single Cooper
pairs through a large gap $\Delta $ for a $Z_{2}$ charge in the bulk of the
array, and a large energy barrier for a $\pi $ vortex crossing: the former
requires rather large quantum phase fluctuations  (thus large $E_{\mathrm{C}%
}/E_{\mathrm{J}}$) and large chain length $N$ whereas the latter needs a
large phase stiffness so that $E_{\mathrm{C}}/E_{\mathrm{J}}$ and the chain
length cannot be too large. To resolve this conflict we shall use the
hierarchical construction which uses iteratively series and parallel
compositions of arrays whose building blocks are rhombi threaded by half a
flux quantum.

Note that whereas it is possible to suppress completely single Cooper pair
tunneling in ideal conditions in a single chain by having exactly half a
flux quantum and a perfect mirror symmetry along the diagonal for each
rhombus, reduction of the $\pi $ vortex tunneling amplitude is only achieved
by choosing a low enough value of $E_{\mathrm{C}}/E_{\mathrm{J}}$ and
remains non-zero for any non-zero value of this parameter. In other words,
for these arrays there is no analogue of the local $Z_{2}$ symmetry that
would prevent single $\pi $ vortex tunneling; these arrays do not implement
exactly the self dual model proposed by Kitaev in his seminal paper on
topological quantum computation.~\cite{Kitaev97} Nevertheless, as we shall
demonstrate, the proper choice of the array geometry that involves parallel
chains of rhombi and a proper choice of the $E_{\mathrm{C}}/E_{\mathrm{J}}$
ratio allow for the exponentially high protection against any local noise.

The protection provided by the system of parallel chains of $\pi $- periodic
elements can be understood as the Hamiltonian version of the error
correction based on the repetition code. In this code one uses $N\times M$
faulty qubits to encode one logical qubit. At the first stage of this
construction one forms $2$ classical states $|C0>$ in each of the $N$
qubits: $|C0>=|00...0>$ and $|C1>=|11...1>$. This classical repetition
protects against $X$ errors that change $0$ into $1$ or backwards. At the
second stage one forms quantum superposition of these two states $|\pm
>=(|C0>\pm |C1>)$ in each of these $M$ rows and takes the direct product of
the resulting states: 
\begin{equation}
\left\vert \pm \right\rangle =\left( 
\begin{array}{c}
\frac{1}{\sqrt{2}}(\left\vert 0000\right\rangle _{1}\pm \left\vert
1111\right\rangle _{1})\otimes  \\ 
\frac{1}{\sqrt{2}}(\left\vert 0000\right\rangle _{2}\pm \left\vert
1111\right\rangle _{2})\otimes  \\ 
\frac{1}{\sqrt{2}}(\left\vert 0000\right\rangle _{3}\pm \left\vert
1111\right\rangle _{3})\otimes  \\ 
\frac{1}{\sqrt{2}}(\left\vert 0000\right\rangle _{4}\pm \left\vert
1111\right\rangle _{4})%
\end{array}%
\right)   \label{RepetitionCode}
\end{equation}%
This step protects against $Z$-errors that change the relative phase of
individual $|0>$ and $|1>$ states and, thus, flip $|\pm >$ states.

In order to establish the correspondence between the array of Josephson
junctions and the array of qubits we identify the rhombus states $\left\vert
\pm \pi /2\right\rangle $ with the states of spin $1/2$ with $S_{z}=\pm 1/2$%
, the low energy states of the rhombi array are equivalent to the states of
the spin array of the same size. The actual Josephson circuit has also high
energy modes with $E_{p}\sim \sqrt{8E_{J}E_{\mathrm{C}}}$, which are separated from
the low energy states by a large gap and thus can be mostly ignored.
Identifying $\left\vert 0\right\rangle $ ($\left\vert 1\right\rangle $) of
the faulty qubit with the even (odd) states of individual rhombus ( $%
\left\vert e\right\rangle =(\left\vert +\pi /2\right\rangle +\left\vert -\pi
/2\right\rangle )/\sqrt{2}$, $\left\vert o\right\rangle =(\left\vert +\pi
/2\right\rangle -\left\vert -\pi /2\right\rangle )/\sqrt{2}$), we see that $%
(\left\vert 0000\right\rangle _{1}\pm \left\vert 1111\right\rangle _{1})/2$
state of a qubit chain corresponds to state rhombi chain in which all rhombi
are in even or odd states with the phase difference across the chain $\pm
\pi /2$ respectively.

The idea of the physical array construction is therefore to ensure that its
ground state corresponds to the 'ferromagnetic' state of one chain in the
charge (odd/even) basis and to the equal phase difference across different
the chains. The former is ensured by the Hamiltonian $H=\sum_{(i,i^{\prime
})}\sigma _{i}^{x}$ $\sigma _{i^{\prime }}^{x}$ where the sum is over all
pairs $(i,i^{\prime })$ belonging to the same chain. The low energy states of
each chain correspond to the codewords in which all spins point in the same
direction along x-axis, i.e. to $\prod_{i}\left\vert e\right\rangle _{i} $
and $\prod_{i}\left\vert o\right\rangle _{i}$ states in each chain. The
parallel connection of the chains ensures that all signs ($+$ or $-$) in (%
\ref{RepetitionCode}) are the same, thereby protecting against $Z-$errors.

Similarly to a repetition code that can determine up to $(N-1)$ $Z$-errors
and $(M-1)$ $X$-errors in the system of $M$ chains of length $N$, the
Hamiltonian version that we discuss in this paper is protected to the order
of $N$ against the dephasing and to the order of $M$ against the decay. The
main difference between the physical systems that we discuss here and the
mathematical error correction models is due to the existence of other low
energy modes that are ignored in their mapping to the chain of $1/2$ spins.

Although this is not directly relevant for engineering protected qubits, we
briefly mention here that these arrays characterized by the two types of
elementary excitations ($Z_2$ charges and $\pi$ vortices) are quite
interesting systems from the viewpoint of statistical physics, because the
unbroken local $Z_2$ symmetry leads to the destruction of the usual Cooper
pair condensate: on any island $j$ connecting two rhombi in the bulk, the
values $\phi_j$ and $\phi_j+\pi$ for the local phase are equiprobable and
therefore the superconducting order parameter $\langle e^{i\phi_j}\rangle$
vanishes. On the other hand, $\langle e^{i2\phi_j}\rangle$ should be
non-zero in order to prevent the proliferation of $\pi$ vortices. These
arrays therefore realize condensates of pairs of Cooper pairs without the
condensation of single Cooper pairs, a regime that may be qualified as
superconducting nematic.~\cite{Park01,Doucot02,Ioffe02}

\section{Model and its local $Z_2$ symmetry}

\label{models}

\subsection{Model}

We shall consider Josephson junction arrays described by the standard
Hamiltonian: 
\begin{equation}
H=4E_{\mathrm{C}}\sum_{ij}\hat{n}_{i}(\mathcal{C}^{-1})_{ij}\hat{n}_{j}-E_{%
\mathrm{J}}\sum_{\langle ij\rangle }\mathcal{J}_{ij}\cos (\phi _{i}-\phi
_{j})  \label{basic_Hamiltonian}
\end{equation}%
where $\mathcal{C}_{ij}$ is the dimensionless capacitance matrix measured in
units of the capacitance $C$ of individual junctions, and the charging
energy $E_{\mathrm{C}}$ is defined by $E_{\mathrm{C}}=e^{2}/(2C)$
(corresponding to single electron charging energy). The dimensionless matrix 
$\mathcal{J}_{ij}$ encodes the Josephson couplings: $\mathcal{J}_{ij}=\pm 1$
when $i$ and $j$ refer to nearest neighbor sites, and $\mathcal{J}_{ij}=0$
otherwise. The signs in the $\mathcal{J}_{ij}$ variables are required in
order to implement the fully frustrated array with one half of a flux
quantum per elementary plaquette. In this case, the product of $\mathcal{J}%
_{ij}$'s over each such plaquette has to be equal to -1. $E_{\mathrm{J}}$ is
the Josephson coupling energy of individual junctions. The number operators $%
\hat{q}_{j}$ are conjugate to the local phases $\phi _{j}$, namely: $\hat{q}%
_{j}=\frac{1}{i}\frac{\partial }{\partial \phi _{j}}$. It will be also
convenient to write down the corresponding Lagrangian, because it involves
the capacitance matrix that is a much more local object than its inverse: 
\begin{equation}
L=\sum_{\langle ij\rangle }\frac{1}{16E_{\mathrm{C}}}(\dot{\phi}_{i}-\dot{%
\phi}_{j})^{2}+E_{\mathrm{J}}\sum_{\langle ij\rangle }\mathcal{J}_{ij}\cos
(\phi _{i}-\phi _{j})
\end{equation}%
To simplify the notation, we have assumed $\hbar =1$ in the expression of
the kinetic term.

For the computation of large arrays and for the qualitative discussion of
their properties it is convenient to model each fully frustrated rhombus by
an effective junction that transfers coherently pairs of Cooper pairs. The
Josephson energy of such an element (that will be called a $Z_2$ junction)
may be written as $-E_{2}\cos2(\phi_i-\phi_j)$ and its charging energy is
incorporated in the Lagrangian by the same term $\frac{1}{16E_{\mathrm{C}}^{%
\mathrm{eff}} }(\dot{\phi}_i-\dot{\phi_j})^{2}$ as before. Note that in the
classical limit, $E_{\mathrm{C}}^{\mathrm{eff}} =E_{\mathrm{C}} $ because a
rhombus involves two chains of two elementary junctions in parallel.

\subsection{Local $Z_2$ symmetry}

\label{local_Z2}

%%%%%%%%%%%%%%%%%%
\begin{figure}[ht]
\includegraphics[width=2.5in]{TwoRhombi}
\caption{Schematics of the simplest protected array (lower pane) and the
microphotogaph of its implementation in the Al nano-structure (upper pane).
The basic element of all protected arrays the rhombus formed by four
junctions (crosses) and four islands (dots), e.g. islands $0,1,2,3$ or $%
0,4,5,6$ which is frustrated by magnetic flux equal to one half of the flux
quantum. In the implementation the islands are horizontal or vertical lines,
the junctions are formed at each intersection. The line style in the
schematics represents the choice of vector potentials used to encode the
magnetic flux: thick lines denote zero vector potential while dashed lines
denote vector potential that shifts by $\protect\pi$ the relation between
current and phase difference accross such junctions. As explained in the
text, the unique properties of the arrays are due to the local $Z_{2}$
transformations, such as the one centered around site $0$ that flips the
chiralities of these two rhombi without changing the phases at the two outer
sites labelled by $1$ and $6$. }
\label{TwoRhombi}
\end{figure}
%%%%%%%%%%%%%%%%%%

As we have emphasized in the Introduction, we shall construct protected
qubits from networks of superconducting islands connected by fully
frustrated rhombi. In such systems single Cooper pairs are localized in
Aharonov-Bohm cages centered at each node of the network and the parity of
the number of Cooper pairs contained in each cage is conserved.~\cite%
{Doucot02} These conserved quantities can be viewed as the generators of
local discrete $Z_{2}$ symmetries which we now define. Let us for instance
consider the geometry of Fig.~\ref{TwoRhombi} with a pair of two adjacent
rhombi that share the node labelled by 0.


The local $Z_{2}$ transformation centered on this node flips simultaneously
the chiralities of the two rhombi adjacent to it while keeping fixed the
phases of the outer sites 1 and 6. More precisely, it corresponds to the
simultaneous changes: 
\begin{eqnarray*}
\phi _{0} &\rightarrow &\phi _{0}+\pi \\
\phi _{2} &\rightarrow &\phi _{3} \\
\phi _{3} &\rightarrow &\phi _{2} \\
\phi _{1} &\rightarrow &\phi _{1} \\
\phi _{4} &\rightarrow &\phi _{5}+\pi \\
\phi _{5} &\rightarrow &\phi _{4}+\pi \\
\phi _{6} &\rightarrow &\phi _{6}
\end{eqnarray*}
Clearly, such transformation leaves invariant the Josephson energy
associated to this circuit, that is: 
\begin{eqnarray*}
H_{\mathrm{J}} &=&-\cos (\phi _{1}-\phi _{2})-\cos (\phi _{2}-\phi
_{0})-\cos (\phi _{1}-\phi _{3}) \\
&+&\cos (\phi _{3}-\phi _{0})-\cos (\phi _{0}-\phi _{4})-\cos (\phi
_{4}-\phi _{6}) \\
&-&\cos (\phi _{0}-\phi _{5})+\cos (\phi _{5}-\phi _{6})
\end{eqnarray*}%
It also preserves the kinetic terms in the Lagrangian, because the above
transformation simply permutes the time derivatives of phase variables
according to the mirror symmetry along the axis joining the central site $0$
with the two extremal points $1$ and $6$. Note that experimentally, this
requires a rather high reproducibility of the junction parameters across the
array, because we need a perfect symmetry between the upper and lower sides
of any rhombus.

%%%%%%%%%%%%%%%%%%
\begin{figure}[ht]
\includegraphics[width=2.5in]{K3Schematics}
\caption{The array with twelve elementary rhombi (shown as crossed circles),
organized in three chains of four rhombi coupled by a central island in the
middle.}
\label{Example_level3}
\end{figure}
%%%%%%%%%%%%%%%%%%


It is also convenient to express the same symmetry operation in the charge
basis. It sends a wave function $\psi $ of $q_{0}$,...,$q_{6}$ into the new
wave function $\psi ^{\prime }$ defined by: 
\begin{eqnarray}
\lefteqn{\psi \prime (q_{0},q_{1},q_{2},q_{3},q_{4},q_{5},q_{6})=}  \notag \\
&&e^{i\pi (q_{0}+q_{4}+q_{5})}\psi
(q_{0},q_{1},q_{3},q_{2},q_{5},q_{4},q_{6})
\end{eqnarray}%
The presence of the phase factor $e^{i\pi q_{0}}$ suggests that this unitary
operation is directly related to the parity of the number of Cooper pairs on
the node $0$. However, this operation is slightly more complicated, because
it also involves the mirror symmetry along the horizontal diagonals of the
two rhombi and additional $\pi $ translations of phases at sites $4$ and $5$%
. A more precise connection can be obtained in the framework of a Bose
Hubbard model. There, the Josephson coupling is represented by a single
boson hopping process. It is therefore possible to diagonalize this
Hamiltonian (in the presence of the magnetic field of half flux quantum per
rhombus); this yields Wannier orbitals localized in finite clusters of sites
centered on the nodes that share two or more rhombi. In this version, the
local $Z_{2}$ generator centered at a given node counts precisely the parity
of the number of bosons in the localized orbitals centered on this node.~%
\cite{Doucot02} In this paper, we shall not use the Bose Hubbard model, but
rather the quantum rotor Hamiltonian version Eq.~(\ref{basic_Hamiltonian}).


In this work we shall mostly consider coupled chains of rhombi and some
generalizations such as the one depicted on Fig.~\ref{Example_level3}. For
this array we have six binary flips associated to the middle islands of the
six two rhombi subchains, and one six-fold flip of all the six rhombi
adjacent to the central island. It is easy to generalize the precise
definition of the local $Z_2$ flips given above to the case of a node
connected to an arbitrary number of rhombi. One can then check that two
local $Z_2$ transformations centered on different nodes mutually commute.

%%%%%%%%%%%%%%%%%%%
\begin{figure}[ht]
\includegraphics[width=2.5in]{nonlocaltransf.eps}
\caption{(Color online) An example of a non local $Z_2$ unitary
transformation in an array with three coupled chains of four rhombi. The red
dashed dotted line indicates the set of rhombi which chirality flips under
the transformation. Phases on all islands located on the left (resp. right)
of this line are unchanged (resp. shifted by $\protect\pi$). For sites
located precisely on the line, they exchange their phase with the site
opposite to them on the same rhombus. For example, phases at sites 2 and 3
are permuted.}
\label{nonlocaltransf}
\end{figure}
%%%%%%%%%%%%%%%%%%

As we discuss below, the low energy spectra of large arrays can be
reproduced in a coarse grained picture in which the array is built from 
\emph{effective} $Z_{2}$ junctions, i.e. elements doubly periodic in the
phase variable. In this picture the local symmetries take a simpler form,
because the $Z_{2}$ generator based at node $j$ shifts $\phi _{j}$ by $\pi $
while leaving all the other phases unchanged. In this situation, it is
exactly the operator that counts the parity of the number of Cooper pairs on
this node.


For the discussion of protected qubits, it is also important to mention the
existence of non local $Z_{2}$ symmetries in which the local chirality
variables are flipped for a subset of rhombi located on a line that cuts the
array in two parts. An example is shown on Fig.~\ref{nonlocaltransf}. For a
rhombus on this line, such as the one with sites labelled 0,1,2 and 3 on the
figure, the transformation rules are as before, namely: 
\begin{eqnarray*}
\phi _{0} &\rightarrow &\phi _{0}+\pi  \\
\phi _{2} &\rightarrow &\phi _{3} \\
\phi _{3} &\rightarrow &\phi _{2} \\
\phi _{1} &\rightarrow &\phi _{1}
\end{eqnarray*}%
Phases on sites located on the left of the line are unchanged whereas those
on sites located on the right are shifted by $\pi $. As a result the total
phase difference $\pi $ between the right and the left boundaries is also
shifted by $\pi $. As mentioned in the Introduction, the physical meaning of
such symmetry is the absence of any process that transfers an odd number
number of Cooper pairs from one boundary to the other.

\section{Hierarchical construction}

\subsection{Motivation}

%%%%%%%%%%%%%%%%%%
\begin{figure}[ht]
\includegraphics[width=2.5in]{HierarchyGeneral}
\caption{Hierarchical construction: the element at level $n+1$ is obtained
by connecting in parallel $K$ identical chains made of two elements at level 
$n$}
\label{HierarchyGeneral}
\end{figure}
%%%%%%%%%%%%%%%%%%

As discussed in the Introduction, a large number of parallel chains of
rhombi provides a good protection from the environment provided that quantum
fluctuations of rhombi in each chain are strong but the phase difference
between the ends of the chains does not fluctuate due to a large number of
chains. The last requirement is not easy to satisfy for moderately long
chains because the phase stiffness decreases exponentially with the length.~%
\cite{Matveev02} The alternative solution is to couple chains in the
transverse direction. If this coupling is achieved by rhombi as well, one
ends up with the full blown decorated lattice. Although possible, this
design has a disadvantage that quantum fluctuations need to be strong enough
to result in a simultaneous flip of many rhombi that have a common site and
requires a large number of rhombi, all this is not easy to achieve
practically. The attractive alternative that we explore here is to use the
hierarchical construction that does not involve additional elements in
transverse direction. This construction has two substeps at each level.
First, two Josephson elements formed at nth level are connected in series
forming a weaker Josephson element. Second, K such Josephson elements are
connected in parallel compensating this weakening and forming an effective
element of the (n+1)th level. The idea is to choose the parameters so that
the effective Josephson energy $E_{2}\cos (2\theta )$ and charging energy
keep constant ratio at each level. This implies that quantum fluctuations
remain strong but do not wash out the phase stiffness in a large structure.
Provided that the $E_{1}\cos (\theta )$ term was significantly smaller than $%
E_{2}\cos (2\theta )$ at the first stage of the hierarchy, it will decrease
exponentially with the level.


\subsection{Hierarchical construction}

The recursion relation that relates the energy scales at the next hierarchy
level to the ones at the previous level can be obtained analytically in the
limit of large branching ratio, $K$. In this limit, $E_{2}/E_{\mathrm{C}} $
for the individual junction is small which allows to treat them
perturbatively.

%%%%%%%%%%%%%%%%%%
\begin{figure}[ht]
\includegraphics[width=2.5in]{g2} \includegraphics[width=2.5in]{GapSerial}
\caption{ (a) Amplitude of the second harmonics $g_2(x)$ of the Josephson
energy versus phase difference across two identical $E_2 \cos(2\protect\theta%
)$ junctions connected in series where $x=E_2/E_{\mathrm{C}} $. (b) The
energy of the lowest excited state corresponding to one Cooper pair at the
middle island. }
\label{J2}
\end{figure}
%%%%%%%%%%%%%%%%%%


A general (super) junction can be characterized by three parameters: 
\begin{equation*}
H=-E_{2}\cos 2\theta -E_{1}\cos \theta +4E_{\mathrm{C}} q^{2}
\end{equation*}
where $q$ is the charge conjugate to the phase measured in units of $2e$.
Two such junctions connected in series give a Josephson energy: 
\begin{equation*}
H=-E_{2}^{\prime }\cos 2\theta -E_{1}^{\prime }\cos \theta +4E_{\mathrm{C}}
^{\prime }q^{2}
\end{equation*}%
Generally, one expects that the parameters of this element are some
functions of the energies characterizing the original junctions:

\begin{eqnarray}
E_{2}^{\prime } &=&g_{2}(\frac{E_{2}}{E_{\mathrm{C}} },\frac{E_{1}}{E_{%
\mathrm{C}} })E_{\mathrm{C}} \\
E_{1}^{\prime } &=&g_{1}(\frac{E_{2}}{E_{\mathrm{C}} },\frac{E_{1}}{E_{%
\mathrm{C}} })E_{\mathrm{C}} \\
E_{\mathrm{C}} ^{\prime } &=&h(\frac{E_{2}}{E_{\mathrm{C}} },\frac{E_{1}}{E_{%
\mathrm{C}} })E_{\mathrm{C}}
\end{eqnarray}%
In the following we shall be mostly interested in the limit of very small $%
\frac{E_{1}}{E_{\mathrm{C}} }\ll 1$. In this limit $g_{2}$ becomes a
function of only one variable, $\frac{E_{2}}{E_{\mathrm{C}} }$. Perturbative
expansion gives: 
\begin{eqnarray}
E_{2}^{\prime } &=&\frac{1}{8}\frac{E_{2}^{2}}{E_{\mathrm{C}} }
\label{Perturbation2} \\
E_{1}^{\prime } &=&\frac{1}{2}\frac{E_{1}^{2}}{E_{\mathrm{C}} }
\label{Perturbation1} \\
E_{\mathrm{C}} ^{\prime } &=&2E_{\mathrm{C}}  \label{Perturbation3}
\end{eqnarray}%
for sufficiently large $E_{\mathrm{C}} $ (perturbation theory in $E_{2}/E_{%
\mathrm{C}} $). In this limit of small $E_{2}/E_{\mathrm{C}} $ the charging
energy of the two elements connected in series is twice $E_{\mathrm{C}} $.
%%%%%%%%%%%%%%%%%%
\begin{figure}[ht]
\includegraphics[width=2.5in]{K3Schematics} \includegraphics[width=2.5in]{K3}
\caption{ (a) The stage intermediate between levels 1 and 2 of the hierarchy
with $K=3$ branching made from the $\cos2\protect\theta$ elementary blocks
used for testing the validity of the recursive approximation. (b) Exact
diagonalization results for the effective $E_2$ (in the units of bare
charging energy) of this element is shown as dots. The curves correspond to
different approximate schemes. The full line is the result of the recursive
approximation discussed in the text applied twice. The dotted line is the
result of the diagonalization of the series of four junctions with
additional capacitance $K-1$ in the middle and an overal multipying factor
of $K^2$, according to the the form suggested by the perturbation theory in $%
E_{2}/E_{\mathrm{C}} $. }
\label{K3}
\end{figure}
%%%%%%%%%%%%%%%%%%


Combining $K$ such chains in parallel we get that at the next stage of the
hierarchy: 
\begin{eqnarray}
E_{\mathrm{C}} ^{(n+1)} &=&(2/K)E_{\mathrm{C}} ^{(n)}  \label{Recursion1} \\
\frac{E_{2}^{(n+1)}}{E_{\mathrm{C}} ^{(n+1)}} &=&\frac{K^{2}}{16}\left[ 
\frac{E_{2}^{(n)}}{E_{\mathrm{C}} ^{(n)}}\right] ^{2}  \label{Recursion2}
\end{eqnarray}%
One can treat the block formed by $2K$ Josephson elements as a single
junction characterized by only one degree of freedom provided that the
energy scales at the next step are much smaller than the ones at the
previous step. This is satisfied if $K\gg 1$ because in this limit the
typical energy scales are set by $E_{\mathrm{C}} ^{(n)}$ that decreases with
the level $n$. Thus, in this limit, the equations (\ref{Recursion1}-\ref%
{Recursion2}) become exact and show the (unstable) fixed point: 
\begin{equation}
\frac{E_{2}^{(n)}}{E_{\mathrm{C}} ^{(n)}}=\frac{16}{K^{2}}
\label{FixedPoint}
\end{equation}

At realistic values of the branching number $K=2-4$ one expects two types of
corrections to this simple recursion. First, the perturbation expansion that
gives (\ref{Perturbation2}) and $E_{\mathrm{C}} ^{\prime }=2E_{\mathrm{C}} $
is not exact. Second, characteristic energy scales at the next level are not
much smaller than the ones at the previous level and that makes coarse
graining only approximate. We now discuss the importance of these
corrections in turn. Calculation of the next order in perturbation theory
gives: 
\begin{eqnarray}
E_{\mathrm{C}} ^{\prime } &=&\left[ 1-\frac{1}{16}\left( \frac{E_{2}}{E_{%
\mathrm{C}} }\right) ^{2}\right] 2E_{\mathrm{C}}  \label{Perturbation12} \\
E_{2}^{\prime } &=&\left[ 1-\frac{7}{256}\left( \frac{E_{2}}{E_{\mathrm{C}} }%
\right) ^{2}\right] \frac{1}{8}\frac{E_{2}^{2}}{E_{\mathrm{C}} }
\label{Perturbation22}
\end{eqnarray}%
that shows that the corrections to the perturbative result become
significant for $E_{2}\gtrsim 2E_{\mathrm{C}} $. Thus, the equation for the
fixed point (\ref{FixedPoint}) remains approximately valid down to $K=3$. At
small branching numbers the time scales do not grow fast with the level and
this introduces an additional source of corrections. To check the accuracy
of the recursion relations and the perturbative approximation we performed a
numerical study comparing the effective Josephson energy $E_{2}$ expected
from the recursion relations at the second level of hierarchy with the
numerical diagonalization of the Hamiltonian that describes the system shown
in Fig.~\ref{K3}a. We check separately the accuracy of the perturbative
expansion by comparing the results (\ref{Perturbation2}), (\ref%
{Perturbation22}) with the diagonalization of two contacts connected in
series, and the accuracy of the recursion by comparing the effective
Josephson energy $E_{2}$ found numerically with 
\begin{equation}
E_{2}=g_{2}(\frac{K^{2}g_{2}(E_{2}/E_{\mathrm{C}} )}{h(E_{2}/E_{\mathrm{C}} )%
})\frac{1}{K}h(E_{2}/E_{\mathrm{C}} )E_{\mathrm{C}}
\label{RecursiveSecondLevel}
\end{equation}%
where $g_{2}(E_{2}/E_{\mathrm{C}} )$ was evaluated numerically (see Fig.~\ref%
{J2}a) and we used the perturbative result $h(x)\approx 2-1/8x^{2}$. The
result of the comparison is shown in Fig.~\ref{K3}b, where one observes that
the recursive approximation (\ref{RecursiveSecondLevel}) reproduces the
exact results with a good accuracy even for $K=3$ in the physically
important range of the parameters $0<E_{2}/E_{\mathrm{C}} <3$. We note that
the results deviate lightly at low values, the reason can be understood by
comparing the results obtained in the recursion approximation with the
direct perturbation theory for the second level device.


The lowest order of the perturbation theory correspond to a sequence of
charge tunnelings and charge-hole creations each of which moves charge in
the same direction. The constraint that the process should occur with the
minimal number of steps implies that a given sequence involves only one (out
of $K$) branches of the device at each side. The remaining $K-1$ branches
provide the additional capacitance attached to the middle island. Thus, in
this order, the circuit can be replaced by a series of four basic elements,
with the central island capacitively coupled to the ground, provided that
the final result is multiplied by $K^{2}$. In the limit of large $K\gg 1$
the energy of the charge at the middle island is low in $1/K.$

The condition that process occurs in the minimal number of steps does 
preclude intermediate states that contain more than one charge or hole.
However, in the limit of large $K$ the dominant contribution comes only from
the processes that involve the intermediate state in which the only charge
is located at the middle island because this state  has parametrically
smaller charging energy $(16/K)E_{\mathrm{C}}$. Keeping only these terms
reproduces exactly the recursion result in which $g_{2}$ and $h$ functions
are approximated by their lowest order expressions (\ref{Perturbation2}), (%
\ref{Perturbation3}): 
\begin{equation}
E_{2}^{pert}\approx \frac{K^{3}}{2^{10}}\left( \frac{E_{2}}{E_{\mathrm{C}}}%
\right) ^{4}E_{\mathrm{C}}  \label{E2_pert_recursive}
\end{equation}%
At finite values of $K$, three effects combine to modify this result. First,
the expression for $g_{2}$ at the first step has to be modified to take into
account the effect of the larger structure on the energy of virtual charged
excitations between two elementary blocks. Specifically, the prefactor $1/8$
now becomes $1/(8+4/K)$. Second, some processes do not involve intermediate
states with the charge located only on the middle island. Summing these
effects, we get: 
\begin{equation}
E_{2}^{pert}=\frac{K^{4}(K+2)}{2^{8}(2K+1)^{2}}\left( \frac{E_{2}}{E_{%
\mathrm{C}}}\right) ^{4}E_{\mathrm{C}}  \label{E2_pert}
\end{equation}%
Third, the energy of the intermediate state with one charge on the middle
island differs from $(16/K)E_{\mathrm{C}}$ due to capacitance
renormalization by virtual processes, but this occurs at higher order in $%
E_{2}/E_{\mathrm{C}}$. For physically relevant $K=3$ the exact perturbation
result (\ref{E2_pert}) differs from the recursive approximation (\ref%
{E2_pert_recursive}) by 20\%.

\subsection{Medium size rhombi arrays}


%%%%%%%%%%%%%%%%%%
\begin{figure}[ht]
\includegraphics[width=2.5in]{1Rhombus}
\caption{The exact diagonalization result $E_2$ for the second harmonics of
a single rhombus as a function of $E_{\mathrm{J}} /E_{\mathrm{C}} $. }
\label{1Rhombus}
\end{figure}
%%%%%%%%%%%%%%%%%%
%%%%%%%%%%%%%%%%%%
\begin{figure}[ht]
\includegraphics[width=2.5in]{RhombusCapacitance.eps}
\caption{The effective capacitance of a single rhombus as a function of $E_{%
\mathrm{J}} /E_{\mathrm{C}} $ determined, as explained in the text, from the
energy gap of a chain of two rhombi with its rightmost island connected to
the ground plane by the capacitance $C_{0}$. The three curves correspond,
from bottom up, to $C_{0}=2$, $C_{0}=6$, and $C_{0}=10$. }
\label{Effective_capacitance}
\end{figure}
%%%%%%%%%%%%%%%%%%


For experimental implementations of protected qubits with superconducting
circuits, the simplest way to get a $\pi $ periodic element is to use a
rhombus frustrated by half a flux quantum. In this sub-section, we shall
discuss computational issues raised by the treatment of networks of rhombi.
>From the perspective of numerical diagonalizations, these networks involve
very large Hilbert spaces, even for small sizes and moderate values of $E_{%
\mathrm{J}}/E_{\mathrm{C}}$. For example, a chain of four rhombi has 12
degrees of freedom which implies diagonalization in $10^{10}-10^{12}$
dimensional space for $E_{\mathrm{J}}/E_{\mathrm{C}}\lesssim 8.$ Therefore,
it is essential to develop approximate but reliable computation schemes for
the low energy spectrum that allow for a significant reduction in the size
of the Hilbert space. A natural coarse-graining procedure is to try to
replace each rhombus by an effective $Z_{2}$ junction. As we have seen in
section~\ref{models}, this operation preserves the existence of a local $%
Z_{2}$ symmetry which implies that only pairs of Cooper pairs can propagate
coherently across the network. 



Generally, the dynamics of  $Z_{2}$ junction is fully characterized its
action $S(\theta )$ which is a non-local function of the phase $\theta $
across it. For the slow phase dynamics the action is separated into kinetic
and potential part, for former can be approximated by quadratic form%
\begin{equation*}
\frac{1}{16E_{\mathrm{C}}^{eff}}\int dt\dot{\theta}^{2}
\end{equation*}
where $E_{\mathrm{C}}^{\mathrm{eff}}$ is its effective charging energy. The
potential energy part can be characterized by the coefficients $E_{n}$ in
the expansion of the Josephson energy $V(\theta )=\sum E_{n}\cos n\theta $.
Numerical diagonalization of the single rhombus shows that all $E_{n}$ with $%
n>2$ are very small in the relevant parameter range $E_{\mathrm{J}}/E_{%
\mathrm{C}}\lesssim 8$ and can be safely neglected. So, in this
approximation the effective $Z_{2}$ junction representing the rhombus is
characterized by three parameters: $E_{1},E_{2}$ and $E_{\mathrm{C}}^{%
\mathrm{eff}}$. The accuracy of this approximation will be discussed below.
The effective $E_{1}$ for an ideal rhombus vanishes exactly at all orders in 
$E_{\mathrm{J}}/E_{\mathrm{C}}$. The effective $E_{2}$ is generated starting
from the fourth order in $E_{\mathrm{J}}/E_{\mathrm{C}}$. At this order, we
have: 
\begin{equation}
E_{2}=\frac{7}{64}\frac{E_{\mathrm{J}}^{4}}{E_{\mathrm{C}}^{3}}
\end{equation}%
The exact result for $E_{2}$ as a function of $E_{\mathrm{J}}/E_{\mathrm{C}}$
is shown on Fig.~\ref{1Rhombus}.


%%%%%%%%%%%%%%%%%%
\begin{figure}[ht]
\includegraphics[width=2.5in]{2RhombiJ2.eps} %
\includegraphics[width=2.5in]{3RhombiJ2.eps}
\caption{(a) The effective stiffness $E^{\prime }_2$ of a chain of two
rhombi, where the phase $\protect\phi$ on the rightmost island is taken as a
control parameter. The plots show the results of the numerical
diagonalization (dots) of the two rhombi chain, the results of the mapping
onto a chain of two $Z_{2}$ junctions whose bare parameters $E_{2}$ and $E_{%
\mathrm{C}}^{\mathrm{eff}} $ are determined as explained in the paper, and
the recursive approximation for this chain of two rhombi. (b) The same
curves are shown for the chain of three rhombi. }
\label{2Rhombi}
\end{figure}
%%%%%%%%%%%%%%%%%%

We now turn to the effective charging energy. At very small $E_{\mathrm{J}}$
the tunneling of the charge can be neglected and the effective charging
energy is given by electrostatics: $E_{\mathrm{C}}^{\mathrm{eff}}=E_{\mathrm{%
C}}.$ However, we expect rather large corrections to this result in the
relevant range of  $2\lesssim E_{\mathrm{J}}/E_{\mathrm{C}}\lesssim 8$ due
to the charge delocalization which should lead to a decrease in the value of 
$E_{\mathrm{C}}^{\mathrm{eff}}$. Estimating the effective capacitance first
raises the issue of defining it precisely. In Eq.~(\ref{Perturbation12}),
this quantity has been computed for a series of two $Z_{2}$ junctions, to
second order in $E_{\mathrm{J}}/E_{\mathrm{C}}$. This computation has been
performed by first assuming a slow variation (on the scale $\hbar /E_{%
\mathrm{C}}$) of the phase $\phi $ across the two junction chain, and then
integrating out the internal phase degree of freedom. It is only in this
limit, and to the lowest order in perturbation theory, that the effective
action keeps the local $\dot{\phi}^{2}$ form that allows the straightforward
definition of a renormalized capacitance. In general, we expect the
effective action to become non-local in time, and even non-linear in $\dot{%
\phi}$. Furthermore, such a path integral approach is not as directly
numerically implementable as a Hamiltonian formulation. Therefore, we shall
adopt here the following procedure. We consider a series of two identical
rhombi, with an additional capacitance $C_{0}$ on its rightmost island. This
capacitance is introduced in order to slow down the quantum dynamics of the
phase $\phi $ on this island, it also mimics the effect of the remaining
part of the array. Then the gap between the ground-state and the first
excited state with an odd charge on the intermediate island is computed by
numerical diagonalization, for various values of $E_{\mathrm{J}}/E_{\mathrm{c%
}}$ and $C_{0}$. This gap is also computed for a series of two $Z_{2}$
junctions, for which $E_{2}$ is chosen to be the one previously computed on
Fig.~\ref{1Rhombus} and where $E_{\mathrm{C}}^{\mathrm{eff}}$ is adjusted in
order to match the gaps for the odd intermediate charge sector of both
systems. The values of $C_{\mathrm{eff}}$ ($E_{\mathrm{C}}^{\mathrm{eff}%
}=2e^{2}/C_{\mathrm{eff}}$) as a function of $E_{\mathrm{J}}/E_{\mathrm{C}}$
for $C_{0}=2,6$, and $10$ are shown on Fig.~\ref{Effective_capacitance}.
These curves show clearly a tendency for $C_{\mathrm{eff}}$ to saturate to a
value around $3.5C$ as soon as $E_{\mathrm{J}}/E_{\mathrm{C}}\gtrsim 4$, and
provided that $C_{0}\gtrsim 5E_{\mathrm{C}}$. This stability of $C_{\mathrm{%
eff}}$ is expected from physical considerations, because a large $C_{0}$
turns $\phi $ into a slow parameter, and a large $E_{\mathrm{J}}/E_{\mathrm{C%
}}$ reduces the magnitude of its quantum fluctuations. This suggests that
this way to map a network of rhombi onto the network of $Z_{2}$ junctions is
likely to provide a good estimate for both the local $Z_{2}$ gap associated
to the energy difference between states with even and odd particle numbers
on a given node, and the phase stiffness for the charge $4e$ condensate. As
a consistency check of this procedure, we compare in Fig.~\ref{2Rhombi} the
phase stiffness $E_{2}^{\prime }$ for a chain of a few rhombi, with the same
quantity for a chain of  $Z_{2}$ junctions. The results for $N=2$ are shown
on Fig.~\ref{2Rhombi}a, and those for $N=3$ on Fig.~\ref{2Rhombi}b. The
agreement between the rhombi chains and the appropriate $Z_{2}$ chains is
very good throughout the whole range of $E_{\mathrm{J}}/E_{\mathrm{C}}$
values considered here.

%%%%%%%%%%%%%%%%%%
\begin{figure}[ht]
\includegraphics[width=2.5in]{K3J2H}
\caption{(Color online) The system of twelve rhombi at the same level of
hierarchy as the system of junctions shown in Fig.\protect\ref{K3} contains $%
31$ quantum rotors and so is too large for a direct numerical calculation.
The upper (blue) line shows the result of the recursive approximation
described in the text. The lower (purple) line uses the approximate mapping
on the K=3 array of $Z_2$ junctions. The agreement between these two methods
is very good. }
\label{12Rhombi}
\end{figure}
%%%%%%%%%%%%%%%%%%

Another check is provided by the calculation of the stiffness for the system
of twelve rhombi at the same level of hierarchy as the system of junctions
shown in Fig.\ref{K3}, using two different approximate methods. The first
one is a recursive approximation starting with the initial $4e$ tunneling
amplitude of a series of two rhombi that has just been computed (see Fig.~%
\ref{2Rhombi}), and with an effective capacitance of half the effective
capacitance of a single rhombus shown on Fig.~\ref{Effective_capacitance}.
The second method uses the mapping onto the $K=3$ network of $Z_2$ junctions
represented on Fig.~\ref{K3}a. The results of both calculations are shown on
Fig.~\ref{12Rhombi}, and the small difference between the two curves is
another indication that these various procedures are mutually consistent.


\subsection{Low energy excitation spectrum}

%%%%%%%%%%%%%%%%%%
\begin{figure}[ht]
\includegraphics[width=2.5in]{K3Gaps.eps} %
\includegraphics[width=2.5in]{K3RGaps.eps}
\caption{Effective phase stiffness $E^{\prime }_{2}$ (increasing curve) and
energy gaps of low energy states (decreasing curves) of a $K=3$ system of two
rhombi chains connected to a central island, computed using the mapping onto
a similar array of $Z_2$ junctions. These data are plotted versus the
effective $E_{2}/E_{\mathrm{C2}}$ of these $Z_2$ junctions (a), and versus the 
nominal $E_{\mathrm{J}} $ of the individual junctions which compose the rhombi 
(b). To convert from $E_{2}/E_{\mathrm{C2}}$ to $E_{\mathrm{J}}/E_{\mathrm{C}}$ 
we used the effective Josephson energy of the rhombus shown on 
Fig.~\protect\ref{1Rhombus} and the effective capacitance shown on 
Fig.~\protect\ref{Effective_capacitance}. The lowest gap corresponds to one 
$Z_2$ charge on the central island, i.e. to
the label $(1,0,0)$ defined in text. The second curve corresponds to two
types of states, that have one $Z_2$ charge on one site in one of the six
branches, with or without a $Z_2$ charge on the central island. These
states, whose energies are very close, are labelled by $(0,1,0)$ and $%
(1,1,0) $. Finally, the upper curve is associated to the two nearly
degenerate levels with labels $(0,1,1)$ and $(1,1,1)$.}
\label{K3_pus_spectrum}
\end{figure}
%%%%%%%%%%%%%%%%%%


Here, we focus on the system with twelve rhombi at the same level of the $K=3
$ hierarchy as shown on Fig.~\ref{K3}. We assume that the total phase
difference  $\phi $ across the system is fixed at the value $0$ (or $\pi $)
that minimizes the energy.  In this case, and in the limit of vanishing
charging energy, the phase $\phi _{\mathrm{C}}$ of the central island can
take two values modulo $2\pi $ in any classical ground-state, namely $0$ and 
$\pi $. Once $\phi _{\mathrm{C}}$ is known, for each subchain composed of
two rhombi, we have two possible chirality configurations. Therefore, this
system exhibits $2^{7}$ classical ground-states. In the presence of finite
charging energy, we expect that quantum fluctuations lift the degeneracy of
these classical states, leaving instead one family of low-lying states whose
energies will be computed below. Above these states, that are relatively
isolated, we expect a much more dense spectrum, that results from the
excitation of one or more localized Josephson plasmon mode. Our general
discussion of the local $Z_{2}$ symmetry shows that a useful way to label
the energy eigenstates is to keep track of the values of the local $Z_{2}$
charges on each of the seven islands. Below, we shall consider a few
sectors, with at most three sites carrying a non-zero value of the $Z_{2}$
charge (that is an odd number of Cooper pairs). We shall then use labels of
the form $(a,b,c)$, with $a,b,c\in \{0,1\}$ to refer to such states, where $a
$, $b$, and $c$ indicate the $Z_{2}$ charge on the central island, on one
particular site on its left, and on one particular site on its right,
respectively.

%%%%%%%%%%%%%%%%%%
\begin{figure}[ht]
\includegraphics[width=2.5in]{K2Gaps.eps} %
\includegraphics[width=2.5in]{K2RGaps.eps}
\caption{Energy gaps and effective phase stiffness  for the $K=2$ system
with each chain consisting of two rhombi connected in series. 
The left pane gives the gaps and the stiffness for the effective $Z_2$ junctions,
the upper pane converts these results to the physical properties of the chain made 
of rhombi. All notations
are the same as on Fig.~\protect\ref{K3_pus_spectrum}. }
\label{K2_pus_spectrum}
\end{figure}
%%%%%%%%%%%%%%%%%%
%%%%%%%%%%%%%%%%%%
\begin{figure}[ht]
\includegraphics[width=2.5in]{K23Gaps.eps} %
\includegraphics[width=2.5in]{K23R.eps}
\caption{Energy gaps and effective phase stiffness for the $K=2$ system with
each chain consisting of three rhombi connected in series. 
The left pane gives the gaps and the stiffness for the effective $Z_2$ junctions,
the upper pane converts these results to the physical properties of the chain made 
of rhombi. All notations are the same as on Fig.~\protect\ref{K3_pus_spectrum}. }
\label{inverted_pus_spectrum}
\end{figure}
%%%%%%%%%%%%%%%%%%

Note that such an array is already too large for a numerical diagonalization
of the corresponding quantum circuit Hamiltonian to be possible. Therefore,
we have to rely on approximate methods, such as the mapping onto an array of 
$Z_{2}$ junctions, as previously discussed. Numerical results for the
effective phase stiffness and the energy gaps of various excitations are
presented on Fig.~\ref{K3_pus_spectrum}



As expected, the lowest gap is obtained while putting a single $Z_{2}$
charge on the central island, that has the smallest charging energy. This
gap becomes exponentially small at large $E_{\mathrm{J}}/E_{\mathrm{C}}$, as
can be expected from the semi-classical picture according to which the phase 
$\phi _{\mathrm{C}}$ of the central island tunnels between its two
equivalent classical minima. Note that the energy barrier between these two
classical states increases when the connectivity $K$ increases. This is in
agreement with the result of Fig.~\ref{K2_pus_spectrum}, which shows a
strong sensitivity of the gap of the $(1,0,0)$ state on the actual value of $%
K$.


By contrast, the energy of the other excited states shown here is much less
sensitive to $K$, because these involve $Z_{2}$ charges on the chain
islands. We find that the energy of two $Z_{2}$ charges on opposite sides ($%
(0,1,1)$ state) is a little less than twice the energy of the single $Z_{2}$
charge ($(0,1,0)$ or $(0,0,1)$). Of course, these states are degenerate, by
a factor $K^{2}$ for the former and a factor $2K$ for the later. A more
surprising observation is that adding an extra $Z_{2}$ charge on the central
island, when one is already present elsewhere, hardly changes the energy.
This can be understood in the small $E_{\mathrm{J}}/E_{\mathrm{C}}$ limit by
considering the inverse capacitance matrix for which the diagonal element $%
C_{00}^{-1}$ on the central island is equal to twice the off-diagonal
elements $C_{0j}^{-1}$. This implies that the states which carry the
classical charges (and not $Z_{2}$ charges here) $(0,1,0)$ and $(-1,1,0)$
have the same electrostatic energy, equal to $(2+1/K)E_{\mathrm{C}}$.
Likewise, states with charges $(0,1,-1)$ and $(-1,1,1)$ have the same
electrostatic energy, $4E_{\mathrm{C}}$, that is slightly less than twice
the previous one. Notice that while this degeneracy is a feature of the
small $E_{\mathrm{J}}/E_{\mathrm{C}}$ limit, it is also expected in the
semi-classical limit of large $E_{\mathrm{J}}/E_{\mathrm{C}}$, because of
the exponentially small tunneling gap between to two classical states of the
central island.

We have also considered another $K=2$ system in which all the four chains
connected to the central island are composed of three rhombi. The numerical
results for the effective phase stiffness and for the lowest energy gaps are
presented on Fig.~\ref{inverted_pus_spectrum}. The comparison with Fig.~\ref%
{K2_pus_spectrum} shows that this latter system has a smaller phase
stiffness but larger gaps than the former $K=2$ array with chains of only
two rhombi.

As we have argued in the Introduction, for experimental implementation of a
protected qubit  one needs to maximize both the energy gap (to reduce the
effective noise-induced tunneling amplitude of a single Cooper pair across
the system) and the effective phase stiffness (to restrict the amplitude of
a $\pi $ shift of the phase on the central island). These conflicting
requirements suggest that a good working point is obtained for the values of 
$E_{\mathrm{J}}/E_{\mathrm{C}}$ where the curve for the effective stiffness $%
E_{2}^{\prime }$ and the curve for the lowest gap cross. The optimal array
is then the one that maximizes the common value of $E_{2}^{\prime }$ and the
gap at this working point. Comparing Figs.~\ref{K3_pus_spectrum}, \ref%
{K2_pus_spectrum}, and \ref{inverted_pus_spectrum}, we conclude that the
most promising array seems to be the one with $K=2$ composed of three rhombi
chains.

\section{Perspectives}

%%%%%%%%%%%%%%%%%%
\begin{figure}[ht]
\includegraphics[width=2.5in]{K3RCrCur.eps}
\caption{The effective phase stiffness $E_{2}^{\prime }$ in various excited
states of the $K=3$ array of two rhombi chains shown in Fig. 6. This
quantity is directly related to the $\Phi _{0}/2$ periodic modulation of the
critical current of the array imbedded in a large SQUID loop threaded by a
magnetic flux $\Phi $. The sensitivity of the critical current to the
quantum state of the array provides a way to probe experimentally these low
energy levels.}
\label{Critical_current}
\end{figure}
%%%%%%%%%%%%%%%%%%


Experimentally, the variation of the total energy for various states as a
function of the total phase difference $\phi $ across the array can be
determined by measuring the critical current of the system composed of the
array in parallel with a large Josephson junction. Motivated by the
theoretical analysis for rhombi chains~\cite{Protopopov04,Protopopov06}, the
Grenoble group has shown that, in $N=8$ rhombi chains, the critical current
as a function of $\phi $ becomes $\pi $- periodic when the flux per rhombus
is close to half a flux quantum.~\cite{Pop08} Furthermore, they observed a
clear evidence of quantum phase fluctuations in a sample with $E_{\mathrm{J}%
}/E_{\mathrm{C}}=8$, where the critical current versus phase oscillations
are rounded and reduced in amplitude. A direct experimental evidence for a
finite $Z_{2}$ charge gap in the finite array corresponding to the geometry
of Fig.~\ref{Example_level3} has been obtained at Rutgers, for a several
samples with the ratio $E_{\mathrm{J}}/Ec$ ranging from 2 to 5.


The next important step would be to measure the energy spectrum of low lying
excitations. This can still be done in principle by critical current
measurements, because, as shown on Fig.~\ref{Critical_current}, the global
phase stiffness is sensitive to the actual quantum state of the whole array.
The selective excitation of the system can be achieved by sending microwave
pulses on an electrode that is capacitively coupled to the central island.

Finally, the goal is to operate such an array as a protected qubit. For the
array shown on Fig.~\ref{Example_level3}, the logical qubit can be encoded
in the two possible values of the central island phase. A better protection
would be obtained by combining two such arrays in series, using again the
central island phase for encoding. As shown by Kitaev, such a system allows
for a universal set of robust quantum gates.~\cite{Kitaev06} We hope that
the research reported here will facilitate further efforts to implement this
very promising class of systems.

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%% BACKMATTER
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

\begin{theacknowledgments}
LI is thankful to LPTMS Orsay and LPTHE Jussieu for their hospitality 
through a financial support from CNRS while BD has enjoyed the hospitality
of the Physics Department at Rutgers University. This work was made possible
by financial support from ANR grant 06BLAN-0218-01, NSF ECS-0608842 and
DARPA HR0011-09-1-0009.
  
\end{theacknowledgments}

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%% The following lines show an example how to produce a bibliography
%% without the help of the BibTeX program. This could be used instead
%% of the above.
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

\begin{thebibliography}{99}
\bibitem{Gladchenko08} S. Gladchenko, D. Olaya, E. Dupont-Ferrier, B. Dou%
\c{c}ot,  L. B. Ioffe, and M. E. Gershenson, arXiv:0802.2295.

\bibitem{Kitaev97} A. Yu. Kitaev, arXiv:quant-ph/9707021, and Ann. Phys.
(New-York) \textbf{303}, 2, (2003).

\bibitem{Park01} K. Park, and D. Huse, Phys. Rev. B \textbf{64}, 134522,
(2001).

\bibitem{Doucot02} B. Dou\c{c}ot and J. Vidal, Phys. Rev. Lett. \textbf{88},
227005, (2002).

\bibitem{Ioffe02} L. B. Ioffe, and M. V. Feigelman, Phys. Rev. B \textbf{66}%
, 224503, (2002).

\bibitem{Matveev02} K. A. Matveev, A. I. Larkin, and L. I. Glazman, Phys.
Rev. Lett. \textbf{89}, 096802, (2002).

\bibitem{Protopopov04} I. V. Protopopov, and M. V. Feigelman, Phys. Rev. B 
\textbf{70}, 184519 (2004).

\bibitem{Protopopov06} I. V. Protopopov, and M. V. Feigelman, Phys. Rev. B 
\textbf{74}, 064516 (2006).

\bibitem{Pop08} I. M. Pop, K. Hasselbach, O. Buisson, W. Guichard, B.
Pannetier, and I. Protopopov, Phys. Rev. B \textbf{78}, 104504, (2008).

\bibitem{Kitaev06} A. Kitaev, arXiv:cond-mat/0609441.
\end{thebibliography}

\end{document}
