Supplementary Material: Micromagnet-free operation of electron spin qubits in Si/SiGe vertical double quantum dots

Supplementary material: Micromagnet-free operation of electron spin qubits in Si/Si1−xGexvertical double quantum dots Abhikbrata Sarkar1and Daniel Loss1, 2, 3, 4 1Department of Physics, University of Basel, Klingelbergstrasse 82, 4056 Basel, Switzerland 2Physics Department, King Fahd University of Petroleum and Minerals, 31261, Dhahran, Saudi Arabia 3Center for Advanced Quantum Computing, KFUPM, Dhahran, Saudi Arabia 4RDIA Chair in Quantum Computing (Dated: December 20, 2025) Geometry-dependency: valley splitting, shear strain and electrical tunability in the vertical DQD FIG. S1. Simulated shear strain in the DQD due to the metal gate electrodes atop the Si/Si1−xGexheterostructure. a) The xz-map of the εxy(x, y, z)shear strain tensor component in the simulated DQD device. b) Similarly, the xy-map of εxy(x, y, z). c-d) The 2D maps of εxy(x, y, z) in the DQD under plunger gate P1. The Si/Si1−xGex/Si (x=0.033 here) layers are indicated in c). The DQD centered at point (390,410,0) nm spans the following range: x∈[370,410] nm, y∈[390,430] nm, z∈[−6.25,6.25] nm. 3D finite element method (FEM) simulations of the strain distribution at 0.1K are performed using the Structural Mechanics module in COMSOL Multiphysics. The metal gate electrodes are modeled as stress-inducing layers on top of the Si/Si1−xGexheterostructure, and the resulting strain tensor is evaluated throughout the DQD region (Fig. S1). Standard material parameters for Si and Si1−xGexfrom literature have been used [1,2]. Few of the relevant Si1−xGexparameters are listed in Tab. S1. The dominant contribution to the large valley splitting originates from coherent interference of intervalley coupling at the two Si/SiGe interfaces defining the vertical double well. We first examine the dependence of the valley splitting on the barrier height, set by the Ge paramaters values lattice constant, al(5.43+0.20x+0.03x2)×10−10 m density, ρ(2.3+3.5x−0.5x2)×103g cm−3 Young’s modulus, E(130.2−28.1x)×109Pa thermal expansion coeff.,α(2.60 + 2.55x)×10−6K−1 heat capacity, Cp19.6+2.9x 72.6x+28.1(1−x)×103J Kg−1K−1 band gap, Eg1.12 −0.41x+ 0.008x2eV Poisson’s ratio, ν0.278 −0.005x relative permittivity, εr11.7 + 4.5x TABLE S1. Si1−xGexparameters as a function of Ge fraction x. For x=0, these parameters correspond to Si. Note that the value of thermal conductivity Khas been taken from the plot in Ref. [2] highlighting its nontrivial variation with x. concentration xin the Si1−xGexbarrier. Varying xfrom 1.1% to 5.5%, we find that for Vb=15 meV, corresponding to x=3.3% (Si0.967Ge0.033), the valley splitting is maximized (Fig. S2a). In the absence of shear strain, varying the barrier width abetween 2 and 3 nm yields oscillations in the valley splitting ranging from 260 µeV down to 25 µeV. When realistic gate-induced shear strain is included, these oscillations are substantially modified, with the valley splitting varying between 275 µeV and 90 µeV (Fig. S2b). Furthermore, at finite plunger fields |Fz|>0, the oscillations are suppressed and a uniformly large valley splitting of order 150 µeV emerges (Fig. S2c). This behavior reflects the interference of interface-induced intervalley coupling contributions from the two barrier interfaces. Since 2k0≈9.5 nm−1, the relative phase of these contributions is sensitive to changes in the barrier width on the scale of a monolayer (∼0.14 nm). Consistently, we find that the valley splitting exhibits peak-to-peak oscillations when ais varied by approximately 0.16 nm. Gate-induced strain breaks the symmetry between the two interfaces, modifying the relative weighting of their contributions and thereby elevating the minimum valley splitting to 90 µeV. Importantly, while the valley splitting exhibits oscillations on the monolayer scale at zero plunger field, gate-induced strain in realistic devices and plunger field strongly suppress this sensitivity, yielding a robust lower bound ∼150 µeV that is insensitive to atomic scale variations in barrier thickness. We also examine the effect of varying the barrier height on the electrical tun- 2 FIG. S2. Valley Splitting dependence on barrier height and width. a) 2D map of the valley splitting Evas a function of plunger voltage Fzand Si1−xGexbarrier height Vb.Vb=15 (x=3.3%) exhibits the largest valley splitting. b) Valley splitting Evas a function of the barrier width aat Fz=0. The circles and solid fit line in black denote the results without strain; while triangular markers and solid fit line in purple show results with strain. The oscillation in the valley splitting can be fit to a period of 0.31 nm. c) 2D map of Evas a function of Fzand barrier width ain the presence of gate-induced shear strain. FIG. S3. Electrical tunability with varying barrier height. a) gfactor as a function of plunger voltage Fzfor different Si1−xGexbarrier heights Vb(plot legends given in meV). b) EDSR Rabi frequency fπas a function of Fzat different barrier heights Vb;B∥ˆ x,Eac∥ˆ z. c) fπvs. Fzat different barrier heights Vb;B∥ˆ x,Eac∥ˆ x. The barrier width is a=2.5nm. ability of the electron spin qubit. With increasing barrier height, the slope of the gfactor increases, while Vb=15 meV (x=3.3%) produces the largest variation (Fig. S3a). As anticipated, the EDSR Rabi frequency for microwave drive via the plunger gate enhances with increasing maximum ∂g ∂Fzclose to Fz=0 (Fig. S3b). EDSR with in-plane ac electric field (Eac∥ˆx) is the fastest for Vb=15, close to plunger voltage where ∂g ∂Fz≈0(Fig. S3). The z-double well potential basis The double well potential is given by, V(r)=VDW(z)+mt 2 ℏ mtL2 x2 x2+ℏ mtL2 y2 y2 .(1) The z-basis states are calculated by solving the following equation: −ℏ2∂2 z 2ml+VDWϕn(z) = enϕn(z). The piecewise function describing the double well potential is as follows: VDW(z, a, L, Vb) =              V0L+a 2≤z < ∞ 0a 2< z < L +a 2 Vb−a 2≤z≤a 2 0−L+a 2< z < −a 2 V0− ∞ < z ≤ −L+a 2. (2) Below, we provide the full form of the transcendental equation and basis solutions for even and odd bound (en< Vb) and unbound (en> Vb) states. Note that γ=q2mL2Vb ℏ2,r=a 2L. en<Vb,even: 3 γπren Vb (1+r)=nπ+ tan−1−F 1+ en F√V2 0−e2 n 1−enF √V2 0−e2 n , where, F= cos γqen Vbπrcosh γq1−en Vbπr−qVb en−1 sin γqen Vbπrsinh γq1−en Vbπr sin γqen Vbπrcosh γq1−en Vbπr+qVb en−1 cos γqen Vbπrsinh γq1−en Vbπr,(3) ϕn(z) = N×                                        num(F) cos γqen Vbπ(1+r)+den(F) sin γqen Vbπ(1+r) e −γV −1 b√V2 0−e2 nπ(1+r)eγV −1 b√V2 0−e2 nπz/L,−∞ < z ≤ −L1 + r num(F) cos −γqen Vbπz L+den(F) sin −γqen Vbπz L,−L(1 + r)< z < −Lr cosh γq1−en Vbπz L,−Lr ≤z≤Lr num(F) cos γqen Vbπz L+den(F) sin γqen Vbπz L, Lr < z < L(1 + r) num(F) cosγqen Vbπ(1+r) +den(F) sinγqen Vbπ(1+r) e −γV −1 b√V2 0−e2 nπ(1+r)e−γV −1 b√V2 0−e2 nπz/L,−∞ < z ≤ −L1 + r, (4) en<Vb,odd: γπren Vb (1+r)=nπ+ tan−1−F 1+ en F√V2 0−e2 n 1−enF √V2 0−e2 n , where, F= cos γqen Vbπrsinh γq1−en Vbπr−qVb en−1 sin γqen Vbπrcosh γq1−en Vbπr sin γqen Vbπrsinh γq1−en Vbπr+qVb en−1 cos γqen Vbπrcosh γq1−en Vbπr, ϕn(z) = N×                                        − num(F) cos γqen Vbπ(1+r)+den(F) sin γqen Vbπ(1+r) e −γV −1 b√V2 0−e2 nπ(1+r)eγV −1 b√V2 0−e2 nπz/L,−∞ < z ≤ −L1 + r −num(F) cos −γqen Vbπz L−den(F) sin −γqen Vbπz L,−L(1 + r)< z < −Lr sinh γq1−en Vbπz L,−Lr ≤z≤Lr num(F) cos γqen Vbπz L+den(F) sin γqen Vbπz L, Lr < z < L(1 + r) num(F) cosγqen Vbπ(1+r) +den(F) sinγqen Vbπ(1+r) e −γV −1 b√V2 0−e2 nπ(1+r)e−γV −1 b√V2 0−e2 nπz/L,−∞ < z ≤ −L1 + r, (5) en>Vb,even: 4 γπren Vb (1+r)=nπ+ tan−1−F 1+ en F√V2 0−e2 n 1−enF √V2 0−e2 n , where, F= cos γqen Vbπrcos γqen Vb−1πr+q1−Vb ensin γqen Vbπrsin γqen Vb−1πr sin γqen Vbπrcos γqen Vb−1πr−q1−Vb encos γqen Vbπrsin γqen Vb−1πr,(6) ϕn(z) = N×                                        num(F) cos γqen Vbπ(1+r)+den(F) sin γqen Vbπ(1+r) e −γV −1 b√V2 0−e2 nπ(1+r)eγV −1 b√V2 0−e2 nπz/L,−∞ < z ≤ −L1 + r num(F) cos −γqen Vbπz L+den(F) sin −γqen Vbπz L,−L(1 + r)< z < −Lr cos γqen Vb−1πz L,−Lr ≤z≤Lr num(F) cos γqen Vbπz L+den(F) sin γqen Vbπz L, Lr < z < L(1 + r) num(F) cosγqen Vbπ(1+r) +den(F) sinγqen Vbπ(1+r) e −γV −1 b√V2 0−e2 nπ(1+r)e−γV −1 b√V2 0−e2 nπz/L,−∞ < z ≤ −L1 + r, (7) en>Vb,odd: γπren Vb (1+r)=nπ+ tan−1−F 1+ en F√V2 0−e2 n 1−enF √V2 0−e2 n , where, F= cos γqen Vbπrsin γqen Vb−1πr−q1−Vb ensin γqen Vbπrcos γqen Vb−1πr sin γqen Vbπrsin γqen Vb−1πr+q1−Vb encos γqen Vbπrcos γqen Vb−1πr, ϕn(z) = N×                                        − num(F) cos γqen Vbπ(1+r)+den(F) sin γqen Vbπ(1+r) e −γV −1 b√V2 0−e2 nπ(1+r)eγV −1 b√V2 0−e2 nπz/L,−∞ < z ≤ −L1 + r −num(F) cos −γqen Vbπz L−den(F) sin −γqen Vbπz L,−L(1 + r)< z < −Lr sin γqen Vb−1πz L,−Lr ≤z≤Lr num(F) cos γqen Vbπz L+den(F) sin γqen Vbπz L, Lr < z < L(1 + r) num(F) cosγqen Vbπ(1+r) +den(F) sinγqen Vbπ(1+r) e −γV −1 b√V2 0−e2 nπ(1+r)e−γV −1 b√V2 0−e2 nπz/L,−∞ < z ≤ −L1 + r. (8) The transcendental equation can be solved on a graph by plotting LHS and RHS and identifying the intersecting points to produce ensolutions for each cases. The normalization coefficient Ncan be obtained by using R∞ −∞ |ϕn(z)|2dz = 1. Also, num() and den() denote the numerator and denominator of the expression in parentheses, respectively. For this work, we approximate V0>> en(infinite outer boundary of the double well), in which case 1+ en F√V2 0−e2 n 1−enF √V2 0−e2 n≈1. 5 FIG. S4. Schematic of the pulsing sequence of hopping between top and bottom Si dots. The in-plane principle axis of the gtensor are denoted in blue (red) for the top (bottom) QD, with an angular difference ϕbetween them. Each (τt+τb)sequence of hopping back and forth generates 2ϕnet rotation of the spin. Hence after nsuch sequences, e.g. single qubit Sgate (π/2rotation) can be achieved if 2nϕ =π/2. The x−ybasis is described by the following set of equations: ϕl(x) = 1 p2ll!Lx√πe−x2 2L2 xHl(x/Lx) ϕm(y) = 1 p2mm!Ly√πe−y2 2L2 yHm(y/Ly). (9) Hopping-based spin rotation mechanism An anisotropic gtensor can be harnessed to realize single qubit rotations via a spin hopping (or shuttling) sequence between the top and the bottom dots of the DQD structure. As demonstrated in figure S4, the spin could be initialized in the top QD in the direction of the applied magnetic field B∥x. The spin is allowed to precess from −π/2to π/2about the principal g-tensor axis for τttime, followed by the spin hopping to the bottom QD. The spin then precesses from −π/2to π/2about the principal g-tensor axis of the bottom QD for τbtime. When shuttled back to the top QD, the angular difference ϕbetween the two g-tensor axes would produce a 2ϕrotation of the original spin orientation ∥x. By repeating this pulsing sequence ntimes, a π/2rotation of the spin about the z-axis, for example, can be achieved if 2nϕ =π/2. The gate-time for this rotation is estimated as Tπ/2 g=n fL, where fL=gµBB his the Larmor frequency of the electron spin qubit (fL≃2.8GHz for the electron spin at |B|= 100 mT). [1] F. Schäffler, Semiconductor Science and Technology 12, 1515 (1997). [2] Ioffe Institute, “SiGe: Basic properties and material parameters,” https://www.ioffe.ru/SVA/NSM/Semicond/ SiGe/basic.html.