Experimental certification of millions of genuinely entangled atoms in a solid

Quantum theory predicts that entanglement can also persist in macroscopic physical systems, albeit difficulties to demonstrate it experimentally remain. Recently, significant progress has been achieved and genuine entanglement between up to 2900 atoms was reported. Here, we demonstrate 16 million genuinely entangled atoms in a solid-state quantum memory prepared by the heralded absorption of a single photon. We develop an entanglement witness for quantifying the number of genuinely entangled particles based on the collective effect of directed emission combined with the non-classical nature of the emitted light. The method is applicable to a wide range of physical systems and is effective even in situations with significant losses. Our results clarify the role of multipartite entanglement in ensemble-based quantum memories and demonstrate the accessibility to certain classes of multipartite entanglement with limited experimental control.


Supplementary Note 1 -Fluorescence measurement
The measurement of the incoherent reemission (the uorescence) in the backward k b and the forward k f modes shown on Fig. 1 (b) in the main text required several corrections. First, the coupling eciencies were measured independently to be equal to more than 70% inside the optical bre. A slight imbalance in the couplings was taken into account. Another correction was necessary due to the non-perfect transmission of the narrow bandpass lter (FWHM of 10 nm at 883.2 nm wavelength) that was used in the backward mode. The transmission was measured to be 65% at the wavelength of the heralded single photon. Finally, another correction is attributed to the PBS used in the backward mode. To apply it the polarization state of the incoherent reemission has to be characterized. Due to the use of polarization preserving quantum memory (consisting of two crystals separated by the half-wave plate [1]) and a bowtie conguration, the polarization state of the uorescence emitted in backward and forward is close to a mixed state. To conrm it we performed a polarization state tomography of the emitted uorescence. We obtained the purity of 51% and the delity of 95% with respect to the completely mixed polarization state. Based on this result a correction of 0.5 was used in all measurements.
In a single crystal, the absorption probability for every spatial position inside the crystal is not the same. Hence, the total probabilities to have incoherent reemission (the uorescence) in the forward or backward directions are not equal. Assuming the total optical depth of the crystal is d one can write the ratio J between emission rate between two modes as (1) However, due to the use of the double-pass conguration in our experiment, this asymmetry disappears. More precisely, the amount of the emitted uorescence in forward and backward modes should be equal. This was veried by the direct and simultaneous measurement of the uorescence in both the forward and backward modes. We found a ratio of at most J = 1.03(10) over a time interval of more than 800 µs after the absorption of the optical pulse, as shown in Supplementary Figure 1 (a).
To measure the SNR the strong coherent state pulses were coupled to the QM prepared in the crystal. Generally, the SNR is fundamentally limited by the number of atoms N involved to the collective re-emission process and can be expressed as where |α| 2 is the absorbed mean photon number of the coherent state pulse, η is the rephasing eciency of the QM and δ is the contribution from the intrinsic noise of the detection system. This expression is correct for a single coherent state pulse or in the limit of the low repetition rate R of the experiment (R 1/T 1 , T 1 being the spin relaxation time). In the general case, taking into account accumulation of the population in the excited state which comes from many coherent state pulses one can write SNR as This expression can be used to verify the source of the uorescence and nally to estimate the number of atoms that are contributing to the collective emission. First, for the xed repetition rate of R = 1 MHz, the input intensity of the coherent state pulses was varied to see the inuence of the noise from the detection system δ (Supplementary Figure 1 (b)). The maximum measured SNR value was 72 dB for the high mean photon number at the input where the detector's noise contribution is negligible.
While decreasing |α| 2 its contribution starts to be more dominant which decreases SNR. The model based on Supplementary Equation (3) and the independently measured parameters explains well the obtained data.
Next, the repetition rate R was varied for the xed |α| 2 which is high enough to make the detector's noise contribution negligible (Supplementary Figure1(c)). In this case higher values of the SNR can be obtained due to the lower accumulation of the population in the excited state (such that R 1/T 1 which corresponds to the highest SNR value). The uorescence lifetime T 1 was measured to be 250 µs from separate measurement. The only free parameter to t the data is the number of atoms N , which was estimated to be 106.0(1) dB ≈ 4.0(1) × 10 10 atoms. The maximal SNR for the measurement with the heralded single photon is expected to be higher since the spectral bandwidth of the AFC structure in this case is ve times larger. This means that the number of atoms that contribute to the absorption of a single photon is larger than the measured N . Nevertheless, we do not correct the measured SNR but directly work with N .

Supplementary Note 2 -Simplication of the minimization problem
Equation (9) in the main text is a constrained minimization problem over 3M complex numbers. Here, we show how to reduce the complexity with a few simple arguments.
First, note that state Eq. (6) in the main text could be subnormalized, because we neglect populations in other subspaces. However, it is straightforward to see that subnormalized states give the same (p 1 , p 2 ) values as the renormalized state mixed with the ground state. Hence, it is sucient to consider normalized states only.
Second, we choose a i ≥ 0 without loss of generality. For the phases of the b i and c i , note that the b i dependent term in Eq. (8) in the main text can be written as Let us write b i = e iϕi |b i | and Eq. (8) in the main text reads Clearly, the phases have to be set to ϕ i =φ and ϑ i = π + 2φ in order to minimize p 2 . Without loss of generality, we setφ = 0 implying that b i ∈ R and c i ≤ 0. Signal-to-noise ratio (SNR) measurement using strong coherent states. (a) The uorescence lifetime measurement in forward and backward spatial modes. A lifetime of the excited state of T 1 = 250 µs was obtained from the exponential t. Both curves overlap almost perfectly up to at least 800 µs after the absorption of the strong coherent pulse. The coherent emission in the forward mode at 50 ns is not shown. (b) The measured SNR in backward mode using strong coherent state pulses as a function of mean photon number per pulse |α| 2 . Due to the dominated contribution from the detector's noise δ for low |α| 2 the SNR value goes down. The solid line is a model curve based on Supplementary Equation (3) and independently measured parameters (with N as a free parameter). The repetition rate was xed to 1 MHz. (c) SNR as a function of the repetition rate of the pulses for |α| 2 = 10 6 . The solid line is a t based on the Supplementary Equation (3) where the number of atoms N is the only free parameter. The estimation gives at least N = 106.0(1) dB which agrees well with the value obtained from the doping concentration. All error bars represent one standard deviation of the measured uncertainty.
for so-called subradiant states, that is, states with lower intensity in forward direction than the incoherent emission.
To summarize, we have a i ≥ 0, b i ∈ R and c i = − 1 − a 2 i − b 2 i , thus reducing the problem to 2M real parameters. The simplied formulas read and Supplementary Note 3 -Formulas from the Lagrange multiplier The partial derivatives of Eq. (10) in the main text are and ∂f ∂λ = f 1 − C. (10) which we insert into ∂f /∂a i = 0 and nd Supplementary Equation (11) can also be written as We notice that from ∂f /∂a i = 0 and ∂f /∂b i = 0 we nd Supplementary Equations (12) and (13), where the right hand side for both is independent of i. Hence, it follows that a i c i = a j c j for all pairs (i, j). Let us x j = 1 and write a ≡ a 1 , b ≡ b 1 and c ≡ c 1 . Supplementary Equations (14) and (15) thus give two equations for two unknowns (a i , b i ) (recall that c i is just a function of a i , b i ). (15)

Inserting Supplementary Equation (14) into Supplementary Equation
The solutions for b i and c i follow accordingly.
Supplementary Note 4 -Asymptotic formula for p 1 To argue that for large M the only relevant conguration is the symmetric one, we calculate p 1 for arbitrary congurations and show that all but the symmetric conguration are asymptotically outside the relevant interval I (see Methods). To this end, note the factor A 2 = Π M i=1 a 2 i in the Supplementary Equations (6) and (7).
For large M , this implies that almost all a i have to be very close to one. However, solutions two to four in Supplementary Equation (16) We therefore see that only congurations with δm = 4 j=2 m j = O(1) asymptotically give nite values of A. More explicitly, inserting Supplementary Equations (20) and (21) into p 1 for a xed conguration C with δm = O(1) gives ) .
A simple calculation shows that the maximum of Supplementary Equation (22)  Zoom in the (p 1 , p 2 ) plane for M = 3. The thick black curve is the full numerical minimization of p 2 given p 1 . The two colored, thinner lines are two relevant congurations (red is the symmetric conguration). One clearly identies p lim 1 1 = 1/3 and p lim 2 1 ≈ 0.29. The full search follows the minimum of the two congurations. Further condence is gained by sampling millions of states (blue dots).