Portraying entanglement between molecular qubits with four-dimensional inelastic neutron scattering

Entanglement is a crucial resource for quantum information processing and its detection and quantification is of paramount importance in many areas of current research. Weakly coupled molecular nanomagnets provide an ideal test bed for investigating entanglement between complex spin systems. However, entanglement in these systems has only been experimentally demonstrated rather indirectly by macroscopic techniques or by fitting trial model Hamiltonians to experimental data. Here we show that four-dimensional inelastic neutron scattering enables us to portray entanglement in weakly coupled molecular qubits and to quantify it. We exploit a prototype (Cr7Ni)2 supramolecular dimer as a benchmark to demonstrate the potential of this approach, which allows one to extract the concurrence in eigenstates of a dimer of molecular qubits without diagonalizing its full Hamiltonian.

O ne of the most intriguing aspects of quantum mechanics are entangled states [1][2][3] , which exhibit correlations that cannot be accounted for in classical physics. Entanglement occurs when a composite quantum object is described by a wave function that is not factorized into states of the object's components, making it impossible to describe a part of the system independently of the rest of it.
Recently, the concept of entanglement has been applied to many different areas of quantum many-body physics, including condensed matter 4,5 , high-energy field theory 6 and quantum gravity 7 . In particular, the development of quantum information science 4,8 has largely increased the interest in the study of entanglement, since it represents an essential resource for quantum information processing applications. Thus, how to experimentally detect and quantify entanglement has become a crucial step for the development of quantum information protocols. However, experimental measurements of entanglement in complex systems are very difficult, because quantum state tomography 9 requires demanding resources and a high degree of control 10,11 . Recently, many theoretical and experimental works have focused on the detection of entanglement [12][13][14][15][16][17][18][19] .
Arrays of weakly coupled molecular nanomagnets represent an ideal playground for investigating quantum entanglement between complex spin systems 4,20 . In particular, supramolecular complexes containing linked antiferromagnetic rings have been demonstrated to be excellent candidates for implementing quantum-computation 21 and quantum-simulation 22 algorithms (see also Supplementary Note 1). A large number of supramolecular complexes of molecular nanomagnets are now being synthesized, but so far entanglement between molecular subunits has only been experimentally demonstrated by exploiting susceptibility as a witness in a (Cr 7 Ni) 2 dimer or by fitting trial model Hamiltonians to electron paramagnetic resonance data 23 . The dimer (see Fig. 1a) consists of two Cr 7 Ni antiferromagnetic rings, and hence is a complex system made of 16 interacting spins, that is, one on each metal site. Strong exchange couplings between the eight spins within each ring lock these spins into a correlated S ¼ 1/2 ground state 24 , which is preserved in the presence of sizeable magnetic fields. Importantly, a weaker inter-ring interaction leads to entanglement between the two composite molecular S ¼ 1/2 spins. The four-dimensional inelastic neutron scattering technique (4D-INS) 25 has the potential to provide a precise characterization of entanglement in such a system. Indeed, the crosssection directly reflects dynamical correlations between individual atomic spins in the molecule, and distinguishes between intra-ring correlations, associated with the composite nature of the molecular qubits, and inter-ring spin-spin correlations, which are a signature of entanglement between qubits.
Here we use this prototype (Cr 7 Ni) 2 dimer as a benchmark system to demonstrate the capability of 4D-INS to investigate entanglement between molecular qubits. By assuming the rings maintain the S ¼ 1/2 ground doublets determined from previous measurements on single-rings 26 , the richness of the experimental information enables us to quantify inter-ring entanglement in the eigenstates through the determination of their concurrence. A magnetic field is used to place the system into a factorized ground state of the two rings and test the occurrence of entanglement in a comfortable parameter regime. Indeed, in this way we can univocally attribute the observed inter-ring correlations to the entanglement in the excited states. This is  possible because of the strong intra-ring couplings which preserve the S ¼ 1/2 ground states. This approach can be applied also to dimers of more complex molecular qubits (described as pseudo-spin 1/2) and it does not require any assumption on the form of the interaction between the qubits.

Results
The (Cr 7 Ni) 2 supramolecular dimer. The Cr 7 Ni antiferromagnetic rings constitute the most studied family of molecular qubits 21,[27][28][29] . The magnetic core is formed by seven Cr ions and one Ni ion arranged at the corners of an octagon. The dominant interaction is the nearest-neighbor antiferromagnetic exchange and leads to an isolated ground doublet behaving as a total spin S ¼ 1/2, which can be used to encode a qubit. These molecules display coherence times sufficiently long for spin manipulations 29 and they have been proposed as prototype for implementing quantum gates 28,30 and quantum simulations 22 .
Here we study the [Cr 7 NiF 3 (C 7 H 12 NO 5 )(O 2 CC(CH 3 ) 15 )] 2 (N 2 C 4 H 4 ) supramolecular dimer, where N 2 C 4 H 4 ¼ N-methyl-D-glucamine (Fig. 1a). The synthesis of (Cr 7 Ni) 2 is described in the 'Methods' section. The Cr 7 Ni subunits have been fully characterized by INS and electron paramagnetic resonance spectroscopies and by low-temperature specific heat and magnetometry measurements 26 . This characterization shows, unequivocally, that the ground state is a S ¼ 1/2 doublet and that it is the only state of the ring populated at the temperatures used in the present work. Indeed, the first excited multiplet lies at about 18 K. The Hamiltonian describing each Cr 7 Ni molecule is (assuming the Ni ion on site 8) where s(i) is the spin operator of the the ith ion (s ¼ 3/2 for Cr 3 þ and s ¼ 1 for Ni 2 þ ). The first two terms correspond to the dominant antiferromagnetic isotropic exchange interaction (with J ¼ 20 K and J 0 ¼ 30 K), while the third term describes the axial single-ion zero-field-splitting terms (with d Cr ¼ À 0.34 K, d Ni ¼ À 7.3 K and the z-axis perpendicular to the ring). The last term represents the Zeeman coupling with an external field B (with g Cr ¼ 1.98 and g Ni ¼ 2.1). The two Cr 7 Ni rings are linked through a pyrazine unit, which provides two N-donor atoms binding to Ni centres in different rings. This leads to a weak exchange coupling between the Ni ions: where A and B label the two rings. We have checked that the effects of the anisotropic zero-field-splitting terms on the results presented in this work are very small (see Supplementary Fig. 1), hence hereafter we neglect the third term in equation (1). Figure 1b reports the calculated magnetic field dependence of the lowest-energy levels of (Cr 7 Ni) 2 . In the presence of an antiferromagnetic ring-ring coupling (j40), the supramolecular system is characterized by an entangled singlet ground state in zero applied field and by an excited triplet. Within this subspace each ring behaves as a total spin S ¼ 1/2, because the strong intraring exchange interactions rigidly lock together the individual spins within each ring. This condition is preserved also in sizeable fields, which can then be used to induce (for B4B cr ) a factorized ferromagnetic ground state for the dimer. In this regime, by focusing on a specific low-temperature transition, we can selectively investigate the entanglement between the rings in the corresponding excited state. Indeed, the origin of the observed inter-ring correlations can be univocally attributed to the excited state involved in the transition (see below) because the ground state is factorized. Moreover, the application of a magnetic field significantly improves the experimental working conditions for detecting inter-ring interactions. In fact, it is easier to resolve a very small splitting d between two peaks in the centre of the energy-transfer range (as it occurs in a sizeable field), than to observe a single-peak centred at d (zero-field case), because it would be partially covered by the elastic signal.
Portraying entanglement with 4D INS. The 4D-INS technique 25 has the potential to yield a deep insight into the entanglement between molecular qubits, by measuring the 4D scattering function S(Q, o) in large portions of the reciprocal Q space and in the relevant energy-transfer 'o range. Indeed, the dependencies of the transition intensities on Q yield direct information on the dynamical spin-spin correlation functions 25 . The zero-temperature spin scattering function is ref. 31 where F d (Q) is the magnetic form factor for the dth ion, |0i and |pi are ground and excited eigenstates, respectively, E p are eigenenergies and R dd 0 are the relative positions of the N magnetic ions within the supramolecular dimer. The products of spin matrix elements are the Fourier coefficients of If |0i is a factorized reference state, and |pi is an excited state where the two rings are entangled, the dynamical correlations of equation (4) are non-zero also for pairs of spins where d is in one ring and d 0 in the other one. Conversely, these inter-ring correlations would be zero if the states of the two monomers were factorized also in |pi, because the corresponding products of spin matrix elements would vanish. The spatial structure of these large-distance correlations produces short-Q modulations in S(Q, E p /:) through the exp(iQ Á R dd 0 ) factors in equation (3), see Fig. 2. Hence, the corresponding pattern of maxima and minima in the measured S(Q, E p /:) is a portrayal of the entanglement between molecular qubits in state |pi. We have measured by 4D-INS a collection of (Cr 7 Ni) 2 single crystals, exploiting the position-sensitive-detectors of the coldneutron time-of-flight spectrometer IN5 at the Institute ILL in Grenoble 32 . This kind of experiment is very challenging because of two conflicting requirements: on the one hand, a detailed Q-dependence of the scattering function over a large range of Q is needed. On the other hand, very high resolution is needed to resolve INS transitions within the dimer's lowest-energy manifold, whose splittings are due to the small inter-ring coupling. In addition, the large single crystals of molecular nanomagnets required for INS are very difficult to grow. To guarantee the occurrence of a factorized reference (ground) state we have applied a magnetic field B ¼ 2.5 T, much stronger than any ring-ring interaction. In addition, measurements have been performed at T ¼ 1.2 K to make the populations of excited states negligible. From the model of equations (1) and (2) two peaks are expected, corresponding to excitations from the factorized ground state towards two excited Bell states (marked by red arrows in Fig. 1b). Figure 1c shows that two inelastic peaks are actually observed at 0.24 and 0.28 meV, both in the energy gain and energy-loss sides (Fig. 1c). These findings are reproduced by equations (1) and (2) with j ¼ 1.1 K. Preliminary measurements were also performed on the high-resolution spectrometer LET at ISIS to estimate the ring-ring coupling.
To illustrate the features in S(Q, o) reflecting inter-ring correlations and entanglement, consider the ideal case of a dimer composed by two perfectly regular and parallel rings, lying in planes parallel to the y-z plane. In this case all intra-ring distance vectors R dd 0 are in the y-z plane, whereas inter-ring vectors have a large component along x (axis perpendicular to the rings). Thus, modulations of S(Q, o) as a function of Q x directly reflect interring dynamical correlations and entanglement through the term exp(iQ Á R dd 0 ) in equation (3). Conversely, intra-ring correlations lead to modulations of S(Q, o) as a function of Q y and Q z . This is illustrated in Fig. 2a,c, that report the calculated INS intensity in the Q x -Q z plane for the two transitions in Fig. 1c, for an ideal (Cr 7 Ni) 2 dimer described by equations (1) and (2). S(Q, o) is characterized by several maxima and minima whose pattern depends on the specific transition and reflects the structure of the involved wavefunctions. In particular, these maps clearly display Q x -dependent modulations, due to the presence of entanglement in the excited |pi states of the dimer. On the contrary, if correlations between the two rings are forced to zero in equation (3), Q x -dependent modulations disappear (Fig. 2e). The residual smooth Q x -dependency of the intensity is merely due to single-ion form factors. A similar behaviour is observed if an Ising inter-ring interaction is considered, as it leads to factorized states (see Supplementary Note 3 and Supplementary Fig. 2). Figure 2b,d,f report the same calculations for a real (Cr 7 Ni) 2 dimer, in which the planes of the rings are not parallel (they form an angle of about 30°) and their normals do not coincide with any axis of the laboratory reference frame (see 'Methods' section). Short-Q modulations of the intensity are still clearly visible for the two transitions to the entangled |pi states of the dimer (Fig. 2b,d) and are absent in Fig. 2f. Longer-period modulations are due to intra-ring correlations in the two non-parallel rings.  Figure 3a,b report the measured Q-dependence of the intensity of the two observed transitions. Short-Q modulations of the intensity are evident for both transitions, with the maxima in the low-energy excitation corresponding to minima in the highenergy one and vice versa. The observed pattern is more complex than in Fig. 2 because dimers with different orientations are present in the single crystals (see 'Methods' section). Nevertheless, these results clearly demonstrate the occurrence of entanglement in the excited states of the supramolecular dimer. Indeed, Fig. 3c,d show that the predicted short-Q modulations coincide with the experimental findings. It is worth stressing that magnetic anisotropy plays a negligible role in these results (see Supplementary Note 2), hence the different orientations of the dimers in the crystals do not affect the possibility of demonstrating the occurrence of entanglement.

Quantification of the entanglement between molecular qubits.
Having experimentally shown the occurrence of entanglement, the next important question is whether is it possible to quantify it. In the following we show that it is possible to extract the concurrence C from the portrait reported in Fig. 3. C is one of the most used measures of the entanglement between a pair of qubits 33 , and its value ranges from 0 for factorized states to 1 for maximally entangled ones. For pure two-qubit states |pi ¼ a|00i þ b|10i þ c|01i þ d|11i, the concurrence is It is evident that C ¼ 0 for a factorized state like |01i or |10i and C ¼ 1 for the Bell states (|10i ± |01i)/ ffiffi ffi 2 p . In the case of a molecular qubit, the |0i and |1i states are typically encoded in the two lowest-energy eigenstates (for example, the ground total spin S ¼ 1/2 doublet in Cr 7 Ni). By considering two molecular qubits (A, B) and restricting to the computational basis (that is, the subspace in which both molecular qubits are in the lowest doublet), the scattering function S(Q, E p /:) for the transitions from the factorized |00i ground state to excited |pi states contains two single-qubit contributions and an interference term: By considering molecular qubits characterized by a well isolated S ¼ 1/2 doublet, the explicit expression of the three contributions is: where the second summation is over the sites d and d 0 of ring A or B, S A;B a are effective spin-1/2 operators acting in the ground doublet of each molecular qubit and c d are the corresponding projection coefficients 34 . This formula can be generalized to other types of molecular doublets. The factorized |00i ground state (obtained by applying a sizeable magnetic field) implies that a ¼ 0 because h00|pi ¼ 0.
All the quantities in I AA , I BB and I AB are here calculated by assuming the S ¼ 1/2 doublet deduced from the model interpreting measurements on the single-ring compound 26 . It is important to note that these quantities can be also determined from measurements on single-qubit compounds (for instance, in Cr 7 Ni the c d coefficients can be extracted by measuring local magnetic moments with nuclear magnetic resonance 35 or polarized neutron diffraction 36 ). Hence, the concurrence C ¼ 2|bc| of the two excited states can be extracted from the data by fitting the observed Q-dependence (Fig. 3) with equation (6). It is worth noting that the inter-qubit interference contribution to S(Q, E p /:) (third term in equation (6)) enables one to discriminate between |pi states with real and complex coefficients. For the present experimental configuration, the observed short-Q modulations in S(Q, E p /:), integrated over an asymmetric range of Q y , point to real wavefunctions (see Supplementary Note 4 and Supplementary Fig. 3), as expected for dominant Heisenberg interactions. Figure  ,b and 4 clearly shows that CC1 for both excited states. It is worth noting that this kind of information cannot be extracted by measuring the energy spectrum as in Fig. 1c, because the mere observation of two peaks is compatible with many possible models with Co1. For instance, an Ising dimer model with two slightly different g values would yield a similar INS energy spectrum, but the two excited states would be factorized and C ¼ 0 (see Supplementary Fig. 2). Furthermore, by determining the concurrence of the eigenstates by this approach, it is also possible to extract the concurrence in the thermodynamic equilibrium state as a byproduct 37 .

Discussion
To summarize, by using the (Cr 7 Ni) 2 supramolecular dimer as a benchmark, we have shown that the richness of 4D-INS enables us to portray entanglement in weakly coupled molecular qubits and to quantify it. This possibility opens remarkable perspectives in understanding of entanglement in complex spin systems. Molecular nanomagnets are among the best examples of real spin systems with a finite size, and are an ideal testbed to address this issue: they are also currently attracting increased attention for quantum information processing [38][39][40] . Tailored finite spin systems can also be assembled on a surface by using a scanning tunnelling microscope 41,42 , but molecular nanomagnets have the advantage that a macroscopic number of identical and independent units can be gathered in the form of high-quality crystals. These enable elusive properties like entanglement to be explored by bulk techniques as 4D INS. It is important to underline that the present method works independently of the specific form of the inter-qubit interaction. Indeed, neither the demonstration of the entanglement through the observed short-Q modulations of the intensity, nor its quantification using equation (6), exploit equation (2). Hence, entanglement can be investigated even between molecular qubits where a sound determination of the full spin Hamiltonian is not possible, as might be the case for molecules containing 4f or 5f magnetic ions. The use of larger single crystals and the consequent increase in the quality of the data, together with additional measurements at higher energies, will even allow the determination of two-spin dynamical correlation functions involving individual spins belonging to different molecular qubits. For instance, this would permit the direct determination of the quantum Fisher information 13 , a witness for genuinely multipartite entanglement and the detection of entanglement between subsystems of the complex spin cluster 15 . Finally, the improvement in the flexibility and in the flux of the new generation of spectrometers, like those under development at the new European Spallation Source (https://europeanspallationsource.se), will further expand the possibilities of this kind of technique.  15 ] 2 (C 4 H 4 N 2 ) (3.0 g) was dissolved in refluxing dichloromethane anhydrous (100 ml) in an Erlenmeyer 500-ml flask, while stirring for 10 min and then anhydrous anisole (C 7 H 8 O) (25 ml) added and the solution was filtered. The filtrated was left undisturbed at ambient temperature under nitrogen for 2 weeks, during which time large well-shaped needles crystals alongside with small crystals including good X-ray quality crystals grew. The crystals were identified by X-ray crystallography as [Cr 7 NiF 3 (Meglu)(O 2 C t Bu) 15 ] 2 (C 4 H 4 N 2 ) (CH 2 Cl 2 ) Á 4(C 7 H 8 O) in short (Cr 7 Ni) 2 dimer, and maintained in contact with the mother liquor to prevent degradation of the crystal quality. Elemental analysis for a powder sample obtained from several large (Cr 7 Ni) 2 15 )] 2 (N 2 C 4 H 4 ) were collected at a temperature of 150 K using a using Mo-K radiation on an Agilent Supernova, equipped with an Oxford Cryosystems Cobra nitrogen flow gas system. Data were measured and processed using CrysAlisPro suite of programs. More details on the crystal structure determination, crystal data and refinement parameters are reported in Supplementary Note 5 and Supplementary Table 1.

Methods
Supplementary Data 1 contains the supplementary crystallographic data for this paper.
Neutron scattering experiments. INS experiments were performed on the IN5 time-of-flight inelastic neutron spectrometer 32 at the high-flux reactor of the Institute Laue-Langevin. The IN5 instrument has a 30 m 2 position sensitive detector divided in 10 5 pixels, covering 147°of azimuthal angle and ± 20°out-of plane. Six (Cr 7 Ni) 2 single crystals were aligned on a copper sample holder with the [010] direction vertical. The crystals were placed with the flat faces perpendicular to the [110] or [1][2][3][4][5][6][7][8][9][10] directions lying on the sample holder, giving two set of crystal orientations. The dimensions of the crystals ranged from a minimum of 5 Â 3 Â 2 mm to a maximum of 8 Â 4 Â 2 mm. Measurements were taken by rotating the crystals (in steps of 1°) about the vertical axis, labelled y in the laboratory frame.
An incident neutron wavelength of 7.5 Å was selected to probe the excitations centred at 0.24 and 0.28 meV, with an energy resolution (full-width at half-maximum) of 22 meV at the elastic line. The magnetic field was applied along the y-vertical axis.
Preliminary measurements were also performed on the high-resolution LET 44 spectrometer at ISIS neutron spallation source using an incident energy of 2.5 meV (40 meV energy resolution), at a temperature of 1.8 K and a field of 5 T.
Data analysis and simulations. Measurements for different rotation angles were combined using the HORACE analysis suite 45 . Corrections for attenuation as a function of rotation and scattering angle based on the slab-like geometry of the sample holder were performed using the formula given in ref. 46. To isolate magnetic signals from the background and to separate in energy the two transitions, the data reported in Fig. 3a,b have been obtained by integrating over the full Q y range and by fitting, for each set (Q x , Q z ), the resulting energy dependence with two gaussians plus a background contribution. The four supra-molecules in the unit cell were included in the calculation of the scattering functions in equations (3) and (6). The two sets of crystal orientations were also taken into account in our simulations, to obtain the intensity maps in Figs 3c,d and 4.
Raw data from the INS experiment were generated at the Insitut Laue-Langevin large-scale facility and are available at the ILL Data Portal with the identifier http://doi.ill.fr/10.5291/ILL-DATA.4-04-457. Derived data supporting the findings of this study are available from the corresponding author upon reasonable request.