Geometrical Formalism for Dynamically Corrected Gates in Multiqubit Systems

The ability to perform gates in multiqubit systems that are robust to noise is of crucial importance for the advancement of quantum information technologies. However, finding control pulses that cancel noise while performing a gate is made difficult by the intractability of the time-dependent Schrodinger equation, especially in multiqubit systems. Here, we show that this issue can be sidestepped by using a formalism in which the cumulative error during a gate is represented geometrically as a curve in a multi-dimensional Euclidean space. Cancellation of noise errors to leading order corresponds to closure of the curve, a condition that can be satisfied without solving the Schrodinger equation. We develop and uncover general properties of this geometric formalism, and derive a recursion relation that maps control fields to curvatures for Hamiltonians of arbitrary dimension. We demonstrate examples by using the geometric method to design dynamically corrected gates for a class of two-qubit Hamiltonians that is relevant for both superconducting transmon qubits and semiconductor spin qubits. We propose this geometric formalism as a general technique for pulse-induced error suppression in quantum computing gate operations.


I. INTRODUCTION
Dynamical gate correction is an important topic in the field of quantum information technology because logical error correction schemes require that the individual gate error be below a given threshold value [1]. So, any quantum error correction protocol can only be implemented after individual gate operations are relatively error-free, often necessitating dynamical decoupling of environmental noise by external pulses. Dynamically corrected gates rely on using precise control of the underlying Hamiltonian to perform rotations such that the error introduced by a given noise source at one point in the pulse will exactly cancel the error introduced by the same noise source at other points in the pulse. This idea, adopted from NMR where pulses are used to reduce spin dephasing such as with the Hahn spin echo effect, is used to extend qubit coherence times [2][3][4][5][6][7][8][9][10][11][12][13][14][15][16]. Hahn spin echo and the closely related CPMG type multiple pulsing techniques are essentially the first dynamical decoupling techniques introduced for noise suppression in NMR and ESR experiments. Similarly, pulse sequences have been proposed to perform quantum gate operations robust to various noise sources [17][18][19][20][21][22][23][24][25][26]. Many of these pulses are derived by finding a series of rotations which combine to produce the desired gate while canceling error to first order or higher. Finding such a sequence often involves numerically solving large systems of nonlinear equations to determine the parameters that define each of the smaller rotations. Solutions to these equations are not unique, and finding one error-correcting pulse sequence only amounts to finding a local maximum in fidelity, and thus many of the resulting pulse sequences are far from optimal. Additionally, many techniques rely on using δ function or square pulses, neither of which can be precisely implemented in physical systems, as waveform generators have bandwidth and amplitude limitations. Thus, general methods for finding smooth, optimal dynamical-decoupling pulses are highly desirable.
One particularly powerful method for generating single-qubit error-correcting pulses is by representing these pulses as curves parametrized by the cumulative error at any given point in time. Ref. 27 first showed that for a simple Hamiltonian, pulses that correct against first-order error can be represented as closed curves in a plane, and that the driving field is proportional to the curvature of the curve. Conditions were given for canceling error at higher orders; specifically, for second order, the total signed area enclosed by the curve must equal zero. This formalism was used to derive the fastest possible pulses implementing specific single-qubit gates given constraints on the magnitude of the driving field [28]. In Ref. [29], this formalism was extended to more general single-qubit Hamiltonians with multiple control fields. Here pulses were represented as curves in three dimensions, with the strengths of the driving fields being related to the curvature and torsion of the curve. There have been other proposed techniques for the reverse engineering of the time-dependent Schrodinger equation for dynamical decoupling, but these ideas turned out to be impractical for actual implementation in quantum computing gate operations [30][31][32].
In this paper, we extend the geometric formalism to systems of multiple qubits subject to quasi-static noise. We show that many of the features of the single-qubit formalism carry over the multiqubit case. Pulses can be represented as curves in a higher-dimensional space, and those which cancel first-order error correspond to closed curves. The length of each pulse is given by the length of the corresponding curve, and the amplitudes of the driving fields are related to the curvature coefficients at each point along the curve. As a demonstration, we derive error-correcting pulses for a class of two-qubit Hamiltonians using curves in six dimensions. We focus on Hamiltonians that describe Ising-type interactions commonly arising in both superconducting transmon qubits and semiconductor spin qubits. We test the performance of these pulses with numerical simulations and confirm that the leading-order infidelity in the presence of noise is canceled. This paper is organized as follows: in Sec. II, we derive the geometric formalism for a general Hamiltonian. In Sec. III, we give an example of a specific two-qubit system and derive a pulse that implements a single-qubit gate on one qubit while suppressing noise on both qubits. We conclude in Sec. IV.

II. MULTIQUBIT GEOMETRICAL FORMALISM
We consider a generic Hamiltonian H 0 (t) that contains time-dependent driving fields. Suppose that there is some error term δH which is small compared to the scale of H 0 (t). In many qubit platforms including superconducting qubits and spin qubits, the time required to perform gates is generally much shorter than the time scale over which δH varies [33][34][35][36][37]. In this case, even dynamic noise can be well approximated as quasi-static, meaning that the noise fluctuation δH is treated as constant during a single pulse, but it can vary from one pulse to the next [20]. This is the situation where dynamical decoupling is most effective, and we employ this quasi-static approximation throughout this work. The full Hamiltonian is then H(t) = H 0 (t)+δH. Our goal is to determine driving fields in H 0 (t) that perform a desired gate while canceling the effects of noise to leading order in δH. We consider the evolution operator in the interaction picture. Let U (t) = T e −i t 0 H(τ )dτ be the time-dependent evolution operator of the full noisy Hamiltonian H, and let R(t) = T e −i t 0 H0(τ )dτ be the ideal noiseless evolution operator, where T represents the time-ordering operator. This means that the time derivative of R is Then the evolution operator in the interaction picture U I (t) satisfies U = RU I , and the time-dependent Schrödinger equation becomes The solution to the Schrödinger equation to first order in δH is given by (3) Note that the integral of H I is a measure of the total accumulated error at any time t, as the gate U actually performed in the presence of noise is equal to the desired gate R only if U I = 1.
In order to help derive noise-canceling gates, we want to represent this traceless part of U I (t) as a curve in some d-dimensional space. To this end, we define S to be the vector space on which H acts. For example, if H describes a n-qubit system, then S is the Hilbert space of n qubits. The set of all traceless Hermitian matrices which act on S themselves form another vector space over R, which we will call V (e.g. H will be an element of V). In the case of n qubits, V is the space of traceless nqubit Pauli strings, which has dimension dim V = 4 n − 1. Define an inner product on this vector space as follows: where dim S is the dimension of the vector space S, and therefore equal to the size of the matrices in V (e.g. for the case of a system of n qubits, dim S = 2 n ). Here, although V = V , we use the notation V , because in what follows, it helps to imagine that we expand V in terms of a basis of Hermitian matrices (e.g., n-qubit Pauli strings) with real coefficients. If these coefficients are time-dependent, then they trace out a curve in a higherdimensional Euclidean space. In what follows, this basis expansion will be kept implicit for the sake of generality and ease of notation. The inner product in Eq. (4) is invariant under unitary frame transformations, which will become important later.
The possible error terms that can be generated by any specific pulse from a Hamiltonian H will span some subspace of V, which we will call V err . The subspace V err is nonuniversal depending on details, and may or may not be equal to V depending on the degrees of freedom present in H 0 and what rotations are possible. For example, for a single-qubit system, V is the 3-dimensional space consisting of the X, Y , and Z Pauli matrices. In Ref. 27, it was shown that for a Hamiltonian with a single Pauli matrix, such as H 0 (t) = Ω(t)X, the error space V err is a 2-dimensional subspace of V. In comparison, Ref. 29 examined a Hamiltonian with two different driving terms such as H 0 (t) = Ω X (t)X + Ω Y (t)Y , and found that V err in fact spans the entire space of V in this case. In general, V err is the space spanned by δH and all vectors to which δH can be rotated by the evolution operator R corresponding to a particular choice of H 0 (t). For example, consider a Hamiltonian of the form where V i are products of Pauli matrices, and k is some arbitrary integer. Then V err will include linear combinations of δH and all vectors of the form Additionally, if any terms in H 0 do not commute with each other, then the time ordering of R will allow sequences of rotations, and thus nested commutators such as −[V j , [V i , δH]] will also be included in V err . Let d be the dimension of V err , and define a ddimensional vector G(t) as follows: This gives the accumulated error present in Eq. (3) as a function of time, divided by |δH| so that G(t) is independent of the actual strength of the noise. Here, the magnitude is given by | A| = A · A for any matrix A. We can picture G(t) as tracing out a curve in V err , with G(0) = 0. If G(t) returns to zero at some later time t pulse , then the rotation performed will be unaffected by error to first order in |δH|, and thus dynamically corrected gates correspond to closed d-dimensional curves by construction. Now consider the time derivative of G(t): Since d G/dt always has unit length, the distance along the d-dimensional curve parametrized by G(t) corresponds to the time t at that point. Thus, the tangent vector to the curve at any point is given by d G/dt, and the normal vector is proportional to the second time derivative of G. Using Eq. (1), we can express this in terms of R, H 0 , and δH as follows: The curvature κ 1 at any point along the curve can be calculated from this, and is given by In this representation of the cumulative error G(t) as a d-dimensional curve, the ideal evolution operator R(t) is a rotation which takes the initial tangent vector δH to the tangent vector d G/dt at any time t along the curve. Because R is a time-ordered exponential, it can be difficult to work with, since it cannot be analytically calculated for most choices of H 0 (t). However, as Eq. (9) shows, this does not prevent us from calculating the curvature, since the curvature at any point on a curve does not depend on the absolute orientation of the curve, but only on points in its own local neighborhood. This is reflected in Eq. (9) by the fact that κ 1 does not depend directly on R. Thus, for a given δH, the curvature provides a direct relationship between G(t) and the Hamiltonian H 0 , without the need to calculate R. If d > 2, the same reasoning demonstrates that the higher-order curvature coefficients are independent of R as well. These can be calculated using the Frenet-Serret equations: Here, the vectors e n (more precisely their normalized counterpartsê n ) define the Frenet-Serret frame at each point along the curve defined by G. The number of vectors is therefore equal to d, the dimension of V err . In the single-qubit case for example where d = 3, e 1 , e 2 and e 3 are the tangent, normal and binormal vectors, respectively. Knowing how these vectors evolve in time is equivalent to knowing how the curve moves through space. These vectors evolve in an interdependent way that is governed by the generalized curvatures: In the single-qubit case, κ 1 is the curvature of the threedimensional space curve, while κ 2 is its torsion. More generally, d − 1 curvatures are needed to characterize a curve in d dimensions. The higher-order derivatives of G(t) can be calculated in a similar fashion to Eq. (8) above, which will generate nested commutators. There will, however, be additional terms due to the time dependence of H 0 . Define the operator C which acts on a matrix V as follows: Then the derivatives of G are given by It is important to note that C does not commute with ∂/∂t, and thus there will be increasingly more terms in this expression for higher values of n. Eqs. (10)- (13) relate the generalized curvatures to the driving fields. While Eq. (13) can be difficult to work with in general, there are specific conditions under which it can be simplified. Specifically, suppose that the following anticommutation relation holds: {H 0 , δH} = 0. Note that if δH is proportional to a Pauli string, then one can always transform to a basis in which this relation holds [29]. For convenience, define Q = δH/|δH|, and assume that Q 2 = 1. Define the operator A n such that the Frenet-Serret basis vectors can be written as: Differentiating Eq. (14) yields the following: On the other hand, the Frenet-Serret equations read: We now show that the two terms in Eq. (15) directly correspond to the two terms in the Frenet-Serret equations. To do this we will show that [A n , Q] = 0 for n mod 4 = 0 or 1, and {A n , Q} = 0 for n mod 4 = 2 or 3. We will make use of the property that for any operators O 1 and O 2 such that O 1 commutes with Q and O 2 anticommutes with Q, the following holds: Consider the first Frenet-Serret vector (the tangent vector to the curve):ê Here A 1 = 1, and thus [A 1 , Q] = 0 holds trivially. Differentiatingê 1 , we obtain the second vector: from which we find that κ 1 A 2 = i{H 0 , A 1 }, and the relationship {A 2 , Q} = 0 follows. Similarly, we can consider the derivative of the nth vector given by Eq. (15). If A n commutes with Q, thenȦ n will also commute with Q, and {H 0 , A n } will anticommute with Q. Similarly, if A n anticommutes with Q, thenȦ n will anticommute with Q, and {H 0 , A n } will commute with Q. In either case, Eq. (17) implies that the two terms in Eq. (15) will be orthogonal to one another. Additionally, one of these two terms will be orthogonal toê n−1 . This means that terms in Eq. (15) can be matched to the terms in Eq. (16), and the curvatures and operators A n can be obtained. Specifically, the following recursion relation holds: The curvatures are obtained by taking the magnitude of this expression, since A n has unit magnitude. Eq. (20) thus provides a general mapping between control fields in the Hamiltonian and curvature coefficients of the curve.

III. DERIVING EXAMPLE PULSES
In this section, we apply the general geometrical formalism derived in Sec. II to a specific two-qubit Hamiltonian, and use it to derive pulses that implement dynamically corrected gates in the presence of quasi-static noise.
Consider the case in which the qubits are coupled by an Ising-type interaction, which creates an energy splitting between the |00 ↔ |01 and |10 ↔ |11 transitions. This type of interaction is common in solid-state qubit systems. For example, it arises in the context of superconducting transmon qubits from both direct capacitive coupling and also from resonator-mediated interactions [38,39]. It also applies to semiconductor spin qubits, both for exchange-based [40,41] and capacitive [42][43][44] inter-qubit couplings. Single-qubit or two-qubit gates can then be performed by driving only one qubit (see e.g., Refs. 25,39,45). The noiseless Hamiltonian in this case is where X i and Z i are the Pauli matrices on qubit i. Suppose that the leading source of noise δH is proportional to Z 2 , corresponding to slow fluctuations in the energy splittings of both qubits, as may arise, for example, from charge and field noise. Then using our geometrical formalism, we can represent the cumulative error as a curve in the 6-dimensional space V err spanned by we calculate the five curvature coefficients as follows: Every pulse Ω(t) corresponds to some 6-dimensional curve given by these curvature coefficients. It is important to realize that the converse is not true: every 6dimensional curve does not necessarily correspond to a pulse that can be created by the Hamiltonian given by Eq. (21). This is because there are in general five degrees of freedom corresponding to a 6-dimensional curve (corresponding to the five curvature coefficients), but in Eq. (21) we have constrained the form of the Hamiltonian such that it has only one independent driving field. Because the Hamiltonian given by Eq. (21) is block diagonal, the system can alternatively be treated as two independent 2 × 2 subsystems corresponding to the two blocks, with the constraint that the driving field Ω must be the same in each block. These subsystems can each be treated using the single-qubit formalism given in Ref. 29. This means that a closed 6-dimensional curve which satisfies Eq. (22) can be mapped to two separate 3dimensional closed curves with equal lengths and curvatures to each other, but different torsions given by E 1 and E 2 .
Consider the case where |E 1 − E 2 | t −1 pulse , where t pulse is a time scale corresponding to the total length of a pulse. Then since κ 3 is proportional to |E 1 − E 2 |, it will be small, and thus the entire curve will lie close to a 3-dimensional subspace of the 6-dimensional space. Then using the 3D formalism inside this subspace will correct curves to leading order in |E 1 − E 2 |t pulse . This happens because the 3D formalism is designed to correct against small Z errors, so if E 1 is close to E 2 , then H 0 can be treated as two copies of the same 2 × 2 Hamiltonian, and the Z 1 Z 2 term can be included with the error in each case. This means that no real 2-qubit gates can be performed using this method, since these pulses correct against the Z 1 Z 2 term used to create the interaction needed for 2-qubit gates in the first place. Thus, we will focus instead on the case where E 1 and E 2 differ significantly.
To demonstrate the utility of this formalism, we numerically derive a dynamically corrected gate that implements a Z rotation on one qubit while canceling the effects of noise on both qubits. In the representation of the error as two 3D curves, this will require generating two curves with differing constant torsion, with the equal curvatures as functions of arclength. There has been much work deriving methods of generating closed curves of constant torsion [46][47][48][49]; however, the analytic methods often used are not applicable here due to the additional constraints of generating two curves with equal curvatures. Thus, we approach the problem by numerically piecing together sets of helices to make two smooth closed curves, which corresponds to a series of square pulses canceling errors similar to Supcode [20]. Alternatively, in the 6-dimensional picture, higher dimensional generalizations of helices (curves with constant curvature coefficients) can be used. These will be five dimensional, since κ 5 is zero, due to dΩ/dt vanishing for a square pulse. Care must be taken for the step function transition between different values of Ω, as dΩ/dt becomes a delta function at these points. At these transitions, Ω(t) can be treated as a steep constant slope over a short interval of time , in the limit where → 0. In this case κ 1 through κ 4 will be finite, and thus will have no effect on the curve as → 0. κ 5 will have a delta function contribution, resulting in a rotation betweenê 5 andê 6 at the point along the curve corresponding to the transition. The exact angle of rotation φ corresponding to a step from Ω 1 to Ω 2 will be given by: We choose to consider curves which are n-fold rotationally symmetric, as these curves are easier to visualize and require fewer parameters to represent. To this end, we consider one curve segment where the Frenet-Serret frame at one endpoint is equal to the Frenet-Serret frame at the other endpoint after having undergone a relative 2kπ/n rotation, for any integer k coprime to n. Any curve segment such as this will produce a closed curve when repeated n times, provided that the displacement vector between endpoints lies within the plane of rotation. For our numerics, we consider 3-fold symmetric curves, which amounts to choosing a periodic pulse with period equal to one third of the total pulse length. Then we numerically adjust the parameters corresponding to the legnths and curvatures of the helices until the conditions on the displacement vector and the Frenet-Serret frames at the endpoints are met. Fig. 1 shows an example of a square pulse derived in this way, with E 2 = 2E 1 . In experiments, square pulses cannot be exactly performed, since pulse generators have bandwidth and pulse rise-time limitations. However, we can numerically search for a smooth pulse similar in shape to the square pulse already obtained. We do this by choosing a pulse shape of the following form: .
(24) The form is meant to approximate a Lorentzian pulse, except that it is periodic with period t p . The curves can then be numerically generated, and the parameters in Eq. (24) adjusted until closed curves are obtained. Because we use more parameters than the dimension of the space, solutions which produce closed curves are not necessarily unique. Uniqueness is not essential here since all we need is a solution providing dynamical decoupling, and the fact that there may be other solutions is not a problem. If finding solutions inside a desired parameter regime is difficult, a third Lorentzian can be added to the pulse to allow for more parameters. Adding more parameters should simplify the process of finding a possible solution as long as we do not insist on unique solutions; however the shape of the resulting pulse will be more complicated as more parameters are used. Using this method, we find the pulse and curves shown in Fig. 2. In Fig. 3, we plot the error curve in the 6-dimensional representation of the same pulse by projecting it into two 3D subspaces. We see that it retains its 3-fold rotational symmetry. In two of these dimensions, e 1 (0) and e 4 (0), the curve covers much more distance than in the other four.
By looking at each 2 × 2 block in the block-diagonal matrix individually, we can use the same method as with the single-qubit case to determine the gates performed. For symmetric pulses like the ones we generated, this will result in e 1 (t pulse ) = e 1 (0), and similarly for e 2 and e 3 , where these vectors belong to the 3-dimensional spaces corresponding to each block. Thus, these symmetric pulses can only perform identity operations in each of the 2×2 blocks. These identity operations can have a different relative sign, resulting in a Z 1 gate, as is the case with the pulse shown in Fig. 2. Removing the symmetry requirement will generate other gates. In order to demonstrate that these pulses do properly correct against error, we numerically calculate the infidelity of the gate using the full noisy Hamiltonian for noise strengths across several orders of magnitude. Here, we consider a single quasistatic Z 2 noise source, and vary its strength, plotting the resulting infidelity versus the noise strength in Fig. 4. We define infidelity as |U − R|, where U is the noisy pulse, and R is the ideal pulse (in this case Z 1 ). We see that the infidelity scales as the square of the noise strength, indicating that the gate cancels error to first order in |δH|. Thus, our proposed geometric scheme does correct for noise in the Ising 2-qubit gate operations.

IV. CONCLUSION
We have shown how to generalize the geometric formalism for producing dynamically corrected single-qubit gates to multiqubit systems. Like the single-qubit case, the cumulative first-order error can be represented as a curve in Euclidean space. Distance along the curve corresponds to the elapsed time from the beginning of the pulse, and the strength of the driving fields can be related to the curvature coefficients at each point on the curve. Critically, using these curvature coefficients circumvents the need to evaluate the time-ordered exponential of the Hamiltonian, which in general cannot be done except in very special cases. We presented equations which show how to calculate these curvature coefficients in terms of the time derivatives and commutators of the Hamiltonian H 0 and noise source δH. We used this formalism to numerically derive dynamically corrected gates for a two-qubit Hamiltonian that is relevant for Ising-coupled qubits, as relevant for superconducting transmons and singlet-triplet spin qubits, and demonstrated that the pulses cancel first-order error by computing the infidelity of the gate as a function of the noise strength.
While this formalism provides a good starting point for deriving dynamically corrected gates for multiqubit systems, there are still several challenges. In particular, it is usually the case that only a few terms in the Hamiltonian can be controlled dynamically. In this case, we need to restrict to the subset of curves that produce Hamiltonians of the desired form. This issue arises when the number of control fields is less than the number of generalized curvatures, in which case we need to find curves for which some of the generalized curvatures remain con-stant. Note that this is an issue even in the single-qubit case if only one control field is present. In this circumstance, one needs to restrict to closed curves of constant torsion. How to do this in general for multiple curvatures is an important open question for future work.