Topological two-dimensional Floquet lattice on a single superconducting qubit

Previous theoretical and experimental research has shown that current NISQ devices constitute powerful platforms for analogue quantum simulation. With the exquisite level of control offered by state-of-the-art quantum computers, we show that one can go further and implement a wide class of Floquet Hamiltonians, or timedependent Hamiltonians in general. We then implement a single-qubit version of these models in the IBM Quantum Experience and experimentally realize a temporal version of the Bernevig-Hughes-Zhang Chern insulator. From our data we can infer the presence of a topological transition, thus realizing an earlier proposal of topological frequency conversion by Martin, Refael, and Halperin. Our study highlights promises and limitations when studying many-body systems through multi-frequency driving of quantum computers.


I. INTRODUCTION
Noisy intermediate-scale quantum (NISQ) computers may not yet offer fully fault-tolerant quantum computing facilities, but they nevertheless constitute a versatile experimental platform with the potential for fundamental research, small-scale computation or quantum simulation [1]. The typical model of a quantum computer is that of a quantum circuit, which is a sequence of gates applied to the qubits. A small set of gates is sufficient to in principle perform any computation on the quantum computer [2]. This allows digital simulation of many-body quantum mechanics, with important applications in quantum chemistry [3,4]. In principle, the time-evolution of any many-body quantum systems can be simulated by applying a Trotterization, which turns the continuous time evolution into a local quantum circuit [5]. This approach to quantum dynamics has now been benchmarked for a range of different models on existing quantum computers [6][7][8][9].
The abstraction into quantum circuits hides the complexity of the underlying platform. In superconducting circuits, the currently leading technology, the set of available gates correspond to a set of carefully calibrated microwave pulses that are applied to its input ports, and a quantum circuit is the concatenation of these pulses [10]. In the absence of driving, the continuous evolution of a superconducting-qubit quantum computer is similar to a disordered Bose-Hubbard model [11], a fact that suggests that these devices may also offer exciting direction in analogue quantum simulation. This perspective has been explored in a series of theoretical and experimental works [12][13][14][15]. If the intrinsic many-body nature of quantum computers is combined with the capacity to apply essentially arbitrary drives, they may serve also as powerful analogue quantum simulators for very large classes of time-dependent Hamiltonians.
The evolution under time-dependent Hamiltonians is incredibly rich and exhibits many novel phenomena, even at the level of individual qubits. A particular example is the temporal topological transition that occurs in the presence of * Both authors contributed equally. quasiperiodic driving, theoretically predicted by Martin, Refael, and Halperin in 2017 [16]. Using a Floquet treatment of the driven qubit, the dynamics is related to the properties of a two-dimensional lattice model, the Bernevig-Hughes-Zhang (BHZ) Chern insulator [17]. As a function of time, the driven qubit explores the whole Brillouin zone, which causes the work done by the two drives to be quantized and proportional to the integer Chern number, which is determined by the parameters in the drives. Recently, temporal topology has been classified, analogously to the classification of topological insulators [18], and extension to larger systems, such as spin-resonator systems [19], or two-qubit systems with interactions [20] have been proposed. While quasiperiodic driving typically maps to systems without boundary, one can in principle also introduce boundaries through quantum feedback [21].
In this work we experimentally demonstrate this temporal topological behaviour of a single-qubit system on an existing quantum device, using continuous driving, implemented with the fine-grained access offered by QISKIT pulse [22]. Choosing a specific driving with two incommensurate frequencies, we observe a topological transition in the temporal dynamics of the qubit, finding good agreement with simulations.
Motivated by this experiment, we theoretically derive the class of Hamiltonians that can readily be implemented on state-of-the-art quantum computers. As one concrete example, this offers an exciting perspective to study strongly interacting Floquet systems [23][24][25][26][27] with an exquisite level of control. The high-level site-specific control as well as highfidelity single-site readout confers quantum computers certain advantage over other quantum simulators based on light [28] or ultracold atoms [29], making them ideally suited for the analogue quantum simulation of generic many-body timedependent Hamiltonians.

II. THEORETICAL DESCRIPTION
In quantum computers based on superconducting transmonlike qubits [30,31], the Hamiltonian describing a single qubit can be cast in the form of a Duffing oscillator with driving (1) arXiv:2012.01459v1 [quant-ph] 2 Dec 2020 The qubit frequency is denoted ω 0 and U quantifies the anharmonicity that separates the lowest two levels that define the qubit from the higher levels of the superconducting circuit. The ideal drive signal is parametrized as [22] where Ω max denotes the maximum Rabi frequency attainable in the system, and we choose the detuning ∆ c = 0. Anticipating our later interpretation as a spin-1/2 particle in a timedependent field, we parametrize the dimensionless drive shape d(t) ∈ C, |d(t)| ≤ 1 in terms of the dimensionless magnetic fieldh with time-dependent phase In order to treat the system as a qubit, the maximum drive strength needs to be much weaker than the anharmonicity Ω max U. Assuming this is fulfilled, applying a rotatingwave approximation and passing into a frame rotating with the qubit frequency, we rewrite Eq. (1) as where we have now introduced the dimensionful magnetic field h(t) = Ω max h (t) and the spin lowering operator σ − . If we pass (again) into a rotating frame with respect to the timedependent Hamiltonian H z (t) = −h z (t)σ z , the Hamiltonian becomes that of a spin-1/2 particle in a time-dependent magnetic field Note the subtlety when passing from Eq. (5) to Eq. (6). Instead of changing the qubit frequency in the lab frame, we re-parameterize time by passing into a rotating frame with respect to H z (t). This constitutes a continuous version of virtual Z gates [32].

III. TOPOLOGICAL FREQUENCY CONVERSION
Given Eq. (6), we can now realize the proposal by Martin et al. [16]. We apply the time-dependent magnetic field where ω 2 /ω 1 = (1 + √ 5)/2. We consider the model in the strong drive limit, i.e., η ω 1 , ω 2 . In this limit we can make a direct connection to the two-dimensional BHZ Chern insulator [17]. We extract the non-trivial Chern number by calculating the work done by each drive.
As a function of M the system shows a crossover from a regime in which the average work done is zero (no topology) to a regime where energy is pumped from one drive to the other at a constant rate Fig. 1. Work done Work done by the two incommensurate drives, calculated using Eq. (12). The experimental data for M = 1, ω 1 = 0.125 measured at 800 points is compared with a numerical simulation of the same setup. The experimental data is fitted with a line (coloured dashed lines) using least-squares regression, for which the 95% confidence interval for the slope is shown as the coloured region. The expected slope is shown as black dashed lines.

A. Floquet lattice
The connection to the BHZ Chern insulator is elucidated by considering the Floquet construction which relates our timedependent zero-dimensional system with a two-dimensional lattice Hamiltonian. This is done by considering the following the time-dependent states where α = 0, 1 label the two basis states for the qubit, and n = (n 1 , n 2 ), with n i ∈ Z and ω = (ω 1 , ω 2 ). The vector n is associated with the sites on an infinite two-dimensional square lattice. The time-dependent Schrödinger equation i∂ t |Ψ(t) = H(t)|Ψ(t) , can then be written as where ψ k = n e −i k· n ψ n , and k = (k 1 , k 2 ) with k i ∈ [0, 2π). The first three terms correspond to the BHZ Chern insulator, discussed in more detail in the Sup. Mat.. The final term is in mixed form and corresponds to a linear potential along the direction ω, which can be interpreted as an electric field. In the strong-drive limit, η ω 1 , ω 2 , this induces a time dependence in the momentum k → k 0 + ωt [16]. Making the substitution for the time-dependent momentum we recover Eq. (7). The dynamics corresponds to the adiabatic evolution of the ground state of Eq. (7) in the Brillouin zone, so long as the gap controlled by M is sufficiently large compared with ω 1 .

B. Chern number
The topological properties of the time-dependent Hamiltonian Eq. (7) manifest as the quantised average energy pumping between the two external drives with frequencies ω 1 and ω 2 . This average pumping rate is set by the Chern number of the corresponding BHZ phase. Note that for this quantization to hold we need that the drives are incommensurate and we are in the strong-drive limit.
To experimentally extract this pumping rate we measure the work done by each of the drives. If we first split the Hamiltonian into the two contributions from each drive, then the work done by each drive over a period T is given by where |Ψ(t) is the state of the qubit at time t under the time- We can therefore extract the Chern number by measuring the dh i /dt as a function of time and computing by integral in Eq. (12).

C. Experimental details
For our simulation we initialize the qubit in the instantaneous eigenstate of Hamiltonian Eq. (7) at t = 0 and then drive the system to realize the evolution under the time-dependent Hamiltonian Eq. (7). In practice, this means that we drive the qubit according to the Hamiltonian in the lab frame Eq. (5), perform state tomography at 800 points in time and then obtain the observables in the rotating frame [corresponding to evolution under Eq. (7)] by applying the rotation where φ(t) is given in Eq. (4). When postprocessing the data, we also mitigate some of the errors due to decoherence and measurement by projecting the measured qubit density matrix onto the Bloch sphere. We find that this method of error mitigation improves quantitative agreement with numerical simulations in all cases we considered. We fix ω 2 = ϕω 1 , where ϕ = (1 + √ 5)/2 is the golden ratio. Ideally we would be like to set the ratio ω 1 /η as small as possible. However, due to the finite coherence time τ ≈ 100µs of the IBM Q device, we must choose ω −1 1 τ. We also choose η = 0.9Ω max in order to avoid driving transitions to higher excited states. Using numerical simulations we found the best compromise was to choose a total simulation time of 20µs and set ω 1 = 0.125η, which corresponds to ω −1 1 ≈ 240ns. In order to improve the fidelity of our results we start the drive with a linear ramp of the drive frequencies ω 1 , ω 2 over a period of 444ns. This reduces transient effects and reduces highfrequency Rabi oscillations seen in the simulation results. The single qubit IBM device we used, codenamed armonk, had qubit frequency ω 0 ≈ 4.97 GHz and maximum Rabi frequency of Ω max ≈ 36.9 MHz.

D. Results
To obtain the work done by each drive [Eq. (12)], we first compute dh i (t)/dt = d h i (t)/dt · σ from the data. We then perform the integration in Eq. (12) numerically. Figure 1 shows the experimental results for the work done, for the case of M = 1. This corresponds to a phase of the BHZ model with Chern number C = −1, resulting in an average linear decrease (increase) of W 1 (W 2 ). The experimental results are in good quantitative agreement with numerical simulations. Furthermore, by fitting a linear curve to this data using least-squares regression, the slope is in close agreement to that predicted by Eq. (13). Fig. 2 shows the extracted Chern number for a range of values of M, compared against numerical simulations of the same setup as well as the ideal result in the strong-drive limit. We find reasonable quantitative agreement with the simulations. Furthermore, the transition between different phases with different Chern number is clearly visible in the experimental data. Beyond the extracted value of the Chern number, the qualitative difference between the phases with |M| < 2 and |M| > 2 is clearly seen by considering the portion of the Bloch sphere covered under time evolution, shown inset Fig. 2. For |M| < 2 we observe that the state of the qubit explores the full Bloch sphere under the dynamics of Eq. (7) reflecting that the surface traced by the time dependent state has non-trivial winding around the origin and hence non-zero Chern number. When |M| > 2 the state of the qubit is instead restricted to either the north or south hemisphere of Bloch sphere and does not wind around the origin and the Chern number is zero.
While overall, the agreement is good, some outliers appear in the experimental data. By plotting the equivalent of Fig. 1 for those points [cf. SI], we observe that at these points the evolution of the qubit departs from the ideal evolution predicted by the numerics. We thus attributed those points to effects that have not been taken into account in our model, most notably excitations to doubly excited states and decoherence, but it is difficult to tell which effect exactly is responsible.

IV. EXTENSION TO QUBIT ARRAYS
In view of the vast research effort in Floquet matter, the question is pertinent whether the presented scheme generalizes to several qubits. To answer this question we must be aware that the coupling between two qubits can be engineered in many ways, as one can choose between direct capacitive or inductive coupling, or indirect coupling via a resonator [33][34][35][36][37]. Indeed, alternative approaches have been taken by IBM (see, e.g., Ref. [38]) and Google [39], with each implementation having their own (dis-)advantages.
In Appendix A we consider two cases in particular. The first is when both qubit frequencies and their coupling is fixed, which is relevant for the devices developed by IBM. In this case, the coupling has to be engineered with time-dependent driving of the qubits in order to implement resonant processes to second order in the Hamiltonian [40][41][42]. This technique produces tuneable XZ, YZ, ZX, and ZY interactions, which are used to engineer CR gates [41,43]. Implementing a similar single-qubit driving as before (Eq. (5)), and passing into a time-dependent reference frame, an effective Hamiltonian of the form can be derived (Appendix A). On bipartite lattices, a (virtual) rotation of every second spin maps this to either an Ising interaction or (in general anisotropic) XY interactions, making this technique very versatile. A drawback is that the Hamiltonian (15) is obtained to second order in the original Hamiltonian and by neglecting quickly rotating terms. It is thus an approximation and care has to be taken that all the steps in its derivation are valid. Nevertheless, these conditions can usually be fulfilled through careful choice of the driving parameters. It may further be possible to actively counteract unwanted effects from this approximation, analogous to the DRAG scheme used to improve the fidelity of digital quantum gates [10]. If tuneable interactions are available and the qubits can be brought into resonance [39], an approach pursued (for example) by the Google team, the driving scheme is simplified and one readily arrives at the general spin Hamiltonian with a magnetic field along z that takes the form where ω i are the qubit frequencies and ω 0 is some arbitrarily chosen reference frequency. We note that tuneable qubit frequencies feature in many implementations [44][45][46][47], but not in all [48]. In the latter case-of fixed qubit frequencies but tuneable interactions-one option might be to drive each qubit with a far off-resonant drive to induce an drive-strengthdependent AC Stark shift. When taking into account the second excited state of each qubit, the qubit array Hamiltonian (16) maps to a Bose-Hubbard model with time-dependent hopping and freely tuneable site-dependent drive and disorder. We note that the interpretation as a (time-independent) Bose-Hubbard model has enabled the experimental measurement of microscopic features of the many-body localized phase [14]. Moving to periodically varying hoppings could allow one to study manybody Floquet models and implement quench and ramp experiments from carefully prepared initial states.

V. DISCUSSION
In our initial discussion, we have emphasized the vast promise that time-dependent analogue quantum simulation on quantum computers holds. However, our experiment has also highlighted some of the difficulties that need to be overcome on the way. Concretely, the experimental results shown in Fig. 2 would improve with greater coherence time of the qubit, as this would allow us to reduce the modulation frequencies ω 1 and ω 2 in the Hamiltonian, which in turn improves the strong-driving and adiabatic modulation approximation. This would make the blue simulated dots in Fig. 2 agree better with the step function. Furthermore, the actual results from the quantum computer differ from the simulation, which highlights the need to improve our model of the quantum computer.
Nevertheless, as we have demonstrated, the temporal topological nature of the qubit time evolution can clearly be extracted despite experimental shortcomings. With the current progress in superconducting-qubit technology, we expect that our understanding and the fidelity of available quantum computers to increase rapidly, allowing for more complex experiments, in particular with more qubits.
Going beyond (quasi)periodic driving, the capability to engineer time-dependent many-body Hamiltonians offers many exciting perspectives to investigate non-equilibrium physics.
For example, ramping through a quantum phase transition might allow one to study Kibble-Zurek scaling or in general the dynamics of phase transitions such as the superfluid to Mott insulator transition [49]. Slow variation could also be used to explore adiabatic algorithms and departures from them. Many-body non-equilibrium physics is notoriously difficult to study with classical computers and this is therefore a prime area of applicability for quantum simulators. Noisy intermediate-scale quantum computers offer a versatile combination of single-site control and readout, and large enough system sizes, and thus promise to support and complement analytical and computational approaches to understand manybody non-equilibrium physics. As we have demonstrated with the single-qubit experiment, this is a realistic outlook.

DATA AVAILABILITY
All of the data presented in this paper, along with the python code used to generate this data and the figures in this paper, can be found in a public GitHub repository [50].
Appendix A: Qubit array Hamiltonian

Fixed qubit frequencies and fixed coupling
If neither the qubit frequencies nor the coupling between them is tuneable, the interaction Hamiltonian is timeindependent. Restricting directly to the qubit subspace and assuming capacitive coupling for simplicity, we have If the coupling J i j is weak compared to the relative detuning, the qubits are effectively decoupled. Coupling between the i th and j th qubit can be turned on by driving either one at the frequency of the detuning between them [40,41]. To see this, let us pass to the frame rotating with the qubit frequencies, and neglect terms rotating at the qubit frequency, which yields where we parametrize the drive in two steps where d (0) i (t), g i j (t) are slow functions (compared to the coupling rate), we have defined the qubit-qubit detuning ∆ i j = ω i − ω j and the sum over j is taken to run over coupling qubits for which J i j 0.
We can now separate the Hamiltonian into a slow and a fast part. The slow Hamiltonian corresponds to the single-qubit Hamiltonian derived in the main text Eq. (5), except for many qubits The fast part is off-resonant and has little effect, except when there are resonant processes. For example, to second order, the hopping becomes resonant again, leading to an Ising term of order |J i j | 2 /∆ i j . Since our original assumption is that the coupling is weak compared to the detuning, J i j ∆ i j , this is negligible. However, the other term that appears at second order, mixing the off-resonant drive g i j with the hopping J i j is not negligible and in fact gives rise to the interaction Combining Eqs. (A4) and (A5) allows one to engineer a wide range of time-dependent Hamiltonians. Note that due to drive enabling a resonant interaction between coupled qubits, the disorder has disappeared. The disorder can be restored by detuning the "coupling drive" g i j (t) = exp i t h ( j) z (t )dt g i j (t) and performing the same transformation on the single-qubit drives d (0) The induced coupling is rich, as we have independently controllable XZ, ZX, YZ, and ZY interactions. If the qubit coupling graph is bipartite, i.e., if it can be separated into two subgraphs A and B such that nodes from A are only connect to B and vice versa, as for example in a square lattice, we can make a unitary transformation on every qubit one of the sublattices (B, say) that rotates Z into X and X into −Z. This allows us to engineer two particularly important interactions (i) We only have coupling XZ, with X on the A lattice and Z on the B lattice. After the unitary transformation, the coupling is XX. Similarly, starting from the coupling ZX, one obtains ZZ.
(ii) We employ both XZ and ZX. The unitary transformation yields g 1 XX + g 2 ZZ, or equivalently g 1 XX + g 2 YY.
Finally, we note that the biggest limitation for the controllability of this setup an effect similar to frequency crowding. As this approach requires one to apply a range of different control tones at different frequencies, unwanted resonances can occur, which may affect the final effective Hamiltonian for the system. This problem particularly concerns lattices with high connectivity, but can typically be avoided by careful consideration and design of the applied tones.

Tuneable frequencies, but fixed interaction
If the qubit frequencies are tuneable, all the techniques from the previous section can be employed, but one gains the additional capacity to preclude frequency crowding.
In addition, one can make use of the direct interaction of the qubits by bringing them close to resonance. This realizes an XY interaction with tuneable disorder.

Tuneable interactions
For a set of qubits that interact via a tuneable XY interaction, the full Hamiltonian reads Parameterizing the drive as before in terms of an effective magnetic field h i (t), passing to a rotating frame with respect to the Hamiltonian H 0 = i (ω 0 − h (i) z (t))a † i a i and restricting to the qubit subspace, we find Choosing h (i) z (t) = h z (t) for all i removes the phase on the interaction terms. The resulting Hamiltonian is This is the Hamiltonian quoted in the main text. The disorder in the energies is inherited from the intrinsic disorder of the platform.
(B1) This Hamiltonian has the formĤ( k) = h( k) · σ = a h a ( k)σ a , and spectrum given by λ ± = ±| h|, or more explicitly, which is gapped except at M = −2, 0, 2. Depending on the value of M this model realizes three different phases with these gap closings marking the phase transitions between them. For |M| > 2 we have a trivial insulator. For −2 < M < 0, we have Chern insulator with Chern number C = +1, and 0 < M < 2 describes a Chern insulator with C = −1.
The eigenstates for the upper and lower band have the particularly simple form which correspond to vectors on the Bloch sphere (anti-)aligned with the field direction h. This can easily be seen by considering the vector of expectation values σ = ± h/h.

Chern number
To compute the Chern number we require the Berry connection A i ( k) and/or the Berry curvature F i j ( k), defined by The Chern number C ∈ Z is then computed by integrating the Berry curvature over the Brillouin zone: For the two-band model, the Berry curvature has a particularly intuitive form, namely The vector h defines a map from the torus (B.Z.) to a twodimensional surface M embedded in three dimensions. The Berry curvature then corresponds to the dot product of a field due to a unit charge at the origin and the normal of the surface S . The Chern number therefore measures the flux through the surface S , that is which measures the degree of map defined by h, or equivalently, the winding of the surface M around the origin. The surface prescribed by h( k) is closed self-intersecting surface with two lobes. the upper (along z-axis) lobe has normal facing outward, while lower has inward facing normal. The parameter M controls the position of this surface along the z-axis. The three phases therefore correspond to when the origin falls outside of the surface, or within one of the two lobes. This is shown schematically using a two-dimensional cross section in Fig. 3. Given this interpretation we can see that when |M| > 2 the Chern number must be zero since h z does not change sign.
Appendix C: Measuring the work done by each drive The computation of the work done by each drive proceeds by three main steps: i) state tomography, to obtain the spin expectation values σ at each time; ii) compute the expectation value dh i /dt ; iii) numerically integrate Eq. (12) to get the work done by each drive.
To perform state tomography we run the experiment three times (with 8192 shots each) and measure the qubit in the X, Y and Z basis. We take the computational or natural basis for the qubit to be the Z basis. Each measurement gives an outcome of 0 or 1 corresponding to the eigenstates of the Pauli-Z operator with positive (+1) and negative (−1) eigenvalues, respectively. By averaging over the 8192 shots we obtain the average expectation value Z , with standard error of the mean 2/ √ 8192 ∼ 2%. To measure in the X basis we apply a Hadamard gate before measuring the qubit, and for the Y basis we apply an S † gate and a Hadamard gate before measuring. We then move to the rotating frame using Eq. (14). Fig. 4 shows an example of the measured Pauli expectation values corresponding to the data shown in Fig. 1 (i.e., M = 1) compared with the exact numerical simulation. The measured expectation values will not in general be consistent with a pure state due to a range of factors including statistical errors, measurement errors, and decoherence. We therefore normalize the vector of expectation values σ , which corresponds to projecting onto the Bloch sphere.
From this projected tomography data we can then compute the expectation values dh where, e.g.x · σ = σ x . In Fig. 5 we show the expectation values computed from the Pauli expectation values shown in Fig. 4, corresponding to Fig. 1 of the main text. The final step is to numerically integrate these expectation values to get the work done by each drive, as given in Eq. (12). We use a simple trapezoid rule using data sampled at 800 points in time to 20 µs.