Introduction

Atom qubits in silicon1 rely on using the potential well naturally formed by the donor atom nucleus to bind the electron spin at cryogenic temperatures allowing for reproducible qubits with a low gate density2. To date, electron spin qubits in natSi have demonstrated \({T}_{2}^{* }=55\) ns3 increasing to \({T}_{2}^{* }=268\)μs in isotopically pure 28Si4 by addressing one of the phosphorus atoms statistically present in an ion implant region5,6. Atomic precision STM lithography allows us to pattern both individual single donor atoms, but also multi-donor molecules in silicon devices. Whilst single donors can trap up to two electrons, donor molecules can host multiple electrons, depending on number of donors present7. This ability to change the number of donors creates an enhanced tunability of the exchange interaction required for two-qubit gates8 and is known to reduce the impact of randomly placed charge impurities on both single and exchange coupled qubits9. Donor molecules have also been shown to produce longer spin relaxation times10 with an inbuilt addressability of spin resonance transitions allowing for lower error rates in multi-qubit control11.

In this paper we present coherent control of an electron spin on a multi-electron 1P–2P phosphorus donor molecule patterned by STM lithography in natSi with a Rabi frequency of 1.2 MHz and a \({T}_{2}^{* }\) dephasing time of 295 ns. Utilising the fifth electron spin on this molecule we observe a single narrow (7 MHz full width half maximum (FWHM)) electron spin resonance (ESR) transition, which we attribute to the reduced hyperfine interaction in our multi-electron system and potential nuclear spin pumping effects. Our results confirm that electron spin coherence times in natSi can be maintained in atomically engineered devices in silicon12, patterned close to heavily doped phosphorus leads and with multiple phosphorus nuclear spins in the qubit. These results are promising for precision engineered phosphorus in silicon architectures, where the transition to isotopically pure 28Si13 promises significant increases in coherence times4.

Results and Discussion

The donor qubit is fabricated using STM hydrogen lithography, where we employ an STM tip to selectively desorb hydrogen atoms from the silicon surface with subsequent dosing with phosphine gas and a high temperature incorporation anneal14. The central device structure after lithography is shown in Fig. 1a. Control structures with a metallic doping density are shown in the bright regions. Bias voltages VL, VC, VR, VB applied to in-plane gates provide electrostatic control over the charge state, and a single-electron transistor (SET), which is biased at VSET ≈ 300 μV, permits single-shot charge detection15 by monitoring the drain current ISET. The SET also acts as an electron reservoir for the donor sites, patterned at a distance of 18.5 nm from the SET.

Fig. 1: Charge stability diagram of a single P donor (L) and a 1P–2P molecule (Rα − Rβ).
figure 1

a Overview STM micrograph of the device during STM lithography. Bright regions are exposed silicon which are subsequently dosed with phosphine to form conductive leads to the central qubit region (dotted white rectangle). b The central device region before dosing, with three sites (L, Rα, Rβ) circled. c The same area following phosphine dosing and imaging. Here additional STM imaging of the surface in the fast (horizontal) scan direction has resulted in the tip un-intentionally creating additional dangling bonds along a scan line (position marked by the red arrow). In the final dosed image surface features containing phosphorus are marked whilst all unmarked bright spots are single dangling bonds incapable of incorporating P atoms. d Measured charge stability diagram of the device, showing ISET as a function of VL and VR (VC = 750 mV, VB = 300 mV). Discontinuities represent charge state transitions on the three sites (blue-L, green-Rα or pink-Rβ dashed lines). Inter-site charge transitions between Rα and Rβ are highlighted in yellow. White numbering gives the electron occupation numbers [L, Rα, Rβ]. Spin readout of Zeeman doublet states is only possible within the purple shaded areas. The operating point for the spin control experiments is given by the triangular red marker.

The dotted rectangular central device region is enlarged in Fig. 1b, c, which show the donor incorporation sites before and after phosphine dosing, respectively. By examining the height profiles we can identify all features in the STM micrographs (see Supplementary Note 1 for a full analysis). To the left of the SET, we observe one PH fragment in site L after dosing. To the right of the SET we observe two incorporation sites. The first site, Rα incorporated a single P atom, indicated by the P–Si heterodimer labelled surrounded by seven dangling bond sites. From previous experiments we know that a single dangling bond cannot absorb and incorporate P into the surface16. The second site, Rβ absorbed one PH2 and two PH species, giving rise to two incorporated P atoms. Collectively, these two sites Rα and Rβ, separated by 8 nm, form a 1P–2P molecular state which can bind up to five electrons. Extrapolating from a previously obtained tunnel coupling of 4.3 GHz for a 2P–3P(3e) system, we expect the tunnel coupling for the fifth electron case in our molecular system to be in the hundreds of GHz taking into account the reduced separation between the quantum dots (16 nm → 8 nm), the additional electrons loaded onto the system (3e → 5e), and the lower number of P atoms comprising the quantum dots (2P–3P → 1P–2P)7.

Following lithography, dosing and annealing, 50 ± 5 nm of epitaxial silicon was grown at 250 C and electron beam lithography used to pattern surface ohmic contacts and a microwave antenna for ESR control11. Using etched registration markers, we align the antenna at a lateral distance of 300 ± 50 nm from the buried donors to produce an oscillating magnetic field B1 perpendicular to the substrate at the donors’ position. Measurements were performed in a dilution refrigerator at an electron temperature of ~250 mK and under the application of a static magnetic field of B0 = 1.45 T along the [110] crystallographic axis, within the plane as indicated in Fig. 1a.

With a positive potential on the two centre gates (VC = 750 mV, VB = 300 mV), we vary gate voltages VL and VR to obtain the charge stability diagram of the SET and donor system, shown in Fig. 1d. From the charge stability diagram we observe breaks in the SET current due to charge transitions from the three donor sites (L, Rα and Rβ) to the SET with three distinct slopes (blue, green and pink). The difference in slope of these charge transition lines in the stability diagram equates to a unique relative lever arm, α of the left and right gate to each dot, confirming the presence of three different charge sites where phosphorus has been incorporated. We identify several accessible charge transitions within the achievable gate range, giving up to a total of six electrons in the system in the upper right corner. At more negative gate voltages (VL < 0 and VR < 0) the system is fully ionised, see Supplementary Note 2. Increasing the gate voltages lowers the electrochemical potential of the bound states, progressively adding electrons in the order Rβ, Rβ, Rα, Rα, Rβ, L (along the VR = VL line).

From charge stability diagrams using different combinations of gates we extract the ratio of the lever arms of each gate set to the donor dots. We then use these ratios, in combination with electrostatic finite element modelling, to triangulate2 and thus confirm the positions of the three donor dot sites L, Rα and Rβ (see Supplementary Note 3). We also calculate the addition energies from the stability diagram, and confirm that the incorporation sites RαRβ indeed consist of a 1P–2P molecular configuration, see Supplementary Note 4. We find that whilst each site has well-defined charge states there is strong capacitive coupling between them with a mutual charging energy \({E}_{{R}_{\alpha }-{R}_{\beta }} \sim 40\) meV. Moreover, we find that the spin parity is determined by the total number of electrons across both RαRβ and not their individual electron number. For example, if Rβ were an isolated quantum dot, spin readout17,18 should be possible when Rβ is hosting 1 (odd) electron, e.g. in the [0,1,1] charge region. However, we observe that spin readout is only possible in the purple shaded regions in Fig. 1d where the total spin count of both Rα and Rβ is odd. This result strongly suggests that the unpaired spin states on Rα and Rβ form molecular singlet states, arising from the strong exchange coupling J between them. For spin readout we used a magnetic field of B0 = 2.5 T, allowing us to estimate that J > gμB0 ~ 70 GHz.

The ESR results we present next are obtained for the [L; Rα; Rβ] = [0; 2; 3] five-electron charge state, at the position of the red triangle marker in Fig. 1d, where the first four electrons in the molecule occupy magnetically inactive singlet states. We note that we obtained ESR spectra for different electron occupations of the donor molecule with hyperfine values that are in agreement with tight binding calculations of the 1P–2P molecular state and this work will be published elsewhere.

Spin-lattice relaxation in donor based qubits is known to be limited by acoustic-phonon mediated valley repopulation19, with a dependence of the relaxation rate T1 on the magnetic field B0 of the form \({T}_{1}^{-1}={K}_{5}{B}_{0}^{5}\). As shown in Fig. 2a, we find a spin-lattice relaxation time of T1 ≈ 16 s (at B0 = 1.45 T), and a proportional constant K5 = 0.0072 T−5 Hz for the 1P–2P molecule. We can compare these results to earlier measurements for a single P atom where T1 ≈ 6 s17,18 and T1 ≈ 15 s for a 2P quantum dot at B0 = 1.5 T10. Whilst the values are very similar we note that these earlier T1 measurements were recorded for the first electron on either the single donor or the 2P quantum dot. This is in contrast to the fifth electron on our molecular 1P–2P system. Since the fifth electron is most likely to be delocalised across the molecule, our results represent the first spin relaxation measurements of a molecular 1P–2P system and as such are distinct from previous results.

Fig. 2: Spin relaxation and resonance spectrum of the 1P–2P molecule’s five-electron state.
figure 2

a Spin relaxation rate measurement of the five-electron state, as a function of magnetic field. Error bars indicate 95% confidence intervals. The known T1 time of a single electron bound to a single P donor is shown in grey. b Energy level diagrams for each stage of the pulse sequence used in the spin control experiments, showing the potentials of the Zeeman states relative to the SET Fermi energy. c Example readout trace of the SET current as a function of time, for a representative four stage pulse sequence incorporating an ESR control pulse. d Electron spin resonance spectrum of the five-electron state of the 1P–2P molecule at B0 = 1.45 T. The fitted peak width is 7 MHz. The inset displays a narrower range of frequency settings and marks the two frequencies fmw,off, fmw,on used in the coherent driving experiments.

Since the long spin lifetime does not allow for fast initialisation by relaxation, we employ a four-level pulse sequence based on controlled tunnelling to and from the SET, which we show in Fig. 2b. The spin state is first emptied, leaving the molecule in a spin 0, even parity charge state (four electrons total). To initialise the spin, an odd parity electron is loaded deterministically spin-down by aligning the SET Fermi energy between the two Zeeman split spin \(\frac{1}{2}\) states. This state is then isolated from the SET by electrostatically plunging far below the Fermi energy during the application of an ESR control pulse via the microwave antenna. Finally, the Zeeman split state is brought near the system’s Fermi energy, producing a detectable change in ISET depending on the spin state. Each full sequence takes 24 ms, limited by the electron tunnel rate (~3 kHz) to and from the SET, with a single example time trace shown in Fig. 2c. We obtain statistics of the fraction of spin-up outcomes p by repeating this pulse sequence.

The molecule’s ESR spectrum for the five-electron state is shown in Fig. 2d. The figure plots p as a function of microwave frequency fmw for an external field of B0 = 1.45 T, microwave power of PMW ~18 dBm (before ~70 dB attenuation in the transmission line) and microwave pulse duration of 330 ns using a single tone at the given frequency. As shown in more detail in the inset which focuses on the centre region (the shaded area marks the frequency range of the inset), we only observe one resonance at a frequency of fmw = 40.394 GHz. Superimposed is a fit to a single gaussian peak (red line) with a FWHM of δfFWHM = 7 ± 1 MHz. The linewidth can be attributed to Overhauser fluctuations of the background spin bath of 29Si atoms in natural silicon20.

The existence of a single peak in the ESR spectrum of Fig. 2d has two potential explanations. Firstly, that the hyperfine coupling between the P donor nuclear spins and the bound electron spin is small. It has been shown that the hyperfine coupling of the third electron on a 2P quantum dot is reduced by a factor of ≈10 compared to a single electron on a single P atom21. Here the outer electron spin is screened from the positively charged donor cores by the spin-paired electrons of the inner orbitals and results in a lateral spread of the outer electron wavefunction and a reduced charge density and hyperfine interaction at the central nuclear sites11. For the fifth electron on our 1P–2P molecular state an even larger reduction in hyperfine coupling might be expected due to the presence of four electrons shielding the inner core. Depending on the exact geometrical donor configuration21, the hyperfine coupling in donor molecules with one electron typically range between tens to hundreds of MHz. With five electrons on our 1P–2P molecule the hyperfine coupling is likely to be no more than a few MHz. This value is comparable to the Overhauser field broadened ESR linewidth in natSi and prevents us from resolving any splitting of the peak experimentally.

Whilst the observation of a single ESR line can be explained by the reduced hyperfine splitting from the increased lateral extent of the wavefunction of the fifth electron, we cannot rule out the alternate explanation that only one nuclear spin configuration is observed in our measurement. This could arise from a nuclear spin pumping effect driven by the repeated loading and excitation sequence22, effectively forcing the phosphorus nuclei into a fixed spin configuration during the measurement.

We measure the electron g-factor by varying the external magnetic field B0 and observing the resonance condition hfmw = gμBB0 (with h Planck’s constant and μB = 57.9 μeV/T the Bohr Magneton). We find g = 1.99 ± 0.02 in agreement with the accepted value for phosphorus donors in silicon23, with an uncertainty limited by the calibration of our superconducting magnet, see Supplementary Note 5.

In the coherent driving experiments, we employ an interleaved measurement method where we constantly toggle the driving frequency between the two frequencies marked in Fig. 2d, fmw,on = 40.394 GHz (on resonance) and fmw,off = fmw,on − 20 MHz (off resonance). In this way we constantly measure a background signal when exciting off-resonance, which allows us to display the results as Δp = p(fmw,on) − p(fmw,off), independent of classical effects such as heating, charge noise or bias voltage drift, which can alter the background spin-up fraction over time. Both p(fmw,on) and p(fmw,off) are obtained using a single-shot readout protocol shown in Fig. 2c.

Coherent control is achieved by varying the pulse duration for fixed microwave power and frequency and extracting the spin-up fraction using the single-shot spin readout method outlined above. The observed Rabi oscillations in Fig. 3a–d for different settings of microwave power, exhibit non-exponential decay envelopes due to the fluctuating Overhauser field of 29Si nuclear spins. The instantaneous Overhauser field component parallel to B0 contributes an unknown detuning relative to the expected electron resonance condition, causing the spin to accumulate phase at an unknown rate. This detuning also modifies the rate of driven spin rotation during an applied ESR pulse, as greater detuning tilts the effective rotation axis out of the plane perpendicular to B0.

Fig. 3: Coherent spin control of the molecular electronic state.
figure 3

ad Spin-up difference signals Δp for four different mircowave power values as a function of pulse duration, showing clear Rabi oscillations. We note that Δp does not reach 1 since the echo sequences are performed by driving a single frequency tone while the fluctuating Overhauser field in natural silicon detunes the instantaneous resonace from the driven frequency as has been observed previously20. Error bars indicate ±1σ, and black dashed lines indicate fits to a power-law decay. For the highest driving strength in d we obtain a Rabi frequency of fR = 1.18 MHz. e The extracted Rabi frequencies as a function of microwave amplitude \({B}_{1} \sim \sqrt{{P}_{{\mathrm{MW}}}}\) show a linear dependence. f Phase control, demonstrated by a spin echo sequence in which the initial and final \(\frac{\pi }{2}\) rotations differ in phase by ϕ, producing a sinusoidal variation on the recovered spin state. Experimental parameters are PMW = 14 dBm, tπ = 840 ns and τ = 50 ns, and the pulse sequence is indicated in the inset.

For this experiment we are in the so-called weak-driving regime, where the Rabi frequency (≤1.2 MHz) is less than the width of the Overhauser field distribution (7 MHz)24,25. In this regime, the decay profile for times beyond the first oscillation period is well described by a power law envelope together with a phase shift of \(\frac{\pi}{4}\). From fits to the formula:

$${{\Delta }}{p}_{\uparrow }({t}_{p})=\frac{a}{\sqrt{{t}_{p}}}\cos \left(2\pi {f}_{\text{R}}{t}_{p}+\frac{\pi}{4}\right)+c,$$
(1)

with a, c, fR as free parameters, indicated by dashed black lines in Fig. 3a–d, we derive the ensemble average Rabi frequency fR for each power setting. Parameters a and c describe the visibility of the oscillations, and reflect the variation in Overhauser field values across the time ensemble of single-shot measurements together with state preparation and readout errors.

In total ≈11,400,000 single-shot traces are represented in these graphs, from which we extract the expected linear dependence between the amplitude of the driving field (\(\sqrt{{P}_{{\mathrm{MW}}}}\)) and the Rabi frequency as shown in Fig. 3e. At the highest power setting we see a Rabi frequency of fR = 1.18 ± 0.02 MHz, and on average an initial π rotation within about 376 ns.

The Overhauser field variations can be partially overcome using spin echo techniques. In Fig. 3f we demonstrate a spin echo sequence (as indicated in the inset) consisting of a \(\frac{\pi }{2}\) rotation from the initial spin-down state, a short delay of time τ = 50 ns, a refocusing π pulse which serves to correct for any unknown phase accumulated throughout the sequence, and following a second delay of τ, a final \(\frac{\pi }{2}\) rotation of varying phase ϕ relative to the first two pulses. For ϕ = 0, the three pulses together result in a total rotation of angle \(\frac{\pi }{2}+\pi +\frac{\pi }{2}=2\pi\) about the X-axis and the final spin outcome is spin-down. For ϕ = π, the final pulse rotates the state in the opposite direction, leaving the final state as spin-up. For intermediate values of the phase angle the final rotation leaves the spin state as a superposition. The data is fit to a sinusoid:

$${{\Delta }}{p}_{\uparrow }({{\Phi }})=c-a\cos ({{\Phi }}),$$
(2)

superimposed as a black dashed line. The result in Fig. 3f demonstrates two axis control over the state of the electron spin, enabling any single qubit unitary transformation by a combination of X (ϕ = 0) and Y (\(\phi =\frac{\pi }{2}\)) axis rotations.

We also study the coherence properties of our spin system by applying compound pulse sequences. We perform a Hahn echo experiment, in which we initialise the electron spin-down, then apply \(\frac{\pi}{2}, {\pi}\) and \(\frac{\pi}{2}\) pulses about the X-axis, separated by time τ as shown in the inset to Fig. 4a before spin readout. Operating at PMW = 18 dBm, we set the π time in the Hahn echo experiment to tπ = 376 ns and vary the free precession time 2τ in Fig. 4a. For short times, we recover the spin-down state, and for very long evolution times, the spin state becomes completely randomised due to complete dephasing caused by uncontrolled fluctuations in the spin bath. The transition between the two regimes is described by a compressed exponential26 of the form:

$${{\Delta }}{p}_{\uparrow }(2\tau)=c-a\exp \left(-{\left(\frac{2\tau}{{T}_{2}}\right)}^{n}\right),$$
(3)

shown as a black dashed line. While fit parameters c and a define the visibility, the coherence properties are captured by an exponent n = 2.6 ± 1 and the Hahn echo decoherence time T2 = 298 ± 30 μs which are both remarkably close to values measured for phosphorus in bulk natSi26. This result indicates that decoherence is mainly due to a slowly fluctuating bath of 29Si nuclear spins27.

Fig. 4: Decoherence and dephasing times of the molecular five-electron spin state.
figure 4

a Hahn spin echo measurement of the five-electron state, plotting the measured spin-up fraction Δp as a function of total free precession time 2τ. b Spin echo envelope measurement with a Ramsey-type pulse sequence. The initial free precession time τ = 3 μs is fixed, and we plot Δp as a function of the relative delay of the final projection pulse δτ. Error bars indicate ±1σ on the spin-up fraction.

Changing the timing of the last projection pulse by introducing a deviation of δτ from a fixed free evolution time τ = 3 μs (shown in the inset to Fig. 4b results in a Ramsey-type spin echo envelope experiment. Phase coherence with the Larmor precession due to B0 is lost if the instantaneous Overhauser detuning is comparable to \(\frac{1}{\delta \tau }\). We sample this spin dephasing across the time ensemble of single shot measurements to determine the pure dephasing time \({T}_{2}^{* }\). Fitting the data in Fig. 4b to a gaussian distribution with standard deviation σ, we obtain \({T}_{2}^{* }=\sqrt{2}\sigma =295\pm 23\) ns. This value is about five times longer than observed in an ion implemented device measuring an individual P phosphorus atom20 in natural silicon. We note however that electron spin dephasing can be influenced by interaction with the nuclear spin bath, which depends on the time dynamics of the hyperfine interaction under a specific experimental pulse sequence28. Accurately modelling the decoherence time would require a full simulation of the multi-electron few-donor wavefunction29, along with a method such as correlated cluster expansion30 to account for non-classical memory effects in the interacting 29Si atoms. However, we can qualitatively attribute the extended dephasing time \({T}_{2}^{* }\) to the reduced electrostatic confinement and thereby more extended wave function of our five-electron molecular state. This results in a weaker coupling to a larger ensemble of nuclear spins, producing a narrower distribution of Overhauser field values. At the same time, our molecular confining potential does not significantly impact the concentration of 29Si nuclear spins or the timescale of their flip-flop dynamics, leaving the T2 decoherence time observed via Hahn echo comparable to that seen for a single phosphorus in bulk natSi26.

In conclusion, we have demonstrated coherent control of an electron spin bound to a 1P–2P donor molecule in isotopically natural silicon. We found that the inhomogeneous dephasing time \({T}_{2}^{* }\) is five times longer than that observed for a similar experiment performed on a single P donor. Along with the intrinsic addressability of donor molecules, and their advantages for achieving controllable inter-qubit exchange coupling, we confirm that multi-electron states bound to donor molecules provide a long-lived coherent spin resource for quantum computing. This work confirms that the presence of numerous dopant atoms within the planar conductive leads of a donor-defined readout SET, along with multiple nuclei within the donor qubit, does not limit its coherence.

Methods

Device fabrication

The STM hydrogen lithography was done at a pressure below 1e−11 mbar with an Omicron Variable Temperature instrument. A chemically cleaned Si(001) wafer was passivated in a beam of atomic hydrogen after a 1100 C reconstruction anneal. The hydrogen mask was selectively removed by scanning with a tip voltage of 3–6 V and current setpoint of 1–10 nA. After lithography, phosphine dosing and imaging, the wafer was heated to 350 C before a 55 nm layer of epitaxial silicon was grown at a rate of 0.15 nm/min to encapsulate the incorporated donors. The donor layer was electrically contacted by depositing aluminium over contact vias produced by reactive ion etching. The contact structures, as well as the microwave antenna were all defined by electron beam lithography using a PMMA mask.

Cryogenic measurements

Spin resonance measurements were performed at 50 mK in a dilution refrigerator. A superconducting solenoid magnet provided the external magnetic field. DC voltage offsets were generated by Yokogawa 7651 and Stanford Research Systems SIM928 voltage sources. Time-varying voltage pulses were generated by a National Instruments USB6363 DAC/ADC device and added to the static offsets via resistive voltage dividers. The combined gate control signals were filtered by two-stage lumped element RC filters inside the dilution fridge with a low-pass cutoff of 150 kHz. The microwave signals were supplied to the on-chip antenna from a Keysight E8267D vector signal generator (with phase and pulse modulation signals supplied by a Tektronix 5014C arbitrary waveform generator) via a lossy stainless steel coaxial cable running from room temperature to 50 mK. The SET current readout signal was amplified by a Femto DLPCA200 transimpedance amplifier and then electrically decoupled and filtered by a Stanford Research Systems SIM910 JFET isolation amplifier and SIM965 Bessel filter before being digitised by the National Instruments USB6363 DAC/ADC.