Coherent Multispin Exchange Coupling in a Quantum-Dot Spin Chain

Haifeng Qiao, Yadav P. Kandel, Kuangyin Deng , Saeed Fallahi, Geoffrey C. Gardner , Michael J. Manfra, Edwin Barnes, and John M. Nichol 1,† Department of Physics and Astronomy, University of Rochester, Rochester, New York 14627, USA Department of Physics, Virginia Tech, Blacksburg, Virginia 24061, USA Department of Physics and Astronomy, Purdue University, West Lafayette, Indiana 47907, USA Birck Nanotechnology Center, Purdue University, West Lafayette, Indiana 47907, USA School of Materials Engineering, Purdue University, West Lafayette, Indiana 47907, USA School of Electrical and Computer Engineering, Purdue University, West Lafayette, Indiana 47907, USA


I. INTRODUCTION
A unique and valuable feature of electron-spin qubits in quantum dots [1,2] is the voltage-controlled nearestneighbor Heisenberg exchange coupling. Heisenberg exchange coupling results from the interplay of the electronic confinement potential, the Coulomb interaction, and the antisymmetric nature of the electronic wave functions under particle exchange. On a basic level, exchange coupling enables two-qubit [1,[3][4][5] and three-qubit gates [6] for single-spin qubits. Exchange coupling also allows rapid and high-fidelity initialization and readout of pairs of spins. As a result, exchange coupling underlies the operation of electron-spin qubits consisting of two [7,8], three [9][10][11][12][13][14], or more [15][16][17] electrons. Superexchange [18,19] in multielectron systems and extended exchange-coupled spin chains can enable new forms of quantum information transfer [20,21] and the generation of many-body entangled states [22]. Recent experiments exploiting pulsed exchange coupling in spin chains point to the feasibility of these proposals [23]. Heisenberg spin chains are also predicted to generate nonequilibrium quantum phenomena [24,25].
Partly because of these exciting possibilities, independent and automated control of interdot tunnel couplings has been the focus of intense research in quantum-dot arrays [26][27][28][29][30][31]. However, generating multiple independent, nonzero, and coherent exchange couplings in quantum-dot arrays is challenging for many reasons. First, the standard procedure to measure tunnel couplings involves detuning pairs of dots away from the symmetric idling point [32], making it difficult to calibrate exchange couplings under actual experimental conditions. Second, as discussed further below, the nonlinear and nonlocal dependence of the exchange couplings on the confinement gate voltages poses a significant challenge. Third, multiple nonzero exchange couplings generate complicated spectra that do not permit easy measurement and iterative tuning of individual exchange frequencies.
In this work, we demonstrate coherent multispin exchange coupling in a GaAs quadruple quantum dot. We show that the nonlinear and nonlocal dependence of exchange couplings on confinement gate voltages results, in large part, from electronic wave-function shifts during exchange pulses. We model our data using the Heitler-London (HL) expression for exchange coupling between two spins [33], assuming that the barrier-gate pulses used to induce exchange coupling primarily shift the locations of the electrons. The model parameters we use change slightly depending on the number of spins involved, suggesting that additional effects beyond wave-function shifts, including perhaps the quantum-dot potential depths and widths, are also important.
The parameters we extract by fitting our data to the Heitler-London model agree well with electrostatic simulations of the confinement potential of our device. We also show that a simpler, exponential model also fits our data well and can be used to predict gate voltages for independent control of exchange couplings. We demonstrate two-, three-, and four-spin exchange coupling in our four-dot device. These results are applicable to Si qubits, which feature reduced hyperfine coupling and longer electronspin coherence compared to GaAs spin qubits. Our results are also applicable to longer arrays of spin qubits, an encouraging prospect for quantum information processing and the exploration of Heisenberg spin chain physics.

II. DEVICE
We use a quadruple quantum dot in a GaAs/AlGaAs heterostructure with overlapping gates [ Fig. 1(a)] [34][35][36]. Two extra quantum dots placed above the main-dot array serve as charge sensors and are configured for rf reflectometry [37,38]. We tune the confinement potential using "virtual gates" [26,30,[39][40][41] such that each dot contains only one electron. We define virtual-plunger-gate voltages P 1 , P 2 , P 3 , and P 4 as linear combinations of the physicalplunger-gate voltages (p 1 , p 2 , p 3 , p 4 ) such that changes to P i are proportional to changes in the electrochemical potential of dot i. We also define virtual-barrier-gate voltages B 1 , B 2 , and B 3 as the voltage applied to the corresponding physical barrier (b 1 , b 2 , and b 3 ) together with a linear combination of physical plunger voltages chosen such that the chemical potentials of the dots are unchanged by the barrier pulse. In particular, barrier-gate pulses therefore involve voltages applied to physical barrier gates in addition to compensation pulses applied to all four physical plunger gates. The virtual gates are related to the actual voltages via the measured capacitance matrix A through G ¼ Ag, where G ¼ ½P 1 ; P 2 ; P 3 ; P 4 ; B 1 ; B 2 ; B 3 T and g ¼ ½p 1 ; p 2 ; p 3 ; p 4 ; b 1 ; b 2 ; b 3 T . In the following, we use the term "virtual-gate voltages" to mean pulses away from the idling tuning of the device, which is the symmetric operating point [42,43] of the four-dot array with one electron in each dot.
For initialization and readout, we configure the four-spin array into two pairs. We refer to spins 1 and 2 as the left pair, and spins 3 and 4 as the right pair. We initialize the array in the product state j↓↑↓↑i via adiabatic separation of singlets in the hyperfine gradient [7,23]. Here, the arrows indicate the spin states of all four spins. Alternatively, we can prepare a polarized triplet state j↑↑i in either pair by exchanging electrons with the reservoirs while each dot contains a single electron [44]. We measure both pairs via Pauli spin blockade [7,23,45], together with a shelving mechanism [46], to enhance the readout fidelity.
The spin-state Hamiltonian of the quadruple dot is Here, J i is the exchange coupling strength (with units of frequency) between dots i and i þ 1, σ i ¼ ½σ x i ; σ y i ; σ z i is the Pauli vector describing the components of spin i, and h is Planck's constant. Note that B z i is the z component of the magnetic field experienced by each spin, and it includes both a large 0.5-T external magnetic field and the smaller hyperfine field. The quantization axis (z direction) is defined by the external magnetic field direction. The x and y components of the hyperfine field are neglected in this Hamiltonian since their sizes are negligible compared to the external magnetic field. In addition, B z i has units of frequency.

III. EFFECT OF POSITION SHIFTS ON EXCHANGE COUPLING
A single nonzero exchange coupling J i is easily tuned by adjusting the voltage on the relevant barrier B i [42,43]. When we extend the interaction to more than two spins by pulsing another virtual barrier gate B j , however, the original exchange coupling J i is strongly affected. For example, a large pulse to B 3 nominally induces a nonzero J 3 . But adding an additional pulse to B 2 during the evolution rapidly and nonlinearly reduces J 3 before eventually turning on a J 2 [ Fig. 1(b)]. In fact, J 3 reduces to nearly zero before J 2 turns on.
We now show that this effect primarily results from lateral shifts of the quantum dots during a barrier-gate pulse. We have self-consistently calculated the electron density and potential of our device in COMSOL using the Thomas-Fermi approximation [47]. Our simulation replicates the behavior of the physical device with high fidelity. For example, the predicted incremental voltage on p 4 required to transition from one to two electrons in dot 4 is about 12.5 mV. This value corresponds quite closely to the behavior of the device. The same transition in the physical device requires a change of 12 mV on p 4 . Figure 2(a) shows the potential associated with dots 2 and 3 tuned to single occupancy as a function of the barrier voltage B 2 . Our simulations include compensation pulses on the plunger gates to match our use of virtual gates in the actual experiment. The dots clearly move toward each other as B 2 increases.
To obtain more detailed information about this process, we fit our two-dimensional simulated potentials to an equation of the form [33] Vðx; yÞ Here, x and y are coordinates in the plane of the twodimensional electron gas, and V 0 , a, l x , and l y characterize the potential wells of the dots. Double Gaussians of this type are commonly used to model double dots, but usually a separate barrier term is included, as in Ref. [33]. Our simulated potential is shallow enough that a separate barrier term is not required to reproduce the potential we simulate. By fitting our simulated potentials to Eq. (2), we find how the parameters a, l x , and V 0 vary with B 2 [ Fig. 2(b)]. The distances between dots 2 and 3 (2a 2 ) and dots 1 and 2 (2a 1 ) change approximately linearly during the barrier pulse, but in opposite directions, because dot 2 moves closer to dot 3 but farther from dot 1. Based on the simulations, we calculate that a 2 changes by about −0.3 μmV −1 , and a 1 changes by about 0.19 μmV −1 . Other parameters of the confinement potential change as well during the barrier pulse.
Reference [33] computes the exchange coupling between two quantum dots in a potential of the form Eq. (2) in the HL framework. At zero magnetic field, the result is . Note that m is the electron effective mass, ϵ 0 is the permittivity of free space, ϵ is the dielectric constant of the material, and I 0 is the zeroth-order modified Bessel function. In writing this equation, we have assumed that the minima of the double-dot potential occur at x ¼ AEa. We have also ignored the magnetic-field-dependent terms because, for the magnetic field used here (0.5 T), the effective magnetic confinement is still weaker than the electrostatic confinement.
To determine if lateral position shifts can indeed explain our data, we experimentally measure how B j affects J 1 in our device (Fig. 3). We then fit our data to Eq. (3). To parametrize the effect of the B j , we allow for where the α j are fit parameters. We also fit for V 0 , but we constrain l x ¼ l y ¼ 100 nm, which is approximately the value we obtain from our simulations. The fitted values are V 0 ¼ 11.4 meV,  This level of agreement is reasonable, considering that the simulation calculates the semiclassical self-consistent potential associated with the electron density in the quantum dots. Moreover, the potential of each dot is not exactly Gaussian, making it challenging to extract an exact potential depth. Our simulated values of the quantum-dot position shifts (−0.3 μmV −1 and 0.19 μmV −1 ) described above agree well with the fitted values of α 1 and α 2 . This agreement supports our hypothesis that lateral position shifts cause the observed trends in our data.

IV. MODELING THE DEPENDENCE OF EXCHANGE COUPLING ON ALL BARRIER GATES
The nonlinear and nonlocal dependence of exchange couplings on the barrier-gate voltages, which results from position shifts of the quantum dots, poses a challenge to implementing simultaneous exchange coupling between all dots in an extended array. Previous work has investigated how to adjust multiple interdot tunnel couplings iteratively [26,30]. Here, we discuss two different models that allow us to determine the virtual-gate voltages given a set of target exchange couplings. In contrast to previous iterative approaches, our approach generates a predictive model. As discussed further below, predictive tuning of exchange couplings in extended spin chains is especially helpful when multiple exchange couplings are present because the observed spin oscillation frequencies do not correspond to the bare two-spin Heisenberg couplings. We also use our model to control coherent exchange coupling instead of incoherent electron tunneling. Finally, our approach has the advantage that it enables calibrating exchange couplings at the symmetric operating point, where tunnel couplings cannot be easily determined.
Our general approach is to measure how all of the J i depend on the B j and then to fit the parameters of a nonlinear model to the data. Using these fit parameters, we create a model that allows us to generate a set of virtualgate voltages GðjÞ, for a set of target exchange coupling values j ¼ ½j 1 ; j 2 ; j 3 , where the j i are the desired exchange coupling values. We compute the actual gate voltages using the transformation described above. We validate this model by inducing exchange coupling between two, three, and four spins and compare our observations with simulations; we find good agreement.
To calibrate the models, we begin by inducing one strong exchange coupling J i ≫ J j≠i and measuring the effect of the B j on that exchange coupling. For example, to measure how J 1 depends on B j , we initialize the array, as discussed above, with B j ¼ 0. Then, we pulse B 1 from 25 to 65 mV, and we record exchange oscillations at each pulse height. Setting B 1 ¼ 60 mV, which yields a large but stillmeasurable J 1 ∼ 200 MHz, and B 3 ¼ 0 mV, we then pulse B 2 from 0 to 35 mV, and we record exchange oscillations. Setting B 1 ¼ 60 mV and B 2 ¼ 0 mV, we sweep B 3 from 0 to 40 mV and again record exchange oscillations. The pulses on B 2 are not sufficient to induce substantial J 2 due to the large pulse height on B 1 . The pulses on B 3 induce substantial J 3 , but J 1 is not affected by the next-nearestneighbor exchange coupling. We extract the oscillation frequencies through a fast Fourier transform (FFT) of the data (Fig. 3). We repeat this process for the other J i (see Supplemental Material [48]).
We fit the resulting J i vs B j data to a set of equations related to Eq. (3), of the form J i ¼ J HL ðV i 0 ; a i Þ, where and Here, the V i 0 and the α ij are fit parameters. As discussed above, we constrain l x and l y to be the values found from simulations. These equations model our data quite well (Fig. 3), and the parameters we extract from the fits agree reasonably well with our simulations. Values of V i 0 range from 6.4 to 11.4 meV, and values of α ii range from −0.40 to −0.43 μmV −1 . Once the model is calibrated and the parameters found, we choose target exchange coupling values j. We then numerically solve the set of equations The black data points in each panel are obtained from the fast Fourier transform of a data set similar to Fig. 1(b). In (a)-(c), the dark blue line is the fit to the exponential model, and the light blue line is the fit to the HL model. Panels (d)-(f) show the difference between the fits and the data for the two models.
j i ¼ JðV i 0 ; a i Þ for the interdot separations a i , and we invert Eqs. (4)-(6) to find the desired barrier-gate voltages.
While this model (the "HL model") originates from a microscopic theory, the exchange coupling is a highly nonlinear function of the potential parameters, and some a priori knowledge of the quantum-dot confinement potential is desirable. Using the HL model, in practice, also involves numerically solving nonlinear equations, which can be susceptible to errors. An alternative, more robust, model for the dependence of the J i on the B j is motivated by the realization that the part of the expression for J i in Eq. (3) that is most sensitive to the interdot separation is the factor as before, we expect the J i to depend approximately exponentially on the B j . Thus, we introduce the following "exponential model" for our data: Here, δ ij , β i , γ i are fit parameters. Empirically, the fit parameters γ i are required for two reasons. As discussed above, the exchange couplings are not pure exponential functions of the barrier gates. Second, the hyperfine gradient can increase the measured oscillation frequency above the bare exchange frequency. Including γ i in the fit allows us to accommodate these deviations from pure exponential behavior. The exponential model also matches our data quite well (Fig. 3). Typical values of β i and γ i are on the order of 10 MHz, and values of δ ii range from 93 to 114 V −1 . From Eq. (7), we expect that δ ij ≈ −4aα ij =l 2 0 .
Using a ¼ 100 nm, l 0 ¼ 32 nm, and taking a typical value of α ii ¼ −0.4 μmV −1 , we expect δ ii ≈ 156 V −1 , which agrees reasonably well with our fitted values. We have also conducted measurements to confirm that the δ ij do not depend significantly on the barrier-gate voltages, supporting the form of the exponential model (see Supplemental Material [48]). With the model parameters in hand, we choose a set of target exchange coupling values j. Setting J ¼ j, we invert Eqs.
For both the exponential and HL models, we require that the virtual-plunger gates remain fixed at the symmetric operating point, and we transform the virtual-gate voltages to physical gate voltages using the capacitance matrix through g ¼ A −1 G, as discussed above.

V. MODEL VALIDATION
In practice, we prefer to use the exponential model because it features a robust inversion process and fits our data as well as the HL model (Fig. 3). (See the Supplemental Material [48] for a comparison of the gate voltages generated by these two models.) We first validate our approach by sweeping j 3 , the target exchange coupling between spins 3 and 4, linearly from 0 to 200 MHz [ Fig. 4(a)]. The observed oscillation frequency matches our expectation. We also perform the same test on qubits 1 and 2, and qubits 2 and 3, and observe good agreement between the measured and target frequencies.
In our device, we estimate that the residual exchange coupling at the idling point with no barrier-gate pulse is not zero but rather a few megahertz for each pair. In generating the data of Fig. 4(a), which features a swept j 3 , we did not simultaneously require that j 1 ¼ j 2 ¼ 0. If we had required j 1 ¼ j 2 ¼ 0, our model would generate large negative values of B 1 and B 2 . Instead, for this data set, we fixed B 1 ¼ B 2 ¼ 0 since the weak residual exchange coupling does not significantly affect these data.
Next, we induce three-spin exchange oscillations between spins 2-4. The array is initialized in the j↑↑↓↑i state. Preparing j↑↑i on the left side ensures that spins 2-4 remain in the s z ¼ þ 1 2 subspace, regardless of the sign of the local magnetic field gradient. We fix B 1 ¼ 0, j 2 ¼ 70 MHz, and we sweep j 3 linearly from 10 to 150 MHz; we measure oscillations on the right pair. We compare our data to simulated predictions (see Supplemental Material [48] for details on the simulation). The observed frequencies closely match our expectation, confirming that we can correctly set the target frequencies [ Figs. 4(b)-4(c)]. Note the presence of three distinct frequencies in the measured spectrum of Fig. 4(b). These frequencies are not the bare Heisenberg couplings. Instead, they result from the different energy splittings between the singletlike and tripletlike states of three spins [9]. The theoretically predicted [9] low-lying energy splittings of three exchange-coupled spins are shown in red in Fig. 4(b), and they overlap nicely with our measurements [9]. This theoretical prediction assumes zero magnetic gradient between the dots. The presence of a hyperfine gradient in our device causes the experimental data to slightly deviate from the theoretical predictions [ Fig. 4(c)]. However, this deviation is relatively insignificant for exchange strengths above 10 MHz.
The model parameters we extract from the fits alone do not suffice to accurately generate the target three-spin exchange frequencies for both the exponential and HL models. We empirically find that the voltages GðjÞ generate actual exchange frequencies J < j when two or more of the j i are nonzero (see Supplemental Material [48]). To overcome this challenge, we make minor adjustments to the fitted model parameters and compare the observed threespin exchange frequencies with simulations. We repeat this process for a few iterations until the experimental data match closely the simulated predictions, usually within about 10 MHz. In the exponential model, we normally need to modify the values of β i and γ i , while the values of δ ij can remain mostly unchanged. For the HL model, we usually need to slightly increase the confinement strength V 0 . See the Supplemental Material [48] for a comparison of the fitted and adjusted parameters and gate voltages.
The need for a modified parameter set for three-spin exchange coupling may originate because of the following. First, we calibrate our model when only one of the J i is large. However, three-spin exchange coupling requires multiple large J i , for which we need several large, simultaneous barrier-gate pulses. Second, our assumption that only linear position shifts induce exchange coupling likely breaks down at large gate voltages. Indeed, Fig. 2(b) shows that both the characteristic size and confinement energy change during barrier-gate pulses. Third, generating three-spin exchange coupling involves simultaneous large voltage pulses on several plunger and barrier gates, and any errors in our measured capacitance matrix will cause errors in the exchange couplings.
Despite the need for an additional set of parameters for three-spin exchange coupling, the data in Fig. 4 show that our model can still be used in this regime. Indeed, this additional set of parameters also suffices to induce fourspin exchange coupling because the exchange coupling has vanishing dependence on the next-nearest-neighbor barrier gates [ Figs. 3(c)]. To calibrate the two sets of three-spin parameters needed for four-spin coupling, we first tune the model for exchange coupling between spins 2-4, which yields precise values of δ 22 , δ 23 , δ 32 , and δ 33 , together with β 2 , β 3 , γ 2 , and γ 3 . To calibrate the model for exchange coupling between spins 1-3, we leave these parameters fixed and tune δ 11 , δ 12 , δ 21 , β 1 , and γ 1 . With the parameters tuned in this way, we induce simultaneous exchange coupling between all four spins in the array by initializing the array in the j↓↑↓↑i state, and we sweep all three target exchange frequencies linearly from 10 to 40 MHz. The oscillations are measured on both the left and the right pairs

VI. DISCUSSION
A significant source of error in our model is the fluctuating nuclear hyperfine field, which is of order 10 MHz in our device (see Supplemental Material [48]). Empirically, our model can generate exchange couplings that are accurate to about 10 MHz, suggesting that the hyperfine fields are a significant source of error. Our model becomes increasingly inaccurate when exchange frequencies approach the hyperfine field values.
Although different sets of parameters are needed for twoand three-spin exchange coupling, the range of validity of the model is sufficient to accurately predict voltages for a wide range of exchange coupling values in either of these cases. To improve on this model, additional effects beyond lateral position shifts probably need to be considered. As discussed above, it is likely that other parameters of the dots are changing or that the capacitance matrix may require further refinement. More sophisticated modeling of the electrostatic potential, using the actual computed potential, as opposed to a Gaussian approximation, and a more accurate calculation of exchange couplings beyond the two-spin HL framework [49] will likely also help to improve the model.
In the future, it seems likely that fabricating devices with extremely narrow barrier gates may help to reduce the voltages needed for simultaneous exchange coupling between multiple electrons in quantum-dot spin chains, although large pulses may still be required to adjust the exchange coupling from very small to very large values. It is less clear that the challenge of position shifts can be avoided with narrow barrier gates. In any case, we expect that electrostatic simulations can be used effectively to guide the design of future quantum-dot spin chains. We also expect that the use of electrostatic simulations to guide quantum-dot spin-qubit tuning and operation will become a valuable tool as quantum-dot devices increase in complexity.
We expect that our technique can apply to different types of quantum dots. A similar exponential dependence of exchange has been observed in Si/SiGe quantum dots [5,42] and GaAs devices [43]. The form of Eq. (3) suggests that lateral wave-function shifts may be a primary cause of this behavior. A similar exponential model has also been used to describe the behavior of tunnel couplings in a GaAs device with a nonoverlapping, open-gate architecture [31]. We also expect that this model could be applied to Si=SiO 2 devices, where voltage-controlled tunnel coupling between quantum dots has now been demonstrated [50].
The possibility of creating quantum-dot Heisenberg spin chains opens up a wide array of new phenomena to explore. As we have shown in previous work, pulsed exchange couplings between spins in quantum-dot spin chains can enable new forms of quantum information transfer [23,51]. When multiple exchange couplings can be controlled simultaneously, coherent information transfer via antiferromagnetic spin chains becomes a possibility [20][21][22]. Various techniques related to adiabatic quantum information processing, including adiabatic state transfer [52] and adiabatic state and gate teleportation [53], seem within reach. Single-pulse three-qubit gates [6] and the creation of long-range entangled states also become possible. Beyond quantum computing applications, Heisenberg spin chains are essential elements of models for quantum magnetism, and they underlie important nonequilibrium many-body phenomena of great interest, such as effects related to many-body localization [24,54] and time crystals [25].
As an example of the utility of our approach, we return to the experiment described in Fig. 4, where we generate coherent exchange coupling between three spins. These data represent an experiment in which an initial state j↑↑↓↑i evolves under simultaneous exchange coupling between dots 2 and 3 and dots 3 and 4. When J 2 ¼ J 3 , evolution for a total time ð2=9J 2 Þ yields the final state jψi ¼ ð1= ffiffi ffi 3 p Þj↑i ⊗ ðj↑↑↓i þ e i2π=3 j↑↓↑i þ j↓↑↑iÞ. Here, jψi is equivalent to a W state on qubits 2-4, up to singlequbit rotations. For the case shown in Fig. 4(d), MHz. Thus, in a total evolution time of approximately 3.2 ns, we expect to generate a three-qubit entangled state. Based on Monte Carlo simulations, which integrate the Schrödinger equation including initialization errors, charge noise, and hyperfine noise [ Fig. 4(d)] (see Supplemental Material [48]), we expect that we can prepare qubits 2-4 in the state ð1= ffiffi ffi 3 p Þðj↑↑↓i þ e i2π=3 j↑↓↑iþ j↓↑↑iÞ with a fidelity of about 0.60. The dominant source of infidelity is state preparation error associated with thermal population of excited states and the fluctuating hyperfine gradient. Neglecting preparation errors, we expect that the fidelity can exceed 0.9, largely because the total evolution time is much shorter than the relevant dephasing times in this system. The creation of extended W states in larger Heisenberg spin chains is discussed further in Ref. [22].
In addition to generating a many-body entangled state, the case where J 2 ¼ J 3 can also generate a remote entangling operation, as we now discuss. An initial state j↑↑↑↓i evolves after a time ð2=3J 2 Þ to jψi ¼ ð ffiffi ffi 3 p =2Þj↑i ⊗ ðj↓↑↑i þ ðe iπ=2 = ffiffi ffi 3 p Þj↑↑↓iÞ. Here, jψi features entanglement only between qubits 2 and 4. Moreover, the initial state of qubit 3 is unchanged. In fact, it is easily verified that, in the absence of magnetic gradients, this operation is equivalent to an ðS 24 Þ 2=3 operation, where S 24 indicates a SWAP gate between qubits 2 and 4. Compare this to the common entangling operation ðS ij Þ 1=2 between two qubits i and j. Based on our simulations, we expect that with perfect state preparation, we could generate this threequbit state with a fidelity of about 0.86 in a total time of about 9.5 ns. The generation of remote exchange in extended exchange-coupled spin chains is further discussed in Ref. [55].

VII. CONCLUSION
In summary, we have demonstrated simultaneous coherent exchange coupling between two, three, and four spins in a quadruple quantum dot. We have also shown that lateral position shifts of the quantum dots during barrier pulses present a significant hurdle to implementing simultaneous exchange coupling between multiple electron spins. Using a phenomenological model based on a microscopic theory, we can predict the virtual-gate voltages required to generate a set of target exchange couplings. Our method is also scalable and applicable to other types of quantum-dot spin qubits, such as silicon qubits, which offer the possibility of isotopic purification and substantially reduced nuclear hyperfine fields, leading to longer electronspin coherence times. This method enables us to generate a four-site Heinseberg spin chain, which is an exciting prospect for the exploration of the physics associated with interacting spin chains.
Recently, we became aware of a related work demonstrating independent control of tunnel couplings in a quantum-dot array [31].