As the international scientific community seeks consensus on a new standard for redefinition of the SI (International System of Units) second (ref. ²), it is essential that the uncertainty budgets of optical standards are both rigorously evaluated and experimentally tested through comparison. Although several state-of-the-art optical standards now report evaluated fractional uncertainties below 10⁻¹⁸ (refs. ^12,13,14), same-species comparisons to experimentally test these claims are relatively few^6,7,15,16. Meanwhile, inter-species comparisons at the highest precision have shown significant inconsistencies relative to reported uncertainty budgets^17,18,19,20. An evaluated uncertainty budget is a quantitatively testable prediction that two frequency standards of the same species, corrected for their evaluated systematic shifts, should agree to within their combined uncertainties. In this work, we test that prediction directly by adhering to three criteria for a rigorous assessment of our optical frequency standard²¹: (1) a clear set of experimental measurements of atomic properties and environmental factors determining the systematic uncertainty budget; (2) a same-species comparison to demonstrate validity of the uncertainty budget to the degree the measurement precision allows; and (3) a stress test on any systematic that is substantially larger than the measurement precision of the comparison, that is, a frequency comparison in which the systematic is deliberately shifted and observed to give the predicted shift in the clock frequency. Applying these criteria, we report one of the first comprehensive characterization of all known sources of systematic uncertainty for two independent ¹⁷⁶Lu⁺ single-ion optical references, each with evaluated systematic uncertainty near 1 × 10⁻¹⁹, a fourfold improvement on the lowest previously reported^12,14. Their agreement is verified at the level of 5.7 × 10⁻¹⁹, limited by the comparison precision, after 200 h of averaging.

The exceptionally low systematic uncertainty is possible because of the relative insensitivity of the ¹⁷⁶Lu⁺¹S₀ ↔ ³D₁ transition to perturbations compared with other leading contenders. In particular, this transition has the lowest sensitivity to blackbody radiation (BBR) and magnetic fields of any established clock system, and the large atomic mass makes it less susceptible to motional shifts than lighter atomic species. The complete uncertainty assessment, summarized in Table 1, builds on several recent advances, including an improved method for micromotion compensation²², evaluation of the quadrupole shift by ¹⁷⁶Lu⁺ microwave spectroscopy²³ and evaluation of background gas collision effects²⁴.

Frequency comparison at these extreme levels of precision on practical timescales requires exceptionally high stability, a figure of merit for which single trapped ion optical standards have historically lagged relative to neutral-atom optical lattice clocks^25,26,27,28. In this work, we use correlation spectroscopy^29,30,31, which rejects common-mode phase noise from a shared clock laser and enables interrogation times limited only by the atoms. The measured agreement therefore characterizes the atomic references themselves, independent of the performance of any particular clock laser. Here, we demonstrate interrogation times up to 10 s and a comparison measurement instability of 4.8 × 10⁻¹⁶ (τ/s)^−1/2, nearly a factor of 3 improvement over previous work¹⁶ and within a factor of three of the best reported ratio measurement instabilities involving an ionic frequency standard realized by correlation spectroscopy³¹, differential spectroscopy³² or clock comparisons¹⁸ with state-of-the-art laser stabilization³³.

For ¹⁷⁶Lu⁺, a same-species comparison of independent systems provides exceptionally strong verification because all evaluated systematic effects—except for the readily controlled quadratic Zeeman shift—are small compared with the statistical measurement precision. Consequently, even if these known effects were perfectly common mode, they could not mask an error larger than the statistical uncertainty. To undermine the validity of the claimed performance, a hypothetical systematic shift would have to be both of unknown physical origin and common mode across the two independent experimental setups.

All optical clock systems use some form of averaging, in which the frequencies of two or more transitions are combined such that the average frequency is insensitive to perturbations from external fields. Most commonly, this is an average over one or more pairs of Zeeman lines, which cancels the linear Zeeman shift³⁴. In the case of ¹⁷⁶Lu⁺, we use averaging over the transitions from to for F = 6, 7 and 8 at the optical frequencies f_F. Conceptually, this hyperfine averaging technique³⁵ leads to an effective J = 0 level for which the quadratic Zeeman and quadrupole shifts of the (J = 1) states are largely eliminated in the hyperfine average (HA) frequency, . Rather than sequentially interrogate multiple optical transitions, the average can be conveniently realized by Ramsey spectroscopy on a single optical transition by using a hybrid microwave and optical pulse sequence^16,36. As shown in Fig. 1b,c, the interrogation sequence consists of a Ramsey interrogation on the |g⟩ to |8⟩ transition and microwave transitions within the interrogation time to transfer population between the hyperfine states. Timing is such that the effective time in each of the hyperfine states is the same. Steering the laser to the central Ramsey fringe ensures the linear combination is the HA frequency f. This interrogation scheme has several advantages: exclusive use of first-order field-insensitive m = 0 magnetic substates, the HA is realized every interrogation cycle, and it is readily extended to incorporate the hyper-Ramsey (HR)^37,38 technique for suppression of the laser-induced a.c. Stark shift by inclusion of an additional phase-shifted optical π pulse. The timing, pulse areas and phases of the full hyperfine-average–hyper-Ramsey (HA–HR) sequence used in this work are shown in Fig. 1d with further details given in the Methods.

Full size image

a, Simplified experimental scheme for correlation spectroscopy used to compare Lu-1 and Lu-2. Acousto-optic modulators (AOMs) labelled 1a and 2a control the intensity, frequency and relative phase of clock laser pulses delivered to the respective trapped ions, whereas AOM 2b actively compensates for differential path length fluctuations referenced to retroreflecting mirrors near the viewports of the respective vacuum chambers. Microwaves are delivered by horns external to the vacuum chambers. b, Atomic-level structure of ¹⁷⁶Lu⁺ showing the wavelengths of transitions used. c, The 848-nm optical clock transition and ³D₁ microwave clock transitions addressed in the clock interrogation sequence. Ω_α and f_α denote the coupling strengths and frequencies for the fields driving the indicated transitions, respectively. d, Clock interrogation sequence for hyperfine-averaged HR spectroscopy. e, Parity signal observed for a scan of the optical frequency difference with a Ramsey time T_R = 5 s. Each point represents the mean of 400 measurements, with error bars given by quantum projection noise.

Source data

We compare single ¹⁷⁶Lu⁺ ions in two independent ion traps, Lu-1 and Lu-2, using correlation spectroscopy^29,30, which enables comparison of the references at high stability limited only by the differential atomic coherence rather than the laser coherence³¹. The experimental setup is shown schematically in Fig. 1a and is similar to that described in previous work¹⁶. After initialization of both ions in the |g⟩ state, the ions are simultaneously interrogated with the HA–HR sequence followed by detection of population in the ³D₁ state. For total Ramsey time T_R much longer than the laser coherence time, the expectation value of the parity, Π, is given by³¹

where δf ≡ f⁽¹⁾ − f⁽²⁾ is the difference in the HA-equivalent linear combination of laser and microwave field frequencies used within the interrogation sequence for Lu-k, ϕ is a phase step applied only to the final optical π/2 pulse for Lu-2, and is the parity fringe contrast. A typical frequency scan demonstrating minimal loss of parity fringe contrast for a Ramsey time T_R = 5 s is shown in Fig. 1e.

At the typical operating magnetic field of B₀ = 0.1 mT in both experiment chambers, we observe the contrast loss as a function of interrogation time shown in Fig. 2a. In contrast to previous work¹⁶, we are no longer limited by thermal dephasing due to ion heating: the dashed line in Fig. 2a represents the simulated contrast loss due to thermal effects as evaluated for the current system. Although the HA frequency is extremely insensitive to magnetic fields, the individual component transitions are weakly susceptible to magnetic fields through the quadratic Zeeman shifts of the states, for which the |6⟩ and |8⟩ states have linear magnetic sensitivities of about ±4 kHz mT⁻²B₀ ~±400 Hz mT⁻¹. Independent measurements on magnetic sensitive states indicate the ambient magnetic field has flicker instability of 7 nT on the timescales of interest for clock interrogation (about 1−100 s). The observed loss of contrast is consistent with simulations of uncorrelated magnetic-field noise at this level for each ion, shown in Fig. 2a (solid line). Coherence times on the microwave transitions of individual systems, their dependence on the d.c. magnetic field and the observed contrast for HA–HR correlation spectroscopy are all consistent with uncorrelated magnetic field noise at this level.

Full size image

a, Parity contrast inferred from the servo instability (points) compared with the simulated contrast loss due to ion heating (dashed black line) and due to uncorrelated 7 nT magnetic field flicker noise in both chambers (solid black line). b, Allan deviation for the longest continuous measurement with T_R = 5 s. The dashed line is the projection noise limit for full parity contrast (C = 1) and the solid line is observed instability 4.8 × 10⁻¹⁶ (τ/s)^−1/2 corresponding to C = 0.9. c, Frequency difference (black points) from 11 measurements with a combined duration of 8.3 days over a 12.4-day interval (67% uptime). Error bars represent the statistical uncertainties that are assumed to be projection noise limited. The shaded vertical bars show the measurement intervals, with the colour indicating Ramsey time for T_R = 5 s (blue), 7.5 s (orange) and 10 s (green). The weighted mean frequency difference is [−0.1 ± (5.7)_stat ± (1.0)_sys] × 10⁻¹⁹ with reduced chi-squared \({\chi }_{\nu }^{2}=0.75\) for 10 degrees of freedom.

Source data

The results of 11 comparison measurements totalling 200 h are shown in Fig. 2c. With the exception of two initial measurements taken with Ramsey times of 7.5 s and 10 s, a Ramsey time of T_R = 5 s was used. In all cases, the observed instability is consistent with the quantum projection noise limit after accounting for the fringe contrast. Figure 2b shows the Allan deviation of the longest continuous run (37 hours, T_R = 5 s), demonstrating an instability of 4.8 × 10⁻¹⁶ (τ/s)^−1/2. The weighted mean frequency difference of all 11 measurements is [−0.1 ± (5.7)_stat ± (1.0)_sys]× 10⁻¹⁹.

We evaluate a total systematic uncertainty of 1.2 × 10⁻¹⁹ for Lu-1, 1.3 × 10⁻¹⁹ for Lu-2 and 1.0 × 10⁻¹⁹ for the difference, with individual contributions summarized in Table 1. The only systematic shifts with magnitude larger than 10⁻¹⁸ are the quadratic Zeeman shift and the BBR shift, which is a testament to the insensitivity of ¹⁷⁶Lu⁺ to perturbation. We proceed with a brief discussion on the evaluation of systematic effects, leaving the technical details to the Methods and Extended Data.

The quadratic Zeeman shift is, by a factor of 100, the largest systematic and the only systematic more than an order of magnitude larger than our measurement precision. It is characterized by the quadratic Zeeman coefficient for which a value of α_z = −4.892 64(88) Hz mT⁻² has been previously reported¹⁶. Measurement of the quadratic Zeeman coefficient constitutes a stress test of the reference in accordance with our third criterion, as it calibrates changes in the reference as a function of applied magnetic field. In the interest of scientific rigour and because it is the largest systematic shift, α_z was remeasured in both the years 2024 and 2025, with both measurements being subsequent to a major lab renovation and replacement of the Lu-1 trap. The later measurement data are shown in Extended Data Fig. 1a,b. There is no statistically significant deviation from a quadratic dependence, which highlights the effectiveness of hyperfine averaging operating over a wide range of magnetic fields. The three measurements of α_z are shown in Extended Data Fig. 1c. We take the weighted mean value of α_z = −4.893 27(51) Hz mT⁻², where the statistical uncertainty has been scaled by the square root of the reduced chi-square statistic \({\chi }_{\nu }^{2}=1.85\). It is noted that the probability \(P[{\chi }_{\nu }^{2} > 1.85]=0.16\) for two degrees of freedom and so the three measurements of α_z are statistically consistent. For high-accuracy comparison at 0.1 mT, uncertainty in α_z contributes 1.4 × 10⁻²⁰ uncertainty to the individual optical references, but this contribution is correlated in the comparison. The uncertainty in the difference of quadratic Zeeman shifts is determined from the magnetic field stability and is <10⁻²¹ owing to the tracking and active steering of the magnetic fields in both chambers.

The blackbody radiation shift contributes the largest systematic uncertainty for the individual clocks, at 1 × 10⁻¹⁹, which is limited by the uncertainty in the d.c. polarizability^39,40. The differential BBR uncertainty is limited by the assessment of the environmental temperature, which is detailed in the Supplementary Information. Briefly, the assessment bounds the temperature between that of the ambient environment and the maximum temperature of the local environment of the ion due to radiofrequency (RF) heating of the ion trap. The contribution of RF heating is evaluated by temperature measurements on a test trap of identical construction combined with modelling and numerical simulations of the heat transfer. We estimate the temperatures T₁ = 301.5(1.7) K for Lu-1 and T₂ = 300.8(1.2) K for Lu-2, which capture, within 1σ uncertainty, the minimum and maximum temperature bounds. We note that it would require an increase in the upper bound temperature by 2° to change the BBR uncertainties of the individual systems even at the lowest significant digit reported, 1 × 10⁻²⁰.

Some of the shifts in Table 1 are suppressed by the use of long interrogation times: the a.c. Zeeman shift due to the microwave interrogation pulses, the a.c. Stark shift due to the 848 nm clock laser and HA timing errors due to microwave coupling uncertainties. By contrast, the second-order Doppler shift (SODS) due to thermal motion becomes more significant as the interrogation time is increased because of ion heating. Heating rates for both traps were characterized separately for radial and axial principal axes. At T_R = 5 s, the SODS uncertainty is still at the low 10⁻²⁰ level, but would become the largest source of uncertainty for significantly longer interrogation times.

Both the a.c. magnetic field at the RF trap frequency and the d.c. electric field gradient are characterized by auxiliary measurements on Lu⁺, in contrast to previous work that used Ba⁺ (ref. ⁴¹). The a.c. magnetic field is inferred by measuring an Autler–Townes splitting on the to transition with microwave spectroscopy. The electric field gradient is inferred from the quadrupole shift measured on the ³D₁ microwave transitions, made possible by recent precision measurements of the unperturbed frequencies of the microwave clock transitions²³. Systematic uncertainty in the quadrupole shift on the HA optical frequency is limited by the uncertainty in the residual quadrupole moment^23,42.

At the level of accuracy and precision demonstrated in this work, the shift arising from collisions between the ion and the background gas is a potentially important systematic. We have investigated this through both a classical analysis similar to that in ref. ⁴³ and through a quantum treatment given in ref. ⁴⁴, for which we get agreement from both descriptions. We have shown²⁴ that the shift can be bounded by the classical Langevin collision rate and a factor that characterizes the decoupling of the ion from the clock laser due to the recoil motion. Based on collision rates inferred from these decoupling measurements, we infer background pressures of 3.5 nPa and 5.8 nPa and estimate collision shift bounds of 3.9 × 10⁻²⁰ and 3.2 × 10⁻²⁰ in the respective chambers (Methods).

An ideal optical reference transition has three essential properties⁴⁵: (1) a high frequency and long-lived upper state, enabling good short-term stability; (2) insensitivity to external perturbation, enabling high accuracy; and (3) practical operation with accessible laser wavelengths. The ¹S₀ ↔ ³D₁ transition in ¹⁷⁶Lu⁺ satisfies all of these criteria to a high degree simultaneously, perhaps more than any other reference. The week-long excited-state lifetime removes any practical limit to interrogation time. The cooling, repumping and clock wavelengths are all accessible with commercial solid-state and diode laser sources. The insensitivity to external perturbation has made possible the record-low evaluated uncertainty reported here, although these experiments are operated at room temperature in relatively basic linear Paul traps without magnetic shielding, which may be readily developed into transportable optical clocks⁴⁶ without compromising performance. Moreover, an improved evaluation of the d.c. polarizability as proposed in ref. ⁴⁷ would immediately reduce the total uncertainty of these existing systems to the level of mid-10⁻²⁰, although our measurement precision must advance to take advantage of this ability. These properties in combination make ¹⁷⁶Lu⁺ a compelling candidate as the international community assesses systems for the redefinition of the SI second.

Frequency comparisons are the means by which clock performance is established in practice. The international roadmap towards redefinition of the SI second² sets two thresholds: an evaluated systematic uncertainty of ≲2 × 10⁻¹⁸ and validation through frequency comparisons between institutes demonstrating agreement to an overall uncertainty of ≲5 × 10⁻¹⁸. Between-institute comparisons at this level remain an outstanding challenge, with only a limited number of within-institute same-species comparisons having reached this threshold. Table 2 summarizes all reports, to the best of our knowledge, that meet either threshold. Looking beyond these thresholds, differential gravitational redshift between remote sites becomes an unavoidable limitation. Over long baselines, geopotential differences are presently known to no better than a few centimetres in equivalent height, corresponding to a relative frequency uncertainty ≳2 × 10⁻¹⁸ (refs. ^4,48), well above both the evaluated uncertainty and the verification level of this work. Below 10⁻¹⁸ testing of uncertainty budgets thus requires comparison of references at the same location, with inter-institute validation at this level only possible through transportable references⁴⁶. Of the standards listed in Table 2, fewer than half are accompanied by a same-species comparison, and before this work none had total uncertainty below 10⁻¹⁸.

Accuracy reflects reproducibility, for which we are limited in this report by the measurement precision. But reproducibility alone is insufficient if large systematics are incorrectly evaluated, which is why we have stress-tested the quadratic Zeeman shift multiple times. By contrast, the BBR shift is still too small for practical temperature differences to influence the observed agreement.

At the level of precision demonstrated here, a comparison between remote references is necessarily a direct measure of geopotential difference rather than a test of the optical references. This sensitivity resolves height differences at the millimetre scale, beyond the reach of conventional geodetic techniques^3,4. Likewise, for tests of fundamental physics⁵, an unexplained frequency difference is evidence of new physics only if the references themselves are above suspicion. Fulfilment of this potential therefore requires rigorous, measurement-based assessment and demonstrated reproducibility of the standards. Below 10⁻¹⁸, lutetium stands alone in this regard.

Experimental sequence

Each interrogation beginsonly after both ions have been subjected to Doppler cooling and prepared in the |g⟩ state with >99% probability using the conditional state preparation sequence described in ref. ⁵³. Both ions are interrogated with the HA–HR sequence shown in Fig. 1d, with the parameters τ_L ≈ 4 ms, τ₁ = τ₂ ≈ 12 ms, and a total Ramsey time given by T_R = 3(T + τ₁ + τ₂). Comparisons of the two independent frequency references are carried out by measuring the parity alternately for and, after N interrogation cycles, steering \(\varPi \left(\frac{{\rm{\pi }}}{2}\right)-\varPi \left(-\frac{{\rm{\pi }}}{2}\right)\) to zero by updating δf. For the comparisons at the same 0.1 mT magnetic field, as is the case for all results shown in Fig. 2, the microwave drive frequencies f₁ and f₂ are identical for both ions and fixed throughout. For comparisons at different magnetic fields, as when evaluating the quadratic Zeeman coefficient α_z, the microwave frequencies are necessarily different to compensate quadratic Zeeman shifts. The frequency difference δf is set by an acousto-optic modulator (AOM), labelled AOM 2a in Fig. 1a. Two auxiliary measurements are interleaved with the comparison servo: Rabi spectroscopy of the optical transition to ensure the 848 nm clock laser is near resonance and measurement of the Zeeman splitting for each ion using microwave spectroscopy. The Zeeman splitting is used to infer the magnetic field, and feedback is applied to shim coils to compensate for any slow field drift, with further details in the Supplementary Information. HA–HR interrogation time is typically 87% of the total duty cycle, including the auxiliary measurements.

Magnetic field stability and servo

The short-term stability of the magnetic field is shown in Extended Data Fig. 2 (blue points). This was evaluated from a separate experiment measuring a Zeeman splitting with fast servo attack time. Over the range of Ramsey times of interest in this work (about 1–50 s), the magnetic field noise shows flicker instability of approximately 7 nT. During comparison measurements, the magnetic field is steered to the set point of B₀ = 0.1 mT by a compensation coil with a typical servo attack time of t_ser ≈ 80 s using interleaved measurements of the Zeeman splitting. The orange points in Extended Data Fig. 2 show the field instability inferred from the in-loop measurements of the Zeeman splitting, which averages down ∝(t_ser/τ) for reasons explained in the Supplementary Information. The true instability of the magnetic field is limited by both projection noise and deadtime , which was verified in a test run using interleaved out-of-loop measurements of the Zeeman splitting (Extended Data Fig. 2, green points). This out-of-loop stability is taken as an upper bound on the magnetic-field instability with compensation engaged. The magnetic field contributes uncertainty through the quadratic Zeeman shift of \(\frac{\delta {\nu }_{\mathrm{QZ}}}{{\nu }_{0}}=2{\alpha }_{z}{B}_{0}\delta B(\tau )\approx 2\times 1{0}^{-19}\times {(\tau /{\rm{s}})}^{-1/2}\) beyond the servo attack time.

Detection and background gas collisions

To detect background collisions as much as possible and ensure the ions are sufficiently re-cooled before the next experiment cycle, we use the detection sequence shown in Extended Data Fig. 3a at the end of every Ramsey experiment. The intervals d_i represent adaptive Bayesian state detection⁵⁴ by 646 nm fluorescence. If the ion is detected dark in the initial detection d₀ after the Ramsey sequence, which is about 50% of events, then three attempts are made to shelve on the 848 nm clock transition and detect. The variable 848 nm clock pulse areas are tailored to maximize population transfer even for higher thermally occupied vibrational n states. The probability of at least one successful shelving is >99.5% for the expected thermal distribution accounting for ion heating during the Ramsey dark time T_R. If the ion is not detected bright on any of d_0...3, then it is assumed that a collision has occurred with sufficient energy transfer to significantly reduce either (1) the coupling on the 848 nm transition or (2) the 646 nm fluorescence rate. During clock comparison servo operation, if a collision is detected in this way on either Lu-1 or Lu-2, that interrogation cycle is considered invalid and repeated.

If a collision is detected (d_0...3 all dark), the additional ending sequence shown in Extended Data Fig. 3a is applied. First, a repump pulse (350 nm, 622 nm and 895 nm) and another detection d₄ is attempted. If bright, this is attributed to a lower-energy-transfer collision of type 1: sufficient to reduce efficiency of shelving on the clock transition but not detection. If dark, then the collision is of type 2: energetic enough to disrupt detection. In either case, a cycle of cooling and fluorescence detection with a high threshold, b_i, is repeated until the ion is confirmed bright to ensure the ion is effectively cooled for the next experiment cycle. Extended Data Fig. 3b shows the probability of outcomes (1) and (2) as a function of T_R from which we extract rates Γ⁽¹⁾ and Γ⁽²⁾. To estimate the detectable collision rate Γ, we assume and Γ⁽²⁾ = P_dΓ, where P_d is the probability a collision interferes with bright detection, resulting in outcome (2). Of the lower energy transfer collisions at rate (1 − P_d)Γ, we assume half are not detected because the ion was in ³D₁ at the end of the Ramsey experiment and half are detected as outcome (1). We infer the collision rates Γ = 1.9 × 10⁻³ s⁻¹ for Lu-1 and 5.9 × 10⁻³ s⁻¹ for Lu-2.

After the comparison measurement campaign, additional experiments were performed to estimate the collision rate without the complications introduced by Ramsey spectroscopy. The ions were prepared in the ¹S₀ ground state, and after a delay of 5 s, an attempt was made to reshelve and detect using a similar sequence as shown in Extended Data Fig. 3a. From the probability that all of multiple shelving attempts on the clock transition fail, we directly infer the detectable collision rates Γ = 2.9(3) × 10⁻³ s⁻¹ for Lu-1 and 6.5(6) × 10⁻³ s⁻¹ for Lu-2, in reasonable agreement with the rates inferred from the measurement campaign data.

Our extensive analysis of collision shifts in ion-based optical clocks is given in ref. ²⁴. The collision shift is evaluated with respect to a Langevin collision rate Γ_L, defined as the rate of collisions below a critical impact parameter that result in an inward-spiralling trajectory. We calculate the overestimate of the Langevin collision rate using eq. 31 of ref. ²⁴ and with the cutoff velocity given by the 50% threshold of the Ramsey suppression factor as defined in section 2 of ref. ²⁴. Using ω ≈ 1 MHz for Lu-1 and ω ≈ 500 kHz for Lu-2 yields Γ_L = 1.3 × 10⁻³ s⁻¹ and 2.1 × 10⁻³ s⁻¹ for the respective Langevin rates, corresponding to 3.5 nPa and 5.8 nPa background pressures of molecular hydrogen at 300 K and consistent with pressure gauge readings. An inverted magnetron gauge on Lu-1 reads 4 nPa, and the ion pumps on both chambers read ‘low pressure’ (<13 nPa) at the limit of the sensitivity of the ion pump controllers. With κ = 1 in eq. 60 of ref. ²⁴, we estimate collision shift bounds of 3.9 × 10⁻²⁰ and 3.2 × 10⁻²⁰ for Lu-1 and Lu-2, respectively.

Blackbody radiation

The differential dynamic polarizability, Δα(ω), for the 848-nm clock transition in ¹⁷⁶Lu⁺ has been well characterized^39,40. The BBR shift is given by

where and T₀ = 300 K. Over the practical temperature range of 270–330 K, the uncertainty contribution from Δα(ω) is well approximated by \(9.8\times 1{0}^{-20}{\bar{T}}^{4}\) (ref. ⁴⁰). A detailed temperature assessment is given in the Supplementary Information in which we bound the minimum and maximum temperatures to [299.8, 303.2] K for Lu-1 and [299.6, 301.9] K for Lu-2.

Trap secular frequencies

The Lu-1 trap used in previous work¹⁶ had an unusually high heating rate, which has been resolved after replacing the ion trap with one of the same design⁵³. Moreover, the helical resonators for both Lu-1 and Lu-2 were changed to reduce the RF-drive frequency Ω_RF and obtain higher trap confinement without increased RF power. For 0.25 W of RF power at Ω_RF = 2π × 9.4 MHz and 11.2 MHz, we obtain secular trapping frequencies of (207, 1,063, 1,134) kHz and (138, 492, 543) kHz for Lu-1 and Lu-2, respectively, in which the weakest confinement corresponds to the axial direction. The combined effect of reduced heating rates and the smaller Lamb–Dicke parameters resulting from the increased radial confinement means ion heating is no longer a significant limitation to the interrogation time. The simulated contrast loss shown in Fig. 2a is evaluated by averaging the final optical pulses over a thermal distribution of radial mode occupation accounting for the heating during the Ramsey time.

Thermal second-order Doppler

Thermal motion gives rise to an SODS given by

where T_i are the average temperatures during the Ramsey interrogation for the axial (i = 1) and radial (i = 2, 3) principal axes, respectively, and k_B is Boltzmann’s constant. The geometry factor κ₂ = κ₃ ≈ 2 accounts for equal contributions from the secular motion and intrinsic micromotion for the two radial principal axes, which is valid when the confinement is predominantly ponderomotive⁵⁵. Axial confinement is predominantly due to the static potential, and so only secular motion contributes (κ₁ ≈ 1).

The initial temperatures and heating rates are determined by measuring the temperature immediately after Doppler cooling and after a fixed delay time. The axial temperature is measured by fitting the thermal dephasing of Rabi flopping on the 804 nm E2 transition. Although the 804 nm laser wave vector is at 45° with respect to the ion trap axis and has a projection onto all principal axes, the radial modes have sufficiently low thermal occupation and a small Lamb–Dicke parameter η that they contribute negligibly to the thermal dephasing. From the weighted mean of three measurements over the course of the campaign, we find the initial temperatures T_1,0 = 190.5(5.9) μK for Lu-1 and 189.6(4.7) μK for Lu-2 and heating rates for Lu-1 and 119(17) μK s⁻¹ for Lu-2 (Extended Data Fig. 5). The average temperature during the Ramsey experiment for each mode is given by .

The radial temperatures are measured by spectroscopy on the secular-motion sidebands using the 804 nm E2 transition. The average thermal occupation of one of the radial modes is extracted from the ratio of population transferred on the red and blue sidebands⁵⁶:

For both traps, the initial measurements on both radial modes are consistent with the temperature 150(30) μK, which is also consistent with the initial axial temperatures to within statistical uncertainty and corresponds to approximately three times the Doppler cooling limit. The sideband ratio method is most sensitive for low , so for measuring the radial heating rates, we first apply Zeeman-degenerate Raman sideband cooling⁵⁷ to prepare in the motional ground state and then measure the sideband ratio after some delay. Extended Data Fig. 6 shows the results of radial heating measurements in both chambers, which yield heating rates \(\left[\frac{{\rm{d}}{T}_{2}}{{\rm{d}}t},\frac{{\rm{d}}{T}_{3}}{{\rm{d}}t}\right]=[110(15),29.1(2.9)]\,\mathrm{\mu K}\,{{\rm{s}}}^{-1}\) for Lu-1 and [85.5(5.5), 70.4(5.5)] μK s⁻¹ for Lu-2.

Equation (3) is the SODS for a thermal state in the high-temperature limit, , and does not account for the zero-point fluctuations for the quantum ground state for both secular and micromotion^58,59, which contribute additional SODS of . This term contributes only a 6% correction to the total SODS for Lu-1 but has been included for completeness.

We estimate the total SODS for a Ramsey time T_R as follows:

Excess micromotion

The micromotion shift has two components: an SODS due to oscillatory motion and the a.c. Stark shift due to the RF electric field. Owing to the very low differential polarizability of the clock transition, the a.c. Stark contribution is more than two orders of magnitude smaller than the SODS component and therefore negligible. The micromotion amplitude is evaluated using phase-modulated sideband spectroscopy on the 804-nm E2 clock transition²². The excess micromotion (EMM) SODS is given by

where k₈₀₄ is the wavenumber of the 804-nm clock transition, and β is the modulation depth. As phase-modulated sideband spectroscopy is sensitive to the phase of the micromotion, we distinguish between two quadrature components β = β_m + iβ_p, where β_m is the excess micromotion due to displacement by stray fields, and β_p is due to a phase shift between RF electrodes, which we are not able to compensate. Intrinsic micromotion is accounted for separately with the thermal SODS.

Following the procedure described in ref. ²², we measure the modulation depth in three orthogonal directions to evaluate . The results of all measurements of β_m,i taken for both ion traps are shown in Extended Data Fig. 4a,b, in which the open circles are measured before compensating for the d.c. stray field and closed circles immediately after. To estimate the EMM shift for each comparison measurement, we linearly interpolate the β_m,i measurements, assuming linear growth of the d.c. stray field between measurements. Extended Data Fig. 4c shows the EMM shifts evaluated from the time average of over each comparison interval. We estimate an average relative shift of −1.6(1.3) × 10⁻²⁰ for Lu-1 and −1.4(0.4) × 10⁻²⁰ for Lu-2.

As reported in ref. ²², β_p is negligibly small for Lu-1. For Lu-2, the three orthogonal components of β_p are measured to be [2.22(13), 1.43(12), 1.92(13)] × 10⁻², corresponding to an EMM shift of −2.42(19) × 10⁻¹⁹. β_p was remeasured at the end of the campaign and found to be in statistical agreement with the original measurement.

a.c. Zeeman, RF

Time-varying magnetic fields give rise to an a.c. Zeeman shift, with the dominant contributions coming from currents in the electrodes driven by the RF-trapping potential. The contribution depends on the component of the RF magnetic field perpendicular to the applied d.c. field⁶⁰. The clock shift given by

where is the root-mean-square RF magnetic field amplitude perpendicular to the d.c. field and α_⊥ = 0.20 mHz μT⁻² is the sensitivity coefficient after ³D₁ hyperfine averaging⁶⁰.

Although in previous work B_⊥ was measured using an Autler–Townes splitting on the Ba⁺ clock transition⁴¹, here we use an Autler–Townes splitting on the Lu⁺³D₂ |6, 0⟩ to |5, 0⟩ microwave transition. The ³D₂F = 5 Landé g-factor, g₅ = −0.38575750(19) (ref. ⁵³), is the largest of all ³D₁ and ³D₂F manifolds, and for Ω_RF ≈ 2π × 10 MHz, the linear Zeeman splitting can be brought into resonance with Ω_RF at experimentally accessible magnetic fields of B ≲ 2 mT.

The experimental measurement procedure is as follows. Neodymium permanent magnets are used to bias the magnetic field to the required 1.7 mT for Lu-1 and 2.0 mT for Lu-2. Three pairs of coils are used to fine-tune the magnetic field amplitude and optimize the orientation for the state preparation by aligning to a π-polarized 646-nm laser. This ensures the d.c. magnetic field is aligned in the same direction as during the comparison measurements. From , we transfer sequentially to the state using a microwave π pulse and then the state using an 804-nm–848-nm Raman pair. We observe the Autler–Townes splitting by interrogating the to |5, 0⟩ microwave transition with a π pulse for a range of detunings, after which the remaining population is reshelved to by Raman transfer for detection.

Following the treatment in ref. ⁶⁰, the effective Hamiltonian in the rotating wave approximation describing the microwave interrogation is

where Ω_M is the microwave coupling, Δ_M is the microwave detuning, \({\varOmega }_{B}=\frac{1}{\hbar }\frac{\sqrt{15}}{2}{g}_{5}{\mu }_{B}{B}_{\perp }\) is the coupling between and |5, ±1⟩ states due to the transverse oscillating magnetic field at frequency Ω_RF, and

is the detuning of the RF field from the |5, 0⟩ to |5, ±1⟩ Zeeman splittings, ω_±1. Degeneracy of ω_±1 is lifted by the differential quadratic Zeeman shift for which Δα ≈ 1.23 kHz mT⁻².

Extended Data Fig. 7a shows the simulated spectrum of the Autler–Townes splittings for the operating parameters of Lu-1 and Ω_B =2π × 3 kHz. For both Lu-1 and Lu-2, several spectra were taken at a fixed magnetic field in the vicinity of the Autler–Townes splitting, of which three from Lu-1 are shown in Extended Data Fig. 7b–d. We fit the measured data to simulated spectra to obtain Ω_B. Taking the mean of the Ω_B fit values, with the standard deviation of fit values as the uncertainty, we find to be 0.220(5) μT and 0.303(4) μT for Lu-1 and Lu-2, respectively. The RF-drive voltage to the trap is monitored over the campaign and for both Lu-1 and Lu-2 varied by less than 0.5%.

a.c. Zeeman, microwave

When applying the microwave fields during the Ramsey sequence for hyperfine averaging, there is a probe-induced shift because of the σ_± polarization components off-resonantly coupling to m = ±1 Zeeman states. Evaluation of this shift is discussed in detail in the supplementary material of ref. ¹⁶. The total microwave a.c. Zeeman shift is given by

where Δ_k,F is the shift of |F, m = 0⟩ when microwave coupling Ω_k is on. The microwave field polarizations are characterized by measuring the π times τ_kq for the field k and polarization q, at fixed microwave power, from which the shifts Δ_k,F are evaluated as

where ω_F > 0 is the Zeeman splitting for the hyperfine level F.

This shift is suppressed by ensuring τ_k ≪ T_R, balancing the circular polarization components (τ_k+ ≈ τ_k−), and maximizing the π coupling (τ_k0 ≪ τ_k±). All microwave horns are mounted on rotation mounts and are tuned, through a combination of translation and rotation, to the condition τ_k+ ≈ τ_k− as much as possible at the start of the campaign. All τ_kq were characterized at three points in the campaign, and the evaluated clock shifts were consistently below 1 × 10⁻²⁰.

Quadrupole shift

The hyperfine-averaged residual quadrupole moment has been reported as \(\widetilde{\varTheta }=-2.54(0.25)\times 1{0}^{-4}\,e{a}_{0}^{2}\) (ref. ⁴²), which has been independently confirmed by an alternate method²³. The residual quadrupole shift, , is evaluated by measuring the quadrupole shift, δν_Q,1, on the ³D₁ |7, 0⟩ to |8, 0⟩ microwave clock transition, which are related by

where \(\varTheta {(}^{3}{D}_{1})=0.63862(74)\,e{a}_{0}^{2}\) (ref. ⁶¹). We measure the microwave transition frequencies for both traps by microwave Ramsey spectroscopy with a 10 s interrogation time. The quadrupole shifts δν_Q,1 are inferred using the unperturbed microwave frequencies²³ and accounting for the second-order Zeeman shifts. The residual quadrupole shifts on Lu-1 and Lu-2 are evaluated to be 3.3(3) × 10⁻²⁰ and 1.35(14) × 10⁻¹⁹, respectively. Given the accuracy of the measured unperturbed microwave frequencies²³, the quadrupole shift may be easily suppressed further by setting the magnetic field angle so as to null the microwave quadrupole shifts as much as is required.

Microwave coupling errors

Incorrect microwave pulse areas due to uncertainty in the respective couplings give rise to timing errors in HA–HR spectroscopy. As derived in the supplementary material of ref. ¹⁶, the shifts due to microwave coupling errors are bound by

where Δ_k is the microwave detuning, τ_k is the microwave π pulse duration, and q_k is the fractional error in the microwave coupling for the k = 1, 2 microwave transitions.

These shifts are suppressed for Ramsey times T_R ≫ τ_L, τ_k, as is the case for the comparison experiments reported here. The magnetic field was actively steered to a fixed value in both chambers so the microwave detunings Δ_k are stable and known precisely from the quadrupole shift assessment. The microwave couplings were measured at the beginning, middle and end of the campaign with about 0.2% measurement precision and found to deviate by less than 0.5%, with the exception of Ω₂ on Lu-2, which drifted by 1.5% from the beginning to the end of the campaign. Extended Data Fig. 8, for example, shows this shift evaluated for the worst case of q₂ = 0.015 as a function of the detuning Δ₂. Although the Δ_k are precisely known at the <10 mHz level, which, in principle, allows for a more accurate estimation of this shift, we take a conservative approach and use the full width of the bounding envelope, given by the second term in equation (11), as the uncertainty.

a.c. Stark

For an optical π time of \({\tau }_{{\rm{L}}}=\frac{{\rm{\pi }}}{{\varOmega }_{{\rm{L}}}}=4\,\mathrm{ms}\), the a.c. Stark shift on the |g⟩ to |8⟩ optical transition when interrogating with the 848-nm laser is approximately Δ_S = 2π × 25 Hz. In the absence of additional effects, both the HR³⁷ sequence, Extended Data Fig. 9a, which was used in previous work¹⁶, and the time-reversed HR sequence, Extended Data Fig. 9b, which is used in this work, equivalently suppress the resulting clock shift to \(\frac{2}{{\rm{\pi }}{T}_{{\rm{R}}}}{\left(\frac{\varDelta }{{\varOmega }_{{\rm{L}}}}\right)}^{3}\), in Hz, where Δ = Δ_SP − Δ_S is the error in the frequency step, Δ_SP, applied during the optical interrogation pulses. To set the laser frequency step, Δ_SP, the a.c. Stark shift was measured to better than 1% for both Lu-1 and Lu-2 at the start of the campaign by interleaved self-comparison of Rabi and HR spectroscopy. The value of the Δ_SP was fixed throughout the campaign, and we bound the uncertainty in Δ to 2% of |Δ_S| for Lu-1 and 5% of |Δ_S| for Lu-2 based on the maximum observed variation in the optical couplings, which were measured at the beginning, middle and end of the campaign with measurement precision of about 0.2%.

When including the effects of ion heating, HR schemes present a weak linear dependence on Δ (ref. ⁶²). This linear dependence is evaluated by simulation with the results for the typical operating parameters of Lu-2 shown in Extended Data Fig. 9c. It is noted that the linear slope has the opposite sign for HR and HR reverse schemes, and that a hybrid error signal constructed by averaging the two schemes cancels the thermal effect, restoring the cubic dependence on Δ. The idea of combining sequences to generate a more robust error signal is similar to modified HR⁶³ and generalized HR⁶⁴, but unlike those schemes for which ion heating results in either an offset or increased sensitivity compared with the basic HR scheme⁶², this simple hybrid combination of HR and HR reverse heavily suppresses the linear dependence on ion heating, at least in the regime τ_L ≪ T_R, in which heating during the optical pulse is negligible compared with during the dark time T_R.

For the experiments reported here, the shift was already sufficiently small that only HR reverse was used. From the measured radial heating rates, uncertainty in Δ and the linear dependence determined from simulation, we estimate a total uncertainty from the a.c. Stark effect to be 3.8 × 10⁻²¹ for Lu-1 and 2.5 × 10⁻²⁰ for Lu-2 for T_R = 5 s.

a.c. Quadrupole, RF

The oscillating RF quadrupole field couples off-resonantly to hyperfine transitions giving rise to an a.c. quadrupole shift that is not cancelled by hyperfine averaging⁶⁵. We estimate this shift by assuming an ideal linear Paul trap with RF potential of the form Φ(x, y, z) = ϵ(x² − y²) in the principal-axis frame, where the electric field gradient \({\epsilon }=\frac{m{\varOmega }_{\mathrm{RF}}{\omega }_{{\rm{r}}}}{e\sqrt{2}}\) is determined by the radial pseudo-potential confinement frequency ω_r ≈ 2π × 1,100 kHz for Lu-1 and 510 kHz for Lu-2. The magnetic field in both chambers is aligned to approximately 33(3)° with respect to the ion trap axis (z). Under these conditions, we evaluate⁶⁵ a relative clock shift of the HA frequency of −2.0 × 10⁻²¹ for Lu-1 and −5.7 × 10⁻²² for Lu-2.

As shown schematically in Fig. 1a, the differential path length is actively stabilized with reference to retroreflecting mirrors below the respective experimental chambers that are near the table surface and approximately 20 cm below the trapped ions. The differential phase stability is characterized by out-of-loop measurement of the optical phase with the clock beams directed to a common beam splitter instead of the ions, requiring approximately 2 m of additional unstabilized optical path length. When simultaneously switching on AOMs 1a and 2a, as for the Ramsey pulses in the comparison interrogation sequence, a differential phase chirp is induced by the lock circuitry as shown in Extended Data Fig. 10a. This is well modelled by a damped harmonic oscillation with the fit parameters shown in Extended Data Fig. 10a. By straightforward extension of the analysis given in ref. ⁶⁶, we evaluate a −4.8(8) × 10⁻²¹ systematic shift to the difference frequency between Lu-1 and Lu-2 for T_R = 5 s.

We consider two sources of first-order Doppler shift (FODS): drift of the ion with respect to the trap along the clock probe (vertical) direction and differential drifts of the optical phase.

The measured micromotion modulation depth β_m,3 can be directly related to the displacement of the ion from the trap RF null along the clock interrogation direction²². For small modulation depth, the relationship is linear with scale factors of approximately 0.77 μm/β_m,3 for Lu-1 and 2.0 μm/β_m,3 for Lu-2. Extended Data Fig. 4d shows the evaluated FODS shift for each comparison measurement, where the average velocity is estimated from the inferred displacement between the EMM measurements. The weighted FODS for the entire campaign is estimated to be 2.4(3.5) × 10⁻²³ for Lu-1, −0.3(2.0) × 10⁻²² for Lu-2 and 0.5(2.0) × 10⁻²² for the difference.

Using the same measurement setup as for AOM chirp, we measure the out-of-loop differential phase for long durations for which we observed an instability asymptote of 3.0 × 10⁻¹⁶ (τ/s)⁻¹ as shown in Extended Data Fig. 10b (blue). For typical servo operation at 87% interrogation duty cycle, we evaluate a contribution of 2.0 × 10⁻¹⁷ (τ/s)^−1/2 (orange). This is well below the comparison instability on all timescales. When operating as a single clock, the path length is actively stabilized from the laser to reference mirror instead of differentially between chambers, but can be expected to contribute at the same level or lower because of the shorter unstabilized path.

The RF synthesizers used for the AOMs are based on the AD9912 chip and have approximately 7 μHz resolution. This contributes about 1 × 10⁻²⁰ fractional uncertainty, half of the minimum step size, to the comparison servo. The synthesized microwave frequencies, f₁ and f₂, were identical and generated from the same sources for Lu-1 and Lu-2, set to the nearest 1 mHz. All synthesizers are referenced to a common hydrogen maser, which is calibrated to ≲2 × 10⁻¹⁵ at the start of the measurement campaign⁶⁷. For the comparison measurement, any error in maser accuracy is common mode and contributes no uncertainty to the difference.

However, we note that when realizing a Lu⁺ standard using HA–HR spectroscopy as , the accuracy of the synthesized microwaves f_k is also not a limitation to optical clock accuracy because the servo action steers the laser to compensate any errors in the microwave fields such that this specific linear combination of the fields is the HA frequency. The HA–HR servo is sensitive to the instability of both the laser and microwave fields at the same level, in Hz. In practice, this places a more relaxed requirement on the microwaves as compared with the laser; for example, for a laser with fractional instability , the microwaves need only a fractional stability better than \(\frac{\delta {f}_{k}}{{f}_{k}}\lesssim \frac{\delta {f}_{{\rm{L}}}}{{f}_{k}} \sim 1{0}^{-12}\) to not limit the clock instability.

Gravitational redshift

The differential redshift between the ions is given by , where g ≈ 9.776 m s⁻² is the local gravitational acceleration and δh = h₁ − h₂ is the height difference of the ions. We measure the heights of the ions relative to a laser levelling assembly fixed to the optical table, which is precisely aligned to each ion by the same imaging optics used for state detection. We determine the ion heights h₁ = 7.44(10) mm and h₂ = 11.41(10) mm. Measurements of the table levelling by a precision spirit level with 0.02 mm m⁻¹ resolution are limited by the table surface flatness but bound the table tilt to <0.2 mm m⁻¹. The ion traps are separated by 1.2 m, and table levelling contributes the dominant uncertainty of 240 μm. The height difference of the ions is evaluated to be δh = −3.97(27) mm, corresponding to a −4.31(29) × 10⁻¹⁹ differential shift.

Clock stability

The stability of a single-ion optical clock is quantum-projection-noise limited at the longest achievable interrogation time, which is set by either the optical coherence of the local oscillator or the atomic coherence, whichever is shorter. State-of-the-art cavity-stabilized lasers reach fractional instabilities in the low 10⁻¹⁷ at 1 s (refs. ^28,33) and have enabled demonstrated single-ion clock instabilities as low as 3.5 × 10⁻¹⁶ (τ/s)^−1/2 (ref. ¹²). Although our laboratory lacks a local oscillator of this quality with which to challenge ion clock stability at this level directly, our comparison does probe the limits set by atomic coherence, allowing us to infer the technical limits to interrogation time independent of any clock laser. Comparison of two references by correlation spectroscopy carries a factor of √2 higher instability than a comparison of independent clocks, and a factor of 2 higher than a single standalone clock, for equal interrogation time. Insofar as the present comparison is limited by uncorrelated magnetic field noise, these systems, when combined with state-of-the-art laser stabilization, would support a single-clock instability of about 2.4 × 10⁻¹⁶ (τ/s)^−1/2. Because the ³D₁ lifetime imposes no practical limit, substantially longer interrogation times remain possible, contingent on overcoming magnetic field noise, ion trap heating and background gas collisions. Reducing the d.c. magnetic field lowers the linear sensitivity to field noise at the cost of smaller Zeeman splittings; nevertheless, operation at 25 μT is achievable with modest changes to our laser configuration. By our assessment, this would permit interrogation times up to 30 s in the present system. Even longer coherence times can be anticipated with the use of magnetic shielding and cryogenics, albeit at the cost of increased experimental complexity.

Source data are provided with this paper. Any additional data are available from the corresponding authors upon request.

All codes to analyse the data in this study are available from the corresponding authors upon request.

This project was supported by the National Research Foundation, Singapore, through the National Quantum Office, hosted in A*STAR, under its Quantum Engineering Programme 3.0 Funding Initiative (W25Q3D0007) and under its Centre for Quantum Technologies Funding Initiative (S24Q2d0009).

1. These authors contributed equally: K. J. Arnold, M. D. K. Lee

Authors and Affiliations

1. Centre for Quantum Technologies, National University of Singapore, Singapore, Singapore

   K. J. Arnold, M. D. K. Lee, Qi Zhao, Qichen Qin, Zhao Zhang, N. Jayjong & M. D. Barrett

2. Temasek Laboratories, National University of Singapore, Singapore, Singapore

   K. J. Arnold

3. Department of Physics, National University of Singapore, Singapore, Singapore

   M. D. Barrett

K.J.A. and M.D.K.L. performed the experiments and independently analysed the data. K.J.A. prepared the manuscript. M.D.B. supervised the project. K.J.A., M.D.K.L., Q.Z., Q.Q., Z.Z., N.J. and M.D.B. contributed to the design of the experiment, discussion of results and revision of the manuscript.

Corresponding authors

Correspondence to K. J. Arnold or M. D. Barrett.

Competing interests

The authors declare no competing interests.

Peer review information

Nature thanks the anonymous reviewers for their contribution to the peer review of this work. Peer reviewer reports are available.

Publisher’s note Springer Nature remains neutral with regard to jurisdictional claims in published maps and institutional affiliations.

Open Access This article is licensed under a Creative Commons Attribution-NonCommercial-NoDerivatives 4.0 International License, which permits any non-commercial use, sharing, distribution and reproduction in any medium or format, as long as you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons licence, and indicate if you modified the licensed material. You do not have permission under this licence to share adapted material derived from this article or parts of it. The images or other third party material in this article are included in the article’s Creative Commons licence, unless indicated otherwise in a credit line to the material. If material is not included in the article’s Creative Commons licence and your intended use is not permitted by statutory regulation or exceeds the permitted use, you will need to obtain permission directly from the copyright holder. To view a copy of this licence, visit http://creativecommons.org/licenses/by-nc-nd/4.0/.

Reprints and permissions