Robust quantum-network memory using decoherence-protected subspaces of nuclear spins

The realization of a network of quantum registers is an outstanding challenge in quantum science and technology. We experimentally investigate a network node that consists of a single nitrogen-vacancy (NV) center electronic spin hyperfine-coupled to nearby nuclear spins. We demonstrate individual control and readout of five nuclear spin qubits within one node. We then characterize the storage of quantum superpositions in individual nuclear spins under repeated application of a probabilistic optical inter-node entangling protocol. We find that the storage fidelity is limited by dephasing during the electronic spin reset after failed attempts. By encoding quantum states into a decoherence-protected subspace of two nuclear spins we show that quantum coherence can be maintained for over 1000 repetitions of the remote entangling protocol. These results and insights pave the way towards remote entanglement purification and the realisation of a quantum repeater using NV center quantum network nodes.


INTRODUCTION
Linking multi-qubit nodes into a large-scale quantum network [1][2][3][4] will open up exciting opportunities ranging from fundamental tests [5] and enhanced time-keeping [6] to applications in quantum computing and cryptography [1,[7][8][9]. Pioneering experiments with atomic ensembles [3], single atoms trapped in vacuum [2,4,10,11] and spins in solids [12][13][14] have demonstrated entanglement between two optically connected nodes. Directly extending these schemes to quantum networks involving many nodes and spanning large distances is hindered by unavoidable imperfections, including photon loss and local control errors, which cause the success probability and entanglement fidelity to decay rapidly both with number of nodes and with distance.
These challenges can be overcome via entanglement purification [15] in a repeater-type [16] setting that exploits quantum memories within each node [7,17,18]. Crucially, one needs to control and readout individual qubits within the node as well as create entanglement with remote qubits without inducing decoherence on the other qubits in the node. In principle, each of these tasks can be accomplished probabilistically using detectors and quantum memories for single optical photons [19], but the resulting inefficiency poses a severe challenge for practical quantum network realizations [2,3]. Instead, many efforts are geared towards a layered architecture, as depicted in Fig. 1(a). Here, remote quantum nodes are probabilistically coupled via optical photons, while each node has several deterministically addressable memory qubits that do not interact with and are therefore not disturbed by the optical channel.
A promising candidate for implementing such a quantum network architecture is the nitrogen-vacancy (NV) center in diamond. The NV electronic spin provides an optical interface that can be used to establish entanglement between distant nodes [13,20,21], while nearby nuclear spins can serve as multi-qubit registers [22][23][24][25] with second-long coherence times demonstrated in isolation [26]. Here, we address the open challenge of robust simultaneous execution of these two processesremote entanglement generation and local qubit storage -which is a key prerequisite for entanglement purification and quantum repeaters and therefore a critical task in quantum networks. We implement individual control over five nuclear spin qubits, in which we store quantum states while repeatedly using the electronic spin in a sequence which has previously been used to generate inter-node entanglement [13,20,21]. We study in detail how the fidelity of storage depends on the coupling between electronic and nuclear spins and on the average time the electronic spin is in an unknown quantum state. We then use decoherence-protected subspaces to enhance the robustness of quantum state storage, which enables us to increase the exponential decay constant of the qubit fidelity above 1000 repetitions of the inter-node entangling sequence.

EXPERIMENTAL SETTING
The experiments are performed on a diamond device with a natural abundance of 13  Individual electronic spins (purple spin symbols) are entangled probabilistically over large distances using photons (red curly lines). Each electronic spin is hyperfine coupled to a quantum register of surrounding nuclear spins (yellow) that can be deterministically controlled (green arrows). (b) Electronnuclear coupling. The nuclear spins precess in an external magnetic field ⃗ B. The precession axis and frequency, ω0 or ω−1, (black vectors) depend on the state of the electronic spin, |0⟩ or | − 1⟩, via the hyperfine interaction with parallel component A ∥ and perpendicular component A ⊥ (green vectors). (c) Experimental sequence to generate entanglement between remote NV electronic spins [27]. By optical pumping on the "Reset" transition, the spin is initialized in |0⟩ at time tr. Subsequently, a spin superposition state is created via a microwave π/2 pulse. Spin-photon entanglement is then generated via two optical excitations, separated by a microwave π-pulse that inverts the spin state |0⟩ ↔ | − 1⟩. (d) NV electronic orbital and spin level scheme at cryogenic temperature. The ground states |0⟩ (| ± 1⟩) are optically coupled to the excited states |Ex,y⟩ (|E1,2⟩ and |A1,2⟩, red arrows), respectively. These states decay either directly (red dashed arrows) or via the metastable spin singlet states |S⟩ (blue dashed arrows) to one of the ground states.
to a temperature of 4 K in a Helium bath cryostat and apply a magnetic field of 40 mT along the NV symmetry axis. Before each experimental run we prepare the NV center in the negative charge state and ensure resonance with the lasers [23]. By using spin-selective optical transitions the electronic spin (S = 1, m s = 0 ≡ |0⟩, m s = ±1 ≡ | ± 1⟩) is initialized with a fidelity above 0.99 and read out in single shot with an average fidelity of about 0.94. We employ tailored pulse sequences on the electronic spin [28,29] that yield high-fidelity individual control of five 13 C spin qubits surrounding the NV center studied here (the same as in Ref. [29]). In Table I we list for each nuclear spin qubit the measured hyperfine coupling parameters, the dephasing time T * 2 and the combined fidelity of spin initialization and readout F i,r (see S.I.).
We now focus on nuclear spin coherence during application of the Barrett-Kok inter-node entangling sequence [27] (see Fig. 1c) that was used in recent experiments  (2) TABLE I. Experimentally determined parameters of the five 13 C nuclear spin qubits used in this work: hyperfine couplings A ∥ and A ⊥ , precession frequency difference ∆ω, and combined initialization and readout fidelity Fi,r. The hyperfine parameters are measured via dynamical decoupling spectroscopy [30] and have an uncertainty in the last digit. [13,20,21]. This protocol has two steps. First, entanglement between the electronic spin and the emission time of a single photon is created at both nodes. Subsequently the two photons are measured after interfering on a beamsplitter, probabilistically projecting the electronic spins into a maximally entangled state. Because of unavoidable inefficiencies including photon loss, this sequence has to be repeated many times to generate remote entanglement with a high probability. Each time an entanglement attempt fails the electronic spin of the NV center is projected into an unknown state and has to be reset. This can lead to decoherence of the surrounding nuclear spin quantum memories via the (always-on) hyperfine interaction. The interaction Hamiltonian is in secular approximation: (1) Here, S and I denote the electronic and nuclear spin operators, respectively, and A ∥ (A ⊥ ) denote the parallel (perpendicular) hyperfine coupling strength. If the electronic spin state is |0⟩, the nuclear spin will precess at the Larmor frequency ω 0 = 2π × γ| ⃗ B| around the applied magnetic field ⃗ B (see Fig. 1(b)), where γ = 11 kHz/mT is the nuclear gyromagnetic ratio. If the electronic spin state is | − 1⟩, however, the nuclear spin will precess around a slightly tilted axis and at a different In a sufficiently strong magnetic field, γ| ⃗ B| ≫ A 2 ⊥ + A 2 ∥ , the change in precession axis is quadratically suppressed and nuclear spin decoherence is mainly caused by dephasing due to the linear change in the precession frequency Thus, randomisation of the electronic spin state is expected to lead to dephasing of a nuclear spin on a timescale that is inversely proportional to the parallel hyperfine coupling strength [31,32]. This dephasing can be mitigated by a suited dynamical decoupling sequence [32], which is inherent in the Barrett-Kok entangling sequence ( Fig. 1(c)): for τ = 0, the time interval between the microwave (MW) π/2 pulse and the MW π pulse has the same duration as the time interval between the MW π pulse and the start of the electronic spin reset. Thus, the electronic spin will be in |0⟩ and | − 1⟩ for an equal amount of time, irrespective of the random optical projection. Thus, under the condition that the spin reset is instantaneous and happens at a precisely known time, the dephasing is exactly nullified. However, electronic spin reset by optical pumping is a stochastic process of finite time duration. As the spin state is not known during this process, it induces irreversible dephasing of the nuclear spins. Therefore, it is desirable to use nuclear spins with weak coupling strength and to make the electronic spin reset as short as possible.

ELECTRONIC SPIN RESET
We first investigate the spin reset process and optimize the reset time. We initialize the electronic spin in |0⟩ by applying laser fields that are on resonance either with the transitions | − 1⟩ ↔ |E 1 ⟩ and | + 1⟩ ↔ |E 2 ⟩ or, for comparison, with the transitions | − 1⟩ ↔ |A 1 ⟩ and | + 1⟩ ↔ |A 2 ⟩, see Fig. 1(d). Compared to our previous remote entanglement experiments [13,20,21], the use of higher magnetic fields requires a second laser beam because of the comparably large ground state Zeeman splitting between | − 1⟩ and | + 1⟩ of about 2 GHz. The lasers address different excited states to avoid pumping to a coherent dark state. Both repump laser beams have approximately circular polarization and equal power. The excited states have a lifetime on the order of 10 ns [33]. From the excited states, the NV can decay either back to | ± 1⟩, or to the metastable spin singlet states |S⟩. In addition, spin mixing in the excited states also opens a direct decay channel to |0⟩ [34].
To determine the time it takes to reset the spin, we prepare it in | − 1⟩ and pump it with the reset lasers for a varying duration. After a delay of 2.5 µs to ensure that no population is left in the excited or singlet states [34], we measure the population in |0⟩, see Fig. 2(a). The spin pumping process exhibits a double-exponential decay with a relative amplitude ratio for the fast and slow time scales of around 3:1, which slightly depends on the excited states used. The two timescales of this decay depend on the applied laser power. At high power, the reset timescales saturate (see S.I.) at 29(1) ns and 463 (14) ns when driving transitions to |A 1,2 ⟩ (green) and at 48(1) ns and 432(26) ns when driving transitions to |E 1,2 ⟩ (yellow).
We attribute the slower timescale, which is the same for both configurations within measurement uncertainty, to the decay constant of the metastable singlet states. The fitted value is in the same range as previously reported values [35]. The faster timescale has a contri-bution from both direct spin flip transitions to |0⟩ and transitions into the singlet states which then decay to |0⟩. The difference in fast timescales between the two configurations is explained by different decay rates to the ground states and metastable singlet states from the excited states used [33]. When the laser power is reduced, we observe a gradual increase of both timescales, as shown for the |A 1,2 ⟩ configuration (cyan, blue and black curves).
To obtain additional insight into the spin reset process, we measure the probabilities p i to arrive in the states i = |0⟩, | − 1⟩, and | + 1⟩, again 2.5 µs after applying a repump pulse of varying duration. In Fig. 2(b) we plot 1 − p |0⟩ , p |−1⟩ , and p |+1⟩ . We fit the data to rate equation models (solid lines), one for each repump configuration. These models assume the lifetimes of the individual states as measured in Ref. [33], a lifetime of the metastable singlet states of 440 ns (taken from the above fits), a relative singlet decay ratio R i to the states i of (R |0⟩ : R |+1⟩ : R |−1⟩ ) = (2 : 1 : 1) [34], and full decay of the singlet and excited states before the ground state population is measured. The model uses equal decay of both excited states |A 1,2 ⟩ to | ± 1⟩, while the states |E 1,2 ⟩ decay either to | + 1⟩ or to | − 1⟩ [34]. The decay rates of the excited states to |0⟩ and the rate of excitation and stimulated emission caused by the repump laser are free parameters in the model. The quantitative agreement between data and the model strengthens the hypothesis that the slow timescale of the repump process is dominated by the lifetime of the metastable singlet states. The population of the metastable singlet states (before decay to the ground states) predicted by the models is shown as the green dashed curve.

SINGLE-NUCLEAR-SPIN QUANTUM MEMORY
We now turn to the dephasing of single nuclear spin qubits induced by the stochastic trajectory of the electronic spin state during reset. Using recently developed techniques [29], we first initialize one of the five controlled nuclear spins by a projective measurement (see S.I.). Ideally, this prepares the nuclear spin superposition state 1 √ 2 (| ↓⟩ + | ↑⟩). We then perform N repetitions of the remote entanglement sequence. We omit the short optical π pulses, as they are expected to induce negligible nuclear spin dephasing since they preserve the electronic spin eigenstate and can be timed such that the detrimental effect of occasional spin flips (p F lip < 1 %) [23] is mitigated by the dynamical decoupling sequence. In addition, the fast optical pulses only lead once per sequence to a population of the excited state, whose spin projection is zero and whose 12 ns lifetime is short compared to the reset procedure. We track the dephasing of the nuclear spins by measuring the length of their Bloch vector projection onto the XY plane. We do not include the Z Repump time (ns) projection as it remains small. In addition, we discard the small fraction of experimental runs in which the NV electronic spin is ionized (S.I.). We first investigate and optimize the timings t and τ of the dynamical decoupling sequence shown in Fig. 1(c). We find that the dephasing rate shows a clear minimum when t equals the Larmor period of the nuclear spins (see S.I.), in which case the effect of entanglement between electronic and nuclear spins caused by the tilted rotation axis, as shown in Fig. 1(b), is minimized. We therefore set t = 2π ω0 ≃ 2.3 µs. We then repeat the entanglement sequence N = 200 times and vary the time τ . Assum-ing an exponential distribution of the repumping time t r , one expects to obtain the best possible decoupling when τ is equal to ⟨t r ⟩ [32]. As can be seen in Fig. 3(a), we observe a broad Gaussian distribution centered around an optimal value of τ ≈ 0.44 µs for all four measured nuclear spins, in good agreement with the slow timescale of the initialization process shown in Fig. 2(a). As mentioned, we attribute this timescale to the decay of the metastable singlet states. At first sight, it is surprising that a singlet state which has zero spin projection and thus no hyperfine coupling would induce dephasing. A possible explanation is that the orbital angular momentum of an E-symmetric singlet state induces a magnetic moment that is comparable to that of the electronic spin ground states.
Additional dephasing can result from experimental imperfections. To prevent errors caused by imperfect spin initialization, e.g. when the laser power drifts over time, we apply the repump laser longer than required for the initialization curves to saturate below 0.01, which is 2 µs (1.5 µs) for the |E 1,2 ⟩(|A 1,2 ⟩) repump configuration at 2000 nW, respectively. To prevent errors caused by imperfect MW pulses, we employ a Hermite pulse envelope with a broad spectrum in order to drive the |0⟩ ↔ | − 1⟩ transition irrespective of the spin state of the nitrogen nucleus of the NV center. We employ this pulse in a BB1 composite pulse sequence [36], consisting of five pulses of less than 0.2 µs duration each, to further suppress pulse errors. In this way, we are able to exclude experimental imperfections as a relevant source of the observed dephasing (see S.I. for details).
With the optimized timings and pulses, we then vary the number of repetitions N of the remote entanglement sequence. For all five nuclear spins, we observe an exponential decay of the single-qubit coherence with N , see Fig. 3(b). Even for the nuclear spin with the largest parallel hyperfine coupling -for which the dephasing is fastest -more than N = 100 repetitions of the entanglement sequence can be run before the Bloch vector length drops to 1/e. For the nuclear spin with the smallest coupling strength, this number increases to N ≈ 500. Further improvements could be achieved by using nuclear spins with even lower parallel hyperfine coupling, although this would generally come at the price of an increased duration of local control operations.

DECOHERENCE-PROTECTED SUBSPACE QUANTUM MEMORY
Motivated by the observation that the memory dephasing is mainly determined by the parallel hyperfine coupling strength, we investigate a different approach to increasing the maximum number of repetitions before a qubit is dephased. Instead of encoding the qubit in a single nuclear spin, we can encode in a decoherence- protected subspace (DPS) [37] of two or more nuclear spins. In this way, the net parallel hyperfine coupling can be strongly reduced while the speed of the individual quantum gates remains the same. A natural choice for a DPS with reduced dephasing is given by the basis states | ↓ i ↑ j ⟩ and | ↑ i ↓ j ⟩ of nuclear spins i and j with comparable parallel hyperfine coupling strength. An encoded qubit will then to first order dephase at a rate that is determined by the coupling strength difference ∆ω ≃ 2π × (A ∥,i − A ∥,j ), which can be much smaller than the individual coupling strengths. On the other hand, encoding a qubit in the states | ↑ i ↑ j ⟩ and | ↓ i ↓ j ⟩ will lead to increased dephasing rates.
To experimentally verify these expectations, we create the states (| ↓ i ↑ j ⟩ + | ↑ i ↓ j ⟩)/ √ 2 and (| ↑ i ↑ j ⟩ + | ↓ i ↓ j ⟩)/ √ 2 [29] (see S.I.) and measure the qubit state projection onto the XY plane of the Bloch sphere under the remote entangling protocol. Fig. 4(a) shows the results obtained for nuclear spins 2 and 3. When initializing the qubit in a decoherence-protected (decoherence-enhancing) two-spin state, we observe a strong improvement (reduction) of the maximum number of repetitions. In the DPS case, we can perform more than one thousand repetitions before the Bloch vector length drops to 1/e. This decay constant can be fully explained by the intrinsic dephasing time T * 2 of the nuclear spins (see S.I.). This shows that the dephasing induced by the entanglement protocol has become negligible in this DPS quantum memory.
The coherence of a DPS with small effective coupling strength might also be limited by the population decay (T 1 ) of the individual spins induced by the entanglement protocol. For the current sample, we find an exponential decay of population with N with decay constants in the range of 1000 to 10000 repetitions, depending on the individual nuclear spins used (see S.I.). We note that this effect is not limiting the coherence of the two-qubit DPS investigated here, but we expect it to become relevant for the smallest investigated coupling strengths once the intrinsic dephasing (T * 2 ) is canceled by an echo.

SCALING OF THE DEPHASING RATE WITH HYPERFINE COUPLING STRENGTH AND REPUMP DURATION
Finally, we perform an extensive quantitative investigation of the scaling of dephasing with the coupling strength and with the time it takes to reset the electronic spin. Fig. 4(b) shows the number of entanglement sequence repetitions N 1/e for which the state fidelity of a balanced superposition state decays to 1/e of the initial value [i.e. the fitted decay constant in 4(a)]. We investigate all 5 individual nuclear spins (open circles) and all 20 possible two-spin subspaces (filled circles), whose coupling strength is the sum or difference of the individual ones. The depicted five data sets correspond to the different values of the repumping time constants shown in Fig. 2. To ensure that we only investigate the scaling of the dephasing with coupling strength and repump duration, we correct for the effects of T * 2 decay, which becomes dominant for the leftmost three data points (see S.I.) and could be compensated by a suited echo sequence on the nuclear spins [26].
We compare the data to the model of Blok et al. [32], which assumes an exponentially distributed repump timescale and that the NV stays in the ground state |−1⟩ until it is reset. When the value of τ used in the dynamical decoupling sequence is equal to the average repump time ⟨t r ⟩, a nuclear spin has acquired a phase shift of ∆ω × (t r − τ ) until the electronic spin is reset to |0⟩. In the limit of large N , the binomial probability distribution of required electronic spin resets can be approximated by a Gaussian distribution. For a balanced superposition input state, as investigated here, this leads to a predicted qubit fidelity of Thus, the model correctly predicts the observed exponential dephasing of the qubit with increasing number of repetitions N . As expected, the decay constant depends on the hyperfine coupling strength of the nuclear spin with a faster decay for increased coupling strength. However, the prediction of the model when inserting the measured optimal value of τ = 0.44 µs (red dotted line) does not exhibit quantitative agreement with the measured data. A possible explanation is that the model neglects the double-exponential reset time distribution and the time spent in | + 1⟩, in one of the excited states or in the meta-stable singlet states. These assumptions are certainly not justified in the present experiment.
To account for this, we leave τ in Eq. 3 as a free parameter and introduce an offset parameter C to the coupling strength, ∆ω → (∆ω + C). Setting C to ≈ 2π × 15 kHz leads to reasonable agreement of the model (solid curves in Fig. 4(b)) with the data. The observed fit values (τ = 0.43(3) µs for the |E 1,2 ⟩ and τ = 0.46(1) µs for the |A 1,2 ⟩ repumping configuration) agree within error with both the measured slow repumping timescale, see Fig.  2(a), and the optimal value of τ in Fig. 3(a).

CONCLUSION AND OUTLOOK
We have studied a prototype quantum network node consisting of nuclear spin qubits hyperfine coupled to an optically active electronic spin in a diamond with natural isotope abundance. Repeated application of a remote entangling protocol is observed to cause dephasing of the nuclear spin qubits. We found that this dephasing can be mitigated by using nuclear spins with small parallel hyperfine coupling strengths, and, even more effectively, by encoding in decoherence-protected subspaces of multiple spins. For the smallest coupling strengths investigated here, the storage of quantum states is robust to more than 1000 remote entangling attempts. We expect that our experimental findings can be generalized to other physical systems [38] and other quantum protocols in which a repeated reset of an ancilla qubit with always-on coupling to a memory is required.
In the future, the implementation of high quality optical cavities should allow for a further reduction of the time it takes to reset the electronic spin, as the Purcell effect [39] induced by such resonator increases the probability of direct spin-flips without populating the singlet states. In addition, the development of techniques Encoding in a decoherence-protected (decoherence-enhanced) subspace leads to strongly decreased (enhanced) dephasing shown in red (magenta). The initial fidelity in the two-spin case is slightly reduced because encoding and readout require more control operations on the nuclear spins. (b) Number of sequence repetitions that are possible before the nuclear qubit Bloch vector length drops to 1/e of its initial value, for qubits encoded in single nuclear spins (empty circles) and in two-spin states (filled circles) of different effective coupling strengths ∆ω. The four depicted data sets are taken for increasing repump duration, caused by a reduced repump laser power. The solid curves are fits according to the model presented in the text.
to measure the electronic spin state non-destructively or within a decoupling sequence might fully eliminate the need for probabilistic repumping. Finally, the realization of quantum networking protocols that are based on photon absorption [40][41][42] rather than photon emission may reduce the number of required electronic spin resets until a successful entanglement event is heralded. Even in the absence of such future improvements, we anticipate that the current results will enable first demonstrations of the purification of remote entanglement [15,43] and proofof-principle operation of a quantum repeater [16] based on NV centers in diamond [17].

C spin initialization and readout
In this section, we explain the procedures used to initialize and readout the 13 C spins using the NV electronic spin as an ancilla [29]. We first describe how a single nuclear spin is prepared in the state |X⟩ ≡ (| ↓⟩ + | ↑⟩)/ √ 2. We employ the gate sequence shown in Fig. 5(a). The electronic spin is first reset to |0⟩ using a laser pulse on the repump transitions (see main text). Subsequently, it is rotated to (|0⟩ + |1⟩)/ √ 2 by a single π/2 microwave (MW) pulse. Then, a sequence of MW π pulses performs a selective rotation on one of the nuclear spins, which is initially in a mixed state ρ m . The rotation is conditional on the initial electronic spin state, as described in detail in Ref. [28]. In this way, the nuclear and electronic spins are correlated. After another MW π/2 pulse, a measurement of the electronic spin projects the nuclear spin to | ± X⟩, where the sign depends on the measurement result.
In our experiment, the readout fidelity of the ancilla qubit is asymmetric: If a photon is detected, the NV is with a very high probability (> 99%) in the fluorescent state |0⟩. Therefore, we repeat the initialization sequence until |0⟩ is measured, which prepares the nuclear spin in |X⟩. The same procedure can be used to initialize all controlled nuclear spins sequentially.
The procedure to initialize a single nuclear spin in | ↓⟩ is shown in Fig. 5(b). Here, the electronic spin is prepared in |0⟩. We then apply a reduced swap operation, which consists of two MW rotations of the electronic spin and two rotations of the nuclear spin that are conditional on the electronic spin. Ideally, this deterministically transfers the electronic spin |0⟩ state to the nuclear spin state | ↓⟩ [29].
Finally, we describe the pulse sequence used to prepare a single qubit in a two-spin decoherence-protected subspace. Here, the nuclear spins are prepared in | ↓↓⟩ by sequential application of the scheme described above. Then, the gate sequence shown in Fig. 5(c) prepares the nuclear spins either in (| ↓↓⟩ + | ↑↑⟩)/ √ 2 or in (| ↑↓⟩ + | ↓↑⟩)/ √ 2, again depending on the measurement result [29]. To avoid errors caused by the asymmetric electron readout fidelity, we again condition the start of the experiment on the detection of a photon, and we insert a MW π pulse before readout if we want to prepare the other nuclear spin state.
The techniques to do tomography on the nuclear spin consist of the same set of operations as those for initialization. They are described in detail in the supplementary information of Ref. [29].
To characterize the control we achieve over the five individual nuclear spins, we prepare each of them separately in |X⟩ and then perform tomography using the electronic spin as an ancilla. We correct for ancilla readout errors. The combined nuclear spin initialization and readout fidelity F i,r is the overlap of the reconstructed nuclear spin density matrix ρ with the ideally expected state, i.e. F i,r = ⟨X|ρ|X⟩. We emphasize that measuring F i,r involves several two-qubit quantum gates and singlequbit rotations. The individual fidelity of each of these operations is expected to be higher than the combined value.
Optimization of the dynamical decoupling sequence As can be seen in Fig. 1(b) of the main text, the hyperfine coupling of the NV electronic and the 13 C nuclear spins leads to two effects: first, a change of the nuclear spin precession frequency, and second, a tilt of the precession axis. In the limit of large magnetic fields, the latter is quadratically suppressed. At the magnetic field used in the experiments, however, the nuclear spins may become partly entangled with the electronic spin in each repetition of the entanglement sequence. This induces decoherence upon projection of the latter.
To minimize this decoherence, two strategies can be explored in order to minimize the entanglement. First, keeping the entanglement sequence very short by minimizing the time t between the pulses, such that the overall nuclear spin evolution during the sequence remains small. Second, the timing of the sequence can be adjusted such that it is matched to one of the nuclear spin pre- cession frequencies, i.e. t = 2π ω−1 or t = 2π ω0 . In this case, the nuclear spins have approximately undergone two full rotations (one around the magnetic field and one around the tilted axis) and end up at the same projection along the Z axis as they started, which minimizes entanglement between electronic and nuclear spin. Fig. 6 shows a measurement in which this theoretical expectation is observed. The measurement uses the sequence shown in the inset of Fig. 3(a) in the main text. Single 13 C spins are initialized in | ↓⟩ and the entanglement sequence is applied N = 450 times. We then measure the length of the Bloch vector Z projection for several values of the wait time t in the dynamical decoupling sequence. We observe two distinct maxima, one at short t and one around the inverse of the nuclear spin precession frequency. The position of the latter slightly depends on the sign and magnitude of the hyperfine coupling. We use t = 2π ω0 ≃ 2.3 µs (red dashed line) in all experiments presented in the main text.

Influence of electronic spin reinitialization and MW pulse errors
In addition to the sources of dephasing described in the main text, also experimental imperfections such as MW pulse errors and imperfect electronic spin reinitialization in |0⟩ might lead to decoherence.
To investigate the effect of the former, we repeat the measurement of the DPS with the smallest dephasing, formed by nuclear spins 2 and 3, but this time with single MW pulses rather than the BB1 pulses used in the main text. Albeit this leads to a reduced MW pulse fidelity, we observe no decrease of the decay constant. In addition, we perform a measurement in which we replace the MW π/2 pulse of the entanglement sequence by a π pulse, such that the electronic spin ends up in |0⟩ at the end of the sequence and no repumping is required. We observe that the dephasing now follows a Gaussian decay and is strongly reduced compared to the case of an initial π/2 pulse. The observed decay can be fully explained by the independently measured nuclear spin dephasing time T * 2 , which testifies that the electronic spin reset is the dominant dephasing mechanism in our experiment.
Effects from imperfect reinitialization in |0⟩ seem unlikely from the remaining population in | ± 1⟩, p ±1 < 0.2(2) %, measured after optical pumping for 1.5 µs using 2 µW of laser power. This expectation is confirmed by the observation of an unchanged decay constant upon stepwise increasing the duration of the repumping laser pulse up to 5.5 µs. Similarly, an additional delay of several µs after the repumping, which would increase the time spent in a state that causes dephasing whenever the pumping process is not perfect, has no discernible effect.
We therefore conclude that the dephasing caused by experimental imperfections is negligible in the measurements presented in the main text.
Power dependence of the optical pumping timescales and NV deionization probability In the main text, we demonstrate that a fast reset of the electronic spin is required to minimize decoherence of the coupled nuclear spins. The reset is doubly exponential, with timescales that depend on the applied laser power. Fig. 7(a) shows the power dependence of these timescales measured in the |A 1,2 ⟩ repump laser configuration. We observe a saturation at ∼ 500 nW of applied laser power. To ensure the fastest possible reset even in the case of NV position or laser power drift, we use 2 µW, which is well above saturation, in all measurements presented in the main text except where mentioned otherwise.
Such comparably high repump power has the unwanted effect that the NV can be transferred to the neutral charge state, which leads to a detectable (and thus heralded) loss of the nuclear spin qubit. This deionization is likely a two-photon process [34], such that we expect its probability to increase quadratically with applied laser power.
To investigate the effects of varying reset timescales and hyperfine coupling strengths on the nuclear spin dephasing without added noise, we post-select the analysis of the results presented in the main text to experimental runs in which the NV has not been deionized. To this end, we measure the NV charge state after completion of the experimental sequence by applying resonant lasers on both the cycling and repump transitions for 50 µs. If the NV is in the negative charge state, we detect on average more than 10 photons in this time interval. If we detect less than two photons, the NV has likely been deionized.
To characterize the influence of NV deionization on the results presented in the main text, we measure the deionization probability as a function of the number of electron resets when applying 2 µW of laser power in the |A 1,2 ⟩ repump configuration, i.e. we study the highest powers used in the main text. Note that in this experiment the electron is prepared in | − 1⟩ and thus reset in every repetition. As can be seen in Fig. 7(b), we observe an exponential decay of the probability to remain in the negatively charged state, p N V − , with a decay constant of N d = 2.82(8) × 10 3 repetitions. In the entanglement sequence, where repumping is required in 50 % of all trials (on average), we expect the decay to be two times slower. Albeit the resulting decay is larger than all observed dephasing constants, it becomes evident that a further reduction of the induced decoherence (e.g. via the use of decoherence-protected subspaces with even smaller coupling strengths than the ones investigated in the main text) will likely require a reduction of the maximally applied laser power.

Nuclear spin polarization decay and natural dephasing correction
As explained in the main text, we expect that the dominant decoherence mechanism of the nuclear spins when applying the entanglement sequence is dephasing which is induced by the electronic spin reset. In contrast, the nuclear spin population is expected to be much less sensitive. To investigate this, we initialize a qubit in the eigenstates | ↓⟩ for each single nuclear spin and in | ↓↓⟩ and | ↑↓⟩ for each investigated two-spin configuration. We then repeat the entanglement sequence many times and measure the exponential decay of the qubit Z projection. As can be seen in Fig. 8, the population decay (red) is typically much smaller than the dephasing shown in black (data identical to Fig. 4(b) of the main text).
In addition to the decoherence that is induced by repeatedly applying the entanglement sequence, the nuclear spins also exhibit natural dephasing. In Ref. [29], this dephasing was shown to lead to a Gaussian decay of the state fidelity with a timescale T * 2 that is different for each measured 13 C spin, c.f. Table I of the main text. This dephasing can be mitigated by a suited dynamical decoupling sequence on the nuclear spins, which can e.g. be implemented by using radio-frequency pulses that are resonant to the nuclear spin Larmor frequency [26].
In our experiments, we do not implement such nuclear spin decoupling sequence. Instead, before fitting the decay we correct the obtained data points for their natural dephasing. To this end, we divide by the factor ) 2 ), where t N is the absolute time elapsed after doing N repetitions of the entanglement sequence and T * 2 is the effective natural dephasing time. For single nuclear spins, we take the measured value of T * 2 , whereas we extrapolate the value for a DPS that is formed by the spins i and j as T * 2,ij ≈ 1 √ (1/T * 2,i ) 2 +(1/T * 2,j ) 2 . In Fig. 8, we show both the corrected (green) and uncorrected (black) data, where it becomes evident that natural dephasing only becomes relevant for DPSs with very low effective coupling strengths.
We thank Adam Gali for discussions on the NV singlet states. We acknowledge support from the Dutch Organization for Fundamental Research on Matter (FOM), the Netherlands Organization for Scientic Research (NWO) through a VENI grant (THT), and the European Research Council through a Starting Grant (RH). . The open circles are single nuclear spins, while the filled circles are twospin subspaces. The red data shows the population decay of a qubit that is encoded in an eigenstate, i.e. along the axis of the external magnetic field. This decay is typically much smaller than the dephasing of a qubit prepared in a superposition state (black and green). The black and green curve show the same measurements, analyzed either with (green) or without (black) applying a correction for the intrinsic nuclear spin dephasing time T * 2 (see text). A significant deviation is only observed for the leftmost data points.