A natural heavy-hole flopping mode qubit in germanium

Flopping mode qubits in double quantum dots (DQDs) allow for coherent spin-photon hybridization and fast qubit gates when coupled to either an alternating external or a quantized cavity electric field. To achieve this, however, electronic systems rely on synthetic spin-orbit interaction (SOI) by means of a magnetic field gradient as a coupling mechanism. Here we theoretically show that this challenging experimental setup can be avoided in heavy-hole (HH) systems in germanium (Ge) by utilizing the sizeable cubic Rashba SOI. We argue that the resulting natural flopping mode qubit possesses highly tunable spin coupling strengths that allow for one- and two-qubit gate times in the nanosecond range when the system is designed to function in an optimal operation mode which we quantify.


I. INTRODUCTION
The potential of implementing and manipulating physical qubits in semiconductor systems was first recognized in the late twentieth century [1][2][3][4][5] and has since been demonstrated in countless experiments [6][7][8]. In particular, strongly confined HH states in nano-scale quantum systems have been investigated theoretically [9][10][11][12][13][14][15] and realized in recent years in various architectures ranging from one-dimensional hut-and nanowires to twodimensional hole gases in planar heterostructures [16][17][18][19][20][21][22][23]. Long coherence times due to the weak interaction with the host atomic nuclear spin bath and the absence of a valley degree of freedom in the valence band promote the binary pseudo-spin of such HH states to prime candidates for a reliable qubit [24].
In this paper we propose a qubit built from the spin of the bonding HH state in a planar DQD. Compared to single quantum dot (QD) systems [25], such so-called flopping mode qubits possess a large dipole coupling to an applied alternating or quantized cavity electric field, allowing for coupling strengths beyond the decoherence rate [26][27][28]. An alternative qubit-cavity construction uses multi-electron exchange-only qubits [29]. Among the most promising candidates for a platform for floppingmode spin qubits are electrons in Silicon and Carbon DQDs, which have seen detailed studies regarding their performance and decoherence properties [30][31][32][33][34][35][36]. In these systems a magnetic field gradient perpendicular to the QDQ axis is applied to achieve a coupling between bonding and anti-bonding states of different spin. The magnetic field gradient enables the electron spins to distinguish between different positions in space and may therefore be seen as synthetic SOI. While coherent spin rotations via a cavity field could be achieved within this setup, experimental realizations rely on micro-magnets and it is yet unclear how such systems can be scaled to include a large number of qubits. On contrast, we show that HHs in Ge can form a natural flopping mode qubit, The feature that allows us to achieve a natural HH flopping mode qubit is the cubic Rashba SOI. This type of SOI is characteristic for valence band states and stems from the electric field produced by the atomic nuclei experienced by the HHs in a quantum well which introduces structural inversion asymmetry along the growth direction. Starting from the semi-microscopic Hamiltonian of a HH in a DQD subject to the cubic Rashba SOI and an out-of-plane magnetic field, we derive the ground state Hamiltonian describing the one particle states in the left and right dot, and obtain an explicit expression for the tunnel matrix element. Additionally, we find terms describing spin-flip tunneling and intra-dot spin-flips which are induced by the SOI and the effect of excited orbitals. These terms allow for spin couplings in the bonding state when the system is coupled to a classical alternating or a quantized cavity electric field. We investigate oneand two-qubit gate times and find operation times in the nanosecond range. The performance of the device depends crucially on system parameters such as the applied magnetic field, the inter-dot distance, the strength of the Rashba SOI and the dot detuning, and may thus be optimized via quantum engineering.
The remainder of this paper is structured as follows: In Sec. II we introduce the HH system and derive the DQD Hamiltonian in the orbital ground state. Building on these results, we take into account the effect of excited orbitals in a perturbative fashion in Sec. III. In Sec. IV we change into the basis of bonding and anti-bonding states and investigate the coupling of these states via electric fields. We proceed to look at special parameter cases and derive effective spin-photon couplings allowing for spin rotations in the bonding state in Sec. V. Finally, Sec. VI provides a conclusion and an outlook on possible future research concerning natural flopping mode qubits. On top of the potential surface, we show a schematic of the flopping mode qubit, i.e., a single HH spin in the bonding state subject to an out-of-plane magnetic field B.

II. LOW-ENERGY HEAVY-HOLE STATES
We consider a semiconductor heterostructure with a Ge quantum well such that the HHs in the material are subject to strong confinement along the out-of-plane (z) axis. To model an in-plane DQD we introduce a quartic, locally harmonic confining potential V forming a double quantum well with dot separation 2a (Fig. 1). Moreover, we allow for a static out-of-plane magnetic field B and take into account the cubic Rashba SOI such that the system may be described by the Hamiltonian H = H 0 + H R [9][10][11], where π = p + eA is the canonical momentum, m the in-plane HH mass, g > 0 the out-of-plane HH g-factor and σ z the Pauli matrix along the quantization axis. The cubic Rashba term H R features the spin ladder operators σ ± = (σ x ± iσ y )/2 with in-plane Pauli matrices σ x/y , π ± = π x ± iπ y and the spin-orbit parameter λ R = 3γ s α R E z /2m 0 Ξ, where γ s = 5.11 is the Luttinger parameter in spherical approximation, α R the Rashba coefficient, E z the average electric field induced by the structural inversion asymmetry due to out-of-plane confinement, m 0 the bare electron mass and Ξ the HH-light hole splitting. There exists another cubic Rashba term which turns out to be far less sizeable and may be neglected completely for the case of the valence band states in Ge where the spherical approximation applies [37]. The time-independent Schrödinger equation with the Hamiltonian (1a) may be solved exactly close to the QD centres x = ±a, where the quartic potential (1b) becomes locally harmonic. In locally symmetric gauge, A L/R = B(−y, x ± a, 0)/2, we obtain Fock-Darwin states ψ FD nl (Appendix A) shifted to x = −a (x = a) for the left (right) dot with energies E nl = l ω L + (n + 1) ω, where ω = ω 2 0 + ω 2 L and ω L = eB/2m denotes the Larmor frequency. In what follows we work in the most general gauge for confinement and constant magnetic field along z, where r, c x , c y are arbitrary real constants, and we introduce the factor /el B including the magnetic length l B = /eB such that r, c x , c y are dimensionless. When transforming the left and right dot wave functions into this gauge, they acquire a magnetic phase and read Since there is a finite overlap between the single QD ground states, we orthogonalize to obtain the Wannier states |L/R , where the normalization factor is given by N = 1 − 2γS + γ 2 −1 and γ = 1 − √ 1 − S 2 /S [38,39]. Defining the Wannier-spin states |L/R, s ≡ |L/R |s with pseudo-spin quantum number s ∈ {↑, ↓} and the Pauli operators for position (τ ) and spin (σ), we obtain for the HH Hamiltonian H in the orbital ground state basis {|L, ↑ , |L, ↓ , |R, ↑ , |R, ↓ }, where t c = 3N γ 4 The energy detuning between the two dots is applied via gate voltages in experiments and we add it as a phenomenological control parameter. On a more microscopic level there arises a spin-conserving tunnel matrix element t c due to the DQD potential. On the other hand, there is a spin-flip tunneling termη due to the Rashba SOI. This term plays an important role in the ground state spin coupling in the bonding state discussed in the following. As a consequence of the U(1) gauge invariance of quantum electrodynamics all observable quantities obtained are gauge-invariant, i.e., they do not depend on the constants r, c x and c y . Additional constant gauge phases, which only affect the gauge field Λ but do not change the vector potential under the transformation A → A + ∇Λ, can be undone by choosing an appropriate representative from the ray of states in the Hilbert space for the wave function in each dot. Since the dependence of the tunnel coupling t c on the magnetic field exceeds a pure phase shift [40] (Fig. 5), our results go beyond the standard interpretation of Peierls substitution prescription [41], which always ought to be treated with the appropriate care [42][43][44]. Indeed, the validity of the approximation stating that the tunnel elements are only changed by a phase factor upon the application of an external magnetic field relies on strongly localized Wannier functions such that the area A ψ they enclose in the x-y-plane is much smaller than the area of a flux quantum, A B = π/eB. The DQD potential discussed in the present work yields local eigenfunctions that do not satisfy this assumption. To be specific, the mean area of the Fock-Darwin ground state ψ FD 00 is estimated to be A ψ = π ψ FD 00 |x 2 |ψ FD 00 ∼ π /mω, and thus

III. EFFECT OF HIGHER ORBITAL STATES
To include the effects of higher orbital states we consider the SOI as a small perturbation and work with the single-dot ground states within first-order perturbation theory, where ω Z = gµ B B and ω ± = ω ± ω L . Note that the normalization is omitted as it differs from unity only in second order in the small perturbation parameter ξ R = λ R ( mω 0 ) 3/2 /(3 ω − + ω Z ). The perturbative approach is valid as long as |ξ R | 1, which holds for all relevant magnetic field strengths, B 100 T for α R E z = 10 −11 eVm and Ξ = 100 meV. Due to the cubic SOI the altered ground states each contain a component with opposite spin in the third excited orbital level. As in Sec. II the states in the left and right dot are obtained by shifting the position coordinate to the dot centres x → x ± a and gauge transforming back to common gauge, Eq. (2), We find the overlaps S ⇑/⇓ = ψL S ⇓⇑ = −S ⇑⇓ , with S as given in Eq. (4). Neglecting the overlaps containing powers of ξ R , we obtain the Wannierspin states, which are orthogonal within our approximation, Here, the index σ = ⇑, ⇓ labels the constituents of a Kramers pair for each dot, which we refer to as a pseudospin doublet. We proceed to project the HH DQD Hamiltonian H = H 0 + H R onto the Wannier basis {|L, ⇑ , |L, ⇓ , |R, ⇑ , |R, ⇓ }. To facilitate the integrals appearing in the calculations, we use the behaviour of the unperturbed left and right dot states including a magnetic phase in symmetric gauge under (partial) parity transformations, The spin-conserving matrix elements and the out-ofplane Zeeman term receive corrections only at secondorder in ξ R which we neglect. In contrast, non spinconserving terms appear in first order in the SOI, and as a result the Hamiltonian features corrections to the spin-flip tunneling termη as well as an extra intra-dot spin-flip term χ, where η =η + δη 1 + δη 2 + g 2 µ B b,η, and t c as in Eq. (6) and We make use of the abbreviationsN = N (1 − γ 2 ) and h 0 = H 0 − gµ B Bσ z /2, and the wavefunctions ψ d nl are the unperturbed shifted Fock-Darwin states including a magnetic phase from transforming back to symmetric gauge. The matrix elements appearing in δη 1 are straightforward to calculate but lengthy; they are shown in Appendix B. Analysing the form of the novel terms, we find that the quantity χ acts as an intrinsic magnetic field gradient of strength 2χ. Note that χ(−B) = −χ(B) and hence χ(B = 0) = 0. While the term arises in a perturbative approach in the SOI, its magnitude may be sizeable as it does not depend on the overlap S between the left and right dot states. At B = 1 T one has 2χ 0.6 µeV for a HH-light-hole splitting Ξ = 100 meV and a Rashba coefficient α R E z = 10 −11 eVm, which is of the same order of magnitude as the magnetic field gradients applied in electronic systems, ∆B x ∼ 1.6 µeV [30,31]. Finally, the corrections to the spin-flip tunneling matrix elementη are threefold: δη 1 arises due to the DQD potential and b is due to the out-of-plane magnetic field in our approximation of zero overlap between different pseudo-spins. On the other hand, δη 2 is due to a non-vanishing overlap between states of different pseudo-spin, and the only instance where our approximation breaks down. The δη 2 contribution can be treated within a different approximative approach as detailed in Appendix C 2. A comparison of the matrix elements as a function of the magnetic field is shown in Fig. 5. In the figure the analytical results are complemented by numerical data which takes into account all overlaps between the dot states and also second-order terms in the spin-orbit parameter ξ R (Appendix C 1.). We find excellent agreement, a posteriori justifying the approximations made.

A. Bonding-anti-bonding basis
It is useful to diagonalize the effective DQD Hamiltonian H, Eq. (12), at zero magnetic field and vanishing SOI, yielding the eigenbasis {|−, ⇑ , |−, ⇓ , |+, ⇑ , |+, ⇓ } with energies E ± = ± t 2 c + 2 /4 (Appendix D). Since the |− states are lower in energy compared to the |+ states, we refer to them as bonding and anti-bonding, respectively, a nomenclature inspired by similar states in molecular physics. To understand the mixing of these states in the system under consideration, we proceed to transform the total HH Hamiltonian H at non-zero magnetic field and in the presence of the SOI into the bonding-anti-bonding basis, where we define χ s = χ sin θ, χ c = χ cos θ with the DQD orbital mixing angle θ = arctan ( /2t c ) and K ± = t 2 c + 2 /4 ± gµ B B/2. As one can see, the SOI couples orbital states of different pseudospin. Our aim is to complement this coupling with dipole transitions induced by an electric field, i.e., transitions between orbital states of equal spin, to obtain effective spin rotations in the orbital of lowest energy (Fig. 2).
2. Spin coupling mechanism in the bonding state. We show the unperturbed Zeeman-split energies of the bonding (green lines) and anti-bonding (orange lines) states as a function of the detuning. The red and blue double-headed arrows represent transitions due to the dipole coupling and SOI, respectively, which can be combined to obtain an effective spin coupling in the bonding state (purple dashed arrows). On the right we display the probability densities of bonding (bottom) and anti-bonding (top) states according to the applied Wannier formalism at zero detuning. Brighter colors correspond to larger values.

B. Coupling to electric fields
Both classical electric fields and quantum mechanical photons may be coupled to the DQD system. Neglecting diamagnetic contributions, the minimally coupled Hamiltonian reads where the time-dependent vector potential A e may describe a classical alternating electric field or photons within the framework of cavity quantum electrodynamics. In the present case the light must be polarized along the DQD axis x to couple bonding and anti-bonding states, and we find for the dimensionless momentum operator where the operators for position (τ ) and spin (σ) act in the bonding-anti-bonding basis and the dipole elements read with the typical lateral QD size r = /mω 0 . Hence, the electromagnetic fields achieve a coupling d c between bonding and anti-bonding states of equal spin. Additionally, due to the effect of higher orbital states which are admixed due to the SOI, there are non spin-conserving transitions within a dot of strength d s and d s cos θ as well as transitions between the two dots of strength d s sin θ. For all cases the coupling strength is sensitive to the interdot distance and the applied magnetic field. We now turn to the description of the explicit coupling to a classical alternating and a quantized cavity electric field.
Electric dipole spin resonance (EDSR) aims to manipulate the spin degree of freedom with the aid of an alternating electric field and has been shown to be implementable in a variety of systems [45][46][47]. For the DQD system in this study, we assume a field of amplitude E, driving frequency ω d and polarization along the DQD axis x, E(t) = −E cos(ω d t)e x . The Hamiltonian describing this field is given by where we define the coupling strength Ω c = eErω 0 /ω d for later convenience. Finally, we turn to the description of a quantized cavity field. A single cavity mode of frequency ω c is described by the Hamiltonian H c = ω c b † b, where b † and b denote the photon creation and annihilation operator, respectively. The interaction of the cavity photons with the HHs confined in the DQD for linearly polarized light along x is of the form where g c = e ω 0 /2 0 r mV ω c is the single dot charge coupling strength of the cavity containing the volume of the cavity V and the material specific relative permittivity r ( r = 16 in Ge) [48,49].

V. SPIN ROTATIONS IN THE BONDING STATE
The goal of this study is to define a qubit via the spin states in the low-energy bonding orbital. In Sec. IV we showed that electric fields can couple bonding and antibonding states of equal spin and that the SOI perturbs the bonding spin states such that they couple to the antibonding state of opposite spin. It was argued that the combined effects of the cavity field and the SOI should allow for an effective spin coupling in the bonding orbital. In this section, we quantify this claim and show that in the logical qubit space defined by the perturbed bonding states, the total electric field-HH flopping mode qubit Hamiltonian may either be written in EDSR form, when the qubit is driven by a classical field, or in Rabi form, for the case of interactions with single photons. The EDSR-qubit Hamiltonian (20) may be used to perform single qubit operations, while two-qubit gates can be implemented with the cavity-qubit Hamiltonian (21). As one-and two-qubit gates are universal for quantum computation [50], the interactions described are sufficient to operate a complete quantum computer. In the following we derive explicit expressions for the qubit energy separation ∆ and the effective spin couplings Ω s and g s for two special cases. Additionally, we highlight the strong dependence of the coupling on system parameters such as the inter-dot distance, the dot detuning and the applied magnetic field.

A. Weak magnetic fields and spin-flip tunneling
Assuming inter-dot distances 2a ∼ 100 nm, we find |χ/η| 1 for B 0.1 T, allowing us to neglect χ at weak magnetic fields (cf. Fig. 5). In this limit, we may exactly diagonalize the HH Hamiltonian in Eq. (14) to obtain the eigenstates where φ ± = − arctan η are the spin-orbit mixing angles. On the one hand, an electric field couples the bonding and anti-bonding states of equal spin, e.g., |−, ⇑ and |+, ⇑ . On the other hand, the SOI couples bonding and anti-bonding states of opposite spin via spin-flip tunneling, e.g., |−, ⇑ and |+, ⇓ (cf. Eq. (22a)). The combination of both effects therefore couples the logical states |−, ⇑ and |−, ⇓ . This is seen most easily by expressing the interaction Hamiltonians in Eqs. (18) and (19) in the basis defined by Eqs. (22a) -(22d) in which the HH Hamiltonian H is diagonal. At the level of the qubit space spanned by the logical basis {|0 = |−, ⇓ , |1 = |−, ⇑ }, we find effective floppingmode qubit Hamiltonians as given in Eqs. (20) and (21) with the effective spin couplings Ω s = Ω c d and g s = g c d containing the common dimensionless electric dipole matrix element where d c and d s are as defined in Eq. (17) . By definition, d describes the ratio of spin to charge coupling strengths and its absolute value |d| is thus referred to as the relative spin coupling in what follows. It is shown as a function of the detuning and half the inter-dot distance a in Fig. 3. The coupling changes abruptly in a narrow region around the resonance gµ B B = 4t 2 c + 2 , corresponding to the point where the energies of unperturbed bonding and anti-bonding states |−, ⇑ and |+, ⇓ align in energy. Beyond this critical line, i.e., when E +⇑ < E −⇓ , the excited qubit state in Eq. (22a) changes its character to be predominantly anti-bonding. This allows for strong dipole couplings to the qubit ground state, and the spin coupling mediated by the electric field is increased. Remarkably, the regime where the coupling becomes strongest is also the regime where the qubit is nearly ideally isolated from the rest of the Hilbert space, reducing leakage errors and hence defining an optimal operation mode (Appendix F). We find the same form for the spin coupling if the tunnel matrix element is modulated via the in-plane confinement energy ω 0 instead of the inter-dot distance (Appendix G). The energies of the system read where we define σ = 1 (σ = −1) for spin up (down) states, and the qubit energy gap is Hence, for B = 0.1 T we have ∆/ in the GHz range, allowing resonant coupling of the bonding ground state spin to photons confined in superconducting resonators [51].

B. Symmetric dot configuration
As the most pronounced coupling strengths appeared at zero detuning, = 0, in Sec. V A, we restrict our attention to the case of a symmetric DQD in this section, allowing us to consider arbitrary magnetic field strengths. Upon diagonalizing the HH Hamiltonian in Eq. (14), the new eigenenergies of the system read from which the energy separation between the logical states can be easily computed as ∆ = E −⇑ − E −⇓ . The spin states in the same orbital are degenerate at zero field, B = 0 due to χ(B = 0) = 0, in accordance with Kramers' theorem.
Turning to the mixing between the qubit states, one may look at the eigenstates of (14). They have the same form as in Eqs. (22a) -(22d) with the modified spin-orbit mixing angles . (27) Consequently, we obtain the same form for the effective electric field mediated spin couplings Ω s = Ω c d and g s = g c d with d as in Eq. (24) at θ = 0. The relative coupling |d| as a function of half the inter-dot distance and the applied magnetic field is displayed in Fig. 4. As can be seen from the plot, the coupling strength shows a strong dependence on the applied magnetic field and the inter-dot distance. Indeed, there exists a critical line where the coupling changes abruptly. It is described by the resonance condition gµ B B = 2t c , which is explained physically as in Sec. V A by equal energies of unperturbed bonding and anti-bonding states and a change of character in the excited qubit state when the line is crossed (Appendix F). The subsequent decrease in |d| is due to the exponential decay of the dipole elements d c and d s in both a and B (Eq. (17)). This knowledge can be taken into account in QD manufacturing to maximize the relative spin coupling strength |d| by constructing the system such that the inter-dot distance equals 2a max , where a max is the value of a where |d| becomes maximum for a fixed magnetic field. In a DQD where the inter-dot distance is tunable (or, alternatively, the in-plane confinement energy, see Appendix G), the characteristic step form of |d| may be utilized as a spin-orbit switch by only slightly changing a (Fig. 4(b)). Since the dependence of d on the magnetic field shows a similar behaviour (Fig. 4(c)), varying the magnetic field may also be used to control the qubit-field interaction. Such tunable couplings are of great importance in quantum processing units. Large couplings can be used to quickly manipulate the spin qubit state, while low couplings may be used to isolate a given spin state, e.g., for the purpose of qubit read-out or initialization. Due to the magnetic compression of the oscillator states a smaller dot separation is favourable for larger magnetic fields. We find a characteristic step-like increase in the coupling when varying the control parameter, implying a high degree of control over the coupling strength. The remaining parameters are chosen as in Fig. 3.

C. Feasible qubit gate times
Finally, we may estimate one-and two-qubit gate times which may be achieved with a natural HH flopping mode qubit in present-day experiments. We assume the system to be engineered such that the qubit is operated at the optimal point, where we may reach relative spin couplings |d| > 1/4 (cf. Fig. 4).
One-qubit gates are most efficiently driven by an alternating electric field of magnitude E and frequency ω d . For a QD of confinement energy ω 0 = 1 meV (yielding a lateral size of r ∼ 25 nm) and an electric field at the resonant driving frequency ω d = ∆ ∼ 100 GHz at B = 1 T, we find Ω s ∼ (1 − 10) E [V/m] MHz and thus one-qubit gate times τ 1 ∼ 1/Ω s of less than a nanosecond for E ∼ (0.1 − 1) kV/m.
On the other hand, two-qubit gates may be implemented by harnessing the long-range spin-photon coupling introduced by a bosonic cavity to obtain a controlled interaction between distant spins. For instance, the iSWAP gate may be performed in the dispersive regime in time τ = (4k + 1)π|δ|/2g 2 s with spin qubitcavity detuning δ = ω c − ∆/ and k = 0, 1, 2, . . . [34]. The entangling and universal CNOT gate may then be implemented using two iSWAP gates and one-qubit gates [52,53]. As superconducting resonators have typical coupling strengths of g c ω c /2πω 0 1 − 10 MHz for QDs of lateral size r 10 − 100 nm [51], HH flopping mode systems in Ge allow for fast two-qubit logic with typical gate operation times τ 2 ∼ |δ| [kHz] ns.

VI. CONCLUSION
We show that a flopping mode qubit arises naturally for HHs in Ge, i.e., without the need to create synthetic SOI via a magnetic field gradient. This is made possible by the strong cubic Rashba SOI, which leads to intradot spin flips and inter-dot spin-flip tunneling. Utilizing these processes, we derive an effective electric field mediated spin coupling of the bonding state for different parameter regimes. Additionally, we quantify and physically motivate optimal points of operation and argue that the spin coupling is highly tunable and can reach several ten percent of the charge coupling when the system is engineered accordingly. Consequently, both one-and twoqubit gate times in the nanosecond range are feasible in present-day experiments with natural HH flopping mode qubits. Our results highlight yet another possible application in quantum computing of the versatile platform Ge.
Further research may include the generalization of the system geometry to elliptical QDs where the effect of the SOI is expected to become more diverse due to a less restrictive symmetry, i.e., the SOI may couple the ground state of a given spin to more than one excited level. Moreover, typical relaxation times for flopping mode systems realized with HHs in Ge need to be calculated to estimate their performance under realistic conditions and definitively judge their potential as fast and reliable qubits.
where π = p + eA with A = B(−y, x, 0)/2 in the symmetric gauge. Due to the circular symmetry of the system, it is convenient to express the general wave functions in planar polar coordinates (r, ϕ), where we introduce the main and azimuthal quantum numbers n and l ∈ {−n, −n + 2, ..., n − 2, n}, respectively, L |l| n r 2 /b 2 = (−1) |l| ∂ |l| r L n+|l| r 2 /2b 2 denote the generalized Laguerre polynomials and b 2 = /m ω 2 0 + ω 2 L with Larmor frequency ω L = eB/2m [54,55]. In the main text, we only need the ground state (n = 0) and the third excited states with extremal azimuthal quantum number (n = 3, l = ±3). They read in cartesian cooridnates (x, y), Appendix B: Matrix elements in the spin-flip tunneling term In this appendix we display the matrix elements needed for the calculation of the correction to the spin-flip tunneling term δη 1 in Eq. (13) of the main text. Writing h 0 = H 0 − gµ B Bσ z /2, they read

Numerical computation of matrix elements
In our numerical study we incorporate the overlaps between states of different pseudo-spin, resulting in the overlap matrix M S in the basis |ψ L ⇑ , |ψ L ⇓ , |ψ R ⇑ , |ψ R ⇓ defined in Sec. III of the main text, where the off-diagonal elements collected in δM S are much smaller than unity in the present system. Aiming to orthogonalize the states, we introduce the Wannier matrix W connecting the left and right dot eigenstates to the Comparison of the matrix elements in the DQD system. We show the tunneling element tc and the spin-orbit induced spin-flip terms χ and η as given in Eqs. (6) and (13) as a function of the magnetic field B (solid lines) and the results of the numerical study (circles). For the spin-flip tunneling term η we display the constituent terms separately to highlight the high accuracy of the approximative analytical results. We use the system parameters ω0 = 1 meV, a = 50 nm, Ξ = 100 meV, m = 0.1m0 and αR Ez = 10 −11 eVm.
Wannier states, The orthonormalization condition d, σ|d , σ = δ dd δ σσ for d ∈ {L, R} and σ ∈ {⇑, ⇓} can be moulded into the Wannier equation in matrix form, Solving this system of non-linear equations yields the components of W and thus the Wannier states. We choose the physical solution given by the one with the largest diagonal elements. With the Wannier states at hand, one may then compute the matrix elements by numerical integration for different values of the magnetic field. A comparison between the numerical approach and the analytical results given in Eqs. (6) and (13) is shown in Fig. 5.

Approximative solution to the Wannier equation
When comparing the numerical data to the approximative solutions of Sec. III, we find that all terms are in excellent agreement except for a correction to the spin-flip tunneling element due to a non-vanishing overlap between states of different pseudo-spin, where h 0 = H 0 − gµ B B/2 and W nm denote the components of the Wannier matrix W . We may approximate the solution to Eq. (C3) as W = I−δM S /2, resulting in an error only at quadratic order in the small off-diagonal elements, W M S W † = I + O(δM 2 S ). Consequently, the spin-flip correction reads which is the form displayed in Eq. (13) of the main text. The approximation is expected to hold for small overlaps between different dot states, and we naively expect good agreement with the exact numerical results for a 30 nm The approximation is valid in the green region, while it does not accurately describe the system in the red region. For S 0.3 (to the right of the dashed line) the error caused by neglecting second-order terms is less than ten percent. We use the system parameters g = 10, Ξ = 100 meV, m = 0.1m0 and αR Ez = 10 −11 eVm as well as (a) ω0 = 1 meV and (b) a = 50 nm. Comparison between approximated analytical (solid lines) and exact numerical results (circles) for the relative spin coupling strength |d| as a function of the magnetic field B as given by the line cuts in Fig. 4(b). In each plot we display the region of half the inter-dot distance a in which |d| shows the characteristic step form. We find excellent agreement between analytics and numerics even for a < 10 nm. The system parameters are set to ω0 = 1 meV, g = 10, Ξ = 100 meV, m = 0.1m0 and αR Ez = 10 −11 eVm. (Fig. 6). Note, however, that this bound is too pessimistic as we can solve the Wannier equation exactly for S ⇑⇓ = 0 and only have to resort to the approximative solution for the computation of δη 2 , which tends to zero as the dot separation is decreased. As a consequence, even for dot separations smaller than 2a = 60 nm the spin coupling is not changed noticeably as can be seen from Fig. 7.

B
[meV] ℏω 0 anti-bonding state |+, ⇓ becomes the dominant contribution. Hence, it becomes more and more susceptible to the electric dipole coupling with the state |−, ⇓ which is the dominant contribution to the ground qubit state |−, ⇓ .
We remark that beyond the critical line, the qubit space is almost ideally isolated from the rest of the Hilbert space under consideration. Using that d s is negligibly small compared to d c and taking into account the above considerations, the momentum operator in the spin-orbit basis reads (Appendix E), Since ν − ±π/2 beyond the critical line, one has cos ν − 0, and the qubit space is decoupled in good approximation.