text:
# Discrete non-Abelian gauge theories in Josephson-junction arrays and quantum computation

## I Introduction

The ideal quantum computer Ekert1996 ; Steane1998 , if it exists, is a very interesting object from a physical point of view. It is a system which may be prepared in any superposition of $`2^K`$ (with $`K10^410^6`$) basic states during a very long time without decoherence. To give an order of magnitude, the decoherence should be much smaller than $`10^6`$ after a characteristic time scale of one operation Preskill1998 . Furthermore, it should allow a set of controlled unitary operations with an accuracy of $`10^6`$ or better. It is natural to think that in order to satisfy these stringent criteria one needs to have a physical system in which decoherence and errors generated by manipulations are exponentially small with respect to the system size or to any other control parameter. We shall refer to such properties as “ideal” in the following. The requirement of an ideal decoherence can be achieved in physical systems with topological order parameters Wen1990 ; Wen1991 ; Read1991  which can be implemented in Josephson-junction arrays Ioffe2002 ; Doucot2003 . However, the physical arrays proposed so far do not allow the ideal implementation of the full set of operations Ekert1996 ; Steane1998  required for quantum computations. In fact, these arrays allow to perform ideally only the simplest operations: a flip of an individual bit and the change of the wave function phase by $`\pm 1`$. It is crucial that they do not allow precise operations that mix the classical states, e.g., the Hadamard rotation $`U_H`$ that transforms the states $`|0`$, $`|1`$ into their mixtures $`(|1\pm |0)/\sqrt{2}`$ requires an application of the correct external perturbation $`V`$ for a given time $`t`$: $`U_H=\mathrm{exp}(iVt)`$ and thus is subjected to errors in both the values of $`t`$ and $`V`$.

Generally, an ideal operation involves an adiabatic change of some parameters that results in a generalized Berry’s phase. Note that all alternative ways to operate the qubit are similar to a spin rotation in an external field. Such operations are not adiabatic and are controlled by the duration of the field pulse which can never be made precise. In fact, it is very difficult to imagine how one can possibly achieve the accuracy of the transformation that is needed for quantum computation with this technique. In contrast, adiabatic transformations can be made, in principle, very precise. The main difficulty with them is to find a system in which the group formed by these transformation is big enough for a quantum computation. The arrays discussed above are implementing the symmetry group $`\mathrm{𝖹𝖹}_2`$ and thus it is not surprising that for them adiabatic changes of the external parameters allow only a very limited set of operations. One expects, following the suggestion by Kitaev Kitaev1997 , that one would get much larger sets of operations in systems allowing more complicated local symmetry transformations, in particular, in the case of discrete non-Abelian groups Propitius1999 .

The purpose of this paper is to consider lattice gauge theories built on more complicated groups than $`\mathrm{𝖹𝖹}_2`$ and to identify some that allow nontrivial adiabatic transformations. After having described the most promising ones, we discuss how they can be implemented, at least, in principle, in Josephson-junction arrays. All solid-state physics implementations of lattice gauge theories involve a lattice in real space but imply a continuous time. The natural framework to describe such models and their experimental realizations is the Hamiltonian formalism. This calls for an explicitly time-independent dynamics which can be achieved by choosing the scalar part of the gauge field $`A_0=0`$. In this case, the vector part $`A_{ij}`$ describes the quantum state of the physical links. However, the Hamiltonian thus constructed still commutes with the large group of local time-independent gauge transformations. It is then convenient to use a description which keeps track only of the underlying physical degrees of freedom. This can be done by restricting the Hilbert space to the space of gauge-invariant states. Note that several studies have already pointed out the existence of close connections between Josephson-junction arrays, lattice gauge theories and generalized statistics for elementary excitations. For instance, incompressible phases analogous to the fractional quantum Hall fluids (FQHF) of two-dimensional electronic systems have been predicted and studied in such arrays. One of them may be viewed as a FQHF of charges Odintsov  induced by an external magnetic field, whereas the other involves a FQHF of vortices Choi ; Stern  in the presence of a uniform fractional background charge on the superconducting islands. An appealing synthesis of these possible phases in the language of lattice gauge theories has been developed in Ref. Diamantini . As we shall show in detail, the lattice gauge models discussed here have the main advantage over the previous ones that effective interactions between elementary excitations are purely of statistical origin, and do not, for instance, depend on their precise locations. This is, of course, highly desirable in the context of quantum computation where the goal is to implement, with high accuracy, a set of universal unitary transformations.

This paper is organized as follows. In Sec. II, we review basic notions of continuous gauge groups and introduce basic notations. We emphasize the connection between the gauge-invariance constraint imposed on the physical states and the Gauss law constraint. Special features arising for non-Abelian groups are discussed. Section III focusses on the Abelian cyclic $`\mathrm{𝖹𝖹}_n`$ group. In this case, we show that imposing a local constraint on the product of group elements along the links of elementary plaquettes (kinematical Gauss law) is equivalent, after a duality transformation, to the requirement of local gauge invariance. We briefly discuss implementations of these two versions of $`\mathrm{𝖹𝖹}_n`$ models in Josephson-junction arrays and the existence of a topologically protected ground state when arrays exhibit macroscopic holes. Section IV addresses the interesting case of theories associated with a finite non-Abelian group. In this situation, the usual local gauge-invariance constraint and the kinematical Gauss law constraint have to be treated separately, since no exact duality transformation is available to connect them. We show that both models present gapped low-energy excitations called fluxons which may be trapped in macroscopic holes of the arrays. We also discuss local violations of the constraints called pseudocharges, and show the existence of a non-Abelian Aharonov-Bohm effect as a pseudocharge winds around a fluxon (and conversely). Section V is dedicated to implementations of these two models, for the simplest non-Abelian groups such as the dihedral groups $`D_3`$ and $`D_4`$. These proposals are based on the frustrated Josephson-junction arrays, where some basic elements have a $`2\pi /n`$ periodic potential energy as a function of the phase difference across them. The computations of the various energy scales driving the physical properties of these arrays, are presented in four appendixes.

## II Continuous gauge theories

Let us start by a brief review of the properties of lattice gauge theories with continuous groups using the language that will be convenient for the discussion of discrete groups. Therefore, we consider a Hamiltonian $`H`$ and a Hilbert space that are invariant under the local gauge transformations built from a continuous group $`G`$.

More precisely, to any set of group elements $`\{g_i\}`$’s labeled by lattice sites, we associate the transformation $`U_{\{g\}}`$ such that $`|\mathrm{\Psi }U_{\{g\}}|\mathrm{\Psi }=|\mathrm{\Psi }`$, and $`[H,U_{\{g\}}]=0`$. For a given lattice, the gauge field $`A`$ is attached to the links and takes its values in the group $`G`$. Local gauge transformations are then given by

$$A_{ij}g_iA_{ij}g_j^1A_{ij}^{}.$$
(1)
In the familiar example of compact quantum electrodynamics (QED) $`A_{ij}=\mathrm{exp}(ia_{ij})`$, $`g_i=\mathrm{exp}(i\varphi _i)`$ and gauge transformation in the variables ($`a_{ij}`$, $`\varphi _i`$) acquires a well-known form $`a_{ij}a_{ij}+\varphi _i\varphi _j`$. Identifying quantum states $`|\mathrm{\Psi }`$ with functions $`\mathrm{\Psi }`$ of classical field configurations $`\{A_{ij}\}`$, the condition of gauge invariance can be explicitly written as

$$\mathrm{\Psi }\{A_{ij}\}=\mathrm{\Psi }\{A_{ij}^{}\}\text{,}$$
(2)
for any choice of $`g_i`$’s. In particular, we may choose transformations localized on a single site $`i`$, so that $`g_j=e`$ if $`ji`$ and $`g_i=g`$, where $`e`$ denotes the neutral element in the group $`G`$. Furthermore, since $`G`$ is a continuous group, it is possible to consider infinitesimal $`g`$’s. The constraint of local gauge invariance (2) can then be expressed as

$$\underset{j(i)}{}\frac{}{a_{ij}^\alpha }|\mathrm{\Psi }=0\text{,}$$
(3)
where $`a_{ij}^\alpha `$ is the generator of the Lie group corresponding to the direction $`\alpha `$ in the vicinity of the neutral element of $`G`$. In the following, the notation $`j(i)`$ stands for all sites $`j`$ connected to $`i`$. The operator conjugated to $`a_{ij}^\alpha `$ is analogous to an electric field $`e_{ij}^\alpha =i/a_{ij}^\alpha `$. In the basis where the electrical field $`e^\alpha `$ is diagonal, Eq. (3) becomes a lattice equivalent of the usual Gauss law $`e^\alpha =0`$ which also reads

$$\underset{j(i)}{}E_{ij}^\alpha |\mathrm{\Psi }=|\mathrm{\Psi }\text{,}$$
(4)
where we have set $`E_{ij}^\alpha =\mathrm{exp}(ie_{ij}^\alpha t)`$, $`t`$ being an arbitrary real parameter. Notice that each site of the lattice corresponds to plaquette of the dual lattice, and the product over the links emanating from a site becomes a product over the links encircling the corresponding plaquette of the dual lattice. This allows for a natural realization of conditions (2) and (4) in physical systems.

At this stage, it is very important to point out a deep difference between Abelian and non-Abelian groups. In the former case, we may diagonalize simultaneously all the generators $`E_{ij}^\alpha `$ since they mutually commute, which is impossible in the non-Abelian case. Because of this, we shall explore two different routes in this work. A natural possibility is to regard Eq. (2) as the basic equation, and as we shall show, it is easy to use this construction for discrete non-Abelian groups. A second approach is to take Eq. (4) as the most fundamental. In this case, we assume that the physical states are also generated by classical configurations of group elements $`E_{ij}^\alpha `$ attached to the links, and we require these classical configurations to satisfy the kinematical constraint

$$\underset{j(i)}{}E_{ij}^\alpha =e\text{.}$$
(5)
Note that in this new scheme, and if $`G`$ is non-Abelian, we can no longer interpret this modified Gauss law as a consequence of a local gauge symmetry. So, these two approaches yield essentially different theories, except in the Abelian case where their equivalence can be precisely established as explained in the following section.

## III Discrete Abelian gauge theories

Let us start with the usual Abelian case and consider more specifically a simple example where the discrete gauge group is the cyclic group $`\mathrm{𝖹𝖹}_n`$, which has also the major advantage that it can be easily realized in a physical system. Indeed, this can be implemented with Josephson-junction elements for which the Josephson Hamiltonian is $`2\pi /n`$ periodic in the phase difference $`(\varphi _a\varphi _b)`$ where $`\varphi _a`$ denotes the phase of the superconducting order parameter at site $`a`$. We discuss the actual construction of such elements in Sec. V. Due to this $`2\pi /n`$-periodicity each link has $`n`$ degenerate ground states characterized by an integer, $`u_{ab}=0,\mathrm{},n1`$ which is defined by

$$\mathrm{\Delta }\varphi _{ab}=\varphi _a\varphi _b=2\pi u_{ab}/n\text{.}$$
(6)
We would like to stress that this type of system exhibits a local $`\mathrm{𝖹𝖹}_n`$ symmetry defined by $`\varphi _a\varphi _a+2\pi m_a/n`$ with integers $`m_a`$. Since we are dealing with quantum systems, it is also natural to interpret this symmetry in a basis which diagonalizes the local particle number operators conjugated to the phase operators. This symmetry then simply means that the particle number on each lattice site is conserved modulo $`n`$. An explicit realization for the $`n=2`$ case has been recently proposed in Refs. Doucot2002 ; Ioffe2002 .

In real systems, charging energies induce quantum fluctuations of the phase variables. Tunneling processes between the various degenerate classical ground states discussed above appear and usually lift completely this classical degeneracy. Physical properties in the low-energy subspace spanned by these classical ground states are then well described by the following Hamiltonian:

$$H_{sites}=\underset{a\mathrm{\Lambda }}{}\underset{m=0}{\overset{n1}{}}r(m)𝒯_a(m)=\underset{a}{}H_a\text{,}$$
(7)
where $`\mathrm{\Lambda }`$ denotes the physical lattice. In order to give a concrete illustration, we work on the geometry displayed in Fig. 1 where $`\mathrm{\Lambda }`$ is the hexagonal lattice and its dual $`\mathrm{\Lambda }_d`$ is the triangular one. Here, the operator $`𝒯_a(m)`$ is defined by

$$𝒯_a(m)|\varphi _a=|\varphi +2\pi m/n_a\text{,}$$
(8)
and $`r(m)`$ is a tunneling amplitude which is expected to be maximal for $`m\pm 1[n]`$.

In this case, the Hamiltonian $`H_a`$ is easily diagonalizable since its eigenstates are “plane waves”

$$|l_a=\frac{1}{\sqrt{n}}\underset{m=0}{\overset{n1}{}}\mathrm{exp}(i2\pi lm/n)|\varphi +2\pi m/n_a\text{,}$$
(9)
and its eigenenergies are given by $`E(l)=2r\mathrm{cos}(2\pi l/n)`$ for $`l=0,\mathrm{},n1`$. Here, we have chosen tunneling amplitudes $`r(m)`$ equal to $`r`$ for $`m=\pm 1`$ and $`0`$ otherwise. Thus the ground state $`|𝒢=_a|𝒢_a`$ is given by a configuration such that $`l=0`$ for all $`a`$.

Although this model is very easy to solve since it only involves single-site operators which mutually commute, it is very fruitful to map it onto a genuine local gauge theory also based on the $`\mathrm{𝖹𝖹}_n`$ group. Indeed, this mapping provides a convenient way to treat physical external perturbations which are coupled to link variables. This construction is based on a description of physical states in terms of the link variables $`u_{ab}`$. The crucial observation is that these variables are not independent as soon as there are loops on the lattice. More precisely, the sum of the integers $`u_{ab}`$’s along any oriented loop is constrained to vanish modulo $`n`$. A natural basis for our system is given by states of the form $`|\mathrm{\Psi }=_{ab}|u_{ab}`$, where the $`u_{ab}`$’s satisfy the constraint. In the following, we shall refer to the description in terms of the $`\varphi _a`$’s as “site” representation, and to the latter (in terms of the $`u_{ab}`$’s) as the “link” representation. This condition on the $`u_{ab}`$’s can be easily interpreted as the Gauss law for an Abelian gauge theory on the dual lattice. To show this, it is convenient to change the basis via the Fourier transform

$$|p_{ab}=\frac{1}{\sqrt{n}}\underset{u=0}{\overset{n1}{}}\mathrm{exp}(i2\pi pu/n)|u_{ab}\text{,}$$
(10)
where $`p`$ is an integer chosen modulo $`n`$. In this basis, a state $`|u_{ab}`$ is an eigenvector for the translation operators $`𝒯_{ab}^{(p)}(p)`$ with the eigenvalue $`\mathrm{exp}(i2\pi pu/n)`$ which provides a simple implementation of these operators in terms of the original fields. Thereafter, we will use the notation $`𝒯_{ab}^{(x)}(y)`$ for the translation operator defined by

$$𝒯_{ab}^{(x)}(y)|x_{ab}=|x+y_{ab}\text{,}$$
(11)
for $`x=u,p`$. Note that $`[(𝒯_{ab}^{(p)}(p)]^n=𝒯_{ab}^{(p)}(np)=\mathrm{𝟏}`$, so that the group composed of the $`𝒯_{ab}^{(p)}`$’s is also $`\mathrm{𝖹𝖹}_n`$. The constraint

$$\underset{plaq.}{}u_{ab}0[n]\underset{plaq.}{}\mathrm{exp}(i2\pi u_{ab}/n)=1\text{,}$$
(12)
for any elementary plaquette is equivalent to requiring the invariance of physical states under the action of the operator $`_{plaq.}𝒯_{ab}^{(p)}(1)`$ which shifts all the $`p_{ab}`$’s around the plaquette by 1.

The last stage of the mapping onto a $`\mathrm{𝖹𝖹}_n`$ gauge theory consists in expressing everything on the dual lattice. Any link $`ab`$ on $`\mathrm{\Lambda }`$ separates two plaquettes which correspond to sites $`i`$ and $`j`$ of $`\mathrm{\Lambda }_d`$. The oriented links of $`\mathrm{\Lambda }`$ and $`\mathrm{\Lambda }_d`$ are thus in one to one correspondence. It is then possible to transform the $`p_{ab}`$’s into the $`p_{ij}`$’s, and similarly the $`𝒯_{ab}`$’s into the $`𝒯_{ij}`$’s. We then recover the Gauss law written as

$$\underset{j(i)}{}𝒯_{ij}^{(p)}(1)|\mathrm{\Psi }=|\mathrm{\Psi }\text{.}$$
(13)
In other words, invariance of physical states under $`_{plaq.}𝒯_{ab}^{(p)}(1)`$ is equivalent to the requirement of invariance under the transformations $`p_{ij}p_{ij}+l_il_j`$ for any pairs of neighboring sites $`i`$ and $`j`$ belonging to $`\mathrm{\Lambda }_d`$, and any integers $`l_i`$ and $`l_j`$. This is exactly the equivalent of condition (2) for a gauge theory associated with the group $`\mathrm{𝖹𝖹}_n`$ on the dual lattice where the gauge “fields” are the $`p_{ij}`$’s.

The simplest Hamiltonian invariant under these transformations reads

$$H=\underset{a}{}\underset{u=0}{\overset{n1}{}}r(u)\underset{b(a)}{}𝒯_{ab}^{(u)}(u)+\underset{ab}{}F(\widehat{u}_{ab})\text{.}$$
(14)
The first term gives quantum dynamics to the links and is equivalent to $`H_{sites}`$. The second term $`F(\widehat{u}_{ab})`$ describes perturbations which lift the degeneracy between the classical ground states. In the physical implementations in Josephson-junction arrays discussed in more detail below, it is usually due to the unprecise value of the flux through the elementary plaquette. For instance, in the realization Ioffe2002  of $`\mathrm{𝖹𝖹}_2`$ theory it is equal to $`u_{ab}E_J\delta \mathrm{\Phi }/\mathrm{\Phi }_0`$. Generally, this term breaks the invariance $`\varphi _a\varphi _a+2\pi m_a/n`$ whereas it preserves the local gauge symmetry $`p_{ij}p_{ij}+l_il_j`$ since it is diagonal in the $`|u`$ representation.

The Hamiltonian $`H`$ can also be written as

$$H=\underset{ijk}{}\stackrel{~}{r}(\widehat{p}_{ij}+\widehat{p}_{jk}+\widehat{p}_{ki})+\frac{1}{n}\underset{ij}{}\underset{p=0}{\overset{n1}{}}\stackrel{~}{F}(p)𝒯_{ij}^{(p)}(p)\text{,}$$
(15)
where $`\stackrel{~}{r}`$ and $`\stackrel{~}{F}`$ denote the Fourier transform defined by

$$\stackrel{~}{K}(p)=\underset{u=0}{\overset{n1}{}}\mathrm{exp}(i2\pi pu/n)K(u)\text{,}$$
(16)
for any function $`K`$. As usual, for a lattice gauge theory, the first term involves the “magnetic” flux $`(\widehat{p}_{ij}+\widehat{p}_{jk}+\widehat{p}_{ki})`$ on the plaquettes $`ijk`$ of $`\mathrm{\Lambda }_d`$, whereas the second term is off diagonal in the $`|p`$ representation and is associated with “electric” fields.

In the absence of this latter term, the model can be solved and its ground state, in the $`|u`$ representation, reads

$`|𝒢`$ $`=`$ $`\left[{\displaystyle \underset{plaq.}{}}\delta \left({\displaystyle \underset{(ab)plaq.}{}}u_{ab}\right)\right]_{ab}{\displaystyle \underset{u=0}{\overset{n1}{}}}|u_{ab}`$ (17)
$`=`$ $`\left[{\displaystyle \underset{i}{}}\delta \left({\displaystyle \underset{j(i)}{}}u_{ij}\right)\right]_{ij}{\displaystyle \underset{u=0}{\overset{n1}{}}}|u_{ij}\text{,}`$ (18)
where the first factor ensures the constraint (12). The ground state can also be simply written in the $`|p`$ representation on the dual lattice

$$|𝒢=\underset{i}{}𝒫_i_{ij}|0_{ij}\text{,}$$
(19)
where

$$𝒫_i=\underset{p=0}{\overset{n1}{}}\underset{j(i)}{}𝒯_{ij}^{(p)}(p)\text{,}$$
(20)
is the projector on the gauge-invariant subspace defined by Eq. (13), and where $`|0`$ stands for the states $`|p=0`$.

The elementary excitations that we shall call fluxons correspond to plaquettes on $`\mathrm{\Lambda }_d`$ for which the corresponding flux is nonvanishing. Such an excitation is produced by an open string of operators

$$\widehat{f}=\underset{(ab)\gamma }{}\mathrm{exp}(i2\pi \widehat{u}_{ab}/n)=\underset{(ab)\gamma }{}𝒯_{ab}^{(p)}(1)\text{,}$$
(21)
where the contour $`\gamma `$ goes over the links of the physical lattice $`\mathrm{\Lambda }`$ that is, over the centers of plaquettes of $`\mathrm{\Lambda }_d`$ on which the gauge field $`p_{ij}`$ is defined (see Fig. 2). Assuming $`r(m)=r`$ if $`m=\pm 1`$ and $`0`$ otherwise, these elementary excitations have energy $`E=2r[1\mathrm{cos}(2\pi /n)]`$. Physically, in the system of Josephson junctions these excitations carry electrical charge $`2e`$. We emphasize, however, that in terms of the equivalent gauge theory these excitations carry no charge, they describe the fluctuations of the gauge field without matter.

Matter appears in this effective theory if one allows sites of the dual lattice where constraint (13) is violated. In terms of the physical Josephson-junction system, it means that we have plaquettes of the physical lattice $`\mathrm{\Lambda }`$ where the sum of the superconducting phases is different from $`0`$, that is vortices in the superconducting phase. Notice that due to $`2\pi /n`$ periodicity these vortices carry $`1/n`$ flux quanta. Thus, it is natural to introduce sites $`i`$ where $`_{j(i)}u_{ij}1[n]`$. In the gauge field language, it implies that gauge transformation at these sites changes the wave function of the gauge field according to

$$\underset{j(i)}{}𝒯_{ij}^{(p)}(p)|\mathrm{\Psi }=\mathrm{exp}(i2\pi p/n)|\mathrm{\Psi }\text{,}$$
(22)
meaning that, at these sites, there is a pseudocharge $`1/n`$. In order to create such an excitation above the ground state, one needs to apply operators $`𝒯_{ij}^{(u)}(1)`$ along an open string on the dual lattice $`\mathrm{\Lambda }_d`$. Because of the presence of this string, it is in general necessary to consider globally neutral charge configurations. This arises because the product of all the generators $`_{j(i)}𝒯_{ij}^{(p)}`$ over all sites $`i`$ is the identity operator, or equivalently because global gauge transformations do not modify the link variables $`p_{ij}`$. When such a pseudocharge moves around a fluxon, it undergoes an Aharonov-Bohm effect with a dephasing $`\mathrm{exp}(i2\pi /n)`$. This is a natural result for a Josephson-junction array where a fluxon carries a physical charge $`2e`$ and where a pseudocharge is a superconducting vortex with flux quanta $`\mathrm{\Phi }_0/n`$. The same result is obtained if we move a fluxon around a fixed pseudocharge as in the Aharonov-Casher effect. A more detailed discussion of this phenomenon is given in the context of non-Abelian gauge theories in Sec. IV.

The classification of the excited states as fluxons and pseudocharges is meaningful if these charges are rare. This is indeed the case if the energy needed to create a vortex of the superconducting order parameter (which becomes pseudocharge in the gauge theory language) is large compared with the Coulomb energy and also with the energy scale $`r`$ of phase fluctuations. Consider qualitatively the opposite case in which the tunneling amplitude $`r`$ is large compared to the potential energy of vortices and to their kinetic energy $`t`$. In a physical Josephson-junction array, this can be realized by constructing the links of the lattice from rhombi Ioffe2002  connected in series to an additional weak Josephson junction such that its Josephson energy $`ϵ_J`$ and its charging energy $`ϵ_C`$ satisfy $`ϵ_J<r`$ and $`ϵ_Jϵ_C`$ Doucot2003 . In these conditions all energy scales associated with weak junctions ($`ϵ_C`$, $`ϵ_J`$, vortex tunneling amplitude $`t`$) are small compared to $`r`$. Thus, in the leading approximation one can neglect the dynamics of the weak junctions and their contribution to the potential energy. The leading process remains the simultaneous flip of three adjacent rhombi described by the first term of the Hamiltonian (14). The condition $`_{j(i)}u_{ij}0[n]`$ no longer holds because weak links favor the spontaneous excitation of vortices. Instead, it is replaced by the condition that all low-energy states have to minimize the kinetic energy part of the Hamiltonian (14), i.e., to satisfy

$$\underset{b(a)}{}𝒯_{ab}^{(u)}(1)|\psi =|\psi .$$
(23)
This condition can be viewed as a Gauss law imposed on the sites of the physical lattice under the gauge group generated by $`𝒯_{ab}^{(u)}`$ \[cf. Eq. (13)\]. In these conditions, the physical array is expected to be in the insulating regime. It is therefore natural that the discrete version (23) of the QED gauge invariance (3) should be restored. Formally, the states satisfying the gauge-invariance condition (23) are described by the same theory as the states satisfying Eq. (13) but the physical meaning of all quantities is reversed. Fluxons correspond now to the fractional vortices in the superconducting phase, while pseudocharges become the physical charges of the insulating array (modulo $`n`$). Finally, note that the degenerate states appearing in arrays with a nontrivial topology are distinguished by the product of $`𝒯_{ij}^{(u)}(1)`$ operators taken over paths defined on the dual lattice.

To conclude this section, we emphasize that when the gauge group is Abelian, there is a duality between the theory defined on the “real lattice” and the theory defined on its dual. This correspondence is established via the Fourier transform which relates the variables $`u`$ and $`p`$ defined on the real and on the dual lattices, respectively. In addition, when both structures are isomorphic, the gauge theory is self-dual, but it is not the case for an arbitrary system. In this sense, $`\mathrm{𝖹𝖹}_n`$ lattice gauge theories on an hexagonal lattice are not self-dual. However, the important characteristic required for the physical operations discussed in our paper is not self-duality but duality which is always satisfied.

## IV Discrete non-Abelian gauge theories

In this section, we consider the case of a finite non-Abelian gauge group $`G`$. As already discussed in Sec. II, it then becomes necessary to consider models resulting from the gauge-invariance requirement and the kinematical Gauss law constraint separately.

### IV.1 Gauge-invariant model

The requirement of local gauge invariance simply means that the wave function $`|\mathrm{\Psi }`$ satisfies $`U_{\{g\}}|\mathrm{\Psi }=|\mathrm{\Psi }`$ for any configuration of the $`g_i`$’s, where $`U_g`$ is defined by

$$U_g=\underset{i}{}\underset{j(i)}{}𝒯_{ij}(g_i)\text{.}$$
(24)
Here, we have adopted multiplicative notation for the group law and we have introduced the left “translation” operator $`𝒯_{ij}`$ defined by

$$𝒯_{ij}(h)|g_{ij}=|hg_{ij}\text{.}$$
(25)
As previously, link variables are oriented in the sense that $`|g_{ij}=|g^1_{ji}`$. The Hamiltonian invariant under those transformations has the form similar to Eq. (15):

$$H=\underset{ijk}{}\stackrel{~}{r}(\widehat{g}_{ij}\widehat{g}_{jk}\widehat{g}_{ki})+\underset{ij}{}\underset{g}{}\stackrel{~}{F}(g)𝒯_{ij}(g)\text{,}$$
(26)
where $`\stackrel{~}{r}`$ and $`\stackrel{~}{F}`$ are now required to be invariant under inner group transformations $`ghgh^1`$ for any $`h`$ in $`G`$. As before, $`\stackrel{~}{r}`$ is associated with the energy cost of a fluxon-type excitation. Note that for a continuous group, the $`g_{ij}`$’s are naturally represented by matrices and $`\stackrel{~}{r}`$ is often chosen as the trace function. By contrast, for a discrete group, we choose here the simplest case for which $`\stackrel{~}{r}(e)=r`$ and $`0`$ otherwise, a choice that differs from previously where we had $`\stackrel{~}{r}(p)=2r\mathrm{cos}(2\pi p/n)`$. As in the Abelian case, this model can be solved in the limit $`\stackrel{~}{F}0`$, where one can construct the ground state explicitly

$$|𝒢=\underset{i}{}𝒫_i_{ij}|e_{ij}\text{,}$$
(27)
where

$$𝒫_i=\underset{h}{}\underset{j(i)}{}𝒯_{ij}(h)\text{,}$$
(28)
is the projector on the gauge-invariant subspace. Note that this expression of the ground state is similar to the one given in Eq. (19) for the Abelian $`\mathrm{𝖹𝖹}_n`$ group.

As in the Abelian case, the excitations are obtained by changing the value of the $`\stackrel{~}{r}`$ on one plaquette. Nevertheless, the non-Abelian nature of $`G`$ prevents us to write down a compact form for the fluxon creation operator $`\widehat{f}`$ which would generalize expression (21). Instead, we first consider the state $`_{(kl)\gamma }𝒯_{kl}(g)_{ij}|e_{ij}`$ where the product goes over the links that are intersected by contour $`\gamma `$ drawn on $`\mathrm{\Lambda }`$ that begins at the boundary of the lattice and ends at the excited plaquette. Then, one projects it onto the allowed Hilbert space getting

$$|f=\underset{i}{}𝒫_i\underset{(kl)\gamma }{}𝒯_{kl}(g)_{ij}|e_{ij}\text{.}$$
(29)
Since the projection operator $`_i𝒫_i`$ is gauge-invariant, one can insert any gauge transformation operator between it and $`_{(kl)\gamma }𝒯_{kl}(g)`$ without changing $`|f`$. Inserting the global gauge transformation with a group element $`h`$, we transform the string of $`𝒯_{kl}(g)`$ operators into the string of $`𝒯_{kl}(g^{})`$ with $`g^{}=h^1gh`$. Thus, all group elements that belong to the same conjugacy class result in the same excitation with energy $`E=r`$.

It is useful to give a more intrinsic description of these fluxons by considering generalized magnetic fluxes along closed paths. For such a path $`\gamma _c`$, the generalized magnetic flux is defined as $`\mathrm{\Phi }(\gamma _c)=_{(ij)\gamma _c}g_{ij}`$. In the absence of fluxon, one has $`\mathrm{\Phi }(\gamma _c)=e`$ independently of $`\gamma _c`$. In the fluxon state $`|f`$, one has $`\mathrm{\Phi }(\gamma _c)=e`$ only if $`\gamma _c`$ does not encircle the plaquette on which the excitation is localized. If the path encloses this plaquette once, then $`\mathrm{\Phi }(\gamma _c)`$ belongs to the conjugacy class of $`g`$, where $`g`$ is the same group element as in Eq. (29) provided the orientation of $`\gamma _c`$ coincides with the one of the string operator used to build $`|f`$. Let us mention three important properties of these generalized fluxes. For this purpose, we need to consider closed paths $`\gamma _c`$ which begin and end at a fixed arbitrary site $`O`$. In this case, we shall say that $`O`$ is the origin of $`\gamma _c`$ or that $`\gamma _c`$ is based at site $`O`$. We have then the following

1. If $`\gamma _c`$ and $`\gamma _c^{}`$ are both based at site $`O`$, and if one may deform $`\gamma _c`$ into $`\gamma _c^{}`$ without going through any localized fluxon, then $`\mathrm{\Phi }(\gamma _c)=\mathrm{\Phi }(\gamma _c^{})`$.
2. If $`\gamma _c`$ and $`\gamma _c^{}`$ involve the same closed loop but with different origins $`O`$ and $`O^{}`$, then $`\mathrm{\Phi }(\gamma _c)`$ and $`\mathrm{\Phi }(\gamma _c^{})`$ belong to the same conjugacy class.
3. $`\mathrm{\Phi }(\gamma _c)`$ is not gauge-invariant in general, but its conjugacy class is gauge-invariant.

For two or more fluxons, the quantum state of the pure gauge theory is no longer determined by the knowledge of the conjugacy class associated with each individual fluxon. This arises because if $`g_1`$ is conjugated to $`g_1^{}`$ and $`g_2`$ to $`g_2^{}`$, then $`g_1g_2`$ is not always conjugated to $`g_1^{}g_2^{}`$. Although all these conclusions are valid for an arbitrary group $`G`$, we shall explicitly consider here only the simplest non-Abelian group $`D_n`$ $`(n3)`$ which contains $`n`$ rotations $`\zeta ^k`$ ($`k=0,\mathrm{},n1`$) and $`n`$ rotations combined with reflections $`\tau \zeta ^k`$. The two generators $`\zeta `$ and $`\tau `$ satisfy the three relations: $`\zeta ^n=e`$, $`\tau ^2=e`$, and $`\tau \zeta =\zeta ^1\tau `$, where $`e`$ is the unit element of the group. For simplicity, let us choose the group $`D_3`$ which contains six elements and three conjugacy classes: $`\{e\}`$, $`\{\zeta ,\zeta ^2\}`$, and $`\{\tau ,\tau \zeta ,\tau \zeta ^2\}`$. If we take $`g_1=g_1^{}=g_2=\tau `$ and $`g_2^{}=\tau \zeta `$, then $`g_1g_2=e`$ and $`g_1^{}g_2^{}=\zeta `$ belong to different conjugacy classes. A complete description of the two-fluxon state requires then to specify the conjugacy classes of the fluxes associated with three sets of paths on $`\mathrm{\Lambda }_d`$: the paths which encircle each fluxon only once without going around the other, and the paths which encircle both fluxons only once (see Fig. 3).

Here, we have used the properties (1) and (2) of general fluxes given above. The non-Abelian nature of $`G`$ induces then a notion of nonlocality when we combine several fluxons. In this sense, these objects do not behave as ordinary independent quasiparticles when we consider their global quantum statistical properties Propitius1999 .

More precisely, the state of the system is characterized by the values of flux calculated for a set of closed paths that begin and end at a given site $`O`$ of the lattice and avoid the fluxon positions. The value of flux is the same for all topologically equivalent paths (with the same origin), i.e., for the paths that can be smoothly deformed into each other without going over the fluxon positions. The set of topologically nonequivalent paths has a natural group structure (provided by path concatenation), called the fundamental group of the lattice with $`n`$ plaquettes removed, we shall denote it by $`\pi _1(n)`$. As is well known Mermin1979 , $`\pi _1(n)`$ is generated by the equivalence classes of $`n`$ paths, each of which winds exactly once around only one particular fluxon. But for $`n2`$, this group is non-Abelian as illustrated in Fig. 4. Notice that specifying a set of winding numbers for a path with respect to each fluxon does not determine the total flux along this path which depends on the order in which the path encircles the fluxons. Formally, any configuration of link variables $`g_{ij}`$ with $`n`$ fluxons establishes a (non-Abelian) group homomorphism from $`\pi _1(n)`$ to $`G`$, which sends the equivalence class of a path $`\gamma `$ based at $`O`$ into $`\mathrm{\Phi }(\gamma )`$. Clearly, $`\mathrm{\Phi }(\gamma _1\gamma _2)=\mathrm{\Phi }(\gamma _1)\mathrm{\Phi }(\gamma _2)`$. In the simpler case of Abelian $`G`$, this subtle effect due to the non-Abelian nature of $`\pi _1(n)`$ is washed out and fluxes do depend only on the set of winding numbers and not on their actual sequence along a path. Finally, we wish to characterize gauge-invariant states with $`n`$ fluxons. As already noted, a gauge transformation $`g_{ij}^{}=h_ig_{ij}h_j^1`$ transforms fluxes $`\mathrm{\Phi }(\gamma )`$ into $`\mathrm{\Phi }^{}(\gamma )=h_0\mathrm{\Phi }(\gamma )h_0^1`$. Thus, equivalence classes of group homomorphisms from $`\pi _1(n)`$ to $`G`$ are in one to one correspondence with the gauge-invariant states for $`n`$ fluxons.

Let us illustrate this for $`n=2`$ and $`G=D_3`$. If fluxons were described by group elements, we would expect $`6^2=36`$ states. But this clearly does not take the gauge-invariance condition into account. If we naively associate conjugacy classes to fluxons, this yields $`3^2=9`$ states for two fluxons. It turns out that the correct counting is $`11`$. We list below the corresponding equivalence classes of group homomorphisms from $`\pi _1(2)`$ to $`D_3`$. Since $`\pi _1(2)`$ has two independent generators $`\gamma _1`$ and $`\gamma _2`$ (see Fig. 4), this homomorphism is completely determined by the knowledge of $`\mathrm{\Phi }(\gamma _1)`$ and $`\mathrm{\Phi }(\gamma _2)`$. In these notations, the eleven possible states of two fluxons are:

$$\begin{array}{c}\mathrm{\Phi }_{00}=\{(e,e)\},\hfill \\ \mathrm{\Phi }_{01}=\{(e,\zeta ),(e,\zeta ^2)\},\hfill \\ \mathrm{\Phi }_{10}=\{(\zeta ,e),(\zeta ^2,e)\},\hfill \\ \mathrm{\Phi }_{02}=\{(e,\tau ),(e,\tau \zeta ),(e,\tau \zeta ^2)\},\hfill \\ \mathrm{\Phi }_{20}=\{(\tau ,e),(\tau \zeta ,e),(\tau \zeta ^2,e)\},\hfill \\ \mathrm{\Phi }_{11}=\{(\zeta ,\zeta ),(\zeta ^2,\zeta ^2)\},\hfill \\ \mathrm{\Phi }_{11}^{}=\{(\zeta ,\zeta ^2),(\zeta ^2,\zeta )\},\hfill \\ \mathrm{\Phi }_{22}=\{(\tau ,\tau ),(\tau \zeta ,\tau \zeta ),(\tau \zeta ^2,\tau \zeta ^2)\},\hfill \\ \mathrm{\Phi }_{22}^{}=\{(\tau ,\tau \zeta ),(\tau ,\tau \zeta ^2),(\tau \zeta ,\tau ),(\tau \zeta ,\tau \zeta ^2),(\tau \zeta ^2,\tau ),\hfill \\ (\tau \zeta ^2,\tau \zeta )\},\hfill \\ \mathrm{\Phi }_{12}=\{(\zeta ,\tau ),(\zeta ,\tau \zeta ),(\zeta ,\tau \zeta ^2),(\zeta ^2,\tau ),(\zeta ^2,\tau \zeta ),(\zeta ^2,\tau \zeta ^2)\},\hfill \\ \mathrm{\Phi }_{21}=\{(\tau ,\zeta ),(\tau \zeta ,\zeta ),(\tau \zeta ^2,\zeta ),(\tau ,\zeta ^2),(\tau \zeta ,\zeta ^2),(\tau \zeta ^2,\zeta ^2)\}.\hfill \end{array}$$
(30)
Note that in this list there are two pairs of different two-fluxon states, $`\mathrm{\Phi }_{11}`$, $`\mathrm{\Phi }_{11}^{}`$ and $`\mathrm{\Phi }_{22}`$, $`\mathrm{\Phi }_{22}^{}`$ which are produced by the combination of nondistinguishable one-fluxon states, namely, states with fluxes $`\mathrm{\Phi }_1=\{\zeta ,\zeta ^2\}`$ and $`\mathrm{\Phi }_2=\{\tau ,\tau \zeta ,\tau \zeta ^2\}`$.

As usual, including static local charges in such a model requires to consider sites which transform differently under the gauge group. More precisely, these sites can be classified according to the representations of the gauge group. For a non-Abelian group, most irreducible representations are no longer one dimensional. This leads us to attach a finite-dimensional representation space spanned by basis vectors $`|\alpha `$ to each of these sites. The corresponding representation of $`G`$ in this subspace is described as usual by a set of matrices $`D_{\alpha \beta }(g)`$. Now, let us assume that we have only one such site. The global wave function may be written as $`|\mathrm{\Psi }=_\alpha |\alpha |\mathrm{\Psi }_\alpha `$. In this slightly generalized setting, the local gauge-invariance constraint becomes $`|\mathrm{\Psi }_\alpha =_\beta D_{\alpha \beta }(g)U_g|\mathrm{\Psi }_\beta `$, if the group element $`g`$ is attached to a site carrying a nontrivial representation. For the $`D_n`$ group that we consider here, the most interesting cases are obtained from the two-dimensional irreducible representation of the group. The wave function of the ground state satisfying this constraint at one site $`o`$, occupied by a pseudocharge in the $`D(g)`$ representation may be written explicitly as

$$|𝒢=\left(\underset{io}{}𝒫_i\right)\stackrel{~}{𝒫}_o|\alpha _{ij}|e_{ij}\text{,}$$
(31)
where $`|\alpha `$ is any vector in the charge Hilbert space, and $`\stackrel{~}{𝒫}_o=_hD(h)_{j(o)}𝒯_{oj}(h)`$.

As for the Abelian case, it is interesting to discuss the situation where a pseudocharge moves around a fluxon and conversely. The state combining a pseudocharge $`|\alpha `$ located at site $`i`$ and the gauge degrees of freedom in the state $`|\mathrm{\Psi }`$ will be denoted by $`|\alpha _i|\mathrm{\Psi }`$. In a trivial background where $`|\mathrm{\Psi }=_{ij}|e_{ij}`$, it is natural to assume that the charge hopping Hamiltonian $`H_t`$ acts according to

$$H_t(|\alpha _i|\mathrm{\Psi })=t\underset{j(i)}{}|\alpha _j|\mathrm{\Psi }\text{,}$$
(32)
where $`t`$ is a positive hopping amplitude. Requiring invariance of $`H_t`$ under local gauge transformations completely determines the action of $`H_t`$ for a nontrivial gauge field background $`|\mathrm{\Psi }`$:

$$H_t(|\alpha _i|\mathrm{\Psi })=t\underset{j(i)}{}[D(g_{ij}^1)|\alpha ]_j|\mathrm{\Psi }\text{.}$$
(33)
Because of the string involved to build the fluxon $`|f`$, it is clear that a pseudocharge moving around it experiences a generalized Aharonov-Bohm effect which is, in general, nondiagonal in the Hilbert space associated with the internal degrees of freedom of the charge. Conversely, we may consider the case of a fluxon associated with the conjugacy class of $`g`$ and moving around a fixed charge. When the fluxon moves from one plaquette to an adjacent one, we obtain the final state by modifying only the field variable $`g_{ij}`$ on the link $`(ij)`$ that has been crossed by the fluxon. This defines an operator $`_{ij}`$ which is equal to $`𝒯_{ij}(g^{})`$ where $`g^{}`$ is a group element belonging to the same conjugacy class as $`g`$. Note that the choice of this representative element inside the class of $`g`$ depends on the actual field configuration. Consider, for instance, the example of two plaquettes shown in Fig. 5.

Thanks to this construction, $`_{ij}`$ commutes with all the local gauge transformations operators $`U_{\{h\}}`$ acting on the gauge degrees of freedom. Let us take an initial state $`𝒫(|\alpha _i|\mathrm{\Psi })`$, where $`𝒫=\left(_{io}𝒫_i\right)\stackrel{~}{𝒫}_o`$. We move the fluxon along a closed path $`\gamma _c`$ which encircles the origin $`o`$ where the pseudocharge is located. This leads to the state:

$`\left({\displaystyle \underset{(ij)\gamma _c}{}}_{ij}\right)𝒫(|\alpha _o|\mathrm{\Psi })`$ $`=`$ (34)
$`𝒫(|\alpha _o`$ $``$ $`\left({\displaystyle \underset{(ij)\gamma _c}{}}_{ij}\right)|\mathrm{\Psi })\text{.}`$
Since $`|\mathrm{\Psi }`$ and $`_{ij}|\mathrm{\Psi }`$ both describe the same fluxon, there is a gauge transformation $`U_{\{h\}}`$ connecting these two states. We choose $`\{h\}`$ such that $`h_i=e`$ for $`i`$ outside the area enclosed by the fluxon trajectory. In this case, $`h_o`$ is conjugated to the flux $`\mathrm{\Phi }=g`$ associated with the fluxon. Our final state is then $`𝒫[D(h_o^1)|\alpha _o|\mathrm{\Psi }]`$. We thus recover the same non-Abelian transformation after a closed cycle as for the moving charge situation.

Furthermore, non-Abelian Aharonov-Bohm effects also appear when a fluxon goes around another fluxon Bais1980 . This situation is shown in Fig. 6, where fluxon 1 carrying initially the flux $`\mathrm{\Phi }_1=g_1`$ goes around fluxon 2 (with initial flux $`\mathrm{\Phi }_2=g_2`$) once counterclockwise. These fluxes are defined as a product of group variable around contours $`\gamma _1`$ and $`\gamma _2`$ starting and ending at point $`O`$ or any other contours topologically equivalent to them. The state with these fluxes can be represented by the configuration of $`g_{ij}`$ that contains two strings against the “flat”, $`g_{ij}=e`$, background. As the first fluxon moves around the second one, its string eventually crosses the string of the static fluxon. After the crossing the string of the first fluxon changes into the string of $`G_1=g_2^1g_1g_2`$ elements. The new fluxes $`\mathrm{\Phi }_1`$ and $`\mathrm{\Phi }_2`$ become $`\mathrm{\Phi }_1=g_2^1g_1g_2`$ and $`\mathrm{\Phi }_2=g_2^1g_1^1g_2g_1g_2`$. Since $`g_2^1g_1g_2=g_2^1g_1^1g_1g_1g_2`$, this state can be reduced to the initial state by a global gauge transformation with the group element $`g=g_1g_2`$, if there are only two fluxons in the system. But if some other fluxons are present somewhere else in the system, besides the two already discussed, they might transform nontrivially under the gauge transformation with $`g`$: $`\mathrm{\Phi }g^1\mathrm{\Phi }g`$ and in a general case the final state cannot be reduced to the initial state by the gauge transformation. This provides another example of nondiagonal Aharonov-Bohm effect, even in the absence of external charges. These considerations have been generalized to composite objects of fluxons and pseudocharges usually called dyons Propitius1999 . The appropriate mathematical framework has been developped in terms of representation theory of the so-called quantum double algebra associated with the group $`G`$ Dijkgraaf1990 ; Bais1992 . We shall not elaborate further on these interesting topics here.

As for Abelian models, the systems with one or several holes have a degenerate ground state. In the limit of zero off-diagonal term in Hamiltonian (26) the fluxons are trapped forever in the holes and have zero energy. The degeneracy of the ground state with one hole is thus equal to the number of different fluxons, i.e., to the number of conjugacy classes. The operator that distinguishes these ground states is $`_Cg`$ where contour $`C`$ encirles the hole. A small but nonzero off-diagonal term in the Hamiltonian leads to an exponentially small amplitude for the fluxon escape from the hole. The operator leading to this escape process is the product of $`𝒯_{ab}(g)`$ along the contour $`\gamma `$ on lattice $`\mathrm{\Lambda }`$ that connects inner and outer boundaries as shown in Fig. 2. This is the same operator that moves fluxons around the system, as discussed above. The gap for the fluxons implies that this process is only virtual and its amplitude is exponentially small. Therefore, the resulting degeneracy splitting is also exponentially small. Note that (as discussed above) the choice of the group element $`g`$ in fluxon motion operator $`𝒯_{ab}(g)`$ depends on the link variables. The contour $`\gamma `$ crosses once any contour $`C`$ and thus changes the flux trapped inside the hole. Since these fluxons undergo non-Abelian Aharonov-Bohm effects, they may be used to construct a complete set of elementary gates for quantum computing, provided the group $`G`$ is “sufficiently” non-Abelian. This has been recently analyzed in great detail by Mochon, first restricting the discussion to states involving only fluxons Mochon1 , and then by extending it to the case where pseudocharges are also involved Mochon2 . This latter construction is specially interesting for our purpose, since it shows that the small $`D_3`$ group (whose practical implementations are discussed in Sec. V) is already capable of universal quantum computation.

### IV.2 Kinematical Gauss law model

As explained at the end of Sec. II, a different model may be constructed by assuming a modified Gauss law in which the constraint is diagonal with respect to the variables $`g_{ij}`$. It reads

$$\underset{j(i)}{}g_{ij}=e\text{,}$$
(35)
and can be directly compared to Eq. (5). This constraint allows for a simple interpretation which is motivated by the physical systems discussed in Secs. III and V. First, we notice that a constraint associated with the sites of $`\mathrm{\Lambda }_d`$ is also naturally attached to the plaquettes on $`\mathrm{\Lambda }`$. We may therefore also write it as $`_{plaq.}g_{ab}=e`$. Second, the constraint that the product of group elements around the plaquette is equal to $`e`$ implies that one can define group elements on the sites such that $`g_{ab}=g_ag_b^1`$. Further, these group elements are defined up to a global transformation $`g_ag_ah`$. For the Abelian systems discussed in Sec. III the group elements on the sites of $`\mathrm{\Lambda }`$ had a transparent physical meaning: they corresponded to the absolute values of the phase of the superconducting islands.

The simplest Hamiltonian preserving the constraint (35) is

$$H=\underset{a}{}\underset{g}{}r(g)\underset{b(a)}{}𝒯_{ab}(g)+\underset{ab}{}F(\widehat{g}_{ab})\text{.}$$
(36)
Note the close similarity with the Abelian model defined in (14). However, for a general non-Abelian group, we do not know how to reexpress the Hamiltonian (36) in terms of a local gauge theory on $`\mathrm{\Lambda }_d`$. Indeed, the key tool in this mapping for the $`\mathrm{𝖹𝖹}_n`$ group was the existence of the Fourier transform for functions defined on $`\mathrm{𝖹𝖹}_n`$. Here, we do not see any obvious way to generalize this trick to functions defined on a non-Abelian group. Nevertheless, if $`F=0`$, the Hamiltonian (36) does exhibit a local symmetry based on the sites of $`\mathrm{\Lambda }`$: $`g_{ab}h_ag_{ab}h_b^1`$ provided $`r(g)`$ depends only on the conjugacy class of $`g`$. In this case, the model can be solved for $`F0`$ because, in this limit, we may interpret the Hamiltonian as a sum of mutually commuting terms attached to the sites of $`\mathrm{\Lambda }`$

$$H_{sites}=\underset{a}{}\underset{g}{}r(g)𝒯_a(g)=\underset{a}{}H_a.$$
(37)
Solution for each site depends somewhat on the structure of the coefficients $`r(g)`$. For the simplest case of positive $`r(g)`$, the ground state is simply $`|𝒢_a=_g|g_a`$. On the whole lattice, the ground state is thus given by: $`|𝒢=_a|𝒢_a`$. This ground state can be also represented in terms of the original link variables:

$`|𝒢`$ $`=`$ $`\left[{\displaystyle \underset{plaq.}{}}\delta \left({\displaystyle \underset{(ab)plaq.}{}}g_{ab}e\right)\right]_{ab}{\displaystyle \underset{g}{}}|g_{ab}`$ (38)
$`=`$ $`\left[{\displaystyle \underset{i}{}}\delta \left({\displaystyle \underset{j(i)}{}}g_{ij}e\right)\right]_{ij}{\displaystyle \underset{g}{}}|g_{ij}\text{.}`$ (39)
We underline that these expressions are the generalizations of (17) and (18) for an arbitrary non-Abelian group.

This model exhibits local excitations attached to the sites of $`\mathrm{\Lambda }`$ or equivalently to the plaquettes of $`\mathrm{\Lambda }_d`$. They are analogous to the Abelian fluxons constructed in Eq. (21) where the site Hamiltonian $`H_a`$ in a system of Josephson junctions involves charging energies. The Hilbert space attached to the site $`a`$ forms the so-called regular representation of $`G`$ where the group element $`h`$ is represented by the left multiplication operator $`𝒯_a(h)`$. As already pointed out, these operators commute with $`H_a`$ so that the eigenstates of $`H_a`$ belong to stable subspaces for the regular representation of $`G`$. For example, if we choose $`G=D_3`$, its regular representation is six dimensional and decomposes into two different one-dimensional representations and two isomorphic copies of the two-dimensional representation of $`D_3`$. At this stage, we should notice that $`H_a`$ also commutes with the representation of $`G`$ by right multiplications (the element $`h`$ being represented by $`(h)`$ so that $`(h)|g=|gh^1)`$. Because of this extra symmetry, the spectrum of $`H_a`$ is composed of two singlets and one quadruplet.

In fact, it is possible to generate the whole set of eigenstates of $`H`$ when $`F=0`$. Let us choose a reference state $`|\mathrm{\Omega }`$ in which the constraint $`_{plaq.}g_{ab}=e`$ is satisfied for all plaquettes. Now, for each site, we consider an eigenfunction $`\varphi _a(g)`$ of $`H_a`$ with energy $`E_a`$. More precisely, we require

$$E_a\varphi _a(g)=\underset{h}{}r(h)\varphi _a(h^1g)\text{.}$$
(40)
We can then construct an eigenstate $`|\mathrm{\Psi }`$ of $`H`$ with the total energy $`E=_aE_a`$ as follows:

$$|\mathrm{\Psi }=\underset{a}{}\underset{g}{}\varphi _a(g)\underset{b(a)}{}𝒯_{ab}(g)|\mathrm{\Omega }\text{.}$$
(41)
In the case where $`\varphi _a(g)`$ is independent of $`g`$, we recover the ground state $`|𝒢`$ given by Eq. (38). Introducing the gauge transformation operators $`U_{\{g\}}`$ as before (where the $`g`$’s now live on the sites of $`\mathrm{\Lambda }`$), we may also write

$$|\mathrm{\Psi }=\underset{a}{}\underset{\{g\}}{}\varphi _a(g_a)U_{\{g\}}|\mathrm{\Omega }\text{.}$$
(42)
This shows that using a gauge transformed reference state $`U_{\{h\}}|\mathrm{\Omega }`$ instead of $`|\mathrm{\Omega }`$ is equivalent to replace the site wave functions $`\varphi _a(g)`$ by $`\varphi _a(gh_a^1)`$. As already discussed, acting on the $`\varphi _a`$’s by these right multiplications produces eigenstates with the same energy. Only sites $`a`$ where $`\varphi _a`$ belongs to a multidimensional representation of $`G`$ generates new states when a gauge transformation with $`h_ae`$ acts on $`|\mathrm{\Omega }`$.

This language is specially useful to describe a different type of excitation, analogous to the pseudocharge in the gauge-invariant model. They are produced in a physical system in which the constraint is modified on a given plaquette so that $`_{plaq.}g_{ab}=u`$ on this plaquette. Note that as it stands, this constraint is not preserved by the action of the Hamiltonian (36), but the conjugacy class of $`_{plaq.}g_{ab}`$ is preserved for any plaquette. This modified constraint may be incorporated in the general construction presented above provided we change the reference state. For instance, we may consider a string $`\gamma `$ on $`\mathrm{\Lambda }_d`$ which runs from the center of the excited plaquette to the boundary of the system. Our previous state $`|\mathrm{\Omega }`$ is then replaced by $`_{(ab)\gamma }𝒯_{ab}(u)|\mathrm{\Omega }`$. Note that in the context of physical implementations, as already mentioned in Sec. III, it is more appropriate to call these excitations vortices.

As in the gauge-invariant model, it is interesting to study the appearance of nondiagonal adiabatic transport of a degenerate fluxon excitation around a pseudocharge. As before, let us describe the situation according to both viewpoints where either one excitation is fixed or the other. For a mobile pseudocharge, we have to transform adiabatically the reference state. Starting from

$$|\mathrm{\Psi }_{in}=\underset{a}{}\underset{\{g\}}{}\varphi _a(g_a)U_{\{g\}}\left(\underset{(ab)\gamma }{}𝒯_{ab}(u)|\mathrm{\Omega }\right)\text{,}$$
(43)
we obtain after one closed path $`\beta `$ of the pseudocharge around the fluxon:

$$|\mathrm{\Psi }_{fin}=\underset{a}{}\underset{\{g\}}{}\varphi _a(g_a)U_{\{g\}}\left(\underset{(ab)\gamma +\beta }{}𝒯_{ab}(u)|\mathrm{\Omega }\right)\text{.}$$
(44)
The projection operator $`_a_{\{g\}}\varphi _a(g_a)U_{\{g\}}`$ ensures that the resulting state has exactly the same energy as the initial one, so that the motion of pseudocharge does not excite fluxons. To make a closer analogy with the gauge-invariant model of Sec. IV.1 we note that this wave function can be also written as

$$|\mathrm{\Psi }_{fin}=\left(\underset{(ab)\beta }{}_{ab}\right)|\mathrm{\Psi }_{in}\text{,}$$
(45)
where operators $`_{ab}`$’s that move the pseudocharge are defined on the links of $`\mathrm{\Lambda }`$ exactly as the $`_{ij}`$’s were defined on the links of $`\mathrm{\Lambda }_d`$ in Sec. IV.1, these operators have the same effect as $`𝒯_{ab}(u)`$ on the flat ($`|\mathrm{\Omega }`$) state and they commute with $`_{b(a)}𝒯_{ab}(g)`$ for any site $`a`$ and any group element $`g`$. Using the same reasonings as for the gauge-invariant model, we see that $`_{(ab)\beta }_{ab}`$ acts in the same way as a gauge transformation $`U_{\{h\}}`$ on the state $`|\mathrm{\Psi }_{in}`$. Since $`h`$ has to be nontrivial (more precisely conjugated to $`u`$) at the site where the fluxon is located, $`|\mathrm{\Psi }_{fin}`$ is degenerate with $`|\mathrm{\Psi }_{in}`$, but these two states may be linearly independent if the fluxon is created in a multidimensional representation of $`G`$.

Alternatively, consider a motion of the fluxon which is due to a small hopping term, $`\widehat{t}(g_{ab})`$ \[arising, e.g. from $`F`$ in Eq. (36)\]. Let us denote by $`\varphi _f(g)`$, the single-site wave function associated with the fluxon. In the absence of a pseudocharge, the fluxon state at site $`a`$ corresponds to

$$|\mathrm{\Psi }_a=\underset{\{g\}}{}\varphi _f(g_a)U_{\{g\}}|\mathrm{\Omega }\text{.}$$
(46)
Acting on this state with the operator $`\widehat{t}(g_{ab})`$ produces the state $`|\mathrm{\Psi }_b`$ provided $`\varphi _f(g_a)\widehat{t}(g_ag_b^1)=\varphi _f(g_b)`$. foot1 . In the presence of a pseudocharge, $`|\mathrm{\Omega }`$ is replaced by $`_{(ab)\gamma }𝒯_{ab}(u)|\mathrm{\Omega }`$. If the link $`(ab)`$ is on the string $`\gamma `$, one has:

$$\widehat{t}(g_{ab})|\mathrm{\Psi }_a=\underset{\{g\}}{}\varphi _f(g_bu_{ab}^1)U_{\{g\}}\left(\underset{(ab)\gamma }{}𝒯_{ab}(u)|\mathrm{\Omega }\right)\text{.}$$
(47)
So, as the fluxon hops from $`a`$ to $`b`$, its internal wave function is modified by a right multiplication associated with an element conjugated to $`u`$ and there is a complete agreement between both descriptions.

Note the close similarity between these results and the ones for the gauge theory model, provided we exchange the roles of pseudocharges and fluxons. For instance, a pseudocharge in the kinematical Gauss law model experiences a nondiagonal Aharonov-Bohm effect when another pseudocharge is moved around it, even in the absence of fluxons. Notice that although the pseudocharges in the Gauge theory are usually associated with independent degrees of freedom, this distinction (with fluxons in kinematical model) is mostly formal. It completely disappears in the alternative (albeit less popular) gauge theory formalism in which one keeps only the dynamical variables $`g_{ij}`$ on the links and creates pseudocharges by replacing the projector $`\stackrel{~}{𝒫}_o`$ in Eq. (31) by $`_h\chi (h)_{j(o)}𝒯_{oj}(h)`$ where $`\chi (h)`$ is the character of the representation of $`G`$ associated with the pseudocharge at charge $`o`$. This approach emphasizes the analogy between the two situations and is more appropriate to describe the implementations discussed in the next section.

In the Abelian case, the system with nontrivial topology, i.e., with holes, acquires a degenerate ground state, since a fluxon located inside the holes has a vanishing energy. Such a fluxon is created by the operator $`_{(ab)\gamma }𝒯_{ab}^{(p)}(1)`$ where the contour $`\gamma `$ in $`\mathrm{\Lambda }`$ connects the inner and the outer boundaries. This operator is also equal to $`_{(ab)\gamma }\mathrm{exp}(i2\pi \widehat{u}_{ab}/n)`$ and is conserved by the Hamiltonian (14). Its eigenstates correspond to linear combinations of various degenerate ground states with fluxons present or absent in the holes. This applies directly to the non-Abelian model of the present section. The degenerate ground states are therefore identified by the value of $`_{(ab)\gamma }g_{ab}`$. In the site representation, it is just the group element $`g_{in}g_{out}^1`$ where $`g_{in}`$, $`g_{out}`$ represent the states of the inner and outer boundaries respectively. These arguments hold provided that each boundary is described by a single group element $`g`$ i.e., that all sites belonging to it are connected by nonfluctuating bonds. In the Abelian case discussed in Sec. III the boundary (denoted by the dashed line in Fig. 2) is just a superconducting wire.

The same arguments as the ones developed in the Abelian case Ioffe2002  show that only nonlocal operators in $`g_{ab}`$ result in the decoherence of the topologically degenerate states. This shows that these degeneracies are protected from decoherence provided that no physical operator couples to the state $`g_a`$ of an individual site.

## V Implementations

As explained in Sec. III the array that implements directly (in a natural basis) the gauge-invariance constraint for the Abelian $`\mathrm{𝖹𝖹}_n`$ group is physically an insulator while the superconducting array becomes equivalent to the same gauge theory only after a dual transformation. Alternatively, the superconducting array can be thought of as implementing directly the modified Gauss law discussed in Sec IV B. In the non-Abelian case these two theories are no more dual to each other, so we discuss these cases separately.

The common element of all arrays discussed here is the Josephson-junction circuit that is characterized by the phase difference $`\mathrm{\Delta }\varphi `$ which can take (in a classical limit) $`n2`$ discrete values $`\mathrm{\Delta }\varphi =2\pi m/n`$; quantum processes leading to transitions between these states. The actual construction of this element for $`n=3`$ (that we shall call masu foot2  below) is shown in Fig. 10. It is built from $`n`$ more basic elements consisting of two Josephson junctions connected in series. The Josephson energy of this element (as a function of the phase difference $`\delta \varphi `$ between the ends) contains all higher harmonics: $`V_0(\delta \varphi )=2E_J|\mathrm{cos}\delta \varphi /2|`$. Now, connecting these basic elements in parallel as shown in Fig. 10 and applying the magnetic field that shifts the relative phases of neighboring circuits by $`2\pi /n`$ gives the desired $`2\pi /n`$ periodic potential. The nonzero charging energy of the junctions leads to transitions between these classical states whose amplitudes are calculated in Appendix A.

### V.1 Gauge-invariance implementation

Here we limit ourselves to the gauge theory with the simplest non-Abelian group $`D_3`$ which is also the permutation group of three elements $`S_3`$. The array implementing this model is shown in Fig. 7. We begin by explaining the qualitative reasoning that leads to this complicated construction. In the $`S_3`$ gauge theory the state of each link has to be described by the permutation of three elements. To implement this in a Josephson-junction array, we construct each bond from two elements (masus) described by discrete phases differences $`\mathrm{\Delta }\varphi _\alpha =2\pi n/3`$ ($`\alpha =1,2`$). Such a bond can be viewed as a discrete transformation from the state of one site characterized by two discrete $`2\pi n/3`$ phases to the state of another site; it is described by one of the nine elements of the $`\mathrm{𝖹𝖹}_3\times \mathrm{𝖹𝖹}_3`$ group. In order to eliminate three extra states and to get the permutation group, we need to introduce a constraint requiring that two phases of each site be different. This is achieved if each site of the lattice contains an element that excludes the configurations with equal phases. The simplest of such an element is a triangle made of three Josephson junctions frustrated by a magnetic flux $`\mathrm{\Phi }_0/2`$. The ground state of this triangle corresponds to the phases $`\varphi _\alpha =2\pi n_\alpha /3`$ with different $`n_\alpha `$. Thus, the state of two masus attached on both sides to the frustrated triangles \[as shown in Fig. 7(b)\] is in a one to one correspondence with the elements of $`S_3`$ group.

The next step in the construction of the array is to ensure that the low-energy states satisfy the gauge-invariance constraint. This is true if the Hamiltonian is dominated by a large “kinetic” term which exchanges two discrete phases on each triangle because such a process changes the states of all three bonds connected to the triangle, further, it is equivalent to the action of a transposition on the states of these three bonds

$$H_T=r_0\underset{a}{}\underset{k=0}{\overset{2}{}}\underset{b(a)}{}𝒯_{ab}(\tau \zeta ^k)\text{,}$$
(48)
where $`r_0`$ is the tunneling amplitude of this process:

$$r_0E_J^{3/4}E_C^{1/4}\mathrm{exp}\left(S_0\sqrt{\frac{E_J}{E_C}}\right)\text{,}$$
(49)
where $`S_0=6.10`$ if junctions in triangles have the same Josephson and charging energies as the ones in masus (see Appendix B). Note that this numerical value is sensitive to the number of nearest neighbors of the effective lattice and is minimized for the hexagonal lattice. As will be clear from the following, in order to be described by the gauge theory, $`r_0`$ should be sufficiently large which is a case easier to achieve for the hexagonal lattice.

The wave functions minimizing the kinetic energy (48) are the ones that satisfy the gauge-invariance constraint

$$\underset{b(a)}{}𝒯_{ab}(\tau \zeta ^k)|\mathrm{\Psi }=|\mathrm{\Psi }\text{,}$$
(50)
since local transpositions generate the whole $`S_3`$ gauge group. The array is faithfully described by the gauge theory if $`r_0`$ sets the largest energy scale in this problem. This is not satisfied if the bonds of the array are just the ones shown in Fig. 7(b), because in this case the product of permutations around the elementary hexagon is unity and it is the energy that is needed to violate this constraint that sets the largest energy scale in the array.

In order to decrease the energy scale associated with the $`_{hex.}g_{ab}1`$ configurations, we need to introduce weak links (shown as loops containing two Josephson junctions in Fig. 7) characterized by small Josephson and charging energies $`ϵ_Jϵ_CE_CE_J`$. In order to avoid a mixture of Josephson junctions with different parameters $`(E_J,E_C)`$ in the same array, we assume that these weak links are constructed from the two usual Josephson junctions with a large Josephson energy $`\stackrel{~}{E}_JE_J`$ and correspondingly small charging energy $`\stackrel{~}{E}_CE_C`$ connected in parallel to form a loop frustrated by the magnetic flux. The effective Josephson energy of such a loop is $`ϵ_J=\stackrel{~}{E}_J(\mathrm{\Phi }/\mathrm{\Phi }_01/2)`$ and its charging energy is $`ϵ_C=\stackrel{~}{E}_C/2`$. The effect of these small energies on the tunneling process that corresponds to transposition is weak and can be ignored. The presence of weak links decreases the energy gap between the six states of the bond associated with the permutation group and the remaining three states that were eliminated from the low-energy spectrum. However, as we show below, we can choose the energies of the weak links, $`ϵ_J`$, $`ϵ_C`$, so that this gap remains large, but the amplitude of the tunneling process, $`r_0`$, remains larger than other scales in the low-energy theory so that the leading term in the effective Hamiltonian is indeed given by Eq. (48).

The Hilbert space satisfying this constraint is still huge: in the absence of other contributions to the Hamiltonian (generated by weak links) all states satisfying the gauge-invariance constraint are degenerate. This degeneracy is lifted when these subdominant terms are taken into account. Weak links give two types of terms: they lead to the dynamics of the discrete variables describing the state of each link and to diagonal potential-energy terms. Let us start with the latter. The main contribution comes in the second order in $`ϵ_J/ϵ_C`$, it is given by

$$U(\{u_i\})=\frac{ϵ_J^2}{2ϵ_C}\mathrm{cos}\left(\frac{2\pi }{3}u_i\right)\text{,}$$
(51)
where $`u_i`$ are discrete variables describing the phase differences on elementary masus and along the sides of the triangles. This energy is minimized if this sum is $`0`$ modulo $`3`$, i.e., if the state of the links is described by the permutation group “acting” on the state of triangles. We assume that the scale of this interaction is relatively large and comparable to $`r_0`$, i.e., that $`\frac{ϵ_J^2}{2ϵ_C}r_0`$. This interaction does not lift the degeneracy of the Hilbert space of the discrete variables that are described by elements of the permutation group. This degeneracy is only lifted by weaker processes which involve longer loops of weak links around each hexagon. Such a loop contains six weak links, it results in the contribution to the potential energy similar to the first term in Eq. (26):

$$V=\underset{abcdef}{}𝒱(g_{ab}g_{bc}g_{cd}g_{de}g_{ef}g_{fa})\text{.}$$
(52)
To calculate the value of the interaction $`𝒱(g)`$, we note that physically this energy is due to the energy of the Cooper pair that starts from one corner of the triangle and moves around the hexagon across the weak links. Furthermore, if $`g`$ is a transposition, the sum of the discrete phases along contours beginning at two corners of the triangle is equal to $`\pm 2\pi /3`$ and one such a sum is $`0`$, while if $`g`$ is a cyclic permutation, all these sums are equal to $`\pm 2\pi /3`$, so $`𝒱(e)=0`$, $`𝒱(\tau \zeta ^k)=2v_0/3`$, $`𝒱(\zeta ^k)=v_0`$, where coefficient $`v_0`$ is the strength of the effective Josephson contact formed by two parallel loops of six weak junctions with a phase difference $`2\pi /3`$ across. To calculate the value of $`v_0`$ we note that two parallel loops are equivalent to one loop with effective Josephson and charging energies $`ϵ_J^{eff}=2ϵ_J`$, $`ϵ_C^{eff}=ϵ_C/2`$ and use the result of the Appendix C for the energy of one loop

$$v_0=\frac{1}{20}\left(\frac{3ϵ_J}{ϵ_C}\right)^5ϵ_J.$$
(53)
We now turn to the dynamical processes which involve phase slips of weak links. These processes change the state of each link individually, similar to the second term in Eq. (26). The potential energy involved in such a process is dominated by the energy of each masu while the charging energy is due to the weak junction (because its charging energy is smaller, it dominates the dynamic term in the action). Thus, the amplitude of such a process is

$$\stackrel{~}{F}(g)E_J^{3/4}ϵ_C^{1/4}\mathrm{exp}(\stackrel{~}{S}),\stackrel{~}{S}=2.44\sqrt{\frac{E_J}{ϵ_C}}\text{.}$$
(54)
For small $`ϵ_C`$ this amplitude becomes exponentially small and so these terms can be neglected in the effective Hamiltonian (26).

Thus, we conclude that the array shown in Fig. 7 provides a realization of the non-Abelian gauge theory discussed in Sec. IV, if $`ϵ_Jϵ_CE_CE_J`$ and

$$\frac{ϵ_J^2}{2ϵ_C}r_0\frac{1}{20}\left(\frac{3ϵ_J}{4ϵ_C}\right)^5ϵ_J.$$
(55)
### V.2 Kinematical Gauss law implementation

As we argue below, the straightforward implementation of the kinematical non-Abelian Gauss law is not possible because it does not protect the system from the coupling to the environment. In this respect, the non-Abelian theories are crucially different from the Abelian ones. In order to explain this difference, we first recollect the physical arguments that show that Abelian arrays (similar to the ones suggested in Ref. Ioffe2002 ) are optimally protected from the environment and then show why non-Abelian ones are different.

For Abelian theories, the Gauss law (4) and (35) is even easier to implement in the physical system than a gauge theory. As discussed in Sec. II, in order to implement an Abelian group $`\mathrm{𝖹𝖹}_n`$ on the triangular lattice shown in Fig. 1 one might assign elements (masus discussed in the preceding subsection) that are $`2\pi /n`$ periodic in superconducting phase to the links of the dual hexagonal lattice. In this system the ’site’ representation acquires a transparent physical meaning – the state of each island is described by the superconducting phase of each site of the hexagonal lattice. Different topological sectors of the system with hole(s) correspond to different phases of the inner boundary. The existence of the protected subspace is related to the $`U(1)`$ physical invariance: superconducting phase of each island (and of the inner boundary) has no meaning by itself, so only the difference of the phases of the neighboring sites might be coupled to physical operators.

The operator conjugated to the superconducting phase is the physical charge of the whole island and this presents some danger since it can be coupled to the external electric fields, especially if it results in a direct coupling of the state of the inner boundary to the electrical fields. Such a coupling, however, becomes exponentially small in the system with a long range order in superconducting phase, more precisely in the systems where $`\mathrm{exp}(in\varphi )`$ is close to $`1`$. In these systems the total charge of the boundary fluctuates strongly, while the conjugate operators appearing in the discrete gauge theory correspond to the total charge modulo $`n`$ and thus are decoupled from the physical electric fields. We compute their weak residual effect on the degeneracy of the protected states in Appendix D and find that it quickly becomes exponentially small \[see Eq. (93)\] with the number of links $`Z`$ connecting the boundary to the inside of the array, i.e., in the thermodynamic limit.

As mentioned in the preceding section, the most straightforward way to generalize this $`\mathrm{𝖹𝖹}_n`$ construction to the simplest non-Abelian group $`S_3`$ would be to construct the bonds characterized by the elements of the permutation group and the sites that can be viewed as results of these permutations. This is achieved in the array shown in Fig. 8 where bonds are constructed from pairs of masus and the sites from triangles frustrated by the magnetic field. This array differs from the one shown in Fig. 7 by the absence of the weak links. All classical states of this array are uniquely specified by the configurations of permutation group elements which satisfy the constraint $`_{hex.}g_{ab}=e`$. A small but nonzero charging energy of each contact would lead to the quantum transitions between different classical states as described by the Hamiltonian (36). In this physical system, these transitions correspond to the tunneling of the phases on each triangle from one configuration to another.

Thus, in ideal conditions, this system is described by the kinematic Gauss law model of preceding section. The situation is not so good for a realistic system because in such implementation, one can not exclude the appearance of physical operators coupled directly to the state of each triangle. In particular, this could distinguish between states with different chiralities (i.e., in which the phase winds clockwise or counterclockwise around the triangle). This is bad news for the effective gauge field model since the phase differences on distant triangles are related by

$$\mathrm{\Delta }\varphi _1=\mathrm{\Delta }\varphi _2+\frac{2\pi }{3}\text{sign}\left(\underset{C}{}g_{ab}\right),$$
(56)
where sign denotes the sign of the permutation and indices $`1`$ and $`2`$ denote two distant triangles. Thus, the terms in the energy that depend on the phase difference $`\mathrm{\Delta }\varphi `$ results in nonlocal terms in the effective Hamiltonian when expressed as a function of the gauge fields $`g_{ab}`$. Such terms ruin the topological protection, especially if these triangles belong to inner and outer borders. To be more specific we shall discuss below the effect of

$$V_{dis.}=V_d\mathrm{cos}(\mathrm{\Delta }\varphi _1\mathrm{\Delta }\varphi _2)\text{,}$$
(57)
which is actually the most dangerous term of this type since it changes the relative energies of different topological sectors of the array shown in Fig. 8 by an amount proportional to $`V_d`$.

The effect of $`V_{dis.}`$ on the array with weak links (discussed in Sec. V.1) is exponentially weak because in this system the relation between phases also contains a contribution from the continuous phases across weak junctions which, in the absence of the superconducting long range order, fluctuate strongly and eliminate the correlation between $`\mathrm{\Delta }\varphi _1\mathrm{\Delta }\varphi _2`$ and $`_Cg_{ab}`$. Indeed, in this case, the relation between $`\mathrm{\Delta }\varphi _1`$ and $`\mathrm{\Delta }\varphi _2`$ acquires an additional term which is due the continuous part of the phase $`\stackrel{~}{\varphi }`$ across the weak links.

$$\mathrm{\Delta }\varphi _1\mathrm{\Delta }\varphi _2=\frac{2\pi }{3}\text{sign}(\underset{C}{}g_{ab})+\underset{C}{}\delta \stackrel{~}{\varphi }\text{,}$$
(58)
where $`\delta \stackrel{~}{\varphi }`$ is the difference of continuous phases across two weak links constituting one bond of the effective lattice. In the insulating limit of this model ($`ϵ_Cϵ_J`$) the continuous phases fluctuate practically independently on each bond with $`\mathrm{cos}(\delta \stackrel{~}{\varphi })(ϵ_J/ϵ_C)`$, suppressing the effect of the noise (57)

$$V_{dis.}^{eff}V_d\left(\frac{ϵ_J}{ϵ_C}\right)^L\mathrm{cos}\left[\frac{2\pi }{3}\text{sign}\left(\underset{C}{}g_{ab}\right)\right]\text{,}$$
(59)
where $`L`$ is the length of the contour connecting triangles $`1`$ and $`2`$.

Alternatively, one can argue that the mapping of the physical system onto the effective model with variables defined on the bonds of the lattice does not produce nonlocal (dangerous) operators if the states of the site are indistinguishable. In the insulating array, this is achieved by making the phase of each triangle (effective site) fluctuate strongly. In the Abelian $`\mathrm{𝖹𝖹}_n`$ array, this is true because absolute values of the phase (which distinguish between different states of a given site) are not measurable, so generally a system in which each site is characterized only by a single phase is a good candidate for an array that can be mapped to a bond model. The problem is that the straightforward implementations of this idea lead to a $`\mathrm{𝖹𝖹}_n`$ group.

In order to add a reflection to the $`\mathrm{𝖹𝖹}_n`$ group, one needs states that are characterized by both the discrete value of the phase and the sign of the conjugate variable, $`q`$. This is possible to achieve if one starts for the system where $`n\times m`$ discrete phases are allowed, and include a large interaction that forms a degenerate doublet of chiral states from the subset of $`m`$ discrete states. To be more specific, consider the low-energy states that are described by the wave functions

$$\mathrm{\Psi }_\pm (l)=\underset{k=0}{\overset{k<m}{}}\mathrm{exp}(\pm 2\pi ik/m)|l,k,$$
(60)
for $`l<n`$. Because the reflection changes the chirality of each state, it is not reducible to a simple rotation, so that combined with the phase rotations, it gives the fundamental representation of $`D_n`$ group. We now discuss the construction of the simplest array that has these general features. The minimal value of $`m`$ that allows one to distinguish left and right states is $`m=3`$. The minimal value of $`n`$ that allows a non-Abelian group is $`n=3`$, however, $`n=m`$ is very difficult to implement because it can not be achieved connecting in series the elements with $`n`$ and $`m`$ periodicity (that would give $`n=m`$ periodicity instead of $`n^2`$ one). Thus, the simplest array of this type corresponds to $`n=4`$ and group $`D_4`$. A possible realization is shown in Fig. 9. Here each link consists of two elements which have rather different values of $`E_J/E_C`$: the masu that is $`2\pi /4`$ periodic operates in the phase regime with $`E_{J4}/E_{C4}>1`$ while the masu that is $`2\pi /3`$ periodic has a different relation $`E_{J3}/E_{C3}1`$. More precisely, the amplitude mixing states that differ by $`2\pi /4`$ is

$$rE_{J4}^{3/4}E_{C4}^{1/4}\mathrm{exp}(3S_4),$$
(61)
where $`S_4=\sqrt{E_{J4}/E_{C4}}`$ must be much smaller than the energy gap

$$\stackrel{~}{\mathrm{\Delta }}E_{J3}^{3/4}E_{C3}^{1/4}\mathrm{exp}(3S_3)$$
(62)
between the doublet of the chiral states and the remaining $`p=0`$ state of the $`m=3`$ set ($`S_3=\sqrt{E_{J3}/E_{C3}}`$).

It is crucial for this design that the states with $`p=\pm 1`$ \[here $`p`$ has the same meaning as in Sec. III, e.g., in Eq. (10)\] have equal energies which are less than the energy of the state with $`p=0`$. This happens if the transition amplitude between the classical minima is real and negative. Generally, the amplitude of the transition that increments the superconducting phase by $`2\pi /3`$ contains the phase factor $`2\pi Q/3`$ where $`Q`$ is the average charge of the superconductor whose phase has been changed (here and below we measure the charge in the units of Cooper pairs). We require that this phase is $`(\pi +2\pi N)`$ where $`N`$ is arbitrary integer. This implies that the system should be biased so that the average charge of the island that experiences this transition is $`Q=3/2+3N`$.

We now discuss the effect of the induced charge on the amplitude of the transition between $`l=0,\mathrm{},3`$ discrete states. For simplicity, we shall assume that the islands which experience these transitions coincide with the islands whose phases change by $`2\pi /3`$ and, thus, have the same charge. These transitions increment phase by $`\pi /2`$; so that their amplitude acquires a phase $`(3/2+3N)\pi /2`$; in the notations of Sec. IV.2, $`r_\zeta =\mathrm{exp}(iQ\pi /2)|r|`$. This amplitude has positive real part for $`Q=9/2,15/2,\mathrm{}`$. The transition between states with $`p=\pm 1`$ appears naturally if the energies of the classical minima with different $`k=0,1,2`$ are not exactly equal because, in the momentum representation, this induces a transition between these states. Such a difference is induced by the deviations, $`\delta \mathrm{\Phi }_3=\mathrm{\Phi }_3\mathrm{\Phi }_0/3`$ of the actual flux through the elementary rhombus from its ideal value, $`\mathrm{\Phi }_0/3`$. Thus, in these conditions the low-energy dynamics of such a system is described by exactly (37) with the parameters $`r_\zeta `$, and $`r_\tau =E_{J3}\delta \mathrm{\Phi }_3/\mathrm{\Phi }_0`$. We emphasize again that the low-energy sector has to be separated from higher-energy excitations by a large gap $`\stackrel{~}{\mathrm{\Delta }}r_\zeta ,r_\tau `$.

## VI Conclusion

In this work, we have proposed several possible implementations of lattice models based on discrete non-Abelian local symmetry groups. The first key ingredient in all these constructions is the existence of a large but discrete manifold of equivalent “classical” states. This manifold appears if the system is built from heavily frustrated elementary blocks, such as Josephson-junction loops threaded by a well adjusted magnetic flux. Further, one needs to select an interaction between these blocks that imposes a set of local constraints in the low-energy sector. We have found two inequivalent (but somewhat similar) types of constraints compatible with the non-Abelian symmetry group. The first type of constraint corresponds to the dynamical selection of the states which are invariant under the gauge group transformations. This is most naturally implemented in insulating arrays where it can be made equivalent to the selection of states that minimize the kinetic energy associated with the local tunneling processes between classically equivalent states. The second type of constraint emerges more naturally in superconducting Josephson-junction arrays (or the ones with small phase fluctuations), where group elements associated with the bonds can be viewed as transformations between the states of the sites of the array. Because the state of each site is uniquely defined, the product of group elements around a closed loop is the unit element. We have called this condition a kinematical Gauss law. This constraint becomes equivalent to the local gauge-invariance condition in Abelian theories where a duality transformation maps one onto the other (Sec. III) but for a non-Abelian symmetry it leads to a different theory in which the roles of charges and fluxons have been reversed (Sec. IV).

In both types of models that we considered, the elementary excitations within the constrained subspace are local in space and have a finite-energy gap. We call these excitations “fluxons”. Generally, these fluxons can be trapped inside a larger hole prepared in the array leading to a degenerate ground state in the thermodynamic limit. In a finite system, this degeneracy is lifted by an exponentially weak tunneling escape process. It is also robust with respect to all noises represented by local operators such as local electric fields (generated by stray charges) or stray magnetic fluxes: the matrix elements of all these operators are exponentially suppressed. This remarkable property may open the way to the interesting applications such as storage of quantum information for which a low decoherence rate is a must. We emphasize that the physical nature of these fluxons depends, however, on the details of implementation and the type of the constraint: the canonical gauge invariance is more easily implemented in insulating arrays in which fluxons have the physical meaning of fractional vortices while kinematical Gauss law is more naturally implemented in superconducting arrays in which fluxons carry physical charges.

Besides fluxons, in any physical array, it is possible to build excitations that violate a local constraint (of course, at some energy cost); in the gauge theory language these excitations are pseudocharges. In the canonical gauge models pseudocharges may provide a multidimensional representation of the discrete gauge group while fluxons are labeled by their conjugacy classes. In the kinematical Gauss law models the role of fluxons and pseudocharges is reversed. The existence of these higher-dimensional representations leads to a very interesting statistical interaction between fluxons and pseudocharges which may be utilized for making exact transformations needed in quantum computations.

Many interesting questions remain to be addressed, we shall list just a few here. First, it is not clear whether it is absolutely unavoidable to emulate gauge theory in order to get the topological protection and non-Abelian statistics. What is the most general class of theories that has these properties? Second, is it possible that some simpler overfrustrated arrays “spontaneously” develop the gauge invariance in the low-energy sector? foot3 . Such an array would be more experimentally accessible and easier to control. Coming back to the proposed arrays it would be very interesting to investigate the experimental manifestations of the nonlocal nature of the fluxon quantum numbers which become apparent only when a few of them are present. As discussed extensively in the mathematical literature Dijkgraaf1990 ; Bais1992 , these quantum numbers interpolate between a set of conjugacy classes (one for each fluxon) and a set of group elements (which is not a gauge-invariant object by itself). It would be important to realize how these subtle but beautiful concepts can be revealed in a physical experiment, and to design some setup that could verify the non-Abelian statistics of the excitations in these models. Finally, it would be very interesting to adapt the construction of the Chern-Simons term (which binds fluxons and charges) known for a continuous groups to the discrete finite gauge groups (even Abelian ones) discussed here and to provide a corresponding physical implementation.

###### Acknowledgements.

L.I. is grateful to ENS and LPTHE for their hospitality which allowed this work to be commenced and finished. This work was partially supported by NSF Grant No. DMR 0210575.

## Appendix A Tunneling in a single rhombus

Here, we calculate the transition amplitude between the degenerate classical states appearing in a single rhombus with diagonal that is placed in the magnetic field with the flux $`2\mathrm{\Phi }_0/3`$ through the whole rhombus. In the following we shall refer to this object as masu. As a function of the phase difference $`\psi `$ across the diagonal the total Josephson energy of such masu is:

$`E(\psi )`$ $`=`$ $`2E_J\left[\right|\mathrm{cos}(\psi /2)|+|\mathrm{cos}(\psi /2+\pi /3)|`$ (63)
$`+\left|\mathrm{cos}(\psi /2\pi /3)\right|].`$
It is periodic in $`\psi `$ with the period $`2\pi /3`$ and has three classically degenerate minima: $`\psi =0,\pm 2\pi /3`$.

The amplitude of transitions $`r_m`$ between these minima can be computed in a quasiclassical approximation, as $`r_mE_J^{3/4}E_C^{1/4}\mathrm{exp}(S_0)`$, and is equal for all transitions. To calculate the action $`S_0`$, we choose two independent superconducting phases as shown in Fig. 10 and rescale $`t`$ into $`t/\sqrt{E_JE_C}`$. In terms of these variables the action is

$`S`$ $`=`$ $`{\displaystyle }[{\displaystyle \frac{3}{8}}\left({\displaystyle \frac{d\theta }{dt}}\right)^2+{\displaystyle \frac{1}{8}}\left({\displaystyle \frac{d\varphi }{dt}}\right)^22\sqrt{3}\mathrm{cos}(\theta )`$ (64)
$`2\mathrm{cos}(\varphi )\mathrm{sin}(\theta )]dt\text{.}`$
The minimum of this action corresponds to the quasiclassical equations of motion:

$`{\displaystyle \frac{d^2\theta }{dt^2}}`$ $`=`$ $`{\displaystyle \frac{8}{\sqrt{3}}}\mathrm{sin}(\theta ){\displaystyle \frac{8}{3}}\mathrm{cos}(\theta )\mathrm{cos}(\varphi ),`$ (65)
$`{\displaystyle \frac{d^2\varphi }{dt^2}}`$ $`=`$ $`8\mathrm{sin}(\theta )\mathrm{sin}(\varphi ).`$ (66)
Solving these equations numerically, we get the expression for the transition amplitude in one masu

$$r_mE_J^{3/4}E_C^{1/4}\mathrm{exp}(S_0)S_0=1.22\sqrt{\frac{E_J}{E_C}}.$$
(67)
## Appendix B Tunneling of four rhombi and triangle

Here, we calculate the transition amplitude of four masus connected by a triangle as shown in Fig. 11. If the junctions constituting the triangle are weak, $`e_JE_J`$ this amplitude can be obtained in the perturbation theory as a product of two transitions that happen consecutively: one involving the left pair of masus and another the right one. The amplitude for each transition is related to the amplitude of the transition in a single masu (67) because the quasiclassical equations of motion for the variables describing two masus joined at the vertex are the same as for a single masu. Thus, this tunneling amplitude contains twice the action (67):

$$r_m^{}E_J^{3/4}E_C^{1/4}\mathrm{exp}(S_0)S_0=2.44\sqrt{\frac{E_J}{E_C}}\text{.}$$
(68)
In this limit, the amplitude of the full process is

$$r_0\frac{(r_m^{})^2}{e_J}\text{,}$$
(69)
where again $`e_J`$ denote the weak Josephson couplings inside the triangles.

In a general case of triangle junctions comparable with the junctions constituting masus, we need to add these terms to the action and solve the full set of equations of motion in order to obtain the action corresponding to this process. The solution of these equations describes a simultaneous variation of phases in all four masus, so the most general action has the form

$`S`$ $`=`$ $`4\sqrt{{\displaystyle \frac{E_J}{E_C}}}{\displaystyle }[{\displaystyle \frac{3}{8}}\left({\displaystyle \frac{d\theta }{dt}}\right)^2+{\displaystyle \frac{1}{8}}\left({\displaystyle \frac{d\varphi }{dt}}\right)^22\sqrt{3}\mathrm{cos}(\theta )`$ (70)
$`2\mathrm{cos}(\varphi )\mathrm{sin}(\theta )+{\displaystyle \frac{3E_C}{8e_C}}\left({\displaystyle \frac{d\theta }{dt}}\right)^2`$
$`+{\displaystyle \frac{e_J}{4E_J}}[2\mathrm{cos}(2\theta )+\mathrm{cos}(4\theta )]].`$
In the most natural case the junctions of the triangle have the same Josephson current per area and the same capacitance per area as the junctions of the masu, so $`e_J/E_J=E_C/e_Ce`$. Minimizing the action, we get the equations of motion

$`{\displaystyle \frac{d^2\theta }{dt^2}}`$ $`=`$ $`{\displaystyle \frac{4}{1+e}}[{\displaystyle \frac{2}{\sqrt{3}}}\mathrm{sin}(\theta ){\displaystyle \frac{2}{3}}\mathrm{cos}(\theta )\mathrm{cos}(\varphi )`$ (71)
$`{\displaystyle \frac{e}{3}}[\mathrm{sin}(4\theta )\mathrm{sin}(2\theta )]],`$
$`{\displaystyle \frac{d^2\varphi }{dt^2}}`$ $`=`$ $`8\mathrm{sin}(\theta )\mathrm{sin}(\varphi ).`$ (72)
Solving these equations numerically, we obtain the dimensionless action $`S(e)`$ as a function of $`e`$ shown in Fig. 12 which should be used in

$$r_0E_J^{3/4}E_C^{1/4}\mathrm{exp}\left(S(e)\sqrt{\frac{E_J}{E_C}}\right)\text{,}$$
(73)
to obtain the value of the transition amplitude. We observe that the action $`S(e)`$ varies only by 25% \[from $`S(0)=4.91`$ to $`S(1)=6.10`$\] when $`e`$ is increased from $`0`$ to $`1`$, so the amplitude of the transition is not too sensitive to the strength of Josephson junctions in the triangles.

## Appendix C Effective Josephson coupling of a chain

Here we calculate the effective coupling by a chain of $`(n+1)`$ weak Josephson contacts. Denoting the phase across each contact by $`\phi _i`$, the phase across the last contact is $`\varphi _i\phi _i`$. After a time rescaling the system is described by a dimensionless Lagrangian

$`L`$ $`=`$ $`L_0+L_J,`$ (74)
$`L_0`$ $`=`$ $`{\displaystyle \frac{1}{2}}{\displaystyle \underset{i}{}}\left({\displaystyle \frac{d\phi _i}{dt}}\right)^2+{\displaystyle \frac{1}{2}}\left({\displaystyle \underset{i}{}}{\displaystyle \frac{d\phi _i}{dt}}\right)^2,`$ (75)
$`L_J`$ $`=`$ $`\lambda \left[{\displaystyle \underset{i}{}}\mathrm{cos}\phi _i+\mathrm{cos}\left(\varphi {\displaystyle \underset{i}{}}\phi _i\right)\right],`$ (76)
where $`\lambda =ϵ_J/(8ϵ_C)`$. Let us assume that $`\lambda 1`$ and treat Josephson energy in perturbation theory. The bare Hamiltonian is quadratic in charges, $`q_i`$, conjugated to the phases $`\phi _i`$

$$H_0=\frac{1}{2}\left[\underset{i}{}q_i^2\frac{1}{n+1}\left(\underset{i}{}q_i\right)^2\right].$$
(77)
The ground state of this Hamiltonian corresponds to all charges $`q_i=0`$ while in the low-energy states some $`q_i=\pm 1`$. The energy of the state with $`k`$ charges $`q_i=1`$ is

$$ϵ_k=k(n+1k)/[2(n+1)].$$
(78)
To find the effective coupling provided by this chain we need to calculate

$$W(\varphi )=\lambda \mathrm{cos}\left(\varphi \underset{i}{}\phi _i\right),$$
(79)
in the lowest order of the perturbation theory in

$$V_{int}=\frac{\lambda }{2}\underset{i}{}\mathrm{exp}(i\phi _i)+h.c.$$
(80)
In the charge representation, one has

$$W(\varphi )=\lambda \underset{i}{}q_i^+\mathrm{cos}\varphi $$
(81)
where $`q_i^+`$ is the operator increasing the charge of junction $`i`$ by $`1`$ and

$$V_{int}=\frac{\lambda }{2}\underset{i}{}q_i^++h.c.$$
(82)
In the leading order, we can limit ourselves to the components of the ground-state wave functions that contain charges $`1`$ at some $`k`$ junctions and $`(nk)`$ opposite charges at the remaining junctions. These components have coefficients $`\mathrm{\Psi }_k`$ that can be found from a recursive relation

$$\mathrm{\Psi }_k=\frac{k\lambda }{2ϵ_k}\mathrm{\Psi }_{k1},$$
(83)
so we get

$$\mathrm{\Psi }_k=\frac{(nk)!}{n!}(n+1)^k\lambda ^k.$$
(84)
Expressing the average $`_iq_i^+`$ through these components via

$$\underset{i}{}q_i^+=\underset{k=0}{\overset{n}{}}\frac{n!}{(nk)!k!}\mathrm{\Psi }_k\mathrm{\Psi }_{nk}\text{,}$$
(85)
we get

$$\underset{i}{}q_i^+=\frac{(n+1)^{n+1}}{n!}\lambda ^n.$$
(86)
Converting back to unrescaled variables we conclude that the effective coupling between the ends of the line consisting of $`n`$ junctions is

$$W_n(\varphi )=\frac{n^n}{(n1)!}\left(\frac{ϵ_J}{8ϵ_C}\right)^{n1}ϵ_J\mathrm{cos}\varphi .$$
(87)
## Appendix D Effect of the electric noise on the protected subspace in $`\mathrm{𝖹𝖹}_n`$ arrays

Quantitatively, consider the island coupled to $`Z`$ links by $`2\pi /n`$ periodic elements. Each element is characterized by the potential energy $`V_{el}(\varphi )=_kV_0(\varphi +2\pi k/n)`$ which has minima at $`\varphi _m=2\pi (m+\frac{1}{2})/n`$. Near each minima the potential is approximately Gaussian:

$$V_{el}\frac{c}{2}E_J(\varphi \varphi _m)^2,$$
(88)
where

$$c=\frac{1}{2}\underset{k=0}{\overset{n1}{}}\left|\mathrm{cos}\left[\frac{\pi }{n}\left(k+\frac{1}{2}\right)\right]\right|,$$
(89)
is a numerical coefficient that is equal to $`1/\sqrt{2}`$ for $`n=2`$ and $`n/\pi `$ for large $`n`$. The dynamics of small phase fluctuations is described by the Lagrangian

$$=\underset{(ij)}{}\left(\frac{n}{32E_C}(\dot{\varphi _i}\dot{\varphi _j})^2+V_{el}(\varphi _i\varphi _j)\right)\text{,}$$
(90)
where $`E_C=e^2/(2C)`$ is the charging energy of each junction. In this approximation, each phase difference is an independent Gaussian variable which is described by the wave functions

$$\mathrm{\Psi }(\delta \varphi )=\mathrm{exp}(\frac{1}{2}\sqrt{\frac{ncE_J}{16E_C}}\delta \varphi ^2).$$
(91)
The canonically conjugate charge fluctuations are therefore also Gaussian with width $`\sqrt{(ncE_J)/(16E_C)}/2`$; charge fluctuations on different junctions add up to the total charge fluctuation on the island with the Gaussian width

$$(Q\overline{Q})^2=\frac{Z}{8}\sqrt{\frac{ncE_J}{E_C}}.$$
(92)
If these fluctuations are restricted to a given value modulo $`n`$, they give slightly different results. Consider, for instance, the difference between average square of the charge between the sectors with different $`Q`$ modulo $`n`$

$$\delta Q^2_k\delta Q^2_l\mathrm{exp}\left(\frac{\pi ^2Z}{4}\sqrt{\frac{cE_J}{n^3E_C}}\right).$$
(93)
This difference quickly becomes exponentially small for a large number of links connecting a given island. Because a topologically protected space can be viewed as spanned by different states of the inner boundary of the array connected to the interior by $`Z1`$ links, the electrical fields have exponentially weak effects on the topologically protected subspace.