A qubit model for U(1) lattice gauge theory

TL;DR

This paper introduces a qubit-based model for U(1) lattice gauge theory that simplifies the gauge link representation, enabling real-time simulations of gauge dynamics in two and three dimensions.

Contribution

It presents a minimal qubit model for U(1) lattice gauge theory using Gauss's law, facilitating quantum simulations of gauge dynamics with few qubits per link.

Findings

Real-time gauge dynamics demonstrated in 2D lattices.

Model extends to 3D with initial real-time dynamics.

Efficient qubit encoding with as few as two qubits per link.

Abstract

A conceptually simple model for strongly interacting compact U(1) lattice gauge theory is expressed as operators acting on qubits. The number of independent gauge links is reduced to its minimum through the use of Gauss's law. The model can be implemented with any number of qubits per gauge link, and a choice as small as two is shown to be useful. Real-time propagation and real-time collisions are observed on lattices in two spatial dimensions. The extension to three spatial dimensions is also developed, and a first look at 3-dimensional real-time dynamics is presented.

Tables

Table 1. Table 1: The correspondence between electric field eigenvalues (named ℓ ℓ \ell ) and states of the qubit register for two cases: 2 qubits and 3 qubits. When referenced in the text, qubits within a register are numbered from right to left: …210.

2 qubits		3 qubits
$\ell$	state	$\ell$	state
2	$|00⟩$	4	$|000⟩$
1	$|01⟩$	3	$|001⟩$
0	$|10⟩$	2	$|010⟩$
-1	$|11⟩$	1	$|011⟩$
		0	$|100⟩$
		-1	$|101⟩$
		-2	$|110⟩$
		-3	$|111⟩$

Equations

$$S$$

$$U_P$$

$$H=\sum_{n}\left(\frac{a^{3}}{2}\vec{E}^{2}(n)-\frac{\beta}{a}\sum_{i=2}^{3}\sum_{j=1}^{i-1}\cos\bigg{(}\theta_{i}(n)+\theta_{j}(n+\hat{i})-\theta_{i}(n+\hat{j})-\theta_{j}(n)\bigg{)}\right)$$

$$[{\theta }_j(n),L_k(n^{\prime })]=i{\delta }_{jk}{\delta }_{n{n}^{\prime }}.$$

$$\left[U_j,L_j\right]$$

$$\left[U_j^{†},L_j\right]$$

$$L_j|0⟩=0$$

$$L_j(U_j{)}^{\ell }|0⟩$$

$$L_j(U_j^{†}{)}^{\ell }|0⟩$$

$$H=\frac{1}{2a\beta}\bigg{(}\sum_{j}L_{j}^{2}-\beta^{2}\sum_{P}(U_{P}+U_{P}^{\dagger})\bigg{)}\,.$$

$$H=L_1^2+L_2^2+L_3^2+L_4^2-{\beta }^2(U_1{U}_2{U}_3^{†}{U}_4^{†}+U_1^{†}{U}_2^{†}{U}_3{U}_4).$$

$$H=4L^2-{\beta }^2(L^{+}+L^{-})$$

$$\left(\begin{array}[]{ccccccccc}64&-\beta^{2}&0&0&0&0&0&0&0\\ -\beta^{2}&36&-\beta^{2}&0&0&0&0&0&0\\ 0&-\beta^{2}&16&-\beta^{2}&0&0&0&0&0\\ 0&0&-\beta^{2}&4&-\beta^{2}&0&0&0&0\\ 0&0&0&-\beta^{2}&0&-\beta^{2}&0&0&0\\ 0&0&0&0&-\beta^{2}&4&-\beta^{2}&0&0\\ 0&0&0&0&0&-\beta^{2}&16&-\beta^{2}&0\\ 0&0&0&0&0&0&-\beta^{2}&36&-\beta^{2}\\ 0&0&0&0&0&0&0&-\beta^{2}&64\\ \end{array}\right)\,.$$

$$L^{+}$$

$$L^{-}$$

$$L$$

$$H=6+2{\sigma }_0^z+4{\sigma }_1^z+4{\sigma }_0^z{\sigma }_1^z-{\beta }^2\left({\sigma }_0^x+\frac{1}{2}({\sigma }_0^x{\sigma }_1^x+{\sigma }_0^y{\sigma }_1^y)\right)$$

$$H$$

$$e^{-i(A+B)t}=e^{-iAt/2}{e}^{-iBt}{e}^{-iAt/2}+O(t^3),$$

$$|\mathrm{init}⟩=\frac{1}{\sqrt{2}}|E_0⟩+\frac{1}{\sqrt{2}}|E_1⟩.$$

$$⟨\psi |H|\psi ⟩$$

$$P_i$$

$$P_i^x$$

$$P_i^{xx}$$

$$P_i^{yy}$$

$$|\psi ⟩=0|00⟩+\frac{1}{\sqrt{2}}\sin \theta |01⟩+\cos \theta |10⟩+\frac{1}{\sqrt{2}}\sin \theta |11⟩.$$

$$H$$

Full text

I Motivation

Numerical simulations of lattice gauge theory started about 40 years ago Creutz:1979dw . During this time, lattice gauge theory has become the central tool for first-principles studies of nonperturbative quantum chromodynamics (QCD), and therefore it plays a vital role in the field of particle physics Tanabashi:2018oca . Lattice gauge theory is also being applied vigorously to theories beyond QCD Brower:2019oor . For QCD and beyond, the standard approach has been Monte Carlo calculations in Euclidean spacetime, which do not provide access to real-time dynamics. Avoiding this Monte Carlo approach is difficult because the required Hilbert space is usually far larger than any classical computer can hope to manage Preskill:2018fag .

The availability of a practical quantum computer would offer a different approach to lattice gauge theory. Calculations could be performed for real time, allowing access to the dynamics that are inaccessible to classical computers. Doubling the size of the Hilbert space would only require the addition of a single qubit. Recent studies have explored various questions related to the quantum computation of lattice gauge theories Osterloh:2005ha ; Buchler:2005fq ; Byrnes:2005qx ; Zohar:2012ay ; Zohar:2012xf ; Banerjee:2012pg ; Banerjee:2012xg ; Wiese:2013uua ; Bazavov:2015kka ; Zohar:2016iic ; Bermudez:2017yrq ; Banuls:2017ena ; Kaplan:2018vnj ; Stryker:2018efp ; Kokail:2018eiw ; Raychowdhury:2018osk ; Zhang:2018ufj ; Hackett:2018cel ; Lamm:2019bik .

This paper describes a qubit model for the simplest lattice gauge theory: compact U(1) in the absence of charged particles. Recall that compactness ensures that this theory is not a trivial theory of free photons, but has a strong coupling region where the theory has confinement and particles akin to glueballs Polyakov:1976fu ; Hasenfratz:1980ue ; Irving:1983ms ; Yung:1988jr ; Hamer:1991dy ; Chaara:1994ar ; Panero:2005iu ; Athenodorou:2018sab . Our model is meant to represent compact U(1) in the strong coupling region. The model is obtained by first writing the Hamiltonian in terms of a minimal set of independent gauge links in a manner that maintains Gauss’s law explicitly, then choosing a particular set of basis states so the Hamiltonian is expressed in terms of Pauli operators that would act on a register of qubits, and then truncating that basis set to a particular subset. Through direct calculations, the truncation error is observed to be small and manageable.

Section II describes the model on the smallest lattice: a single plaquette (four sites in the shape of a square). Section III uses the model on a 1-dimensional row of plaquettes to show the dynamics of excitations. We find that real-time propagation and real-time collisions can be observed quite readily in this model. Section IV applies the model to a 2-dimensional lattice of plaquettes. Section V explains a new issue that arises when extending the model to 3-dimensional lattices, and how it can be handled. A brief conclusion is given in Section VI. All of the calculations for this work were performed using simulators running on classical computers.

II The model on one plaquette

The Euclidean lattice action for U(1) gauge theory is well known Kogut:1979wt . For any hypercubic lattice, it is a sum over all plaquettes,

[TABLE]

where $\beta=1/g^{2}$ and $\theta_{\mu}(n)=agA_{\mu}(n)$, with $g$ the bare gauge coupling and $A_{\mu}(n)$ the vector potential. Choosing the temporal gauge (meaning $A_{t}(n)=0$) allows the continuum limit to be taken in the temporal direction, and brings us to the Hamiltonian for compact U(1) gauge theory:

[TABLE]

where $\vec{E}=\partial_{t}\vec{A}$ is the electric field. The conjugate momentum for the angle $\vec{\theta}$ (because it is a kind of angular momentum, let’s name it $\vec{L}$) is proportional to $\vec{E}$, and obeys the canonical commutation relation,

[TABLE]

This implies that link variables, $U_{j}=e^{i\theta_{j}}$, are ladder operators because they obey the corresponding commutation relations,

[TABLE]

Following Byrnes:2005qx , we can define a state $\left|0\right>$ by

[TABLE]

and then identify a sequence of states with integer eigenvalues,

[TABLE]

where $\ell=\ldots,-2,-1,0,1,2,\ldots$ The Hamiltonian can now be written as

[TABLE]

Let’s rescale energies by dropping the overall factor of $(2a\beta)^{-1}$ from now on. Having arrived at a convenient form for the Hamiltonian of compact U(1) lattice gauge theory, we now apply it to the smallest possible lattice.

Consider a 2$\times$2 lattice with Dirichlet boundary conditions. This lattice has four sites, four gauge links and only one plaquette so

[TABLE]

Because there are no charges in a pure gauge theory, Gauss’s law ensures that electric flux is conserved at every lattice site so $L_{i}^{2}$ is the same for all four links. Moreover, a local gauge transformation at any lattice site allows one adjoining link to be rotated to the identity, and this can be done at three sites on our lattice, leaving a minimal Hamiltonian:

[TABLE]

where $L^{\pm}$ are the ladder operators. Evaluating the Hamiltonian in the basis of electric field eigenstates and truncating to $-4\leq\ell\leq 4$ yields the Hamiltonian matrix

[TABLE]

Because the on-diagonal entries grow quadratically while the off-diagonals remain constant as we approach the corners of the matrix, truncation has a small effect on the lowest eigenvalues at strong coupling (i.e. small $\beta^{2}$). To be specific, in Fig. 1 it was necessary to truncate to a 4$\times$4 Hamiltonian in order to display a visible deviation, and even then the ground state deviation is barely seen. Notice that the 4$\times$4 case corresponds to non-symmetric limits because $\ell=-1,0,1,2$.

To run someday on a quantum computer (and today on a classical simulator), the states have to be mapped to qubits, and operators have to be expressed as products of Pauli operators acting on those qubits. A proof of principle was provided in Byrnes:2005qx though the authors acknowledged that it was a qubit-inefficient implementation. A more efficient alternative here will be used here.

An optimal description of the $\ell=-1,0,1,2$ truncation needs only two qubits. Including a third qubit doubles the number of basis states. Our mapping conventions for those two cases are displayed in Tab. 1.

The operators that appear in the Hamiltonian, Eq. (12), can be written as follows:

[TABLE]

where $\sigma_{i}$ acts on the $i$th qubit in the register. All terms containing $\sigma_{2}$ can simply be omitted for the case of a two-qubit register. Conversely, an extension to a register with more than three qubits is straightforward though unnecessary for our purposes.

Translating from the qubit ladder operators $\sigma^{\pm}$ to standard Pauli operators, the two-qubit Hamiltonian is

[TABLE]

and the three-qubit Hamiltonian is

[TABLE]

The time evolution operator in a quantum theory is $e^{-iHt}$ which requires exponentiating the Hamiltonian. The second-order Suzuki-Trotter formula Lloyd ,

[TABLE]

is used to rewrite the time evolution, which is due to sums of non-commuting terms in the Hamiltonian, into products of unitary operators. High level quantum computer programming languages being developed today provide for exponentiated operators such as $e^{i\sigma^{j}t}$ for $j=x,y,z$ .

As an explicit example of time evolution, choose $\beta^{2}=1$ and begin with an initial state (at $t=0$) that is equally weighted between the ground and first excited states,

[TABLE]

As shown in Fig. 2, the probability of remaining in the initial state oscillates between zero and unity. Several oscillations pass before a distinction is noticed between the calculations with two-qubit and three-qubit registers.

For a final calculation on the one-plaquette lattice, we solve for the ground state eigenvalue and eigenstate. Of course these can be found easily by matrix diagonalization (recall Eq. (13) and Fig. 1), but now an implementation of the variational principle for a quantum computer Peruzzo will be used. To begin, write the two-qubit Hamiltonian’s expectation value as a function of probabilities,

[TABLE]

The quantum computer can be repeatedly prepared in state $\left|\psi\right>$, and then measurements of the appropriate spin states will determine the probabilities in Eqs. (22-25). Those, in turn, give the Hamiltonian’s expectation value.

For a simple trial state, we anticipate that the ground state will be dominated by $\ell=0$, that it will have essentially equal contributions from $\ell=\pm 1$, and that it will have a negligible contribution from $\ell=2$. Therefore we use

[TABLE]

A four-step procedure will create this state in the quantum register. Step 1: Begin with $\left|00\right>$, Step 2: Apply X(1) to get $\left|10\right>$, Step 3: Apply $R^{y}(2\theta,0)$ to get $\cos\theta\left|10\right>+\sin\theta\left|11\right>$, Step 4: Apply $CR^{y}(-\tfrac{\pi}{2},0,1)$ to get $\left|\psi\right>$.

A calculation of the Hamiltonian’s expectation value is plotted in Fig. 3. Our simple trial state identifies an optimal value for $\theta$ that gives a good estimate of the exact ground state eigenvalue. As the plot indicates, many measurements are needed before a precise determination is obtained.

III A one-dimensional array of plaquettes

To explore the real-time propagation of excitations, we now extend the lattice from a single plaquette to a row of plaquettes. The Hamiltonian is again Eq. (10) and the factor of $(2a\beta)^{-1}$ can again be dropped. Application of Gauss’s law and gauge invariance reduces the Hamiltonian to one independent gauge link per plaquette. For example, a row of three plaquettes gives

[TABLE]

which is simply three copies of the single-plaquette Hamiltonian plus an extra term for each gauge link that is shared between plaquettes. The operators for each of the three links are given by Eqs. (14-16). The extension to any number of plaquettes in a row is straightforward.

For definiteness, choose a row of seven plaquettes described by two qubits per gauge link. The time evolution of that system can be computed on a 14-qubit register from any chosen initial state. Consider the initial state where the first plaquette begins in its first excited state and all others begin in the ground state. Time evolution will allow us to observe the motion of that excitation in real time. The exact eigenvalues and eigenstates for an individual plaquette can be determined by diagonalizing the appropriate 4$\times$4 truncation of Eq. (13), and that is what will be used to build our initial state.

The upper plot in Fig. 4 shows our results at $\beta^{2}=1$, calculated by expressing Eq. (27) in terms of Pauli operators and then using the second-order Trotter formula, Eq. (19), to translate into unitary operators acting on the qubit register. Most calculations were performed on a homemade simulator written in fortran, with some checks carried out in the wavefunction simulators of pyQuil Smith and ProjectQ Steiger .

The excitation that begins at $L_{1}$ moves across the lattice to $L_{7}$ and then bounces back. Only a modest dependence on the number of qubits is observed. The small jagged effects in the $L_{1}$ and $L_{2}$ curves are reminders that neighbouring plaquettes share a gauge link at their boundary, leading to the mixing terms in the Hamiltonian of Eq. (27). This is no problem, it just means that our chosen initial state reflects the fact that the plaquettes on our lattice do not factorize. Another consequence is that the seven curves in Fig. 4 do not add exactly to unity, because the probabilities defined by these curves are not exactly mutually exclusive. The deviation from unity is roughly 5 percent and, as must be the case, the sum is larger (not smaller) than unity at all times.

The lower plot in Fig. 4 shows the time evolution of a different initial state. In this case, $L_{1}$ and $L_{7}$ (at opposite ends of the lattice) both begin in the single-plaquette first excited state while the others ($L_{2}$ through $L_{6}$) begin in the single-plaquette ground state. Time evolution shows the excitations traveling to meet at the centre of the lattice and then traveling to the lattice walls where they bounce toward one another again. Comparison of the upper and lower plots indicates that the interaction between two traveling excitations has caused them to traverse the lattice in less time than a single excitation, since the large peak for $L_{7}$ has moved from time $\gtrsim$20 (upper plot) to time $\lesssim$19 (lower plot). An attractive interaction between the two excitations would cause this effect, while a repulsive interaction would cause the opposite effect. The classical analogy is two balls rolling on side-by-side tracks, where one track has a valley and the other has a hill. Both balls begin and end at the same height and same speed, but the ball with the valley (attractive potential) takes less time and therefore wins the race.

Figure 5 indicates that the attractive interaction persists on a longer lattice, while Fig. 6 shows a slightly repulsive interaction for a larger value of $\beta^{2}$. Exploration of a few other $\beta^{2}$ values led to Fig. 7 which suggests an attractive region at strong coupling and a repulsive region at weaker coupling. As the gauge coupling vanishes ($\beta^{2}\to\infty$), the interaction vanishes.

IV A two-dimensional lattice

The model presented in this work is directly applicable to any planar lattice. As was true for a one-dimensional array of plaquettes, Gauss’s law ensures that any two-dimensional array of plaquettes has one independent gauge link per plaquette. Consider the example of a 4$\times$4 lattice, which is nine plaquettes in the shape of a square. According to Eq. (10), the Hamiltonian is

[TABLE]

After expressing this Hamiltonian in terms of Pauli operators and applying the Trotter formula, Eq. (19), to a chosen initial state, its time evolution is obtained.

Fig. 8 shows the case of $\beta^{2}=1$ for two qubits per gauge link (giving 18 qubits in total). The upper plot begins with $L_{1}$ in its first excited state and all others in the ground state, so it corresponds approximately to one excitation beginning in a corner of the lattice. We say “approximately” because $L_{1}$ affects the boundaries of neighbouring plaquettes as well. This was an approximation for the row of plaquettes in Sec. III also, but it is more significant now because this square lattice has more boundaries per plaquette than were present in the row of plaquettes. A qualitative interpretation is clear from the graph. The excitation travels diagonally across the lattice and then bounces back to the original corner. The probability of exciting the intermediate corners en route is noticeably smaller than the probability of exciting the center plaquette.

The lower graph in Fig. 8 shows the time evolution of an initial state where $L_{1}$ and $L_{9}$ begin in the first excited state with all others in the ground state. Qualitatively, the graph suggests that the two excitations travel from their opposing corners, meet in the center of the lattice, travel to the corners that have not yet been visited, return to the center, and then return to their original corners. This particular sequence of events depends on the initial condition and on the choice of $\beta^{2}$. The lower graph in Fig. 8 does not have any peaks that reach above 0.5 after the initial peak, so it would be challenging to analyze this time evolution from probabilities measured on a quantum computer.

V A three-dimensional lattice

In contrast to the lower-dimensional cases, a cubic spatial lattice will have less than one independent gauge link per plaquette, and that fact is a significant issue for implementation of the model. The smallest example is the $2^{3}$ lattice, which is a single cube and therefore has six plaquettes, but Gauss’s law says that only five gauge links are independent. To prevent the magnetic part of the Hamiltonian (i.e. the terms proportional to $\beta^{2}$) from ballooning into a huge number of Pauli-operator terms, it is important to choose the independent links in a clever way. Additionally, it is important to choose the independent links in such a way that this cube will fit into a larger lattice. For example, a $2^{2}\times 3$ lattice is two cubes side by side but the number of independent gauge links is 5+4=9, not 5+5=10.

With $2^{2}\times N$ lattices in mind, we use the labeling displayed in Fig. 9. For the first plaquette, notice that $L_{1}^{+}$ is a raising operator for the left side, $L_{1}^{+}L_{2}^{+}$ is a raising operator for the top, $L_{2}^{+}L_{3}^{+}$ is a raising operator for the back, $L_{3}^{+}L_{4}^{+}$ is a raising operator for the front, $L_{4}^{+}L_{5}^{+}$ is a raising operator for the bottom, and $L_{5}^{+}$ is a raising operator for the right side. The second cube has $L_{5}$ already in place for its left side, and adds $L_{6}$, $L_{7}$, $L_{8}$ and $L_{9}$ as its new links. Each cube has the same layout of independent links, so the row of cubes can be made arbitrarily long by attaching identical cubes to either end.

From Eq. (10), the Hamiltonian for a single cube is

[TABLE]

Ensuring that no $\beta^{2}$ term has a product of more than two $L_{i}^{\pm}$ operators was one reason for choosing this particular $L_{i}$ basis.

At $\beta^{2}=0$, the smallest eigenvalue of this Hamiltonian is zero, the first excited state has the value 4 and is 12-fold degenerate, and the next state has value 6. Increasing $\beta^{2}$ breaks the degeneracies. Choosing two qubits per gauge link gives four basis states per gauge link and makes this Hamiltonian a $2^{10}\times 2^{10}$ matrix which, in contrast to Eq. (13), is not tridiagonal. The effect of truncating to two qubits per gauge link is demonstrated by comparing the dashed curves to the solid curves in Fig. 10, which shows the smallest eigenvalues in the strong-coupling region as calculated by matrix diagonalization.

To calculate time evolution, we use the three-cube lattice of Fig. 9. A first thought is to put $L_{1}$ in the single-plaquette first excited state and all other plaquettes in the single-plaquette ground state. However, every gauge link on this lattice is shared among multiple plaquettes, and the time evolution does not display clear probability peaks like those observed for the planar lattices in Sec. III. Another initial condition is to put the entire left cube into the single-cube first excited state, leaving the other cubes in the single-cube ground state. That option has gauge links shared between neighbouring cubes, but the upper panel of Fig. 11 does show a traveling excitation for some time until the signal disperses. The lower panel of Fig. 11 begins with both end cubes in the first excited state and the center cube in the ground state, leading to rapid dispersal of the signal. A longer lattice and a different choice for the initial state could allow the excitation to persist for a longer time, but we will leave this for future work.

Eq. (29) is a convenient choice for a row of cubes, but it would not be convenient for a plane of cubes. The reason is that only two sides of the cube have gauge links ($L_{1}$ and $L_{5}$) that are ready to serve a neighbouring cube. Now we will write down a different labeling for a cube:

[TABLE]

which corresponds to Fig. 12. Cubes of this type can immediately be used side by side to fill a plane, because $L_{1}^{+}$, $L_{2}^{+}$, $L_{4}^{+}$ and $L_{5}^{+}$ are the raising operators for the four sides of the cube and are directly usable by neighbouring cubes. However, Eq. (30) requires significantly more Pauli factors than Eq. (29) due to the triple products of $L_{i}^{\pm}$ operators in the terms proportional to $\beta^{2}$, so there is a price to be paid when moving beyond a row of cubes.

A final step is to extend the lattice into a 3-dimensional volume. This is necessarily more expensive than for a row or plane of cubes, but there is a recommendation that maintains Gauss’s law and requires no more than five $L_{i}^{\pm}$ operators in a single $\beta^{2}$ term. This remains true for any size of the lattice volume. The recommendation is to build cube towers of any height following the template in Fig. 13, where the floor’s raising operator is $L_{1}^{+}$ and every wall plaquette is a product of one or two $L_{i}^{\pm}$ factors.

Towers like this can immediately be put side by side to create any lattice volume. The Hamiltonian for the two-cube example in Fig. 13 is

[TABLE]

To add another cube on top of the tower shown in Fig. 13, we would define $L_{10}$ around the entire 8-link path on the front of the tower, just as $L_{6}$ appears in the entire path on the front of the present tower. Similar definitions apply to $L_{11}$, $L_{12}$ and $L_{13}$ on the other sides of the tower. An ever higher tower requires ever longer paths for the newest gauge links, representing a nonlocality in the Hamiltonian that was not required for a row or plane of plaquettes.

VI Conclusions

In this work, we have presented a model for compact U(1) lattice gauge theory in 2 or 3 spatial dimensions. Gauss’s law is implemented fully, leaving no redundant gauge field degrees of freedom. The model can be implemented with any number of qubits per gauge link. The smallest eigenvalues of the single plaquette and the single cube, plotted as functions of $\beta^{2}$ in Figs. 1 and 10 respectively, indicate that a truncation to just two qubits per gauge link is already a useful model.

This work reports on calculations performed on $2\times N$ lattices, i.e. a row of plaquettes, showing that simple initial conditions produce visible propagation and collisions of excitations. Calculations on nine plaquettes forming a square lattice produce similar results, though excitations disperse more quickly when this extra spatial dimension is available. Extending the model to a 3-dimensional volume requires a non-locality that was not present in the planar Hamiltonian, and an implementation of the 3D model has been presented in this work.

U(1) is the simplest gauge theory, and the qubit model presented here is useful even with an aggressive truncation. Therefore this model provides an opportunity to explore many aspects of quantum computation using only a small number of qubits. An interesting next step could be to design initial states that allow excitations to propagate for longer times on planar or 3D lattices. We envision such studies being performed by simulators running on classical computers, as was done in the present work, but perhaps it will not be too long until the model can be implemented directly on quantum computers as well.

Acknowledgments

RL is grateful to the organizers of Lattice for BSM Physics 2019 for the opportunity to present this work and to the participants for valuable conversations. This work was supported in part by the Natural Sciences and Engineering Research Council (NSERC) of Canada. TRIUMF receives federal funding via a contribution agreement with the National Research Council (NRC) of Canada.