# Searching for time-dependent high-energy neutrino emission from X-ray
binaries with IceCube

The IceCube Collaboration
(a complete list of authors can be found at the end of the proceedings)

###### Abstract

X-ray binaries are long-standing source candidates of Galactic cosmic rays and
neutrinos. The compact object in a binary system can be the site for cosmic-
ray acceleration, while high-energy neutrinos can be produced by the
interactions of cosmic rays in the jet of the compact object, the stellar
wind, or the atmosphere of the companion star. We report a time-dependent
study of high-energy neutrinos from X-ray binaries with IceCube using 7.5
years of muon neutrino data and X-ray observations. In the absence of
significant correlation, we report upper limits on the neutrino fluxes from
these sources and provide a comparison with theoretical predictions.

Corresponding authors: Qinrui Liu1∗, Ali Kheirandish2
1 Wisconsin IceCube Particle Astrophysics Center (WIPAC) and Department of
Physics, University of Wisconsin-Madison, Madison, WI 53706, USA
2 Department of Physics; Department of Astronomy & Astrophysics; Center for
Multimessenger Astrophysics, Institute for Gravitation & the Cosmos, The
Pennsylvania State University, University Park, PA 16802, USA
∗ Presenter

##  1 Introduction

Cosmic rays (CRs) with energies up to several PeV, the "knee" in the CR
spectrum, are believed to be of Galactic origin. However, where and how these
CRs are accelerated remains an open question. Interactions of very high energy
CRs in the Galaxy will lead to the production of pions, which subsequently
decay into gamma rays and neutrinos, with energies reaching hundreds of TeV.
As electromagnetic processes could also contribute to high-energy gamma-ray
emission, only the detection of high-energy neutrinos would be a smoking gun
for such CR interactions (i.e., hadronic interactions) as they are the only
way to produce neutrinos. The sources of the vast majority of high-energy
neutrinos detected by IceCube are yet to be identified. The isotropic
distribution of high-energy neutrino’s arrival direction suggests dominant
contributions from extragalactic sources. The Galactic contribution to the
diffuse neutrino flux is constrained to $\sim$14% above 1 TeV [1]. Studies
have been conducted to identify Galactic point-like sources, extended regions,
and the diffuse emission produced by CRs interacting with the interstellar
medium. Nevertheless, recent searches for correlations do not show remarkable
signals yet [2, 3, 4]. X-ray binaries (XRBs) are binary systems consisting of
a compact object (neutron star (NS) or black hole (BH)) and a non-compact
companion star. These systems are bright in X-rays and sometimes in gamma
rays. XRBs have been proposed as sites of CR acceleration and hadronic
interactions since the 1980s. XRBs with jets, often regarded as a smaller
version of quasars and referred to as microquasars, have been widely discussed
in the context of hadronic processes in jets. Protons can be accelerated in
the jet, and pions are generated through interactions with the external
radiation field of the accretion disk and/or internal synchrotron photons.
Other discussions focus on hadronuclear interactions, e.g., jet-cloud/wind
interactions when the jet is traversing the matter field of the ejected clouds
or stellar wind from the companion star. For other XRBs where there is no
collimated beam present, hadronic interactions can happen in a wider shocked
region. CR acceleration can take place in the magnetosphere of a spinning NS
and CRs can then further interact with matter from either the accretion disk
or the companion star. See e.g. [5, 6, 7, 8] for theories of neutrino
production in XRBs. Some XRBs have been observed at TeV energies, which
illustrates the capability of these sources to accelerate particles to high
energies.

XRBs are known for their outburst and periodic emission. Thus, it is
reasonable to hypothesize that the possible neutrino emission is related to
either the periodicity or the X-ray outburst activity, which might stem from a
change in the power or target material. Time-dependent analyses can be
performed based on such hypotheses, which benefit from the suppression of the
background, which is dominated by the atmospheric neutrino flux. Both time-
integrated and time-dependent analyses searching for high-energy neutrino
emission have been performed by IceCube and ANTARES, e.g., [9, 10], without
significant detection. Here, we present a study focusing on XRBs using the
IceCube muon track data searching for correlation with the X-ray outburst and
persistent emission of the possible neutrino flux from XRBs, covering an ample
list of sources.

##  2 Analysis

Figure 1: 90% sensitivity and 5$\sigma$ discovery potential of a flaring
source V404 Cyg when varying the threshold (Bayesian blocked light curve in
Fig. 2 ) and a comparison to the time-integrated case, which indicates an
improvement in sensitivity. The spectrum shown here is $E^{-2}$.

This search uses an unbinned maximum likelihood method, which follows the one
described in [11, 12], to seek an excess of neutrino events (signal) above the
background. In both the time-dependent and the time-integrated analyses, the
likelihood function describing the signal includes both spatial clustering and
energy information. In the time-dependent analysis, a unique temporal term is
incorporated in the likelihood, which incorporates a correlation between
neutrinos and X-ray light curves. As the majority of the data is expected to
be background events that are uniform in time, the likelihood function of the
background is constructed with time-randomized data for the time-dependent
analysis and right ascension randomized data for the time-integrated analysis.
The test statistic is obtained by maximizing the likelihood function w.r.t. a
set of parameters, which include the number of signal events ($n_{s}$) and the
spectral index ($\gamma$) for both analyses. For the time-dependent analysis,
in addition to $n_{s}$ and $\gamma$, time-related parameters introduced are
the threshold of a light curve $f_{th}$ for picking flares and the time lag
$T_{lag}$ between the X-ray and the neutrino emission.

The time-dependent analysis focuses on searching for a correlation between the
neutrino emission and the X-ray activity of a source. For this purpose, hard
X-ray light curves are used to construct the time probability density function
(PDF). Light curves are obtained from hard X-ray data reported by Swift/BAT in
the energy range 15-50 keV
111https://swift.gsfc.nasa.gov/results/transients/index.html[13] and MAXI in
the energy range 10-20 keV 222http://maxi.riken.jp/top/slist.html [14]. The
X-ray light curve data are binned in days, and a Bayesian block algorithm is
applied to find the optimal segmentation of the data and identify flares [15].
After the light curves are divided into blocks, the value of each block can be
fitted as a constant, taking into account the uncertainty of each data point.
The normalized blocked light curves then act as the temporal PDF. Fig. 1 shows
the sensitivity of the time-dependent analysis compared to the sensitivity of
the time-integrated analysis from the direction of V404 Cyg, where the
expected improvement is demonstrated.

The sources studied are from the Galactic high-mass XRB (HMXB) catalog [16]
and the Galactic low-mass XRB (LMXB) catalog [17], which include 301 sources.
TeV sources from TeVCat 333http://tevcat.uchicago.edu [18] which are not in
the HMXB or LMXB catalog are added. Starting from the initial source list,
sources without available Swift/BAT or MAXI hard X-ray light curves are
removed. As we are only interested in sources active with flaring or variable
behaviors in X-ray, the variability and excess variance of the light curves
are evaluated such that sources with weak emission are taken out. This step is
applied only to the X-ray data in the time frame overlapped by the neutrino
data sample. If both the Swift/BAT and MAXI light curves pass the selection
criteria, the Swift/BAT data is preferred to be used for one source. After
this selection, there are 102 sources in the initial source list left to be
analyzed.

Figure 2: The temporal PDF before normalization (Bayesian blocks) and the
event distribution within $1.5^{\circ}$ around V404 Cyg in the data sample of
2015 at the time indicated by MJD. The Bayesian blocks have been shifted by
the best-fit time lag (-0.5 days) and the dashed gray line indicates the best-
fit threshold (0.011 cm-3 s-1). Vertical lines represent neutrino events. The
color shows the energy proxy while the height shows the weight of each event
in the likelihood function.

We complement the study with a time-integrated search for neutrino signals
from four notable sources: Cyg X-3, LS 5039, LSI +61 303, and SS 433.
Additionally, two time-integrated stacking tests are conducted for
microquasars and TeV sources separately, with the method used in [4] and an
equal weighting scheme when considering the relative contribution of each
source.

For all searches, we use 7.5 years of all-sky muon track data collected
between 2011-05-13 and 2018-10-14, corresponding to a livetime of 2711 days.
The data sample being used consists of high-quality through-going muon track
events from the entire sky, yielding a total of 1502612 events. Details of the
data sample are described in [19].

##  3 Results & Discussion

Analysis | Name | TS | $\hat{n}_{s}$ | $\hat{\gamma}$ |  $p$-value | 90% CL upper limits
---|---|---|---|---|---|---
Flare | V404 Cyg | 8.3 | 5.4 | 4.0 | 0.754 (0.014) | 0.91
Time-integrated | Cyg X-3 | 6.8 | 44.6 | 3.3 | 0.036 (0.009) | 1.51
TeV XRB stacking | - | 0.1 | 7.7 | 3.5 | 0.587 | 1.22
microquasar stacking | - | 0 | 0 | - | 1 | 7.32
Table 1: The most significant source in the flare/time-integrated analysis
with the TS, and the best-fitted $\hat{n}_{s}$ and $\hat{\gamma}$. Both post-
trial and pre-trial (bracketed) p-values are shown. The results of the 2
stacking tests are also listed. The 90% CL upper limits are parameterized as
$dN_{\nu_{\mu}+\bar{\nu}_{\mu}}/dE_{\nu}=\mathcal{F}_{\nu_{\mu}+\bar{\nu_{\mu}}}\left(E_{\nu}/\rm{TeV}\right)^{-2}\,\cdot
10^{-4}\;\rm{TeV^{-1}cm^{-1}}$ for the flare analysis and
$dN_{\nu_{\mu}+\bar{\nu}_{\mu}}/dE_{\nu}=\phi_{\nu_{\mu}+\bar{\nu_{\mu}}}\left(E_{\nu}/\rm{TeV}\right)^{-2}\,\cdot
10^{-12}\;\rm{TeV^{-1}cm^{-1}s^{-1}}$ for time-integrated analyses.  Figure 3:
The relation between the time-integrated flux at 1TeV and the flaring time for
V404 Cyg. The dashed black line is the flaring time converted from the best-
fit threshold and the red triangle shows the 90% CL upper limit. The orange
line is the 5$\sigma$ discovery potential in IceCube. Purple lines illustrate
the estimated sensitivity at 90% CL and 5$\sigma$ discovery potential in
IceCube-Gen2. The shaded regions are the time-integrated neutrino flux
prediction assuming an E-2 spectrum with an energy cutoff at 100 TeV estimated
following the jet model [20]. The uncertainties are from flux densities in
different frequencies in VLA radio measurements during the flaring time in
2015. The two colors correspond to varying the energy fraction of the jet
carried by accelerated protons $\eta_{p}$. Figure 4:  Red and purple lines
indicate a comparison between current upper limits and estimated 10 yr
sensitivity (light) & discovery potential (dark) in IceCube-Gen2 for Cyg X-3.
As the high-energy neutrino events nearby cut at several TeV, an exponential
cutoff at 5 TeV is also applied for computing upper limits. The shaded regions
show predictions of $pp$ [21] and $p\gamma$ [22] scenarios. The inclusion of a
cutoff is also to be compared to the shaded pink region which includes a
cutoff of CR energy at 100 TeV with the spectral index ranging from 2.4-2.7.
3.25 corresponds to the best-fit spectral index. The gray shaded region shows
the uncertainty from the collision radius.

In the search for correlation of high-energy neutrinos and the flaring
activity of XRBs, the lowest p-value is found for the signal events from the
microquasar V404 Cyg, a low-mass BH XRB, with a pre-trial p-value of 0.014.
However, the $p$-value increases to 0.754 after taking into account the trials
for the number of sources in the catalog. V404 Cyg underwent a major X-ray
flaring episode in 2015. There are 5 sub-TeV neutrino events within
$1.5^{\circ}$ of the source during the time of this flare, and the best-fit
threshold indicates a time duration of 11 days, as shown in Fig. 2. This giant
flare was observed with a duration of approximately 13 days by Swift/BAT [23].

In the time-integrated analysis, both the tests on individual sources and
stacked search find no signal with sufficient statistical significance. The
prominent excess in the point source search is found for Cyg X-3, which
exhibits pre-trial p-value of $9\times 10^{-3}$, leading to a post-trial
$p$-value 0.036 after considering the 4 trials. In the flare analysis, Cyg X-3
has a pre-trial $p$-value 0.09, less significant than the time-integrated
results. Within $1^{\circ}$ around the source location, there are 44 events
above 1 TeV, and the most energetic one among them has deposited energy about
5 TeV, leading to a soft best-fit spectrum. Since there is no significant
signal found, we set the 90% confidential level (CL) upper limits to the
neutrino flux from the sources studied. A summary is shown in Table. 1.

For microquasars, relativistic jets are expected to be the CR acceleration
sites. Possible neutrino emission is expected from the beam dump on either
radiation from the compact object itself or gas from the companion star.
Parameters for neutrino flux prediction in [20], based on the photohadronic
model of [5] can be constrained for some microquasars. Nevertheless, the
simplified estimation has large uncertainties. For V404 Cyg, the X-ray flare
in June 2015 was observed in multiple wavelengths, and the jet activity during
that outburst was studied, e.g., in [24, 25]. A simple estimation of the
neutrino flux using the jet model can be performed with the radio jet
information when the source is in an outburst state. The upper limits reported
here are compared to the time-integrated flux estimation in Fig. 3. The
collision region is estimated from the flaring duration. The values of jet
parameters in the estimation are from [24, 25], and the spectrum is assumed to
be a power-law with an index of 2 and an exponential cutoff of 100 TeV.

For Cyg X-3, one of the microquasars identified as a gamma-ray source in early
observations, many predictions have been calculated in the past decades
depending on different models for microquasars. For a comparison to the upper
limits, we take [22] and [21], which discussed the general $p\gamma$ and $pp$
scenarios based on the AGILE observation respectively, shown in Fig. 4. What
needs to be mentioned is that Cyg X-3 lies in the direction of the Cygnus X
region and is close to the Cygnus OB2 association but with a further distance
compared to the Cygnus X region. The possibility of contamination from the
Cygnus X complex cannot be excluded.

The next generation of the IceCube experiment, IceCube-Gen2, will provide a
factor of eight increase in volume [26], leading to an expected $\sim$5-time
increase in the effective area compared to IceCube, corresponding to an
improvement in sensitivity by the same order, which advances the
identification of neutrino sources. Here, we extend the study to IceCube-Gen2
and estimate the sensitivity and discovery potential for V404 Cyg, as an
example of a flaring source and Cyg X-3, for persistently emitting sources.
The estimated improvement can be seen in Fig. 3 and Fig. 4. The effective
areas of muon tracks are computed from the proposed IceCube-Gen2
configuration, and the projection is evaluated similar to that in [26] without
considering a contribution from the existing IceCube detector. In comparison
with theoretical calculations, it demonstrates the power to either identify
those sources or rule out models with IceCube in the future.

##  4 Summary

A Galactic contribution to the high-energy neutrino flux observed by IceCube
is expected. We present a study of neutrino emission from XRBs, long-standing
candidates for the Galactic sources of CRs and neutrinos. We performed a time-
dependent analysis based on the assumption of flaring neutrino emission. In
parallel, a time-integrated search is also performed on 4 notable sources and
2 stacked lists. In the absence of any significant excess, we set upper limits
on the neutrino emission in the scenarios discussed. The results of the most
significant sources in this search are compared to models of neutrino
production in XRBs. Our estimation of the improved detectability by IceCube-
Gen2 due to higher neutrino event statistics demonstrates the potential for
the future detection and presents a promising outlook of identifying Galactic
cosmic-ray accelerators in the upcoming years.

## References

  * [1] IceCube Collaboration, M. G. Aartsen et al. Astrophys. J. 849 no. 1, (2017) 67.
  * [2] ANTARES, IceCube Collaboration, A. Albert et al. Astrophys. J. Lett. 868 no. 2, (2018) L20.
  * [3] IceCube, HAWC Collaboration, A. Kheirandish and J. Wood PoS ICRC2019 (2020) 932.
  * [4] IceCube Collaboration, M. G. Aartsen et al. Astrophys. J. 898 no. 2, (2020) 117.
  * [5] A. Levinson and E. Waxman Phys. Rev. Lett. 87 (2001) 171101.
  * [6] L. A. Anchordoqui, D. F. Torres, T. P. McCauley, G. E. Romero, and F. A. Aharonian Astrophys. J. 589 (2003) 481–486.
  * [7] G. E. Romero, D. F. Torres, M. M. K. Bernado, and I. F. Mirabel Astron. Astrophys. 410 (2003) L1–L4.
  * [8] W. Bednarek Astrophys. J. 631 (2005) 466.
  * [9] IceCube Collaboration, M. G. Aartsen et al. Astrophys. J. 807 no. 1, (2015) 46.
  * [10] A. Albert et al. JCAP 04 (2017) 019.
  * [11] J. Braun, J. Dumm, F. De Palma, C. Finley, A. Karle, and T. Montaruli Astropart. Phys. 29 (2008) 299–305.
  * [12] J. Braun, M. Baker, J. Dumm, C. Finley, A. Karle, and T. Montaruli Astropart. Phys. 33 (2010) 175–181.
  * [13] H. A. Krimm et al. Astrophys. J. Suppl. 209 (2013) 14.
  * [14] M. Matsuoka, K. Kawasaki, S. Ueno, H. Tomida, M. Kohama, M. Suzuki, Y. Adachi, M. Ishikawa, T. Mihara, M. Sugizaki, et al. Publications of the Astronomical Society of Japan 61 no. 5, (2009) 999–1010.
  * [15] J. D. Scargle, J. P. Norris, B. Jackson, and J. Chiang Astrophys. J. 764 (2013) 167.
  * [16] Q. Z. Liu, J. van Paradijs, and E. P. J. v. d. Heuvel Astron. Astrophys. 455 (2006) 1165.
  * [17] Q. Z. Liu, J. van Paradijs, and E. P. J. v. d. Heuvel Astron. Astrophys. 469 (2007) 807.
  * [18] S. P. Wakely and D. Horan, “Tevcat: an online catalog for very high energy gamma-ray astronomy,” in International Cosmic Ray Conference, vol. 3, pp. 1341–1344.  2008\.
  * [19] IceCube Collaboration, M. G. Aartsen et al. Astropart. Phys. 92 (2017) 30–41.
  * [20] C. Distefano, D. Guetta, E. Waxman, and A. Levinson Astrophys. J. 575 (2002) 378–383.
  * [21] N. Sahakyan, G. Piano, and M. Tavani Astrophys. J. 780 (2014) 29.
  * [22] P. Baerwald and D. Guetta Astrophys. J. 773 (2013) 159.
  * [23] A. Segreto, M. Del Santo, A. D’Aí, V. La Parola, G. Cusumano, T. Mineo, and J. Malzac The Astronomer’s Telegram 7755 (2015) 1.
  * [24] J. C. A. Miller-Jones et al. Nature 569 (2019) 374–377.
  * [25] A. J. Tetarenko et al. Mon. Not. Roy. Astron. Soc. 482 no. 3, (2019) 2950–2972.
  * [26] IceCube Gen2 Collaboration, M. G. Aartsen et al.

## Full Author List: IceCube Collaboration

R. Abbasi17, M. Ackermann59, J. Adams18, J. A. Aguilar12, M. Ahlers22, M.
Ahrens50, C. Alispach28, A. A. Alves Jr.31, N. M. Amin42, R. An14, K.
Andeen40, T. Anderson56, G. Anton26, C. Argüelles14, Y. Ashida38, S. Axani15,
X. Bai46, A. Balagopal V.38, A. Barbano28, S. W. Barwick30, B. Bastian59, V.
Basu38, S. Baur12, R. Bay8, J. J. Beatty20, 21, K.-H. Becker58, J. Becker
Tjus11, C. Bellenghi27, S. BenZvi48, D. Berley19, E. Bernardini59, 60, D. Z.
Besson34, 61, G. Binder8, 9, D. Bindig58, E. Blaufuss19, S. Blot59, M.
Boddenberg1, F. Bontempo31, J. Borowka1, S. Böser39, O. Botner57, J.
Böttcher1, E. Bourbeau22, F. Bradascio59, J. Braun38, S. Bron28, J. Brostean-
Kaiser59, S. Browne32, A. Burgman57, R. T. Burley2, R. S. Busse41, M. A.
Campana45, E. G. Carnie-Bronca2, C. Chen6, D. Chirkin38, K. Choi52, B. A.
Clark24, K. Clark33, L. Classen41, A. Coleman42, G. H. Collin15, J. M.
Conrad15, P. Coppin13, P. Correa13, D. F. Cowen55, 56, R. Cross48, C. Dappen1,
P. Dave6, C. De Clercq13, J. J. DeLaunay56, H. Dembinski42, K. Deoskar50, S.
De Ridder29, A. Desai38, P. Desiati38, K. D. de Vries13, G. de Wasseige13, M.
de With10, T. DeYoung24, S. Dharani1, A. Diaz15, J. C. Díaz-Vélez38, M.
Dittmer41, H. Dujmovic31, M. Dunkman56, M. A. DuVernois38, E. Dvorak46, T.
Ehrhardt39, P. Eller27, R. Engel31, 32, H. Erpenbeck1, J. Evans19, P. A.
Evenson42, K. L. Fan19, A. R. Fazely7, S. Fiedlschuster26, A. T. Fienberg56,
K. Filimonov8, C. Finley50, L. Fischer59, D. Fox55, A. Franckowiak11, 59, E.
Friedman19, A. Fritz39, P. Fürst1, T. K. Gaisser42, J. Gallagher37, E.
Ganster1, A. Garcia14, S. Garrappa59, L. Gerhardt9, A. Ghadimi54, C. Glaser57,
T. Glauch27, T. Glüsenkamp26, A. Goldschmidt9, J. G. Gonzalez42, S. Goswami54,
D. Grant24, T. Grégoire56, S. Griswold48, M. Gündüz11, C. Günther1, C.
Haack27, A. Hallgren57, R. Halliday24, L. Halve1, F. Halzen38, M. Ha Minh27,
K. Hanson38, J. Hardin38, A. A. Harnisch24, A. Haungs31, S. Hauser1, D.
Hebecker10, K. Helbing58, F. Henningsen27, E. C. Hettinger24, S. Hickford58,
J. Hignight25, C. Hill16, G. C. Hill2, K. D. Hoffman19, R. Hoffmann58, T.
Hoinka23, B. Hokanson-Fasig38, K. Hoshina38, 62, F. Huang56, M. Huber27, T.
Huber31, K. Hultqvist50, M. Hünnefeld23, R. Hussain38, S. In52, N. Iovine12,
A. Ishihara16, M. Jansson50, G. S. Japaridze5, M. Jeong52, B. J. P. Jones4, D.
Kang31, W. Kang52, X. Kang45, A. Kappes41, D. Kappesser39, T. Karg59, M.
Karl27, A. Karle38, U. Katz26, M. Kauer38, M. Kellermann1, J. L. Kelley38, A.
Kheirandish56, K. Kin16, T. Kintscher59, J. Kiryluk51, S. R. Klein8, 9, R.
Koirala42, H. Kolanoski10, T. Kontrimas27, L. Köpke39, C. Kopper24, S.
Kopper54, D. J. Koskinen22, P. Koundal31, M. Kovacevich45, M. Kowalski10, 59,
T. Kozynets22, E. Kun11, N. Kurahashi45, N. Lad59, C. Lagunas Gualda59, J. L.
Lanfranchi56, M. J. Larson19, F. Lauber58, J. P. Lazar14, 38, J. W. Lee52, K.
Leonard38, A. Leszczyńska32, Y. Li56, M. Lincetto11, Q. R. Liu38, M.
Liubarska25, E. Lohfink39, C. J. Lozano Mariscal41, L. Lu38, F. Lucarelli28,
A. Ludwig24, 35, W. Luszczak38, Y. Lyu8, 9, W. Y. Ma59, J. Madsen38, K. B. M.
Mahn24, Y. Makino38, S. Mancina38, I. C. Mariş12, R. Maruyama43, K. Mase16, T.
McElroy25, F. McNally36, J. V. Mead22, K. Meagher38, A. Medina21, M. Meier16,
S. Meighen-Berger27, J. Micallef24, D. Mockler12, T. Montaruli28, R. W.
Moore25, R. Morse38, M. Moulai15, R. Naab59, R. Nagai16, U. Naumann58, J.
Necker59, L. V. Nguyễn24, H. Niederhausen27, M. U. Nisa24, S. C. Nowicki24, D.
R. Nygren9, A. Obertacke Pollmann58, M. Oehler31, A. Olivas19, E.
O’Sullivan57, H. Pandya42, D. V. Pankova56, N. Park33, G. K. Parker4, E. N.
Paudel42, L. Paul40, C. Pérez de los Heros57, L. Peters1, J. Peterson38, S.
Philippen1, D. Pieloth23, S. Pieper58, M. Pittermann32, A. Pizzuto38, M.
Plum40, Y. Popovych39, A. Porcelli29, M. Prado Rodriguez38, P. B. Price8, B.
Pries24, G. T. Przybylski9, C. Raab12, A. Raissi18, M. Rameez22, K. Rawlins3,
I. C. Rea27, A. Rehman42, P. Reichherzer11, R. Reimann1, G. Renzi12, E.
Resconi27, S. Reusch59, W. Rhode23, M. Richman45, B. Riedel38, E. J. Roberts2,
S. Robertson8, 9, G. Roellinghoff52, M. Rongen39, C. Rott49, 52, T. Ruhe23, D.
Ryckbosch29, D. Rysewyk Cantu24, I. Safa14, 38, J. Saffer32, S. E. Sanchez
Herrera24, A. Sandrock23, J. Sandroos39, M. Santander54, S. Sarkar44, S.
Sarkar25, K. Satalecka59, M. Scharf1, M. Schaufel1, H. Schieler31, S.
Schindler26, P. Schlunder23, T. Schmidt19, A. Schneider38, J. Schneider26, F.
G. Schröder31, 42, L. Schumacher27, G. Schwefer1, S. Sclafani45, D. Seckel42,
S. Seunarine47, A. Sharma57, S. Shefali32, M. Silva38, B. Skrzypek14, B.
Smithers4, R. Snihur38, J. Soedingrekso23, D. Soldin42, C. Spannfellner27, G.
M. Spiczak47, C. Spiering59, 61, J. Stachurska59, M. Stamatikos21, T.
Stanev42, R. Stein59, J. Stettner1, A. Steuer39, T. Stezelberger9, T.
Stürwald58, T. Stuttard22, G. W. Sullivan19, I. Taboada6, F. Tenholt11, S.
Ter-Antonyan7, S. Tilav42, F. Tischbein1, K. Tollefson24, L. Tomankova11, C.
Tönnis53, S. Toscano12, D. Tosi38, A. Trettin59, M. Tselengidou26, C. F.
Tung6, A. Turcati27, R. Turcotte31, C. F. Turley56, J. P. Twagirayezu24, B.
Ty38, M. A. Unland Elorrieta41, N. Valtonen-Mattila57, J. Vandenbroucke38, N.
van Eijndhoven13, D. Vannerom15, J. van Santen59, S. Verpoest29, M. Vraeghe29,
C. Walck50, T. B. Watson4, C. Weaver24, P. Weigel15, A. Weindl31, M. J.
Weiss56, J. Weldert39, C. Wendt38, J. Werthebach23, M. Weyrauch32, N.
Whitehorn24, 35, C. H. Wiebusch1, D. R. Williams54, M. Wolf27, K. Woschnagg8,
G. Wrede26, J. Wulff11, X. W. Xu7, Y. Xu51, J. P. Yanez25, S. Yoshida16, S.
Yu24, T. Yuan38, Z. Zhang51

1 III. Physikalisches Institut, RWTH Aachen University, D-52056 Aachen,
Germany
2 Department of Physics, University of Adelaide, Adelaide, 5005, Australia
3 Dept. of Physics and Astronomy, University of Alaska Anchorage, 3211
Providence Dr., Anchorage, AK 99508, USA
4 Dept. of Physics, University of Texas at Arlington, 502 Yates St., Science
Hall Rm 108, Box 19059, Arlington, TX 76019, USA
5 CTSPS, Clark-Atlanta University, Atlanta, GA 30314, USA
6 School of Physics and Center for Relativistic Astrophysics, Georgia
Institute of Technology, Atlanta, GA 30332, USA
7 Dept. of Physics, Southern University, Baton Rouge, LA 70813, USA
8 Dept. of Physics, University of California, Berkeley, CA 94720, USA
9 Lawrence Berkeley National Laboratory, Berkeley, CA 94720, USA
10 Institut für Physik, Humboldt-Universität zu Berlin, D-12489 Berlin,
Germany
11 Fakultät für Physik & Astronomie, Ruhr-Universität Bochum, D-44780 Bochum,
Germany
12 Université Libre de Bruxelles, Science Faculty CP230, B-1050 Brussels,
Belgium
13 Vrije Universiteit Brussel (VUB), Dienst ELEM, B-1050 Brussels, Belgium
14 Department of Physics and Laboratory for Particle Physics and Cosmology,
Harvard University, Cambridge, MA 02138, USA
15 Dept. of Physics, Massachusetts Institute of Technology, Cambridge, MA
02139, USA
16 Dept. of Physics and Institute for Global Prominent Research, Chiba
University, Chiba 263-8522, Japan
17 Department of Physics, Loyola University Chicago, Chicago, IL 60660, USA
18 Dept. of Physics and Astronomy, University of Canterbury, Private Bag 4800,
Christchurch, New Zealand
19 Dept. of Physics, University of Maryland, College Park, MD 20742, USA
20 Dept. of Astronomy, Ohio State University, Columbus, OH 43210, USA
21 Dept. of Physics and Center for Cosmology and Astro-Particle Physics, Ohio
State University, Columbus, OH 43210, USA
22 Niels Bohr Institute, University of Copenhagen, DK-2100 Copenhagen, Denmark
23 Dept. of Physics, TU Dortmund University, D-44221 Dortmund, Germany
24 Dept. of Physics and Astronomy, Michigan State University, East Lansing, MI
48824, USA
25 Dept. of Physics, University of Alberta, Edmonton, Alberta, Canada T6G 2E1
26 Erlangen Centre for Astroparticle Physics, Friedrich-Alexander-Universität
Erlangen-Nürnberg, D-91058 Erlangen, Germany
27 Physik-department, Technische Universität München, D-85748 Garching,
Germany
28 Département de physique nucléaire et corpusculaire, Université de Genève,
CH-1211 Genève, Switzerland
29 Dept. of Physics and Astronomy, University of Gent, B-9000 Gent, Belgium
30 Dept. of Physics and Astronomy, University of California, Irvine, CA 92697,
USA
31 Karlsruhe Institute of Technology, Institute for Astroparticle Physics,
D-76021 Karlsruhe, Germany
32 Karlsruhe Institute of Technology, Institute of Experimental Particle
Physics, D-76021 Karlsruhe, Germany
33 Dept. of Physics, Engineering Physics, and Astronomy, Queen’s University,
Kingston, ON K7L 3N6, Canada
34 Dept. of Physics and Astronomy, University of Kansas, Lawrence, KS 66045,
USA
35 Department of Physics and Astronomy, UCLA, Los Angeles, CA 90095, USA
36 Department of Physics, Mercer University, Macon, GA 31207-0001, USA
37 Dept. of Astronomy, University of Wisconsin–Madison, Madison, WI 53706, USA
38 Dept. of Physics and Wisconsin IceCube Particle Astrophysics Center,
University of Wisconsin–Madison, Madison, WI 53706, USA
39 Institute of Physics, University of Mainz, Staudinger Weg 7, D-55099 Mainz,
Germany
40 Department of Physics, Marquette University, Milwaukee, WI, 53201, USA
41 Institut für Kernphysik, Westfälische Wilhelms-Universität Münster, D-48149
Münster, Germany
42 Bartol Research Institute and Dept. of Physics and Astronomy, University of
Delaware, Newark, DE 19716, USA
43 Dept. of Physics, Yale University, New Haven, CT 06520, USA
44 Dept. of Physics, University of Oxford, Parks Road, Oxford OX1 3PU, UK
45 Dept. of Physics, Drexel University, 3141 Chestnut Street, Philadelphia, PA
19104, USA
46 Physics Department, South Dakota School of Mines and Technology, Rapid
City, SD 57701, USA
47 Dept. of Physics, University of Wisconsin, River Falls, WI 54022, USA
48 Dept. of Physics and Astronomy, University of Rochester, Rochester, NY
14627, USA
49 Department of Physics and Astronomy, University of Utah, Salt Lake City, UT
84112, USA
50 Oskar Klein Centre and Dept. of Physics, Stockholm University, SE-10691
Stockholm, Sweden
51 Dept. of Physics and Astronomy, Stony Brook University, Stony Brook, NY
11794-3800, USA
52 Dept. of Physics, Sungkyunkwan University, Suwon 16419, Korea
53 Institute of Basic Science, Sungkyunkwan University, Suwon 16419, Korea
54 Dept. of Physics and Astronomy, University of Alabama, Tuscaloosa, AL
35487, USA
55 Dept. of Astronomy and Astrophysics, Pennsylvania State University,
University Park, PA 16802, USA
56 Dept. of Physics, Pennsylvania State University, University Park, PA 16802,
USA
57 Dept. of Physics and Astronomy, Uppsala University, Box 516, S-75120
Uppsala, Sweden
58 Dept. of Physics, University of Wuppertal, D-42119 Wuppertal, Germany
59 DESY, D-15738 Zeuthen, Germany
60 Università di Padova, I-35131 Padova, Italy
61 National Research Nuclear University, Moscow Engineering Physics Institute
(MEPhI), Moscow 115409, Russia
62 Earthquake Research Institute, University of Tokyo, Bunkyo, Tokyo 113-0032,
Japan

### Acknowledgements

USA – U.S. National Science Foundation-Office of Polar Programs, U.S. National
Science Foundation-Physics Division, U.S. National Science Foundation-EPSCoR,
Wisconsin Alumni Research Foundation, Center for High Throughput Computing
(CHTC) at the University of Wisconsin–Madison, Open Science Grid (OSG),
Extreme Science and Engineering Discovery Environment (XSEDE), Frontera
computing project at the Texas Advanced Computing Center, U.S. Department of
Energy-National Energy Research Scientific Computing Center, Particle
astrophysics research computing center at the University of Maryland,
Institute for Cyber-Enabled Research at Michigan State University, and
Astroparticle physics computational facility at Marquette University; Belgium
– Funds for Scientific Research (FRS-FNRS and FWO), FWO Odysseus and Big
Science programmes, and Belgian Federal Science Policy Office (Belspo);
Germany – Bundesministerium für Bildung und Forschung (BMBF), Deutsche
Forschungsgemeinschaft (DFG), Helmholtz Alliance for Astroparticle Physics
(HAP), Initiative and Networking Fund of the Helmholtz Association, Deutsches
Elektronen Synchrotron (DESY), and High Performance Computing cluster of the
RWTH Aachen; Sweden – Swedish Research Council, Swedish Polar Research
Secretariat, Swedish National Infrastructure for Computing (SNIC), and Knut
and Alice Wallenberg Foundation; Australia – Australian Research Council;
Canada – Natural Sciences and Engineering Research Council of Canada, Calcul
Québec, Compute Ontario, Canada Foundation for Innovation, WestGrid, and
Compute Canada; Denmark – Villum Fonden and Carlsberg Foundation; New Zealand
– Marsden Fund; Japan – Japan Society for Promotion of Science (JSPS) and
Institute for Global Prominent Research (IGPR) of Chiba University; Korea –
National Research Foundation of Korea (NRF); Switzerland – Swiss National
Science Foundation (SNSF); United Kingdom – Department of Physics, University
of Oxford.


# Multiparameter quantum metrology using strongly interacting spin systems

Min Jiang  Hefei National Laboratory for Physical Sciences at the Microscale
and Department of Modern Physics, University of Science and Technology of
China, Hefei 230026, China  CAS Key Laboratory of Microscale Magnetic
Resonance, University of Science and Technology of China, Hefei 230026, China
Synergetic Innovation Center of Quantum Information and Quantum Physics,
University of Science and Technology of China, Hefei 230026, China  Yunlan Ji
Hefei National Laboratory for Physical Sciences at the Microscale and
Department of Modern Physics, University of Science and Technology of China,
Hefei 230026, China  CAS Key Laboratory of Microscale Magnetic Resonance,
University of Science and Technology of China, Hefei 230026, China  Synergetic
Innovation Center of Quantum Information and Quantum Physics, University of
Science and Technology of China, Hefei 230026, China  Qing Li  Hefei National
Laboratory for Physical Sciences at the Microscale and Department of Modern
Physics, University of Science and Technology of China, Hefei 230026, China
CAS Key Laboratory of Microscale Magnetic Resonance, University of Science and
Technology of China, Hefei 230026, China  Synergetic Innovation Center of
Quantum Information and Quantum Physics, University of Science and Technology
of China, Hefei 230026, China  Ran Liu  Hefei National Laboratory for Physical
Sciences at the Microscale and Department of Modern Physics, University of
Science and Technology of China, Hefei 230026, China  CAS Key Laboratory of
Microscale Magnetic Resonance, University of Science and Technology of China,
Hefei 230026, China  Synergetic Innovation Center of Quantum Information and
Quantum Physics, University of Science and Technology of China, Hefei 230026,
China  Dieter Suter  Experimentelle Physik III, Universit$\ddot{\textrm{a}}$t
Dortmund, 44221 Dortmund, Germany  Xinhua Peng<EMAIL_ADDRESS>Hefei
National Laboratory for Physical Sciences at the Microscale and Department of
Modern Physics, University of Science and Technology of China, Hefei 230026,
China  CAS Key Laboratory of Microscale Magnetic Resonance, University of
Science and Technology of China, Hefei 230026, China  Synergetic Innovation
Center of Quantum Information and Quantum Physics, University of Science and
Technology of China, Hefei 230026, China

###### Abstract

Interacting quantum systems are attracting increasing interest for developing
precise metrology. In particular, the realisation that quantum-correlated
states and the dynamics of interacting systems can lead to entirely new and
unexpected phenomena have initiated an intense research effort to explore
interaction-based metrology both theoretically and experimentally. However,
the current framework of interaction-based metrology mainly focuses on single-
parameter estimations, a demonstration of multiparameter metrology using
interactions as a resource was heretofore lacking. Here we demonstrate an
interaction-based multiparameter metrology with strongly interacting nuclear
spins. We show that the interacting spins become intrinsically sensitive to
all components of a multidimensional field when their interactions are
significantly larger than their Larmor frequencies. Using liquid-state
molecules containing strongly interacting nuclear spins, we demonstrate the
proof-of-principle estimation of all three components of an unknown magnetic
field and inertial rotation. In contrast to existing approaches, the present
interaction-based multiparameter sensing does not require external reference
fields and opens a path to develop an entirely new class of multiparameter
quantum sensors.

##  I Introduction

Quantum metrology exploits quantum resources of well-controlled quantum
systems to measure small signals with high sensitivity, and has great
potential for both fundamental science and concrete applications [1, 2, 3, 4,
5, 6, 7, 8]. Interacting quantum systems have been recently recognized as an
unprecedented playground to develop precise metrology [9, 10, 11, 12, 13, 14,
15, 16], where the interactions between probe systems can be a valuable
quantum resource to improve the precision and other metrological performance,
such as the achievement of metrological sensitivity beyond the Heisenberg
limit [17, 18, 19, 20], the robustness against thermal noise [13],
criticality-enhanced metrology [16, 21, 22, 23], and the evasion of rapid
thermalization [24, 25, 11]. Moreover, the use of interactions that correlate
probe systems can relax the experimental challenge of preparing entangled
probe states [14, 26].

However, the current framework of interaction-based quantum metrology mainly
focuses on single-parameter estimations [19, 13, 18, 14, 17, 15, 24, 16, 21,
20], neglecting many possible applications under realistic conditions. Indeed,
numerous metrological applications are intrinsically multiparameter, such as
measuring electric, magnetic, or gravitational fields. As a consequence, the
field of multiparameter metrology has been intensively studied both
theoretically and experimentally [27, 28, 29, 30, 31, 32, 33, 34, 35, 36]. The
potential combination of interaction-based metrology and multiparameter
estimation would open up exciting possibilities for developing precise
metrology in applied and fundamental physics, such as quantum imaging [37,
38], sensor networks [39], navigation [40, 41, 42], the detection of vector
fields [32, 33, 43, 34, 35], as well as the searches for spin-gravity coupling
and ultralight dark matter axions [44, 45, 46]. Despite these appealing
features, a demonstration of multiparameter metrology using the resource of
interactions was heretofore lacking.

Figure 1: Experimental schematic of interaction-based multiparameter
metrology. a, Diagram of the experimental set-up. The polar and azimuthal
angles ($\theta$, $\phi$) parameterize the orientation of the unknown magnetic
field $\mathbf{B}(\theta,\phi)$ in the laboratory frame. The information of
$\mathbf{B}(\theta,\phi)$ is encoded into the nuclear magnetization of
molecules comprising of strongly interacting nuclear spins, which are read out
by an atomic vapour sensor. b, Basic procedure of multiparameter metrology,
including probe state preparation, encoding unknown magnetic field, and
nuclear-spin probe readout. After being polarized in a Halbach magnet, the
nuclear spins in liquid-state molecules contained in a 5-mm tube is
pneumatically shuttled to a sensing and readout region (at a distance of 1 mm
above the atomic vapour); During the pneumatic transfer, the guiding solenoid
and three pairs of Helmholtz coils (not shown, at the bottom of the guiding
solenoid) are used to generate a guiding field ($\approx 10^{-4}$ T) along an
axis (${x}$, ${y}$ or ${z}$) for initializing the spin state; The nuclear spin
magnetization of the liquid, which is read out with a 87Rb vapour sensor (c,
see Appendix A1 for more details), evolves under the unknown magnetic field;
Then, the obtained strongly-interacting nuclear-spin spectrum is used to
determine the magnitude and orientation of the magnetic field.

Here, we demonstrate an approach to interaction-based multiparameter
estimation with the use of strongly interacting nuclear spins in molecules. We
find that, as the interactions among nuclear spins become significantly larger
than their Larmor frequencies, a different regime emerges where the strongly
interacting spins can be simultaneously sensitive to all components of a
multidimensional field, such as three-dimensional magnetic field and inertial
rotation. Moreover, we observe that strong interactions between nuclear spins
can increase their quantum coherence times as long as several seconds, leading
to enhanced measurement precision. We experimentally illustrate the power of
interaction-based multiparameter estimation by applying this approach to
simultaneously estimating all three components of a magnetic field. Our work
also provides a recipe for precisely measuring magnetic pseudo-fields
originating from inertial rotation and constituting a new zero-field vector
gyroscope. The present interaction-based metrology should be generic and
promises to develop an entirely new class of multiparameter quantum sensors
using a broad range of strongly interacting quantum systems, and will impact
diverse metrological applications.

We would like to emphasize the difference of multiparameter metrology between
this work and others. Although a variety of multiparameter quantum sensors
have demonstrated the capability of measuring magnetic fields and inertial
rotations [32, 33, 43, 34, 42, 35, 47, 36], they all make use of non-
interacting probe systems, such as electron and nuclear spins. For example,
all nuclear magnetic resonance (NMR) magnetometers tested so far use only a
single spin species (for example, proton rich substances or noble gas 3He)
that yields a single resonance line at the precession frequency [48, 49, 50].
Therefore, these magnetometers measure only the magnitude of the magnetic
field, providing no information on its orientation. In conventional methods,
external reference fields are additionally required to access the orientation
information [32, 33, 34, 35]. This unavoidably produces significant static and
radio frequency fields contamination, which are classical noise and would be a
bottleneck for reaching the fundamental precision beyond the standard quantum
limit. Unlike previous works, our work uses a different scheme in which
strongly interacting spins are intrinsically sensitive to multidimensional
fields. Thus, our work can avoid the classical noise caused by previously-
required reference fields, thus opening a feasible route towards realizing
multiparameter quantum sensors with quantum noise-limited sensitivity.

##  II Interaction-based multiparameter estimation scheme

Nuclear and electron spin systems are progressively becoming a paradigm for
precision metrology. Here nuclear-spin systems are used as a platform to
develop multiparameter metrology. We focus on the estimation of all three
components of an unknown magnetic field. In contrast to existing works based
on non-interacting spins [32, 33, 43, 34, 42, 35], we employ a collection of
molecules, where each molecule contains at least two different types of
nuclear spins that interact with each other by an exchange interaction (a
Heisenberg interaction). When such nuclear spins are imposed in a magnetic
field $\mathbf{B}(\theta,\phi)$, they can be described by the Hamiltonian,

$\begin{cases}\mathcal{H}(\theta,\phi,B)=\mathcal{H}_{\textrm{spins}}+\mathcal{H}_{\textrm{int}},\\\
\mathcal{H}_{\textrm{spins}}=-\sum_{j}2\pi\gamma_{j}\mathbf{I}_{j}\cdot\mathbf{B}(\theta,\phi),\\\
\mathcal{H}_{\textrm{int}}=\sum_{j>i}^{n}2\pi{{J_{ij}}}\mathbf{I}_{i}\cdot\mathbf{I}_{j},\end{cases}$
(1)

where $\gamma_{j}$ denotes the gyromagnetic ratio of the $j$th spin (here we
set $\hbar=1$), $\mathbf{I}_{j}=(\hat{I}_{jx},\hat{I}_{jy},\hat{I}_{jz})$ is
the spin angular momentum operator of the $j$th spin,
$\mathbf{I}_{i}\cdot\mathbf{I}_{j}$ is called a Heisenberg interaction, and
${J_{ij}}$ is the strength of the Heisenberg interaction between the $i$th and
$j$th spins. As shown in Fig. 1a, the unknown magnetic field
$\mathbf{B}(\theta,\phi)=(B\sin\theta\cos\phi,B\sin\theta\sin\phi,B\cos\theta)$
is parameterized with the polar angle $\theta$, azimuthal angle $\phi$ and
magnitude $B$, which are unknown parameters to be estimated.

Our interaction-based multiparameter estimation includes three stages: probe
preparation, encoding with the strongly interacting spins, and probe readout,
as shown in Fig. 1b. The initial probe state of nuclear spins is prepared as
$\rho_{0}$, which is described in Appendix C1 for details. In the encoding
stage, the key ingredient is the use of the quantum dynamics of the strongly
interacting spins to encode the unknown parameters,

$\rho(\theta,\phi,B,t)=e^{-i\mathcal{H}(\theta,\phi,B)t}\rho_{0}e^{i\mathcal{H}(\theta,\phi,B)t},$
(2)

where $\rho(\theta,\phi,B,t)$ is the time-dependent density matrix of nuclear
spins. Then, an atomic vapour sensor as a detector continuously measures the
time-evolution of the nuclear magnetization
$\textrm{Tr}[\rho(\theta,\phi,B,t)\hat{\mathbf{O}}_{z}]$ [51, 52, 53, 54].
Here, $\hat{\mathbf{O}}_{z}=\sum_{j}\gamma_{j}{\hat{I}}_{jz}$ represents the
observable operator of $z$-magnetization. The observation values can be
Fourier transformed into a spectrum. Consequently, the information of
$\mathbf{B}(\theta,\phi)$ in the Hamiltonian $\mathcal{H}(\theta,\phi,B)$ can
be encoded into the spectrum of the interacting spins, where each spectrally
separating line corresponds to the transition between relevant eigenstates
$|\Psi_{i}\rangle$ and $|\Psi_{j}\rangle$ of $\mathcal{H}(\theta,\phi,B)$ and
its amplitude is proportional to a product of
$\langle\Psi_{i}|\rho_{0}|\Psi_{j}\rangle\langle\Psi_{j}|\hat{\mathbf{O}}_{z}|\Psi_{i}\rangle$
[see Eqs. (11) and (12) in Appendix].

Figure 2: Spectra of 13CHn molecules measured at various orientations of the
magnetic field. The initial spin state is prepared by the guiding field along
the $z$ axis. (a, b) Formic acid (H13COOH, where the acidic proton is
negligible due to rapid exchange), (c, d) formaldehyde (13CH2O), (e, f)
acetonitrile (13CH3CN) are detected parallel ($\theta=0$) and perpendicularly
($\theta=\frac{\pi}{2},\phi=0$) to the applied magnetic field. The splitting
patterns match well the transitions shown in the inset (see Appendix B1 for
details). The amplitudes of the experimental spectra (blue lines) are
asymmetric [54] and are in good agreement with those of the simulated spectra
(red lines).

To show the explicit dependence on $\theta$ and $\phi$ of the spectral
amplitudes, we introduce a convenient orthonormal coordinate system, where
$\mathbf{z^{\prime}}=\mathbf{B}/B$,
$\mathbf{x^{\prime}}=\frac{\partial\mathbf{z^{\prime}}}{\partial\theta}$, and
$\mathbf{y^{\prime}}=(\frac{1}{\sin\theta})\frac{\partial\mathbf{z^{\prime}}}{\partial\phi}$
(see Appendix C2 for details), as shown in Fig. 1a. Consequently,
$\langle\Psi_{i}|\hat{\mathbf{O}}|\Psi_{j}\rangle=\mathcal{P}(\theta,\phi)\langle\Psi_{i}|\hat{\mathbf{O}^{\prime}}|\Psi_{j}\rangle$,
where $\mathcal{P}(\theta,\phi)$, depending on $\theta$ and $\phi$, is a
$3\times 3$ transformation matrix [see Eq. (15)] and
$\hat{\mathbf{O}}^{\prime}_{l^{\prime}}=\sum_{j}\gamma_{j}\hat{I}_{jl^{\prime}}$
($l^{\prime}=x^{\prime},y^{\prime},z^{\prime}$). Importantly, we note that the
transition amplitudes in the new defined coordinates
$\langle\Psi_{i}|\hat{\mathbf{O}}^{\prime}|\Psi_{j}\rangle$ are independent of
the magnetic field orientation ($\theta,\phi$) [see Eq. (17)], where
$\langle\Psi_{i}|\hat{\mathbf{O}}^{\prime}_{z^{\prime}}|\Psi_{j}\rangle$ and
$\langle\Psi_{i}|\hat{\mathbf{O}}^{\prime}_{x^{\prime},y^{\prime}}|\Psi_{j}\rangle$
respectively correspond to the observation parallel or perpendicular to the
magnetic field and hence the zero- or single-quantum resonance line. Thus, the
observed amplitudes of zero- and single-quantum lines therefore can be used to
estimate the multiple parameters of $\theta$, $\phi$ and $B$ [see Eq. (20) in
Appendix], and in turn the magnetic field vector $\mathbf{B}(\theta,\phi)$.

We now explain why the strong interaction is a valuable resource for
multiparameter estimation. The situation is considered where the Heisenberg
interactions are much stronger than those Larmor frequencies
($|J|\gg|\gamma_{j}B|$); in this situation, the nuclear spins in each molecule
strongly interact with each other. Strong interactions between spins can
modify the energy level structure of the system, enabling the direct detection
of some forbidden transitions that are pivotal for our multiparameter
estimation. To show in more detail, we focus on two-spin case. The eigenstates
of the two-spin system in an arbitrary magnetic field are
$|\Psi_{1}\rangle=\mid\uparrow\uparrow\rangle$,
$|\Psi_{2}\rangle=\cos\xi\mid\uparrow\downarrow\rangle+\sin\xi\mid\downarrow\uparrow\rangle$,
$|\Psi_{3}\rangle=\cos\xi\mid\uparrow\downarrow\rangle$
$-\sin\xi\mid\downarrow\uparrow\rangle$,
$|\Psi_{4}\rangle=\mid\downarrow\downarrow\rangle$ [55], where
$\mid\uparrow\rangle$ and $\mid\downarrow\rangle$ are the eigenstates of spin
operators $\hat{I}_{jz^{\prime}}$ along the quantization axis $z^{\prime}$
(defined by the direction of magnetic field, see Fig. 1a) and
$\tan(2\xi)=\frac{J}{(\gamma_{1}-\gamma_{2})B}$. The transition between
$|\Psi_{2}\rangle$ and $|\Psi_{3}\rangle$ corresponds to a zero-quantum
transition, and the transition probability is
$p=|\langle\Psi_{2}|\hat{\mathbf{O}}_{z^{\prime}}|\Psi_{3}\rangle|=1/2|(\gamma_{1}-\gamma_{2})\sin(2\xi)|$.
We immediately see that the strong interacting regime
($|J|\gg|\gamma_{1}B|,|\gamma_{2}B|$) is of particular interest, since
$|\Psi_{2}\rangle$ and $|\Psi_{3}\rangle$ are entangled states between two
spins and their transition probability
$p\rightarrow|\gamma_{1}-\gamma_{2}|/2$. For non-interacting ($J=0$) and
weakly interacting ($|J|\ll|\gamma_{1}B|,|\gamma_{2}B|$) cases, $p\rightarrow
0$, so the zero-quantum line should never be detected directly. In contrast,
strongly interacting spins can produce oscillating magnetization along the
magnetic-field direction ($z^{\prime}$). In this situation, zero-quantum
spectral lines are detectable and moreover they are very narrow lines due to
long-lived coherence mechanism [56]. Based on these unforbidden transitions in
strongly interacting spin systems that encode unknown parameters,
multiparameter estimation becomes possible.

##  III Spectra of strongly interacting spin systems

It is important to experimentally measure the transition lines of the strongly
interacting nuclear spins under the measured magnetic field. We consider a
star topology with one central nuclear spin, and $n$-1 equivalent spins in
liquid-state molecules, such as 13CHn molecules, as the testbed. A liquid of
such molecules is contained in a glass tube with a 100-$\mu$L volume and
polarized with a permanent magnet ($\approx 1.3$ T). Then the liquid is
pneumatically shuttled to a sensing region $\sim 1$ mm above the atomic vapour
cell (shown in the inset of Fig. 1a). During the shuttling, a guiding field
($\approx 10^{-4}$ T) is applied on the spins, allowing adiabatic transfer of
the polarized nuclear spins into the sensing region and preserving the spin
polarization [51, 52, 53, 54]. The nuclear-spin coherence times
($\tau_{\textrm{coh}}$) of the spin clusters are typically a few seconds (see
Appendix A2). As a result, the spin polarization is initialized through
controlling the direction of the guiding field [see Eq. (9) in Appendix].
After the guiding field is turned off within 10 $\mu$s, subsequent evolution
of the spin state produces oscillating magnetization under the magnetic field
to be estimated [see Eq. (11)]. The magnetization encodes the information
about the magnetic field and is read out by an atomic vapour sensor (Fig. 1c,
see Appendix for details).

Figure 3: Magnetic field sensing with formic acid molecules. a, The procedure
of determining a magnetic field vector. b, As one example of reconstructing
magnetic field vector, three initial spin states for formic-acid sensor are
independently prepared by applying a guiding field along ${x}$, ${y}$ and
${z}$, in turn. Then the guiding field is immediately turned off. The
corresponding time-domain signals under an estimated $B(\theta,\phi)$ are
transformed into Fourier spectra (blue lines). The magnetic-field-vector
reconstruction yields $\theta\approx 1.289$ rad, $\phi\approx 0.047$ rad and
$B\approx 1.0788\cdot 10^{-7}$ T. The theoretical simulations are shown in red
lines. c, Estimated precision of $\theta$, $\phi$ when spectral amplitude
suffers from Gaussian noise $\sim\mathcal{N}(0,\sigma^{2})$. Here the variance
$\sigma=0.01$ for simulation (corresponding to SNR=100) and the uncertainty
$\sigma_{\theta}\approx 0.009$ rad and $\sigma_{\phi}\approx 0.017$ rad. d,
Histograms of magnetic field strength measurements. The measurement of a
constant applied $2.9\cdot 10^{-8}~{}\textrm{T}$ field is repeated 300 times
(30 s for each measurement). The precision of measuring the magnetic field
magnitude is given by the corresponding standard deviations $10^{-11}$ T.

We first measure different quantum transitions of 13CHn molecules. When one
measures the spins along an axis parallel to the magnetic field $\mathbf{B}$
(i.e., $\theta=0$), the selection rules for allowed transitions are (referred
as zero-quantum transitions): $\Delta f=0,\pm 1$, $\Delta m_{f}=0$. For
example, formic acid yields a single resonance line at $222.2$ Hz (Fig. 2a).
On the other hand, if one measures the spins along an axis perpendicular to
the magnetic field $\mathbf{B}$ (i.e., $\theta=\pi/2$), the selection rules
for allowed transitions are (referred as single-quantum transitions): $\Delta
f=0,\pm 1$, $\Delta m_{f}=\pm 1$. In this case, formic acid yields two
resonance lines centered at 222.2 Hz (Fig. 2b) split by
$(\gamma_{h}+\gamma_{c})B$. Here, gyromagnetic ratios $\gamma_{c}\approx
10.7077$ MHz/T for carbon spins and $\gamma_{h}\approx 42.5775$ MHz/T for
proton spins. Spectra of other 13CHn molecules show similar patterns, as shown
in Fig. 2c-f. The detailed splitting patterns are analysed in Eqs. (4)-(7) in
Appendix B1. Additionally, zero- and single-quantum transitions can be
simultaneously allowed when one measures the spins along an axis neither
parallel nor perpendicular to the magnetic field. This implies that the
information of magnetic-field orientation can be extracted from the observed
spectral lines of 13CHn molecules. Similar phenomena have been explored in the
context of optical spectroscopy of multi-electron atoms [57, 58] but, to our
knowledge, our work represents the first experimental demonstration in NMR
spectroscopy.

Figure 4: Scheme for a vectorial quantum spin gyroscope based on strongly
interacting nuclear spins. a, The proposed vectorial quantum spin gyroscope is
the same as the set-up of Fig. 1a, which operates in zero magnetic field and
experiences an inertial rotation $\bm{\Omega}$. The rotation vector
$\bm{\Omega}$ can be obtained from the rotationally induced spectra. b,
Simulated rotation-induced spectra of acetonitrile (13CH3CN) in an inertial
rotation (with the rotation frequency of 2 Hz) and various conventional
magnetic field strengths along the $x$ axis. c,d,e,f, Measurement of inertial
rotations using formic acid (H13COOH) for different orientations of the
rotation axis in zero magnetic field. The rotation frequency $\Omega=2$ Hz,
and the rotation orientations are shown in the corresponding inset.

##  IV Demonstration of multiparameter metrology

As a proof-of-principle demonstration, we use formic acid to estimate all
three components of a magnetic field. As shown in Fig. 3a, in one single scan
measurement, the amplitudes of three spectral lines at frequency $\nu_{0}$ for
zero-quantum and at frequency $\nu_{\pm 1}$ for single-quantum are
simultaneously recorded by the normalized form of
$\mathbf{A}_{\textrm{exp}}=[\textrm{a}_{1}^{(\textrm{exp})}(\nu_{-1}),\textrm{a}_{2}^{(\textrm{exp})}(\nu_{0}),\textrm{a}_{3}^{(\textrm{exp})}(\nu_{1})]$.
The magnetic field induces the Zeeman splitting $\Delta=|\nu_{1}-\nu_{-1}|$ of
the doublet and, accordingly, the magnitude of the measured field can be
determined in advance through $\Delta=(\gamma_{c}+\gamma_{h})B$ [the
splittings $\Delta$ for arbitrary 13CHn molecules are shown in Eqs. (4)-(7) in
Appendix]. Furthermore, we use a model [see Appendix C] to calculate the
amplitudes of relevant NMR lines, labelled by
$\mathbf{A}_{\textrm{sim}}(\theta,\phi)=[\textrm{a}_{1}^{(\textrm{sim})}(\nu_{-1}),\textrm{a}_{2}^{(\textrm{sim})}(\nu_{0}),\textrm{a}_{3}^{(\textrm{sim})}(\nu_{1})]$
when $\theta$ and $\phi$ are given. Conversely, the parameters of $\theta$,
$\phi$ can be obtained by minimizing the norm
$|\mathbf{A}_{\textrm{sim}}(\theta,\phi)-\mathbf{A}_{\textrm{exp}}|$. However,
multiple pairs of $(\theta,\phi)$ may result in the same spectrum if they
satisfy certain symmetric conditions. To remove this ambiguity, we measure the
nuclear-spin spectra and three sets of $\mathbf{A}_{\textrm{exp}}$ under
various initial states $\rho_{0}$ by controlling the guiding field along
${x}$, ${y}$ and ${z}$ [51], in turn. The corresponding NMR spectra are shown
in Fig. 3b. The $\theta,\phi$ can be unambiguously determined by finding the
simulated $\mathbf{A}_{\textrm{sim}}(\theta,\phi)$ that simultaneously best
matches the measured $\mathbf{A}_{\textrm{exp}}$ under three initial states.
The magnetic field is determined to be such that $\theta\approx 1.289$ rad and
$\phi\approx 0.047$ rad, and $B\approx 1.0788\cdot 10^{-7}$ T. The simulated
Fourier spectra with red lines are in a good agreement with the experimental
spectra.

In order to benchmark the precision of multiparameter estimation in our
protocol, we simulate the performance of $\theta$, $\phi$ deviation under the
influence of spectral amplitude variance. The variances
$\\{\delta\textrm{a}_{1}^{(\textrm{exp})},\delta\textrm{a}_{2}^{(\textrm{exp})},$
$\delta\textrm{a}_{3}^{(\textrm{exp})}\\}$ are assumed to obey Gaussian noise
$\sim\mathcal{N}(0,\sigma^{2})$, which comes mainly from the photon-shot noise
of the probe beam in our experiment. Figure 3c shows the $\theta$, $\phi$
joint statistical distribution with $\sigma=0.01$ (corresponding to our
current experimental signal-to-noise ratio (SNR) of 100). The uncertainties
(standard deviations) of $\theta$, $\phi$ values are respectively estimated to
be $\sigma_{\theta}\approx 0.009$ rad and $\sigma_{\phi}\approx 0.017$ rad, as
shown in their histograms of deviations. One could increase, for example, the
nuclear-spin initial polarization and NMR detector sensitivity, to decrease
$\sigma$ and then achieve better orientation resolution.

We next discuss how precisely we can determine the magnitude of a magnetic
field vector. The magnetic magnitude can be achieved by exploiting the Zeeman
splitting of strongly interacting nuclear spins, for example, using the
splitting $\Delta$ in the formic acid spectrum. For other 13CHn molecules,
there exists similar Zeeman splitting that can be used to measure the magnetic
field magnitude [see Eqs. (4)-(7) in Appendix]. We use a set of such
measurements to quantify the estimation precision. The uncertainties of
$\Delta$ in the frequencies extracted from spectral profile fits are in the
order of $0.3~{}\textrm{mHz}$. As shown in Fig. 3d, the statistical standard
deviation of absolute magnetic field is $10^{-11}~{}\textrm{T}$.

Nuclear spins can also be used to measure inertial rotation. Conventional
nuclear spin gyroscopes [40, 41] usually rely on measuring the frequency shift
of the freely precessing non-interacting spins (such as 3He) in a bias
magnetic field, limiting the gyroscope to detection of rotation along the
magnetic field direction. Here, we show that it could overcome this limitation
by using strongly interacting nuclear spins, enabling to estimate three
dimensional rotations. Similar to magnetic field sensing, strongly interacting
nuclear spins can also be sensitive to magnetic pseudo-fields originating from
inertial rotation. A magnetic pseudo-field
$\mathbf{B}_{\Omega}^{j}=\bm{\Omega}/\gamma_{j}$ is induced on the $j$th spin
[40, 59, 60] when the nuclear spin is probed with an atomic vapour sensor in a
rotating frame, as shown in Fig. 4a. Here $\bm{\Omega}(\theta,\phi)$ is the
rotation vector and parameterized with the polar angle $\theta$ and azimuthal
angle $\phi$ and magnitude $\Omega$. The corresponding Hamiltonian is

$\begin{cases}\mathcal{H}(\theta,\phi,\Omega)=\mathcal{H}_{\Omega}+\mathcal{H}_{\textrm{spins}}+\mathcal{H}_{\textrm{int}},\\\
\mathcal{H}_{\Omega}=-2\pi\bm{\Omega}(\theta,\phi)\cdot\sum_{j}\mathbf{I}_{j},\end{cases}$
(3)

where $\mathcal{H}_{\textrm{spins}}$ and $\mathcal{H}_{\textrm{int}}$ are
described in Eq. (1).

Simulations are performed to investigate the rotational effect on the spectra
of 13CHn molecules. Figure 4b shows the spectra of acetonitrile in different
magnetic fields and a rotation along the $x$ axis. With decreasing applied
magnetic-field strength, the resonance lines gradually change from a multiplet
to a doublet. More generally, under inertial rotation and zero magnetic field,
arbitrary molecules yield fewer resonance lines than those in real magnetic
field [see Eq. (8) in Appendix]. Thus, in contrast to real magnetic field
sensing, a variety of molecular species that have no overlapping lines can
function as rotation sensors in zero field.

The spectra of strongly interacting nuclear spins can be used to the
simultaneous estimation of the rotation speed and the orientation of the
rotation axis. As an illustration, we consider the case of nuclear spins
initialized with the guiding field along $z$. When the rotation axis is along
the sensitive $z$-axis of atomic vapour sensor, the resonance line of formic
acid remains a single line (Fig. 4c, bottom). When the rotation axis is along
$x$, the resonance line splits into a doublet with a splitting of $2\Omega$
(Fig. 4e, bottom). When the rotation axis is along $\theta=\pi/4$, $\phi=0$,
the resonance line splits into a triplet (Fig. 4d, bottom). The above
described spectra can still be explained by the theories [see Eqs. (4)-(7) in
Appendix] for conventional magnetic field, if $\gamma_{h}\mathbf{B}$ and
$\gamma_{c}\mathbf{B}$ are replaced by $\bm{\Omega}$. The measurement of the
rotation vector $\bm{\Omega}$ can be performed with the same protocol as the
one for the vector magnetometry above. As demonstrations, Figure 4c-f show the
simulated spectra of formic acid for four physical rotations with different
orientations of the rotation axis. Clearly, they show significant differences.
Using the same procedure as for conventional field sensing, we can obtain the
rotation frequency and the orientation of the rotation axis. Based on the
experimental precision and orientation resolution of our sensor for the real
magnetic field, the precision of the rotation frequency is about $2.7\cdot
10^{-4}$ Hz (corresponding to the bias stability of gyroscope [40]) and the
rotation orientation resolution is about 0.02 rad. With further optimization
(see conclusion and outlook for details), it should be possible to enhance the
single-scan precision with several orders of magnitude improvement, which is
sufficient for the measurement of Earth’s rotation (about $10^{-5}$ Hz).

##  V Conclusion and Outlook

In this paper, we have introduced a new approach to multiparameter metrology,
using a set of strongly interacting nuclear spins as a quantum resource. Using
this technique, we have demonstrated a new type of sensor capable of detecting
all three components of a magnetic field or the axis of rotation of an overall
sample rotation. This new measurement scheme has promising applications in
geophysical surveys, high-precision navigation, and testing fundamental
physics. Moreover, our technique naturally avoids the field contamination
induced by previously-required reference fields [32, 33, 34, 35], and opens a
feasible route towards realizing quantum noise-limited multiparameter sensors.

A further improvement of the experimental sensitivity can be anticipated. For
example, the probe states used in this work are thermal equilibrium states
with low nuclear-spin polarization. The extension to hyperpolarized states
with near-unity spin polarization through many well-developed techniques, such
as parahydrogen-induced polarization (PHIP) [61, 62], dynamic nuclear
polarization (DNP) [63] and spin-exchange optical pumping (SEOP) [41]
techniques, could yield at least five orders of magnitude improvements over
currently achievable sensitivity. Moreover, it is also possible to improve the
readout sensitivity of the atomic vapour sensor. For example, the application
of magnetic gradiometer to the detection of nuclear-spin magnetization, as
recently demonstrated in Ref. [52], could enhance the multiparameter
estimation sensitivity by a factor of about ten. With recent quantum control
technologies [51, 64] in ultralow-field NMR [65, 54, 53, 66], the NOON of the
nuclear spins in each molecule can be created which offers an enhanced
sensitivity to the multiparameter estimation at a fundamental level [67].

Although demonstrated for strongly-interacting nuclear-spin systems, the
protocol presented here can be directly applied to other strongly interacting
spin systems and design a new class of multiparameter quantum sensors. We
suggest future theoretical and experimental studies of the interaction-based
multiparameter metrology with electron-nuclear [68, 69, 70] or electron-
electron [71, 72] interacting spin systems, taking a fresh look at many
techniques, such as spectroscopy, magnetometry, imaging, navigation, and even
biocompass mechanism. For example, the application to nitrogen-vacancy centers
[73, 74], which offer the capability of realizing nanometer-scale spatial
resolution, could realize a completely new nanometer-scale vector magnetic-
field microscope. It can also be applied to explore the magnetic-field
orientation dependence of electron-electron spin systems, such as radical
pairs in magnetic proteins, which closely relate with the hitherto mysterious
biocompass mechanism [71, 72]. We are therefore convinced that this new type
of interaction-based multiparameter quantum metrology opens a promising new
field of research.

## Acknowledgements

We acknowledge Dmitry Budker, Ren-Bao Liu, Roman P. Frutos and Teng Wu for
valuable discussions, and Jiankun Chen for drawing the set-up diagram. This
work was supported by National Key Research and Development Program of China
(grant no. 2018YFA0306600), National Natural Science Foundation of China
(grants nos. 11661161018, 11927811, 12004371), Anhui Initiative in Quantum
Information Technologies (grant no. AHY050000), USTC Research Funds of the
Double First-Class Initiative (grant no. YD3540002002).

## APPENDIX A: Experimental details

### 1\. Experimental set-up

As depicted schematically in Fig. 1a, the liquid ($\sim 100$ $\mu$L) contained
in a $5$-mm NMR tube is polarized in a Halbach magnet ($B_{p}\approx
1.3~{}\textrm{T}$), and then pneumatically shuttled to a sensing and readout
region ($\sim 1$ mm above the atomic vapour cell) inside a five-layer
$\mu$-metal magnetic shield. During the pneumatic transfer, a guiding solenoid
and an additional set of mutually orthogonal Helmholtz coils are used to
generate a guiding field along an axis (${x}$, ${y}$ or ${z}$) for preparing
the initial spin state $\rho_{0}$, as described in Eq. (9). The magnitude of
the guiding field ($B_{g}\approx 10^{-4}~{}\textrm{T}$) satisfies $B_{g}\gg
J/(\gamma_{h}-\gamma_{c})$. After the shuttling, the guiding field is switched
off within 10 $\mu$s. The nuclear spins of the liquid then evolves under the
external magnetic field and, accordingly, generate an oscillating
magnetization signal detected with a 87Rb vapour sensor [53, 65, 51, 52]. The
87Rb vapour is resistively heated to $180$ ∘C. As shown in Fig. 1c, the 87Rb
atoms in a vapour are pumped with a circularly polarized laser beam
propagating in the $x$ direction. The laser frequency is tuned to the center
of the buffer-gas (N2) broadened and shifted $\textrm{D}1$ line. The magnetic
field is measured via optical rotation of linearly polarized probe laser light
at the $\textrm{D}2$ transition propagating in the $y$ direction. The
sensitivity of the 87Rb vapour sensor along the $z$ axis is optimized to
$\approx 20~{}\textrm{fT}\cdot\textrm{Hz}^{-1/2}$.

### 2\. Strongly-interacting nuclear-spin samples

As used in our experiments, the liquids including formic acid (H13COOH),
formaldehyde (13CH2O), acetonitrile (13CH3CN), and acetic acid (13CH3COOH) are
from Sigma-Aldrich, and are contained in standard $5$-mm NMR tubes. The
samples were flame-sealed under vacuum in NMR tubes after five freeze-pump-
thaw cycles to remove dissolved oxygen, which is otherwise a significant
source of relaxation. Nuclear spin coherence times $\tau_{\textrm{coh}}\approx
10.4$ s and scalar spin-spin coupling strengths $J\approx 222.2$ Hz for formic
acid, $\tau_{\textrm{coh}}\approx 0.8$ s and $J\approx 163.9$ Hz for
formaldehyde, $\tau_{\textrm{coh}}\approx 4.7$ s and $J\approx 136.25$ Hz for
acetonitrile, and $\tau_{\textrm{coh}}\approx 8.8$ s and $J\approx 129.5$ Hz
for acetic acid.

## APPENDIX B: Splitting of strongly interacting spins

### 1\. Zeeman splitting

When the magnetic field is sufficiently smaller than the scalar spin-spin
couplings (i.e., $\gamma_{h}B,\gamma_{c}B$ $\ll|J|$), the nuclear spins in a
13CHn molecule are strongly coupled to each other. In this situation, the good
quantum numbers are the total spin angular momentum $f=f_{c}+f_{h}$ and the
corresponding magnetic quantum number $m_{f}$ with possible values $-f\leq
m_{f}\leq f$. Here, $f_{c}=\pm 1/2$ and $f_{h}=n/2-k$ are respectively the
total 13C and H spin angular momentum, and $k=0,1,...,n/2$ for even $n$ or
$k=0,1,...,(n-1)/2$ for odd $n$. The eigenstates of 13CHn molecules are
denoted by $|{fm_{f};k}\rangle$ [65, 75], which can be approximately those
eigenstates of the zero-field Hamiltonian
$\mathcal{H}=\sum_{j>i}^{n}{{J_{ij}}}\mathbf{I}_{i}\cdot\mathbf{I}_{j}$. The
examples of 13CH, 13CH2, 13CH3 energy levels are presented in Fig. 2.

First, we consider the case of real magnetic fields. When the ensemble of
13CHn molecules are measured parallel to the external magnetic field, the
selection rules for allowed transitions are $\Delta f=0,\pm 1$, $\Delta
m_{f}=0$. Because the protons are equivalent in a 13CHn, there is an
additional selection rule $\Delta k=0$ [65]. Accordingly, one could observe
$n-2k$ resonance lines centered at $\frac{1}{2}J(1+n-2k)$, with the
frequencies

$\begin{array}[]{l}\nu_{f,m_{f};k}^{f^{\prime},m_{f};k}=\frac{1}{2}J(1+n-2k)+\frac{2m_{f}(-\gamma_{h}+\gamma_{c})}{1+n-2k}B,\end{array}$
(4)

where $\nu_{f,m_{f};k}^{f^{\prime},m^{\prime}_{f};k^{\prime}}$ denotes the
transition frequency between the eigenstates $|{fm_{f};k}\rangle$ and
$|{f^{\prime}m^{\prime}_{f};k^{\prime}}\rangle$. The adjacent zero-quantum
Zeeman splitting is

$\begin{array}[]{l}\Delta_{\textrm{ZQ}}=\frac{2}{1+n-2k}(\gamma_{h}-\gamma_{c})B.\end{array}$
(5)

When the ensemble of 13CHn molecules are measured perpendicular to the
external magnetic field, the selection rules for allowed transitions are
$\Delta f=0,\pm 1$, $\Delta m_{f}=\pm 1$, and $\Delta k=0$. Accordingly, one
could observe $2(n-2k)$ resonance lines centered at $\frac{1}{2}J(1+n-2k)$,
with the frequencies

$\begin{array}[]{l}\nu_{f,m_{f};k}^{f^{\prime},m_{f}\pm
1;k}=\frac{1}{2}J(1+n-2k)+[\frac{2m_{f}(-\gamma_{h}+\gamma_{c})}{1+n-2k}\pm\frac{(n-2k)\gamma_{h}+\gamma_{c}}{1+n-2k}]B.\end{array}$
(6)

The Zeeman splitting of the central doublet, which have the highest signal-to-
noise ratio compared to others, is

$\begin{array}[]{l}\Delta_{\textrm{SQ}}=\frac{2}{1+n-2k}[\gamma_{h}+(n-2k)\gamma_{c}]B.\end{array}$
(7)

For formic acid with $n=1$ and $k=0$, the Zeeman splitting of the doublet is
$\Delta=(\gamma_{c}+\gamma_{h})B$, as used in the text to determine the
strength of magnetic field. The above discussions only consider the first-
order Zeeman effects of 13CHn molecules. The higher-order Zeeman effects cause
a systematic error of measuring magnetic-field strengths and can be seen in
Refs. [75, 44].

### 2\. Rotation-induced splitting

We now consider the splitting of 13CHn molecules under an inertial rotation.
The magnetic pseudo-field induced by the inertial rotation is described in the
text. In the absence of real magnetic field, the above Eqs. (4)-(7) can also
be used in the case of magnetic pseudo-fields originating from inertial
rotation. In this case, $\gamma_{h}B$ and $\gamma_{c}B$ should be replaced by
$\Omega$. Remarkably, in contrast to the case of real magnetic field,
$\nu_{f,m_{f};k}^{f^{\prime},m_{f};k}$ and
$\nu_{f,m_{f};k}^{f^{\prime},m_{f}\pm 1;k}$ become independent of the magnetic
quantum number $m_{f}$. As a result, the spectra of molecules under inertial
rotation and zero field are very simple:

$\begin{cases}\nu_{f,m_{f};k}^{f^{\prime},m_{f};k}(\Omega)=\frac{1}{2}J(1+n-2k),\\\
\Delta_{\textrm{ZQ}}(\Omega)=0,\\\ \nu_{f,m_{f};k}^{f^{\prime},m_{f}\pm
1;k}(\Omega)=\frac{1}{2}J(1+n-2k)\pm\Omega,\\\
\Delta_{\textrm{SQ}}(\Omega)=2\Omega.\end{cases}$ (8)

(1) “$\Delta_{\textrm{ZQ}}=0$” implies that zero-quantum transition lines
remain single lines under rotation, whereas zero-quantum lines are usually
split into a multiplet in conventional magnetic fields, such as those shown in
Fig. 2c,e. (2) “$\Delta_{\textrm{SQ}}=2\Omega$” implies that single-quantum
transition lines are doublets with a splitting of $2\Omega$. This is
significantly different from that in real magnetic fields, where the spitting
depends on $n$ and $k$ and the number of splitting lines is usually large
(Fig. 2b,d,f).

## APPENDIX C: Multiparameter estimation process

The multiparameter estimation presented in this work includes three stages:
probe preparation, encoding with the strongly interacting spins, and probe
readout. In the following, we provide the detailed theory of such three
stages.

### 1\. Probe preparation

As described in the part of experimental set-up in Appendix A1, the nuclear
spins are polarized in a permanent magnet and then pneumatically shuttled to a
sensing region. When the guiding field is turned off within 10 $\mu$s, the
spin state remains the high-field equilibrium state. In the high-temperature
approximation, the relevant initial spin state is [51]

$\begin{array}[]{l}\rho_{0}=\mathbbm{1}/2^{n}-\sum_{j}\varepsilon_{j}\mathbf{I}_{j}\cdot\hat{\mathbf{k}}_{g},\end{array}$
(9)

where $\varepsilon_{j}=\gamma_{j}B_{p}/k_{\textrm{B}}T\sim 10^{-6}$,
$k_{\textrm{B}}$ is the Boltzmann constant, $T$ is the temperature,
$\mathbf{I}_{j}=(\hat{I}_{jx},\hat{I}_{jy},\hat{I}_{jz})$ is the spin angular
momentum operator of the $i$th spin, and $\hat{\mathbf{k}}_{g}$ is the unit
vector along the guiding field orientation. For example,
$\hat{\mathbf{k}}_{g}=(0,0,1)$ corresponds to the case of $z$-guiding field.
Thus, the initial state can be controlled by the guiding field direction
$\hat{\mathbf{k}}_{g}$.

### 2\. Encoding unknown parameters

Here we use the encoding of an unknown magnetic field as an example to
describe the basic encoding process. The process is general and suitable for
encoding inertial rotation. In Appendix B, we have described the process to
determine the strength of an unknown magnetic field or inertial rotation.
Thus, we focus on how to estimate their orientation in this part. The
orientation of magnetic field $\mathbf{B}$ is parameterized with the polar
angle $\theta$ and azimuthal angle $\phi$, as shown in the inset of Fig. 1a.
In the measured magnetic field, subsequent evolution of the quantum state of
the spins in each molecule produces oscillating magnetization, which is
detected by the 87Rb vapour sensor. The observable in our experiment is the
total $z$ magnetization,

$S(\rho_{0},\theta,\phi,t)\propto\textrm{Tr}[\rho(\theta,\phi,B,t)\sum_{j}\gamma_{j}\hat{I}_{jz}],$
(10)

where $\rho(\theta,\phi,B,t)$ is the time dependent density matrix, as
described in Eq. (2). Consequently, the component of the sensor signal,
evolving with the frequency of
$\nu_{f,m_{f};k}^{f^{\prime},m_{f}^{\prime};k}$, can be written with

$\begin{array}[]{l}S_{fm_{f};k}^{f^{\prime}m_{f}^{\prime};k^{\prime}}(\rho_{0},\theta,\phi,t)\propto\\\
\Re_{fm_{f};k}^{f^{\prime}m_{f}^{\prime};k^{\prime}}(\rho_{0},\theta,\phi)\cdot\textrm{cos}(\Phi-2\pi\nu_{f,m_{f};k}^{f^{\prime},m_{f}^{\prime};k^{\prime}}t)e^{-t/\tau_{\textrm{coh}}},\end{array}$
(11)

where the amplitude
$\Re_{fm_{f};k}^{f^{\prime}m_{f}^{\prime};k^{\prime}}(\rho_{0},\theta,\phi)$
and phase $\Phi$ are defined by the following relation

$\begin{array}[]{l}\Re_{fm_{f};k}^{f^{\prime}m_{f}^{\prime};k^{\prime}}(\rho_{0},\theta,\phi)e^{i{\Phi}}\equiv\\\
\langle{fm_{f};k}|\rho_{0}|{f^{\prime}m_{f}^{\prime};k^{\prime}}\rangle\langle{f^{\prime}m_{f}^{\prime};k^{\prime}}|\hat{\mathbf{O}}_{z}|{fm_{f};k}\rangle,\end{array}$
(12)

the spin operator $\hat{\mathbf{O}}_{l}\equiv\sum_{j}\gamma_{j}\hat{I}_{jl}$
($l=x,y,z$), and $e^{-t/\tau_{\textrm{coh}}}$ is due to the decoherence
effect. The amplitude
$\Re_{fm_{f};k}^{f^{\prime}m_{f}^{\prime};k^{\prime}}(\rho_{0},\theta,\phi)$
can be determined by the observed NMR spectrum, as seen in “Probe readout”.

To achieve the clear relation between the amplitude
$\Re_{fm_{f};k}^{f^{\prime}m_{f}^{\prime};k^{\prime}}(\rho_{0},\theta,\phi)$
and the unknown parameters $\theta$, $\phi$, we introduce a convenient
orthonormal coordinate system, defined by using only $\theta$ and $\phi$. In
such a coordinate system (see the inset of Fig. 1a), the quantization axis
($\mathbf{z}^{\prime}$) is chosen along the magnetic field, i.e.,

$\begin{cases}\mathbf{z^{\prime}}\equiv\mathbf{B}/B=(\sin\theta\cos\phi,\sin\theta\sin\phi,\cos\theta),\\\
\mathbf{x^{\prime}}\equiv\frac{\partial\mathbf{z^{\prime}}}{\partial\theta}=(\cos\theta\cos\phi,\cos\theta\sin\phi,-\sin\theta),\\\
\mathbf{y^{\prime}}\equiv(\frac{1}{\sin\theta})\frac{\partial\mathbf{z^{\prime}}}{\partial\phi}=(-\sin\phi,\cos\phi,0).\end{cases}$
(13)

Thus the spin operator in the laboratory frame can be transformed to
$\hat{\mathbf{O}}^{\prime}_{l^{\prime}}=\sum_{j}\gamma_{j}\hat{I}_{jl^{\prime}}$
($l^{\prime}=x^{\prime},y^{\prime},z^{\prime}$) by

$\hat{\mathbf{O}}=\mathcal{P}(\theta,\phi)\hat{\mathbf{O}}^{\prime},$ (14)

where
$\hat{\mathbf{O}}=(\hat{\mathbf{O}}_{x},\hat{\mathbf{O}}_{y},\hat{\mathbf{O}}_{z})^{T}$,
$\hat{\mathbf{O}}^{\prime}=(\hat{\mathbf{O}}^{\prime}_{x^{\prime}},\hat{\mathbf{O}}^{\prime}_{y^{\prime}},\hat{\mathbf{O}}^{\prime}_{z^{\prime}})^{T}$,
and $\mathcal{P}(\theta,\phi)$ is a $3\times 3$ transformation matrix

$\mathcal{P}(\theta,\phi)=\begin{pmatrix}\cos\theta\cos\phi&-\sin\phi&\sin\theta\cos\phi\\\
\cos\theta\sin\phi&\cos\phi&\sin\theta\sin\phi\\\
-\sin\theta&0&\cos\theta\end{pmatrix}.$ (15)

Further, the element
$\langle{f^{\prime}m_{f}^{\prime};k}|\hat{\mathbf{O}}_{l}|{fm_{f};k}\rangle$
in the laboratory frame can be expressed in the new defined coordinate system
as

$\langle{f^{\prime}m_{f}^{\prime};k}|\hat{\mathbf{O}}|{fm_{f};k}\rangle=\mathcal{P}(\theta,\phi)\cdot\langle{f^{\prime}m_{f}^{\prime};k}|\hat{\mathbf{O}}^{\prime}|{fm_{f};k}\rangle.$
(16)

Moreover, the
$\langle{f^{\prime}m_{f}^{\prime};k}|\hat{\mathbf{O}}^{\prime}|{fm_{f};k}\rangle$
is independent of the orientation of the magnetic field ($\theta,\phi$) and
thus can be calculated without knowing any information of the magnetic field
$\mathbf{B}$ in ultralow-field regime. The detailed methods to calculate
$\langle{f^{\prime}m_{f}^{\prime};k}|\hat{\mathbf{O}}^{\prime}|{fm_{f};k}\rangle$
can be seen in Ref. [76]. For the typical 13CH molecules, for example, we have

$\begin{cases}\langle{0,0;0}|\hat{\mathbf{O}}^{\prime}|{1,0;0}\rangle=(0,0,\frac{\gamma_{c}-\gamma_{h}}{2})^{T},\\\
\langle{0,0;0}|\hat{\mathbf{O}}^{\prime}|{1,-1;0}\rangle=(\frac{\gamma_{c}-\gamma_{h}}{2\sqrt{2}},-\frac{i(\gamma_{c}-\gamma_{h})}{2\sqrt{2}},0)^{T},\\\
\langle{0,0;0}|\hat{\mathbf{O}}^{\prime}|{1,1;0}\rangle=(\frac{-\gamma_{c}+\gamma_{h}}{2\sqrt{2}},-\frac{i(\gamma_{c}-\gamma_{h})}{2\sqrt{2}},0)^{T}.\\\
\end{cases}$ (17)

Based on Eqs. (16) and (17),
$\langle{f^{\prime}m_{f}^{\prime};k}|\hat{\mathbf{O}}_{z}|{fm_{f};k}\rangle$
can be represented by the only parameters $\theta$ and $\phi$, i.e.,

$\begin{array}[]{l}\langle{f^{\prime}m_{f}^{\prime};k}|\hat{\mathbf{O}}_{z}|{fm_{f};k}\rangle=\hat{\mathbf{z}}\cdot\mathcal{P}(\theta,\phi)\cdot\langle{f^{\prime}m_{f}^{\prime};k}|\hat{\mathbf{O}}^{\prime}|{fm_{f};k}\rangle,\end{array}$
(18)

where the unit vector $\hat{\mathbf{z}}=(0,0,1)$. For the density matrix
$\rho_{0}$, the element
$\langle{fm_{f};k}|\rho_{0}|{f^{\prime}m_{f}^{\prime};k^{\prime}}\rangle$ has
similar result

$\begin{array}[]{l}\langle{f^{\prime}m_{f}^{\prime};k}|\rho_{0}|{fm_{f};k}\rangle=\\\
(-B_{p}/k_{\textrm{B}}T)\hat{\mathbf{k}}_{g}\cdot\mathcal{P}(\theta,\phi)\cdot\langle{f^{\prime}m_{f}^{\prime};k}|\hat{\mathbf{O}}^{\prime}|{fm_{f};k}\rangle.\end{array}$
(19)

According to Eqs. (12), (17), (18), (19), the amplitude
$\Re_{fm_{f};k}^{f^{\prime}m_{f}^{\prime};k^{\prime}}(\rho_{0},\theta,\phi)$
can be finally represented by the two unknown parameters $\theta$ and $\phi$,
i.e.,

$\begin{array}[]{l}\Re_{fm_{f};k}^{f^{\prime}m_{f}^{\prime};k^{\prime}}(\rho_{0},\theta,\phi)=\\\
B_{p}/k_{\textrm{B}}T\cdot|\hat{\mathbf{k}}_{g}\cdot\mathcal{P}(\theta,\phi)\langle{f^{\prime}m_{f}^{\prime};k}|\hat{\mathbf{O}}^{\prime}|{fm_{f};k}\rangle|\times\\\
|\hat{\mathbf{z}}\cdot\mathcal{P}(\theta,\phi)\langle{f^{\prime}m_{f}^{\prime};k}|\hat{\mathbf{O}}^{\prime}|{fm_{f};k}\rangle|.\end{array}$
(20)

Equation (20) is the central result of our interaction-based multiparameter
estimation. We make two important observations. (1)
$\Re_{fm_{f};k}^{f^{\prime}m_{f}^{\prime};k^{\prime}}(\rho_{0},\theta,\phi)$
is represented with two parts, $\theta$, $\phi$-dependent [i.e.,
$\mathcal{P}(\theta,\phi)$] and $\theta$, $\phi$-independent one [e.g.,
$\langle{f^{\prime}m_{f}^{\prime};k}|\hat{\mathbf{O}}^{\prime}|{fm_{f};k}\rangle$].
This clearly shows the relation between the multiple unknown parameters and
$\Re_{fm_{f};k}^{f^{\prime}m_{f}^{\prime};k^{\prime}}(\rho_{0},\theta,\phi)$
that can be experimentally determined in NMR spectrum; see “Probe readout” for
details. (2) The use of strong interactions between nuclear spins permits
nonzero transition probability
$\langle{f^{\prime}m_{f}^{\prime};k}|\hat{\mathbf{O}}^{\prime}_{z^{\prime}}|{fm_{f};k}\rangle$
(corresponding to zero-quantum transition; see Fig. 2a,c, e), offering an
efficient route to encode the unknown parameters.

### 3\. Probe readout

The amplitude
$\Re_{fm_{f};k}^{f^{\prime}m_{f}^{\prime};k^{\prime}}(\rho_{0},\theta,\phi)$
can be experimentally determined with the observed NMR spectrum, detected with
a 87Rb vapour sensor. One detailed example with formic acid molecules is
presented in Fig. 3b. Specifically, the amplitudes of zero- and single-quantum
lines for 13CH molecules, $\Re_{0,0;0}^{1,0;0}(\rho_{0},\theta,\phi)$ and
$\Re_{0,0;0}^{1,\pm 1;0}(\rho_{0},\theta,\phi)$ can be determined from the
triplet (see Fig. 3b), and lastly can be used to extract the parameters
$\theta$ and $\phi$ of the magnetic field $\mathbf{B}(\theta,\phi)$.

## References

  * [1] Giovannetti, V., Lloyd, S. & Maccone, L.  Advances in quantum metrology.  _Nat. Photon._ 5, 222 (2011).
  * [2] Degen, C. L., Reinhard, F. & Cappellaro, P.  Quantum sensing.  _Rev. Mod. Phys._ 89, 035002 (2017).
  * [3] Pezze, L., Smerzi, A., Oberthaler, M. K., Schmied, R. & Treutlein, P.  Quantum metrology with nonclassical states of atomic ensembles.  _Rev. Mod. Phys._ 90, 035005 (2018).
  * [4] Braun, D. _et al._ Quantum-enhanced measurements without entanglement.  _Rev. Mod. Phys._ 90, 035006 (2018).
  * [5] Aasi, J. _et al._ Enhanced sensitivity of the LIGO gravitational wave detector by using squeezed states of light.  _Nat. Photon._ 7, 613–619 (2013).
  * [6] Budker, D. & Romalis, M.  Optical magnetometry.  _Nat. Phys._ 3, 227–234 (2007).
  * [7] Safronova, M. _et al._ Search for new physics with atoms and molecules.  _Rev. Mod. Phys._ 90, 025008 (2018).
  * [8] Lang, J., Liu, R.-B. & Monteiro, T.  Dynamical-decoupling-based quantum sensing: Floquet spectroscopy.  _Phys. Rev. X_ 5, 041016 (2015).
  * [9] Álvarez, G. A., Suter, D. & Kaiser, R.  Localization-delocalization transition in the dynamics of dipolar-coupled nuclear spins.  _Science_ 349, 846–848 (2015).
  * [10] Lucchesi, L. & Chiofalo, M. L.  Many-body entanglement in short-range interacting fermi gases for metrology.  _Phys. Rev. Lett._ 123, 060406 (2019).
  * [11] Kong, J. _et al._ Measurement-induced, spatially-extended entanglement in a hot, strongly-interacting atomic system.  _Nat. Commun._ 11, 1–9 (2020).
  * [12] Peyronel, T. _et al._ Quantum nonlinear optics with single photons enabled by strongly interacting atoms.  _Nature_ 488, 57–60 (2012).
  * [13] Dooley, S., Hanks, M., Nakayama, S., Munro, W. J. & Nemoto, K.  Robust quantum sensing with strongly interacting probe systems.  _NPJ Quant. Inf._ 4, 1–7 (2018).
  * [14] Nolan, S. P., Szigeti, S. S. & Haine, S. A.  Optimal and robust quantum metrology using interaction-based readouts.  _Phys. Rev. Lett._ 119, 193601 (2017).
  * [15] Zhou, H. _et al._ Quantum metrology with strongly interacting spin systems.  _Phys. Rev. X_ 10, 031003 (2020).
  * [16] Frérot, I. & Roscilde, T.  Quantum critical metrology.  _Phys. Rev. Lett._ 121, 020402 (2018).
  * [17] Roy, S. & Braunstein, S. L.  Exponentially enhanced quantum metrology.  _Phys. Rev. Lett._ 100, 220501 (2008).
  * [18] Napolitano, M. _et al._ Interaction-based quantum metrology showing scaling beyond the Heisenberg limit.  _Nature_ 471, 486–489 (2011).
  * [19] Boixo, S., Flammia, S. T., Caves, C. M. & Geremia, J. M.  Generalized limits for single-parameter quantum estimation.  _Phys. Rev. Lett._ 98, 090401 (2007).
  * [20] Zwierz, M., Pérez-Delgado, C. A. & Kok, P.  General optimality of the Heisenberg limit for quantum metrology.  _Phys. Rev. Lett._ 105, 180402 (2010).
  * [21] Liu, R. _et al._ Experimental Adiabatic Quantum Metrology with the Heisenberg scaling.  _arXiv preprint arXiv:2102.07056_ (2021).
  * [22] Chu, Y., Zhang, S., Yu, B. & Cai, J.  Dynamic framework for criticality-enhanced quantum sensing.  _Phys. Rev. Lett._ 126, 010502 (2021).
  * [23] Rams, M. M., Sierant, P., Dutta, O., Horodecki, P. & Zakrzewski, J.  At the limits of criticality-based quantum metrology: Apparent super-Heisenberg scaling revisited.  _Phys. Rev. X_ 8, 021022 (2018).
  * [24] Rovny, J., Blum, R. L. & Barrett, S. E.  Observation of discrete-time-crystal signatures in an ordered dipolar many-body system.  _Phys. Rev. Lett._ 120, 180603 (2018).
  * [25] Kominis, I., Kornack, T., Allred, J. & Romalis, M. V.  A subfemtotesla multichannel atomic magnetometer.  _Nature_ 422, 596–599 (2003).
  * [26] Boixo, S. _et al._ Quantum metrology: dynamics versus entanglement.  _Phys. Rev. Lett._ 101, 040403 (2008).
  * [27] Szczykulska, M., Baumgratz, T. & Datta, A.  Multi-parameter quantum metrology.  _Adv. Phys.: X_ 1, 621–639 (2016).
  * [28] Vidrighin, M. D. _et al._ Joint estimation of phase and phase diffusion for quantum metrology.  _Nat. Commun._ 5, 1–7 (2014).
  * [29] Hou, Z. _et al._ Zero–trade-off multiparameter quantum estimation via simultaneously saturating multiple Heisenberg uncertainty relations.  _Sci. Adv._ 7, eabd2986 (2021).
  * [30] Roccia, E. _et al._ Multiparameter approach to quantum phase estimation with limited visibility.  _Optica_ 5, 1171–1176 (2018).
  * [31] Hou, Z. _et al._ Minimal tradeoff and ultimate precision limit of multiparameter quantum magnetometry under the parallel scheme.  _Phys. Rev. Lett._ 125, 020501 (2020).
  * [32] Seltzer, S. & Romalis, M.  Unshielded three-axis vector operation of a spin-exchange-relaxation-free atomic magnetometer.  _Appl. Phys. Lett._ 85, 4804–4806 (2004).
  * [33] Patton, B., Zhivun, E., Hovde, D. & Budker, D.  All-optical vector atomic magnetometer.  _Phys. Rev. Lett._ 113, 013001 (2014).
  * [34] Thiele, T., Lin, Y., Brown, M. O. & Regal, C. A.  Self-calibrating vector atomic magnetometry through microwave polarization reconstruction.  _Phys. Rev. Lett._ 121, 153202 (2018).
  * [35] Li, R. _et al._ A dual-axis, high-sensitivity atomic magnetometer.  _Chin. Phys. B_ 26, 120702 (2017).
  * [36] Liu, J. & Yuan, H.  Control-enhanced multiparameter quantum estimation.  _Phys. Rev. A_ 96, 042114 (2017).
  * [37] Spagnolo, N. _et al._ Quantum interferometry with three-dimensional geometry.  _Sci. Rep._ 2, 1–6 (2012).
  * [38] Genovese, M.  Real applications of quantum imaging.  _J. Optics_ 18, 073002 (2016).
  * [39] Komar, P. _et al._ A quantum network of clocks.  _Nat. Phys._ 10, 582–587 (2014).
  * [40] Donley, E. A.  Nuclear magnetic resonance gyroscopes.  In _SENSORS, 2010 IEEE_ , 17–22 (IEEE, 2010).
  * [41] Walker, T. G. & Happer, W.  Spin-exchange optical pumping of noble-gas nuclei.  _Rev. Mod. Phys._ 69, 629 (1997).
  * [42] Kornack, T., Ghosh, R. & Romalis, M.  Nuclear spin gyroscope based on an atomic comagnetometer.  _Phys. Rev. Lett._ 95, 230801 (2005).
  * [43] Hurwitz, L. & Nelson, J.  Proton vector magnetometer.  _J. Geophys. Research_ 65, 1759–1765 (1960).
  * [44] Wu, T., Blanchard, J. W., Kimball, D. F. J., Jiang, M. & Budker, D.  Nuclear-spin comagnetometer based on a liquid of identical molecules.  _Phys. Rev. Lett._ 121, 023202 (2018).
  * [45] Garcon, A. _et al._ Constraints on bosonic dark matter from ultralow-field nuclear magnetic resonance.  _Sci. Adv._ 5, eaax4539 (2019).
  * [46] Jiang, M., Su, H., Garcon, A., Peng, X. & Budker, D.  Search for axion-like dark matter with spin-based amplifiers.  _arXiv preprint arXiv:2102.01448_ (2021).
  * [47] Hou, Z. _et al._ “Super-Heisenberg” and Heisenberg Scalings Achieved Simultaneously in the Estimation of a Rotating Field.  _Phys. Rev. Lett._ 126, 070503 (2021).
  * [48] Legchenko, A., Baltassat, J.-M., Beauce, A. & Bernard, J.  Nuclear magnetic resonance as a geophysical tool for hydrogeologists.  _J. Appl. Geophys._ 50, 21–46 (2002).
  * [49] Gross, S. _et al._ Dynamic nuclear magnetic resonance field sensing with part-per-trillion resolution.  _Nat. Commun._ 7, 1–7 (2016).
  * [50] Waters, G.  A measurement of the earth’s magnetic field by nuclear induction.  _Nature_ 176, 691–691 (1955).
  * [51] Jiang, M. _et al._ Experimental benchmarking of quantum control in zero-field nuclear magnetic resonance.  _Sci. Adv._ 4, eaar6327 (2018).
  * [52] Jiang, M. _et al._ Magnetic gradiometer for the detection of zero-to ultralow-field nuclear magnetic resonance.  _Phys. Rev. Appl._ 11, 024005 (2019).
  * [53] Tayler, M. C. _et al._ Invited review article: Instrumentation for nuclear magnetic resonance in zero and ultralow magnetic field.  _Rev. Sci. Instrum._ 88, 091101 (2017).
  * [54] Jiang, M. _et al._ Interference in atomic magnetometry.  _Adv. Quantum Technol._ 3, 2000078 (2020).
  * [55] Béné, G. J.  Nuclear magnetism of liquid systems in the earth field range.  _Phys. Rep._ 58, 213–267 (1980).
  * [56] Sarkar, R., Ahuja, P., Vasos, P. R. & Bodenhausen, G.  Long-lived coherences for homogeneous line narrowing in spectroscopy.  _Phys. Rev. Lett._ 104, 053001 (2010).
  * [57] Zeeman, P.  The effect of magnetisation on the nature of light emitted by a substance (1897).
  * [58] Condon, E. U., Condon, E. & Shortley, G.  _The theory of atomic spectra_ (Cambridge University Press, 1935).
  * [59] Wood, A. _et al._ Magnetic pseudo-fields in a rotating electron–nuclear spin system.  _Nat. Phys._ 13, 1070–1073 (2017).
  * [60] Heims, S. & Jaynes, E.  Theory of gyromagnetic effects and some related magnetic phenomena.  _Rev. Mod. Phys._ 34, 143 (1962).
  * [61] Adams, R. W. _et al._ Reversible interactions with para-hydrogen enhance nmr sensitivity by polarization transfer.  _Science_ 323, 1708–1711 (2009).
  * [62] Theis, T. _et al._ Parahydrogen-enhanced zero-field nuclear magnetic resonance.  _Nat. Phys._ 7, 571–575 (2011).
  * [63] Maly, T. _et al._ Dynamic nuclear polarization at high magnetic fields.  _J. Chem. Phys._ 128, 02B611 (2008).
  * [64] Jiang, M. _et al._ Numerical optimal control of spin systems at zero magnetic field.  _Phys. Rev. A_ 97, 062118 (2018).
  * [65] Ledbetter, M. _et al._ Near-zero-field nuclear magnetic resonance.  _Phys. Rev. Lett._ 107, 107601 (2011).
  * [66] Jiang, M. _et al._ Zero-to ultralow-field nuclear magnetic resonance and its applications.  _Fundamental Research_ 1, 68–84 (2021).
  * [67] Jones, J. A. _et al._ Magnetic field sensing beyond the standard quantum limit using 10-spin noon states.  _Science_ 324, 1166–1168 (2009).
  * [68] Bermudez, A., Jelezko, F., Plenio, M. B. & Retzker, A.  Electron-mediated nuclear-spin interactions between distant nitrogen-vacancy centers.  _Phys. Rev. Lett._ 107, 150503 (2011).
  * [69] Zhao, N., Hu, J.-L., Ho, S.-W., Wan, J. T. & Liu, R.  Atomic-scale magnetometry of distant nuclear spin clusters via nitrogen-vacancy spin in diamond.  _Nat. Nanotechnol._ 6, 242–246 (2011).
  * [70] Schweiger, A. & Jeschke, G.  _Principles of pulse electron paramagnetic resonance_ (Oxford University Press on Demand, 2001).
  * [71] Xiao, D.-W., Hu, W.-H., Cai, Y. & Zhao, N.  Magnetic noise enabled biocompass.  _Phys. Rev. Lett._ 124, 128101 (2020).
  * [72] Qin, S. _et al._ A magnetic protein biocompass.  _Nat. Mater._ 15, 217–226 (2016).
  * [73] Shi, F. _et al._ Sensing and atomic-scale structure analysis of single nuclear-spin clusters in diamond.  _Nat. Phys._ 10, 21–25 (2014).
  * [74] Mamin, H. _et al._ Nanoscale nuclear magnetic resonance with a nitrogen-vacancy spin sensor.  _Science_ 339, 557–560 (2013).
  * [75] Appelt, S. _et al._ Paths from weak to strong coupling in nmr.  _Phys. Rev. A_ 81, 023420 (2010).
  * [76] Sjolander, T. F., Tayler, M. C., King, J. P., Budker, D. & Pines, A.  Transition-selective pulses in zero-field nuclear magnetic resonance.  _J. Phys. Chem. A_ 120, 4343–4348 (2016).