11institutetext: Observatoire de Paris, LERMA, Collège de France, CNRS, PSL
University, Sorbonne University, 75014, Paris — 11email:
<EMAIL_ADDRESS>22institutetext: Inter-University Centre for
Astronomy and Astrophysics, Post Bag 4, Ganeshkhind, Pune 411 007, India
33institutetext: Department of Space, Earth and Environment, Chalmers
University of Technology, Onsala Space Observatory, 43992 Onsala, Sweden
44institutetext: Ioffe Institute, Polyteknicheskaya 26, 194021 Saint-
Petersburg, and HSE University, Saint-Petersburg, Russia 55institutetext:
National Radio Astronomy Observatory, 520 Edgemont Road, Charlottesville, VA
22903, USA 66institutetext: Max-Planck-Institut für Radioastronomie, Auf dem
Hügel 69, D-53121 Bonn, Germany 77institutetext: Astrophysics Research Centre,
University of KwaZulu-Natal, Durban 4041, South Africa 88institutetext: School
of Mathematics, Statistics & Computer Science, University of KwaZulu-Natal,
Durban 4041, South Africa 99institutetext: National Radio Astronomy
Observatory, Socorro, NM 87801, USA 1010institutetext: Institut
d’Astrophysique de Paris, CNRS-SU, UMR 7095, 98bis bd Arago, 75014 Paris,
France 1111institutetext: French-Chilean Laboratory for Astronomy, IRL 3386,
CNRS and U. de Chile, Casilla 36-D, Santiago, Chile 1212institutetext: Green
Bank Observatory, PO Box 2, Green Bank, WV 24944, USA

# PKS 1413$+$135: OH and H i at $z=0.247$ with MeerKAT

F. Combes PKS 1413$+$135: OH and H i at $z=0.247$ with MeerKATPKS 1413$+$135:
OH and H i at $z=0.247$ with MeerKAT N. Gupta PKS 1413$+$135: OH and H i at
$z=0.247$ with MeerKATPKS 1413$+$135: OH and H i at $z=0.247$ with MeerKAT S.
Muller PKS 1413$+$135: OH and H i at $z=0.247$ with MeerKATPKS 1413$+$135: OH
and H i at $z=0.247$ with MeerKAT S. Balashev PKS 1413$+$135: OH and H i at
$z=0.247$ with MeerKATPKS 1413$+$135: OH and H i at $z=0.247$ with MeerKAT P.
P. Deka PKS 1413$+$135: OH and H i at $z=0.247$ with MeerKATPKS 1413$+$135: OH
and H i at $z=0.247$ with MeerKAT K. L. Emig PKS 1413$+$135: OH and H i at
$z=0.247$ with MeerKATPKS 1413$+$135: OH and H i at $z=0.247$ with MeerKAT
H.-R. Klöckner PKS 1413$+$135: OH and H i at $z=0.247$ with MeerKATPKS
1413$+$135: OH and H i at $z=0.247$ with MeerKAT D. Klutse PKS 1413$+$135: OH
and H i at $z=0.247$ with MeerKATPKS 1413$+$135: OH and H i at $z=0.247$ with
MeerKAT K. Knowles PKS 1413$+$135: OH and H i at $z=0.247$ with MeerKATPKS
1413$+$135: OH and H i at $z=0.247$ with MeerKATPKS 1413$+$135: OH and H i at
$z=0.247$ with MeerKATPKS 1413$+$135: OH and H i at $z=0.247$ with MeerKAT A.
Mohapatra PKS 1413$+$135: OH and H i at $z=0.247$ with MeerKATPKS 1413$+$135:
OH and H i at $z=0.247$ with MeerKAT E. Momjian PKS 1413$+$135: OH and H i at
$z=0.247$ with MeerKATPKS 1413$+$135: OH and H i at $z=0.247$ with MeerKAT P.
Noterdaeme PKS 1413$+$135: OH and H i at $z=0.247$ with MeerKATPKS 1413$+$135:
OH and H i at $z=0.247$ with MeerKATPKS 1413$+$135: OH and H i at $z=0.247$
with MeerKATPKS 1413$+$135: OH and H i at $z=0.247$ with MeerKAT P. Petitjean
PKS 1413$+$135: OH and H i at $z=0.247$ with MeerKATPKS 1413$+$135: OH and H i
at $z=0.247$ with MeerKAT P. Salas PKS 1413$+$135: OH and H i at $z=0.247$
with MeerKATPKS 1413$+$135: OH and H i at $z=0.247$ with MeerKAT R. Srianand
PKS 1413$+$135: OH and H i at $z=0.247$ with MeerKATPKS 1413$+$135: OH and H i
at $z=0.247$ with MeerKAT J. D. Wagenveld PKS 1413$+$135: OH and H i at
$z=0.247$ with MeerKATPKS 1413$+$135: OH and H i at $z=0.247$ with MeerKAT

(Received: November 2022; accepted: December 2022)

The BL Lac object PKS 1413$+$135 was observed by the Large Survey Project
MeerKAT Absorption Line Survey (MALS) in the L-band, at 1139 MHz and 1293-1379
MHz, targeting the H i and OH lines in absorption at $z=0.24671$. The radio
continuum might come from the nucleus of the absorbing galaxy or from a
background object at redshift lower than 0.5, as suggested by the absence of
gravitational images. The H i absorption line is detected at a high signal-to-
noise ratio, with a narrow central component, and with a red wing, confirming
previous results. The OH 1720 MHz line is clearly detected in (maser)
emission, peaking at a velocity shifted by -10 to -15 km s-1 with respect to
the H i peak. The 1612 MHz line is lost due to radio frequency interference.
The OH 1667 MHz main line is tentatively detected in absorption, but not the
1665 MHz line. Over 30 years a high variability is observed in optical depths,
due to the rapid changes of the line of sight caused by the superluminal
motions of the radio knots. The H i line has varied by 20% in depth, while the
OH-1720 MHz depth has varied by a factor of $\sim$ 3\. The position of the
central velocity and the widths also varied. The absorbing galaxy is an early-
type spiral (maybe S0) seen edge-on, with a prominent dust lane, covering the
whole disk. Given the measured mass concentration and the radio continuum size
at centimeter wavelengths (100 mas corresponding to 400 pc at $z=0.25$), the
width of the absorption lines from the nuclear regions are expected up to 250
km s-1. The narrowness of the observed lines ($<$ 15 km s-1) suggests that the
absorption comes from an outer gas ring, as frequently observed in S0
galaxies. The millimetric lines are even narrower ($<$ 1 km s-1), which
corresponds to the continuum size restricted to the core. The radio core is
covered by individual 1 pc molecular clouds, whose column density is a few
1022 cm-2, which is compatible with the gas screen detected in X-rays.

###### Key Words.:

galaxies: ISM, quasars: absorption lines, quasars: individual: PKS 1413$+$135

##  1 Introduction

Since its discovery PKS 1413$+$135 has been a mystery, raising questions
regarding the classification schemes of active galactic nuclei (AGN). It is
classified as a BL Lac for its highly compact variable polarized continuum
emission, corresponding to a relativistic jet oriented close to the line of
sight (Beichman et al. 1981). Bregman et al. (1981) found only weak emission
lines in the optical spectrum of the galaxy associated with the radio source,
supporting the BL Lac character. They determined a redshift of $z=0.26$ for
the galaxy from calcium absorption lines. However, the galaxy might not be the
host of the BL Lac object. The vast majority of BL Lac objects and radio-loud
quasars are hosted by elliptical galaxies, while PKS 1413$+$135 appears to be
hosted by an Sb-c spiral galaxy (McHardy et al. 1991). The galaxy appears
almost edge-on, and therefore the jet axis should be almost in the galaxy
plane. The jet should then trigger significant activity, such as a narrow-line
region, but none is detected in the galaxy. Very long baseline interferometry
(VLBI) measurements have shown that PKS 1413$+$135 might be a compact
symmetric object (CSO), a young radio source with a two-sided parsec-scale
jet; more precisely, the jet structure is slightly bent, like a wide-angle
tail (WAT) source (Perlman et al. 1996, 2002). This means that the jet is not
along the line of sight and does not benefit from beaming. However, the jet
orientation might have changed with scale since the core appears beamed. The
source was also classified as a red quasar since the extinction is rather
high, with a column density higher than 2$\times$1022 cm-2 and visual
extinction A${}_{v}>$ 30 (Stocke et al. 1992); this high extinction explains
the deficit of soft X-rays from PKS 1413$+$135\. Although the radio source
coincides with the galaxy center, it was thought that the radio source might
be a background source, given the absence of emission lines and strong thermal
infrared flux.

Carilli et al. (1992) detected H i 21 cm in absorption in the galaxy at
$z=0.24671$, with a high column density ($>$ 1022 cm-2), supporting the high
extinction found by Stocke et al. (1992), and confirming the classification of
the galaxy as a late-type spiral. Molecular clouds were also observed in
absorption by Wiklind & Combes (1994), supporting a high gas column-density.
Monitoring over two years revealed high variation in the molecular absorption
opacity spectrum, attributed to the superluminal character of the BL Lac. This
high gas column density exacerbates the problem of the absence of thermal
emission expected from the gas heated by an internal active nucleus. An
alternative interpretation of the variability of BL Lac objects is that the
background AGN suffers microlensing by the stars of a foreground galaxy
(Ostriker & Vietri 1990). For better efficiency of stellar microlensing, while
avoiding macrolensing, the AGN should be close to the foreground galaxy. The
Very Long Baseline Array (VLBA) data from Perlman et al. (2002) provide a high
resolution map of the HI absorption, which is mainly in front of the eastern
mini-jet 25 mas from the radio core. The molecular gas absorbs in front of the
core, and may be also the OH line. The core happens to be offset by 13 mas
(=50 pc) from the gravity center of the spiral galaxy. This offset might favor
a background location for the AGN, although Perlman et al. (2002) favor a
patchy absorption from the spiral galaxy. If the AGN is behind the spiral, it
must be at z $<$0.3 to avoid arcs and lensing images, which are not observed,
and thus close to the spiral ($z=0.25$) we should see its galaxy host in HST
images.

Recently, Vedantham et al. (2017b, a) discovered symmetric achromatic
variability (SAV) in PKS 1413$+$135 in the forms of a 1yr long U-shaped
perturbation of the light curve, with a symmetric feature 5 yr later. This
nicely corresponds to the signature of gravitational milli-lensing when knots
in the relativistic jet cross caustics created by lenses of mass between 103
and 106 M⊙ within the intervening galaxy. Combining all VLBI monitoring at
several wavelengths, Readhead et al. (2021) conclude that PKS 1413$+$135 is a
background BL Lac quasar, with radio jets aligned along the line of sight and
located at 0.25 $<$ z $<$ 0.5. The intervening galaxy is a nearly edge-on
Seyfert 2 spiral, providing the milli-lensing episodes responsible in part for
the variability, the other part being intrinsic to the BL Lac. They show that
the jet appears quasi-perpendicular to the galaxy major axis, and has an
extension of 110 mas. Since the optical spectrum from Vedantham et al. (2017b)
should reveal signatures from the background galaxy, which are not detected,
the localization of the radio continuum source is not yet settled. It is
possible that the BL Lac is not hosted by a massive elliptical, but by a
galaxy half as luminous as the foreground Seyfert; thanks to the high dust
obscuration observed (by a factor of $\sim$10 in the near-infrared) and the
luminosity dimming at z=0.5, the expected flux ratio between the two galaxies
could be as high as 100.

The size of the radio source at centimeter wavelengths is 0.11” $\sim$ 400pc
(Perlman et al. 1996). This broadens the absorption spectrum because the gas
actually involved in the absorption traces the velocity gradient over that
length in the disk for H i 21cm and OH 18cm. The millimeter lines are much
narrower since the continuum is only significant in the core (Wiklind & Combes
1994).

Kanekar & Chengalur (2002) detect the main OH line at 1667 MHz, while Darling
(2004) detects the conjugate satellite lines (1612 MHz in absorption and 1720
MHz in emission), but none of the main OH lines. The 1720 and 1612 MHz line
rest frequencies have different dependences on the fine structure constant
$\alpha$ and the proton-electron mass ratio $\mu$, and they have been used to
constrain the variation of these fundamental constants (Kanekar & Chengalur
2002; Kanekar et al. 2004, 2018; Darling 2004). Because the lines are
conjugate, and their optical depths are perfectly symmetric, they are a
privileged tool, eliminating the possible biases that the lines are not coming
from the same gas along the line of sight. In the present work we report new H
i and OH observations of PKS 1413$+$135 carried out as part of the MeerKAT
Absorption Line Survey (MALS; Gupta et al. 2017). We quantify the time
variability of the various components, and through comparison with previous
millimeter-wave observations we draw some conclusions on the nature and
location of the absorptions.

This paper is structured as follows. Section 2 presents the details of the
observations and data analysis with ARTIP. The results in terms of optical
depth and column density are quantified in Section 3, as well as the time
variability. We discuss in Sect. 4 the origin of the various absorption
components in comparison with the dense molecular gas absorptions. Section 5
summarizes our conclusions. Throughout the paper the velocity scale is defined
with respect to $z=0.24671$ (barycentric reference frame), which corresponds
to the main atomic (Carilli et al. 1992) and molecular (e.g., CO, HCO+, HCN,
and OH) absorption components detected at millimeter wavelengths (Wiklind &
Combes 1994) or centimeter wavelengths (Kanekar & Chengalur 2002; Darling
2004). To compute distances we adopt a flat $\Lambda$CDM cosmology, with
$\Omega_{m}$=0.29, $\Omega_{\Lambda}$=0.71, and the Hubble constant H0 = 70 km
s-1Mpc-1. At the distance of the PKS 1413$+$135 absorber, 1 arcsec corresponds
to 3.9 kpc in physical units.

##  2 Observations and data analysis

The field centered at PKS 1413$+$135 was observed as part of MALS on October
16, 2020, using 59 antennas of the MeerKAT-64 array and the 32K mode of the
SKA Reconfigurable Application Board (SKARAB) correlator. The total bandwidth
of 856 MHz centered at 1283.9869 MHz was split into 32,768 frequency channels.
The resultant frequency resolution is 26.123 kHz, or 6.1 km s-1, at the center
of the band. The correlator dump time was 8 seconds. The objects 3C286 and
PKS1939–638 were observed for flux density scale, delay, and bandpass
calibrations. The compact radio source J1347+1217 was observed for complex
gain calibration. The total on-source time on PKS 1413$+$135 was 56 mins. The
data were processed using the latest version of the Automated Radio Telescope
Imaging Pipeline ARTIP based on CASA 6.3.0 (Casa22). The details of ARTIP and
the calibration steps are provided in Gupta et al. (2021).

The spectral line processing of the MALS datasets through ARTIP partitions the
wideband into 16 spectral windows (SPWs) labeled SPW-0 to SPW-15 (for details,
see Gupta et al. 2021). In short, the measurement sets for each SPW are
processed for continuum imaging with self-calibration. We performed three
rounds of phase-only with a solution interval of 1 min, and one round of
amplitude-and-phase self-calibration with a solution interval of 5 min. These
self-calibration gain solutions were applied to the visibilities prior to
continuum subtraction and cube imaging. The H i 21cm and the four OH lines for
$z=0.24671$ would fall in SPWs 4, 8, and 9. Here, we mainly focus on the
stokes $I$ radio continuum and spectral line properties of the target obtained
using the robust = 0 weighting of the visibilities for these SPWs, as
implemented in CASA. The frequency ranges covered by these SPWs and the
corresponding spectral rms are summarized in Table 1. In the table we also
provide the synthesized beams and continuum flux densities for the continuum
images generated using the visibilities from these SPWs.

Table 1: Properties of PKS1413+135 observations for relevant SPWs.

SPW-Id. | Frequency range | Spectral rms | Synthesized beam | Total flux densitya
---|---|---|---|---
| (MHz) | (mJy beam-1) |  | (Jy)
SPW-4b | 1083.4 - 1143.6 | 1.0 |  $16\aas@@fstack{\prime\prime}7\times 8\aas@@fstack{\prime\prime}1$, $-21.8^{\circ}$ | 1.3635$\pm$0.0032
SPW-8c | 1297.4 - 1357.6 | 0.9 |  $14\aas@@fstack{\prime\prime}0\times 6\aas@@fstack{\prime\prime}8$, $-22.7^{\circ}$ | 1.1682$\pm$0.0067
SPW-9d | 1350.9 - 1411.1 | 0.9 |  $13\aas@@fstack{\prime\prime}5\times 6\aas@@fstack{\prime\prime}5$, $-22.8^{\circ}$ | 1.1335$\pm$0.0056
111 $a$$a$footnotetext: Obtained using a single-Gaussian component fitted to
the radio emission. $b$$b$footnotetext: Covers the H i 21cm line.
$c$$c$footnotetext: Covers the OH 1612, 1665, and 1667 MHz lines.
$d$$d$footnotetext: Covers the OH 1720 MHz line.

Figure 1:  MeerKAT stokes $I$ H i 21cm spectrum of PKS 1413$+$135\. The
individual Gaussian components fitted to the absorption profile (dotted
lines), the total fit (dashed line), and the residuals (dot-dashed lines) are
also shown. For comparison, we also show the residuals with respect to the GBT
spectrum reproduced using the Gaussian components provided in Darling (2004)
(dot-dashed green line). The individual XX and YY MeerKAT spectra are shown as
points and crosses, respectively.

##  3 Results

###  3.1 H i absorption

Figure 2:  MeerKAT Stokes-$I$ spectrum (blue) of OH 1720 MHz line at
$z=0.24671$ toward PKS 1413$+$135\. The single-Gaussian component (peak =
0.0076$\pm$0.0031 Jy; shift = -10.37$\pm$0.99 km s-1; sigma = 2.8$\pm$1 km
s-1) fitted to the emission line is shown as a dashed line. The emission
corresponds to a weakly masing line, while the 1612 MHz emission is seen in
absorption (Darling 2004). The RFI at this frequency prevented MeerKAT from
observing the emission.

Figure 1 displays MeerKAT H i 21cm absorption in front of PKS 1413$+$135 at
$z=0.24671$. Given the asymmetry of the line, the best fit providing a
satisfactory model of the total spectrum required five Gaussian components.
The characteristics of all the components are listed in Table 2. These
components may not be physical, but are just practical to completely reproduce
the global profile. There is a wide component in the red wing, which was also
seen by Darling (2004), although it was thought at that time that it might
come from the radio frequency interference (RFI). The global integrated H i
21cm optical depth is $\int\tau$dv$=10.86\pm 0.32$ km s-1. Assuming a
homogeneous and optically thin gas, at spin temperature $T_{\rm s}$= 100K,
with a covering factor $f_{c}^{\tiny\mbox{H\,{\sc i}}}$ = 1, and using the
equation

$N{(\mbox{H\,{\sc i}})}=1.823\times 10^{18}\leavevmode\nobreak\ {T_{\rm
s}\over f_{\rm c}^{\tiny\mbox{H\,{\sc i}}}}\int\leavevmode\nobreak\
\tau(v)\leavevmode\nobreak\ {\rm d}v\leavevmode\nobreak\ {\rm cm^{-2}}$ (1)

we obtain the H i column density $N$(H i) = $(1.98\pm 0.06)\times 10^{21}$
cm-2.

Table 2: Multiple Gaussian fits to the H i 21cm absorption. ID | Center | $\sigma$ | $\tau_{p}$
---|---|---|---
| (km s-1) | (km s-1) |
1 | -47.1 $\pm$ 2.1 | 11.1 $\pm$ 2.2 | 0.00434 $\pm$ 0.00057
2 | 32.9 $\pm$ 4.0 | 15.9 $\pm$ 5.2 | 0.0147 $\pm$ 0.0045
3 | 72.1 $\pm$ 9.2 | 24.9 $\pm$ 4.8 | 0.0108 $\pm$ 0.0017
4 | -2.19 $\pm$ 0.024 | 6.35 $\pm$ 0.064 | 0.463 $\pm$ 0.015
5 | -0.322$\pm$ 0.86 | 11.7 $\pm$ 0.89 | 0.0963 $\pm$ 0.013

###  3.2 OH main and satellite lines

Figure 2 shows that the 1720 MHz satellite line is well detected with MeerKAT,
while Fig. 3 displays the spectrum of the OH main lines. There is only a hint
of a detection of the 1667 MHz line, which was already detected by Kanekar &
Chengalur (2002). The 1665 MHz line is not detected. In local thermodynamic
equilibrium (LTE), the relative strengths of the 18cm lines are expected to
follow the ratio 1612:1665:1667:1720 MHz = 1:5:9:1. Since we tentatively
detect the 1667 MHz line, there is no surprise that there is only an upper
limit for the 1665 MHz line.

Figure 3:  MeerKAT Stokes-$I$ spectrum (blue for robust=0 and orange for
robust=1) of the OH main lines absorption at $z=0.24671$ toward PKS
1413$+$135\. The spectral rms noise is 0.9 mJy beam-1 channel-1. Vertical
arrows indicate the expected positions of the 1665 and 1667 MHz lines. There
is a tentative dip in the flux at the location of 1667 MHz line corresponding
to $\int\tau$dv=0.03$\pm$0.01 km s-1.

We detect the 1720 MHz line in emission (see Fig. 2), while the expected
absorption at 1612 MHz is undetected due to RFI. The two satellite lines have
been clearly detected by Darling (2004) and Kanekar et al. (2018) to be
exactly symmetric, with their sum of negative and positive integrated optical
depths cancelling out. The OH satellite lines were also observed to exhibit
conjugate behavior in PKS1830-211 (Combes et al. 2021).

Table 3: Multiple epochs of HI and OH absorptions at $z=0.24671$

Epoch (HI) | S(Jy) | $\tau_{HI}$ | Ref
---|---|---|---
1992 | 1.25$\pm$0.15 | 0.34 $\pm$ 0.04 | (1)
2003 | 1.37$\pm$0.02 | 0.68 $\pm$ 0.03 | (2)
2020 | 1.363$\pm$0.003 | 0.463 $\pm$ 0.015 | (3)
Epoch (OH) | S(Jy) | $\tau_{1720}$ | Ref
2003 | 1.085$\pm$0.002 | 0.0097 $\pm$ 0.0017 | (2)
2003 |  | 0.014 $\pm$ 0.002 | (4)
2010-12 |  | 0.0071 $\pm$ 0.0003 | (5)
2020 | 1.133$\pm$0.006 | 0.0047 $\pm$ 0.0018 | (3)
222 $a$$a$footnotetext:  Ref (1): Carilli et al. (1992) VLA with resolution
1.28 km s-1. Ref (2): Darling (2004) GBT with 0.8 km s-1. Ref (3): this work.
Ref (4): Kanekar et al. (2004) WSRT with 1.1 km s-1. Ref (5): Kanekar et al.
(2018) Arecibo with 0.18 km s-1.

###  3.3 Time variability

The radio continuum emission from the background blazar is found to vary in
time, sometimes with large factors, likely due to milli-lensing (Peirson et
al. 2022). However, the optical depth of the H i and OH absorption lines vary
even more (see Table 2 and Fig. 4).

Figure 1 shows that the difference in the H i optical depth between Darling
(2004) and the present work is detected at 6$\sigma$, while the peak optical
depth in 2003 was twice that of 1992 (see Table 2). The bulk of the H i
absorption does not come from the core of the radio source, but from one knot
in the radio jet, the eastern mini-lobe at 25 mas (=100 pc) from the core, as
shown with the VLBA by Perlman et al. (2002). Since the jet has superluminal
motions, it is expected that the line of sight explored is varying relatively
quickly, especially in a clumpy medium (Wiklind & Combes 1997). The continuum
source is not extended more than a few tens of milliarcseconds, and the
absorbing medium has the size of an interstellar cloud.

We compared all available OH 1720 MHz spectra from 2003 with the present
MeerKAT spectra. Figure 4 shows that the peak optical depth is varying, as are
the central velocity and the width of the spectral features. The absorption
peaks between 5% and 14% of the continuum, with a central velocity from -10 to
-15 km s-1 with respect to the CO and H i velocities. The FWHM varies between
7 and 11 km s-1. Although observations at a smaller number of epochs are
available for the H i absorption, it appears that it is slightly less
variable, as expected for a less clumpy component. As already noted in Combes
et al. (2021) for PKS1830-211, the atomic component is more diffuse, more
extended both in radius and in height in the galaxy plane.

Additional sources of variation might be due to the milli-lensing events,
observed through the symmetric achromatic variability events monitored by
Vedantham et al. (2017b, a). These events usually last more than one year, and
they occurred clearly in 1993–1995, 2009–2010, 2014–2015 (Peirson et al.
2022). They also affected some absorption measurements in 2010 (Kanekar et al.
2018), and also in the molecular component (Wiklind & Combes 1994, 1997).

Figure 4:  MeerKAT OH 1720 MHz spectrum taken in 2020 (blue), compared to the
spectrum observed in 2003 by Kanekar et al. (2004) with WSRT (green), the
spectrum observed in 2003 by Darling (2004) with GBT (red), and the spectrum
observed in 2010–2012 by Kanekar et al. (2018) with Arecibo (black). All
spectra have been normalized to their observed continuum.

##  4 Kinematic interpretation of the spectra

PKS 1413$+$135 has been considered for a long time as a puzzling source. It
shows two sets of characteristics that appear contradictory. On one side it
exhibits a BL Lac character with a radio jet relativistically boosted toward
the observer, with apparent transverse velocity up to v$\sim$ 1.5c measured
through proper motions with VLBA (e.g., Lister et al. 2016). This would mean
that the jet orientation is close to the line of sight; however, it also has a
counterjet, which suggests a re-orientation somewhere along the jet.
Paradoxically, it has also been classified as a CSO, which are typically small
and young objects that have jet axes not aligned close to the line of sight.
CSOs are often thought to have re-starting jet activity, and PKS 1413$+$135 in
2019 revealed violent variations (more than an order of magnitude) in
$\gamma$-rays (Principe et al. 2021). It is the first CSO detected at very
high energies, larger than 100 GeV (Blanch et al. 2022). Gan et al. (2022)
show that the source becomes harder in frequency when brighter, and that the
recent flare of PKS1413+135 might be linked to its re-started jet activity.

Powerful radio jets are not generally associated with spiral galaxies, and in
the case of PKS 1413$+$135 the radio source core corresponds to the nucleus of
an edge-on early-type spiral galaxy (McHardy et al. 1994). Many authors (e.g.,
McHardy et al. 1991; Perlman et al. 2002) have questioned its association with
this galaxy since no infrared activity is observed, contrary to what is
expected for the host of a powerful AGN that has a radio jet aligned along its
plane. The true galaxy host could be a background galaxy, hard to detect since
hidden by the foreground edge-on galaxy. The latter could have some lensing
effect on the BL Lac, but since no lensed images are seen, its redshift may be
limited to the range 0.2467 $<$ z $<$ 0.5 (Readhead et al. 2021). Only milli-
lensing events are possible, and would correspond to the SAV events observed
by Vedantham et al. (2017b).

The various absorption lines coming from the foreground edge-on galaxy and
observed in front of the background radio source provide precious information.
This data could bring new clues to the nature of this puzzling source.

###  4.1 Galaxy model

The H i spectrum is remarkable for its main narrow central component at zero
velocity, and its redshifted wing, extending up to 100 km s-1. It is also
striking that it is not seen toward the radio core, but mainly from a
counterjet knot, 25ms (=100 pc) northeast from the core (Perlman et al. 2002),
on the minor axis of the galaxy. All the other absorption lines come from the
core region, where the X-ray absorption indicates a large column density N(H)
$\sim$ 6 $\times$ 1022 cm-2 (Stocke et al. 1992).

The systemic velocity was set such that the main narrow H i absorption is
centered at zero velocity. All lines peak at this velocity, except the OH 18cm
satellite lines. This might also be the systemic velocity of the foreground
galaxy since most of the absorption is coming from the minor axis of the
rotating spiral. On the contrary, the two OH satellite lines appear at V=-10
or -15 km s-1, according to the observing epoch (see Fig. 4). The fact that
the H i absorption does not come from the core (almost coincident with the
foreground galaxy nucleus) is not surprising since there is a well-known H i
depletion in spiral galaxy centers (e.g., Walter et al. 2008).

Figure 5: Comparison between the CO(1-0) and H i lines. All spectra are
normalized to their respective continuum. The blue curve is the H i
absorption, observed in 1992 with the VLA by Carilli et al. (1992) with
channels 2.57 km s-1, and the red curve is from this work (H i spectrum
obtained in 2020 with 6 km s-1 channels). Superposed on the CO(1-0) spectrum
observed in 1996 (black) from Wiklind & Combes (1997) is a CO(1-0) spectrum
(green) obtained with the Institut de Radioastronomie Millimétrique (IRAM)
interferometer in 2008, showing no variation.

The edge-on spiral galaxy appears to be of very early type, perhaps a
lenticular. This is supported by the large and massive bulge, as
photometrically fitted by McHardy et al. (1994). With the bulge alone the
effective radius is best fit as 11 kpc, but when the bulge and disk are fitted
together the effective bulge radius is found to be 9.7 kpc. The disk is quite
thick, with characteristic height-to-size ratio h/R = 0.24. The mass
concentration is high, and a Sersic index of n=4 is used. With this mass
concentration, for a galaxy of total mass $\sim$ 3 $\times$ 1010 M⊙ the
velocity width at the nucleus is expected to be on the order of 300 km s-1.

Since all absorption lines are quite narrow, with a FWHM lower than 15 km s-1,
except for the weak red wing of the H i spectrum, it is likely that the
absorbing gas is not located in the nuclear region, but is orbiting in an
outer gaseous ring. These gaseous rings are frequently observed in lenticular
galaxies (e.g., van Driel & van Woerden 1991). The ring likely corresponds to
the dust lane observed crossing the disk in the HST image of McHardy et al.
(1994). Because the dust lane covers the whole disk extent, it must be an
outer ring.

###  4.2 Nature of the absorbing gas

The H i 21cm and OH 18cm lines are narrow, but the millimetric lines are even
narrower. As reported by Wiklind & Combes (1997), the CO(1-0) line is composed
of two main components, of FWHM = 0.6-1.9 km s-1, the deepest absorption being
at velocity V=-0.5 km s-1 and the other at +0.6 km s-1. According to the epoch
observed, a third weaker line might appear wider (FWHM 0.6-7 km s-1) in the
red wing located at V=1-2 km s-1. The position, intensity, and widths of these
components varying with time, by factors of $\sim$2 over two years, due to the
change in line of sight toward the continuum core, with the knots moving at
superluminal speeds.

Figure 5 highlights the difference in widths between the centimeter and
millimeter absorption lines. This is easily interpreted due to the difference
in extent of the background continuum at these two frequencies. As shown in
Perlman et al. (1996, 2002), the radio source contains a jet and counterjet,
of size 100 mas $\sim$ 400 pc, but only at centimeter wavelengths. Already at
43 GHz the source is reduced to a point-like core, with the beam size ($<$ 1
mas). The millimeter pencil beam of the continuum may encounter only one
cloud, of size 1 pc, corresponding to the observed width, according to the
size–linewidth relation (Solomon et al. 1987). The medium is certainly clumpy,
and several components may occur on the edge-on gaseous ring, for a kiloparsec
width of the ring. For H i and OH absorption lines the continuum covers a
broader region of the ring, explaining the observed wider absorptions. A
schematic view of the various components is sketched in Fig. 6.

Figure 5 also shows that the CO(1-0) spectrum has not varied between 1996 and
2008. This might indicate that the core component, toward which the millimeter
absorptions are detected, is not actually moving, and that the superluminal
motions only affect the jet. At centimeter wavelengths, the absorptions occur
toward the more extended continuum emission, including the jets, for both OH
and H i , and this explains their variations. Although the H i and OH
absorptions are not centered at the same point, their extensions overlap, and
that is why their variations are somewhat correlated: their optical depths
both decreased between 2003 and 2020.

Figure 6: Schematic view of the atomic and molecular absorptions. The HST
legacy archive image of the absorbing galaxy in H-band (F160W) is plotted on
the left, and the 5 GHz VLBI image from Perlman et al. (1996) showing the jet
(west) and counterjet (east) is on the right. The contours are spaced by a
factor of $\sqrt{2}$ in brightness. The component N is the nucleus or radio
core. The counterjet components are indicated by A, B, and C. The HI
absorption is centered in B, and the molecular components in N. The jet
structure is aligned on the minor axis of the galaxy.

##  5 Conclusions

We report results from the MALS project on the MeerKAT array toward the BL Lac
source PKS 1413$+$135\. A strong H i absorption is detected, decomposed into
five Gaussian components, with high signal-to-noise ratio (S/N$\sim$ 30 in
channels of 6 km s-1 ). The main OH component at 1667 MHz is tentatively
detected, but not the 1665 MHz component. The 1720 MHz satellite line is
clearly detected, while radio frequency interference (RFI) prevented us from
detecting the expected absorption line at 1612 MHz. The 1720 MHz satellite
line is observed in stimulated (maser) emission, due to radiative pumping
through far-infrared radiation. Previous observations (Darling 2004; Kanekar
et al. 2018) have shown that the two satellite lines are perfectly conjugate,
with the 1612 MHz line in absorption. The absorption lines have been observed
since 1992. The H i component reveals a 30% variation in its depth in 2020
with respect to 2003 (see Table 2), but only slightly in velocity position or
width. The 1720 MHz OH satellite line shows more variation: a factor of 3 in
optical depth, in velocity position, and in width (see Fig. 4). These
variations can easily be explained by the change in the line of sight, related
to the superluminal motion of the background radio continuum knots.

Due to discovered SAV in PKS 1413$+$135, attributed to milli-lensing events,
it is possible that the continuum source belongs to a background AGN, at
redshift lower than 0.5 (Peirson et al. 2022). The size of the jet and
counterjet covers 400 pc at the distance of the foreground absorbing galaxy in
the rest frame 18 and 21cm wavelengths. The very narrow (FWHM $<$ 15 km s-1)
OH and H i absorption lines suggest that the absorbing interstellar medium in
the foreground edge-on galaxy is located in an outer gas ring. The foreground
galaxy is an early-type spiral, likely a lenticular given its very massive
bulge. The width sampled by the gas in the 400 pc nuclear region is expected
to be much larger (250 km s-1). In the millimeter domain, the continuum size
is much smaller ($<$ 4 pc), restricted to its core. Along this pencil beam the
line of sight may encounter individual parsec-scale clouds, accounting for the
extremely narrow (FWHM$\sim$ 1 km s-1) absorption in CO, HCO+, or HCN (Wiklind
& Combes 1997). The typical surface density of molecular clouds (a few 1022
cm-2) corresponds to the optical depth found in X-ray data (Stocke et al.
1992).

###### Acknowledgements.

The MeerKAT telescope is operated by the South African Radio Astronomy
Observatory, which is a facility of the National Research Foundation, an
agency of the Department of Science and Innovation. KM acknowledges support
from the National Research Foundation of South Africa. The MeerKAT data were
processed using the MALS computing facility at IUCAA
(https://mals.iucaa.in/releases). The Common Astronomy Software Applications
(CASA) package is developed by an international consortium of scientists based
at the National Radio Astronomical Observatory (NRAO), the European Southern
Observatory (ESO), the National Astronomical Observatory of Japan (NAOJ), the
Academia Sinica Institute of Astronomy and Astrophysics (ASIAA), the CSIRO
division for Astronomy and Space Science (CASS), and the Netherlands Institute
for Radio Astronomy (ASTRON) under the guidance of NRAO. The National Radio
Astronomy Observatory is a facility of the National Science Foundation
operated under cooperative agreement by Associated Universities, Inc. We made
use of the NASA/IPAC Extragalactic Database (NED).

## References

  * Beichman et al. (1981) Beichman, C. A., Neugebauer, G., Soifer, B. T., et al. 1981, Nature, 293, 711
  * Blanch et al. (2022) Blanch, O., Sitarek, J., & Striskovic, J. 2022, The Astronomer’s Telegram, 15161, 1
  * Bregman et al. (1981) Bregman, J. N., Lebofsky, M. J., Aller, M. F., et al. 1981, Nature, 293, 714
  * Carilli et al. (1992) Carilli, C. L., Perlman, E. S., & Stocke, J. T. 1992, ApJ, 400, L13
  * Combes et al. (2021) Combes, F., Gupta, N., Muller, S., et al. 2021, A&A, 648, A116
  * Darling (2004) Darling, J. 2004, ApJ, 612, 58
  * Gan et al. (2022) Gan, Y.-Y., Zhang, J., Yao, S., et al. 2022, arXiv e-prints, arXiv:2206.09367
  * Gupta et al. (2021) Gupta, N., Jagannathan, P., Srianand, R., et al. 2021, ApJ, 907, 11
  * Gupta et al. (2017) Gupta, N., Srianand, R., Baan, W., et al. 2017, ArXiv e-prints [arXiv:1708.07371]
  * Kanekar & Chengalur (2002) Kanekar, N. & Chengalur, J. N. 2002, A&A, 381, L73
  * Kanekar et al. (2004) Kanekar, N., Chengalur, J. N., & Ghosh, T. 2004, Phys. Rev. Lett., 93, 051302
  * Kanekar et al. (2018) Kanekar, N., Ghosh, T., & Chengalur, J. N. 2018, Phys. Rev. Lett., 120, 061302
  * Lister et al. (2016) Lister, M. L., Aller, M. F., Aller, H. D., et al. 2016, AJ, 152, 12
  * McHardy et al. (1991) McHardy, I. M., Abraham, R. G., Crawford, C. S., et al. 1991, MNRAS, 249, 742
  * McHardy et al. (1994) McHardy, I. M., Merrifield, M. R., Abraham, R. G., & Crawford, C. S. 1994, MNRAS, 268, 681
  * Ostriker & Vietri (1990) Ostriker, J. P. & Vietri, M. 1990, Nature, 344, 45
  * Peirson et al. (2022) Peirson, A. L., Liodakis, I., Readhead, A. C. S., et al. 2022, ApJ, 927, 24
  * Perlman et al. (1996) Perlman, E. S., Carilli, C. L., Stocke, J. T., & Conway, J. 1996, AJ, 111, 1839
  * Perlman et al. (2002) Perlman, E. S., Stocke, J. T., Carilli, C. L., et al. 2002, AJ, 124, 2401
  * Principe et al. (2021) Principe, G., Di Venere, L., Orienti, M., et al. 2021, MNRAS, 507, 4564
  * Readhead et al. (2021) Readhead, A. C. S., Ravi, V., Liodakis, I., et al. 2021, ApJ, 907, 61
  * Solomon et al. (1987) Solomon, P. M., Rivolo, A. R., Barrett, J., & Yahil, A. 1987, ApJ, 319, 730
  * Stocke et al. (1992) Stocke, J. T., Wurtz, R., Wang, Q., Elston, R., & Jannuzi, B. T. 1992, ApJ, 400, L17
  * van Driel & van Woerden (1991) van Driel, W. & van Woerden, H. 1991, A&A, 243, 71
  * Vedantham et al. (2017a) Vedantham, H. K., Readhead, A. C. S., Hovatta, T., et al. 2017a, ApJ, 845, 90
  * Vedantham et al. (2017b) Vedantham, H. K., Readhead, A. C. S., Hovatta, T., et al. 2017b, ApJ, 845, 89
  * Walter et al. (2008) Walter, F., Brinks, E., de Blok, W. J. G., et al. 2008, AJ, 136, 2563
  * Wiklind & Combes (1994) Wiklind, T. & Combes, F. 1994, A&A, 286, L9
  * Wiklind & Combes (1997) Wiklind, T. & Combes, F. 1997, A&A, 328, 48


# Criticality-Enhanced Precision in Phase Thermometry

Mei Yu  Naturwissenschaftlich-Technische Fakultät, Universität Siegen, Siegen
57068, Germany  H. Chau Nguyen  Naturwissenschaftlich-Technische Fakultät,
Universität Siegen, Siegen 57068, Germany  Stefan Nimmrichter
Naturwissenschaftlich-Technische Fakultät, Universität Siegen, Siegen 57068,
Germany

###### Abstract

Temperature estimation of interacting quantum many-body systems is both a
challenging task and topic of interest in quantum metrology, given that
critical behavior at phase transitions can boost the metrological sensitivity.
Here we study non-invasive quantum thermometry of a finite, two-dimensional
Ising spin lattice based on measuring the non-Markovian dephasing dynamics of
a spin probe coupled to the lattice. We demonstrate a strong critical
enhancement of the achievable precision in terms of the quantum Fisher
information, which depends on the coupling range and the interrogation time.
Our numerical simulations are compared to instructive analytic results for the
critical scaling of the sensitivity in the Curie-Weiss model of a fully
connected lattice and to the mean-field description in the thermodynamic
limit, both of which fail to describe the critical spin fluctuations on the
lattice the spin probe is sensitive to. Phase metrology could thus help to
investigate the critical behaviour of finite many-body systems beyond the
validity of mean-field models.

††preprint: APS/123-QED

##  I Introduction

Temperature is one of the key thermodynamical quantities that determine the
feasibility and performance of quantum experiments and technological
applications such as quantum sensors [1] and quantum computers [2, 3]. Its
precise measurement is therefore crucial, but also extremely challenging in
the often desired regime of low temperatures [4, 5, 6, 7, 8]. To estimate the
temperature of a quantum system, one brings in a probe, i.e., the quantum
thermometer, couples it to the system and measures it. Repeating this
procedure, either in parallel with independently prepared copies, or
sequentially with the same continuously coupled probe, one infers the
temperature from the statistics of measurement outcomes. One can measure the
probe after it has already achieved thermal equilibrium with the system [9,
10, 11, 6, 7, 12, 13, 14], or one can monitor its non-equilibrium dynamics
[15, 16, 17, 18, 19, 20, 21, 22, 23, 24, 25, 26, 27].

For the assessment and optimization of quantum thermometry protocols, two main
paradigms are pursued in the literature. Either one characterizes and and
maximizes the globally achievable accuracy averaged over a range of possible
temperature values, typically based on Bayesian estimation and data analysis
methods [28, 12, 29, 30, 31]. Or one focuses on the locally achievable
precision of the estimate around fixed temperatures, as described by the
Cramér-Rao bound in the limit of asymptotic data [32, 33, 9]. The main goal is
then to optimize the measurement strategy with respect to the probe’s quantum
Fisher information (QFI) on temperature. This local paradigm is the most
widespread one [10, 11, 5, 7, 34, 21], and we will adopt it here to
investigate temperature estimation near the critical point of an interacting
spin system.

Physical systems close to a quantum or classical phase transition are
extremely sensitive to small fluctuations of the relevant parameters, which
can be exploited for quantum metrology and enhanced sensing. Previous
theoretical studies in critical quantum metrology were focused on quantum
critical systems in the ground state as probes for the estimation of
Hamiltonian parameters [35, 36, 37, 38, 39], e.g., the magnitude of an
external magnetic field, or intra-system coupling strengths. However,
exploiting the phase transition for enhanced temperature estimation is still a
open problem.

In principle, one can infer the temperature of an equilibrium system from
energy measurements, but in a strongly interacting many-body system, this
would require global multi-particle measurements that are difficult to
implement without disturbing the system. In certain types of systems away from
criticality, measurements of extensive observables could be a viable
alternative [5]. Another one is the minimally invasive “local” probing of
small subsystems in order to infer the global temperature [6, 40, 41]. In any
case, an exponentially bad sensitivity at low temperatures is expected both
for local and global strategies, unless the probed system is gapless [8, 42].

In this work, we study phase thermometry: a non-invasive temperature
estimation strategy in which a probe is brought into contact with a localized
part of the sample system and undergoes pure dephasing, such that the sample
temperature can be estimated by measuring the probe coherence over time.
Previous works have considered phase thermometry of large bosonic samples with
Ohmic spectral densities [43, 24] as well as non-interacting Fermi gases [21,
22, 23, 44].

Here, we set out to explore temperature estimation of interacting systems
close to a phase transition and investigate phase thermometry of a two-
dimensional Ising model at zero field [45]. Ising spin lattices may represent
various physical scenarios such as arrays of Rydberg atoms [46, 47, 48, 49],
systems of nitrogen vacancy (NV) centers [50, 51, 52, 53, 54, 55], a spin-
polarized Fermi gas [56, 57], and macroscopic magnetic systems [45].

A single additional spin immersed in the sample and interacting with a
controlled range of surrounding sample spins will experience temperature-
dependent dephasing and thus serve as the phase probe. We will show that its
sensitivity is critically enhanced at the phase transition point, allowing us
to read out a significant portion of the temperature information contained in
the state of those surrounding sample spins. We compare our Monte-Carlo
simulation results to two well-known effective descriptions: the Curie-Weiss
model [45] in which each spin couples to the uniform sample average of all
other spins, and the mean-field approximation. Both capture the critical
sensitivity boost qualitatively, but they fail to predict both the critical
point and the exact enhancement quantitatively.

We structure the manuscript as follows: section II briefly reviews how to
characterise the temperature estimation precision by means of the Fisher
information. In section III, we introduce our model for the phase thermometer
immersed in a thermal Ising spin lattice, the time evolution of which we solve
using a previously developed framework [58]. We then present our main result,
the quantitative assessment of thermometric performance in terms of the Fisher
information around the critical point, in section IV. It starts with the
effective Curie-Weiss model, which is analytically tractable in the
thermodynamic limit, followed by the exact numerical treatment of a finite
zero-field Ising lattice of $20\times 20$ spins and the associated mean-field
approximation. The scaling of the Fisher information is evaluated as a
function of the number of lattice spins the probe interacts with, and we also
consider a high-temperature expansion to approximate the behaviour in the
paramagnetic regime far above the critical temperature. Finally, we conclude
in section V.

##  II Precision of temperature estimates

Quantum estimation protocols generically start from a quantum system
initialized in a fiducial state $\rho_{p}(0)$, which picks up information
about a parameter, here the temperature $T$ of a thermal environment, as it
evolves to a state $\rho_{p}(t,T)$ at time $t$. One then performs a
measurement, described by a positive operator-valued measure (POVM)
$\bigl{\\{}\coprod_{x}\bigl{\\}}$ satisfying $\int dx\coprod_{x}=\openone$,
which results in a random outcome $x$ at the likelihood $p(x|T)={\rm
tr}\\{\rho_{p}(t,T)\coprod_{x}\\}$. The protocol is repeated $M$ times to
gather an appropriately large data sample
$\mathbf{x}=\\{x_{1},x_{2},...,x_{M}\\}$, from which one can construct a
temperature estimator $\tilde{T}$ yielding the estimate
$\tilde{T}(\mathbf{x})$. One speaks of an unbiased estimator when the
expectation value of the estimator over the likelihood of outcomes at a given
true value $T$ matches the true value, $\mathbb{E}[\tilde{T}(\mathbf{x})]=T$
[32, 28]. This may hold merely approximately or asymptotically, in the limit
of large data, as often the case with the commonly used maximum-likelihood
estimator, for example.

In equilibrium thermometry, the system is allowed to exchange energy with the
environment until it equilibrates to a thermal state after a sufficiently long
time, $t\to\infty$. One then measures in the energy basis to infer its thermal
populations and thereby estimate the temperature. In phase thermometry, a
probe system is prepared in a superposition state that is then subjected to a
thermal sample inducing dephasing over time, which one detects by performing
phase-sensitive interferometric measurements on the probe at finite times $t$.
We will focus on this case here, but for comparison, we will also consider
directly measuring the equilibrium state of a number of spins in the thermal
sample.

The statistical uncertainty of the temperature estimate can be defined with
help of a suitable _relative_ error quantifier. This allows us to compare the
achievable precision at different temperature scales. A suitable common choice
is the relative mean-square deviation from the true value,
$\mathbb{E}[(\tilde{T}-T)^{2}]/T^{2}$, which for unbiased estimators can be
obtained from the data as the variance divided by the squared mean,
$\mathrm{Var}[\tilde{T}]/(\mathbb{E}[\tilde{T}])^{2}$. This figure of merit
for the local precision of temperature estimates at a given true temperature
$T$ obeys the (quantum) Cramér-Rao bound [59, 60], The statistical uncertainty
of the temperature estimate can be defined with help of a suitable _relative_
error quantifier. This allows us to compare the achievable precision at
different temperature scales. A suitable common choice is the relative mean-
square deviation from the true value, $\mathbb{E}[(\tilde{T}-T)^{2}]/T^{2}$,
which for unbiased estimators can be obtained from the data as the variance
divided by the squared mean,
$\mathrm{Var}[\tilde{T}]/(\mathbb{E}[\tilde{T}])^{2}$. This figure of merit
for the local precision of temperature estimates at a given true temperature
$T$ obeys the (quantum) Cramér-Rao bound [59, 60],

$\frac{\mathbb{E}\big{[}\tilde{T}(\mathbf{x})-T)^{2}\big{]}}{T^{2}}\geq\frac{1}{MT^{2}\mathcal{F}_{T}}\geq\frac{1}{MT^{2}\mathcal{F}^{Q}_{T}}.$
(1)

Here, the first inequality is specific to the chosen quantum measurement and
given in terms of the Fisher information (FI) of the associated likelihood
with respect to temperature,

$\mathcal{F}_{T}=\int dx\,p(x|T)\left[\frac{\partial\ln p(x|T)}{\partial
T}\right]^{2}.$ (2)

It quantifies the sensitivity of the likelihood to small temperature changes
around $T$.

The more fundamental second inequality is obtained by optimizing the FI over
all possible POVMs,
$\mathcal{F}^{Q}_{T}=\max_{\\{\coprod_{x}\\}}\mathcal{F}_{T}$. The so defined
quantum Fisher information (QFI) is a function of the state $\rho_{p}(t,T)$
and sets the ultimate precision bound, i.e., the maximum local temperature
sensitivity, attainable for any chosen measurement on that state.

The QFI can be given explicitly in terms of the symmetric logarithmic
derivative (SLD) operator $\hat{L}_{T}$, which is the self-adjoint solution to
the Lyapunov equation

$\frac{\partial\rho_{p}(t,T)}{\partial
T}=\frac{1}{2}\left[\rho_{p}(t,T)\hat{L}_{T}+\hat{L}_{T}\rho_{p}(t,T)\right].$
(3)

The solution depends on the probe state, and its eigenbasis gives the
corresponding optimal projective measurement [61]. The QFI then reads as

$\mathcal{F}^{Q}_{T}={\rm
tr}\left[\rho_{p}(t,T)\hat{L}_{T}^{2}\right]=2\sum_{nm}\frac{|\langle\psi_{m}|\partial_{T}\rho_{p}(t,T)|\psi_{n}\rangle|^{2}}{\rho_{n}+\rho_{m}}.$
(4)

Here, the second expression follows by expanding the probe state in its
eigenbasis, $\rho_{p}(t,T)=\sum_{n}\rho_{n}|\psi_{n}\rangle\langle\psi_{n}|$,
which leads to

$\hat{L}_{T}=2\sum_{nm}\frac{\langle\psi_{m}|\partial_{T}\rho_{p}(t,T)|\psi_{n}\rangle}{\rho_{n}+\rho_{m}}|\psi_{m}\rangle\langle\psi_{n}|.$
(5)

In the following, we will conveniently express the temperature in terms of the
inverse temperature scale parameter $\beta=1/k_{B}T$. This reparametrization
leaves the relative Cramér-Rao bounds in (1) invariant,
$T^{2}\mathcal{F}_{T}=\beta^{2}\mathcal{F}_{\beta}$ and thus
$T^{2}\mathcal{F}^{Q}_{T}=\beta^{2}\mathcal{F}^{Q}_{\beta}$, as one can easily
check using (2). These dimensionless bounds will be used to characterise the
performance of our phase thermometer protocol.

##  III Zero-field Ising model with a single phase probe

Figure 1: Two dimensional Ising lattice model. The probe spin is depicted on
top of the lattice layer and interacts with system spins inside of a cluster
$C$ (dark gray circle) with uniform coupling strength $g$. The system spins
interact with each other at coupling strength $J$.

We consider a finite two-dimensional Ising spin lattice with nearest-neighbour
interactions as the thermal sample whose temperature we seek to estimate. It
is described by the Hamiltonian

$\hat{H}_{s}=-J\sum_{\langle
ij\rangle}\hat{\sigma}_{z}^{i}\otimes\hat{\sigma}_{z}^{j}-h\sum_{i=1}^{N}\hat{\sigma}_{z}^{i},$
(6)

where $N$ is the total spin number, $h$ the single-site energy due to an
external field, $J$ the nearest-neighbour coupling energy, and $\langle
ij\rangle$ denotes any possible pair of nearest neighbours assuming periodic
boundary conditions. The dimensionless spin operator $\hat{\sigma}_{z}^{i}$
describes the Z-component of a spin-$1/2$ on the $i$-th lattice site with
eigenvalues $m_{i}=\pm 1$. We will focus on the zero-field case $h=0$ and on
$J>0$ in the following, which exhibits a transition between a ferromagnetic
and a paramagnetic phase. We describe the sample in thermal equilibrium at
temperature $T$ by the Gibbs state
$\rho_{s}=e^{-\beta\hat{H}_{s}}/\mathcal{Z}$, with the partition function
$\mathcal{Z}={\rm tr}(e^{-\beta\hat{H}_{s}})$.

Even though the above Ising spin lattice Hamiltonian is a highly simplified
model for ferromagnetism, it can be representative of a large class of systems
(including realistic ones) and admit a quantitative understanding of these
systems near their critical point [45]. This is because, by virtue of
universality, the behavior of a system in the vicinity of a phase transition
becomes essentially independent of its microscopic details.

As our non-invasive phase probe, we consider a single spin-$1/2$ (or two-level
system) that interacts with some lattice spins in its vicinity through the
Ising Z-Z coupling; see Fig. 1 for a sketch. For simplicity, we restrict the
probe-sample interaction to a finite cluster of $n$ lattice spins within a
certain radius and assume a uniform coupling rate $g$. This will allow us to
compare the temperature information acquired by the probe with the information
one could obtain if one were to directly measure the local configuration of
the $n$ lattice spins. The joint probe-sample Hamiltonian is

$\displaystyle\hat{H}=\hat{H}_{s}+\frac{\hbar\omega_{p}}{2}\hat{\sigma}^{p}_{z}+\frac{\hbar
g}{2}\hat{\sigma}^{p}_{z}\otimes\hat{Z}_{n},$ (7)

with $\hbar\omega_{p}$ the probe energy, $\hat{\sigma}_{z}^{p}$ the Pauli-Z
matrix of the probe (with eigenvalues $\pm 1$), and
$\hat{Z}_{n}=\sum_{i=1}^{n}\hat{\sigma}^{i}_{z}$ the total spin of the
$n$-site cluster. In the extreme case $n=N$, the probe senses the collective
fluctuating magnetization of the whole sample.

The probe-sample coupling commutes with the probe Hamiltonian, which implies
that the diagonal terms of the probe density operator in the energy basis
remain constant. A generic initial probe state that is uncorrelated to the
Gibbs state $\rho_{s}$ of the sample,

$\rho_{p}=\begin{pmatrix}p&c\\\ c^{*}&1-p\end{pmatrix},$ (8)

will undergo pure dephasing due to averaging over the random phase shifts
caused by the thermally occupied spin configurations of the sample. The
reduced probe state evolves as

$\displaystyle\rho_{p}(t)$ $\displaystyle=$ $\displaystyle{\rm
tr}_{s}\left[e^{-i\hat{H}t/\hbar}\rho_{p}\otimes\rho_{s}e^{i\hat{H}t/\hbar}\right]$
(9) $\displaystyle=$
$\displaystyle\begin{pmatrix}p&ce^{-i\omega_{p}t}r(t,\beta)\\\
c^{*}e^{i\omega_{p}t}r^{*}(t,\beta)&1-p\end{pmatrix},$

with the decoherence factor

$r(t,\beta)={\rm tr}[\rho_{s}e^{-ig\hat{Z}_{n}t}]=\frac{1}{\mathcal{Z}}{\rm
tr}[e^{-\beta\hat{H}_{s}}e^{-ig\hat{Z}_{n}t}].$ (10)

This factor has the form of a characteristic function that generates the
moments of the $n$-site cluster’s magnetization $\hat{Z}_{n}$, via
$\langle\hat{Z}_{n}^{m}\rangle=(i/g)^{m}\partial_{t}^{m}r(0,\beta)$. Hence,
one could in principle infer the full spin distribution of the cluster by
sampling the decoherence factor at many different times.

Here, we will be concerned with the maximum information about the sample
temperature the probe can acquire. Starting from a maximally phase-coherent
initial state, $\rho_{p}=|+\rangle\langle+|$ with $p,c=1/2$, the QFI of the
probe state at an interrogation time $t$ with respect to $\beta$ becomes

$\displaystyle\mathcal{F}_{\beta}^{Q}(t)$ $\displaystyle=$
$\displaystyle\frac{4|\partial_{\beta}r|^{2}+\big{[}r^{*}\partial_{\beta}r-r\partial_{\beta}r^{*}\big{]}^{2}}{4(1-|r|^{2})}.$
(11)

For the results presented in the following, we evaluate the maximum of this
QFI with respect to $t$, showing the highest possible phase thermometry
precision in our setting.

Crucially, our treatment is exact at this point; it encompasses effects such
as coherence revivals, which go beyond Markovian dephasing models in the weak-
coupling limit [62]. For an efficient numerical modelling, we make use of a
binary array-based representation of spin networks [58]: Let the binary array
$\underline{b}\in\\{0,1\\}^{\times N}$ and the associated state vector
$|\underline{b}\rangle$ represent a specific $N$-spin configuration of the
sample lattice, such that
$\hat{\sigma}_{z}^{i}|\underline{b}\rangle=(1-2b_{i})|\underline{b}\rangle$
for all $i$. Short- or long-range spin-spin interactions can then be described
in terms of adjacency matrices. In particular, let
$\smash[b]{\underline{\underline{\Gamma}}}$ be the $N\times N$ adjacency
matrix of the lattice, taking the value $\Gamma_{ij}=-J$ for nearest-neighbour
pairs $\langle ij\rangle$ and zero otherwise. The energy of a given spin
configuration is then
$E(\underline{b})=(2h+8J)\sum_{i=1}^{N}b_{i}-2\underline{b}\cdot\smash[b]{\underline{\underline{\Gamma}}}\underline{b}$,
and its thermal occupation probability is $p(\underline{b})=e^{-\beta
E(\underline{b})}/\mathcal{Z}$.

Now, let us split the binary array,
$\underline{b}=(\underline{b}^{\prime},\underline{b}^{\prime\prime})$, into
the configuration $\underline{b}^{\prime}$ of the $n$-spin cluster interacting
with the probe and the configuration $\underline{b}^{\prime\prime}$ of the
remaining $N-n$ spins. Given the $1\times n$ adjacency matrix
$\underline{\Gamma}^{\prime}$ of the probe-sample interaction, which in our
uniform coupling model amounts to simply
$\underline{\Gamma}^{\prime}=g(1,1,\ldots,1)$, the decoherence factor reads as

$\begin{split}&r(t,\beta)=\sum_{\underline{b}^{\prime}}\left(\sum_{\underline{b}^{\prime\prime}}p(\underline{b}^{\prime},\underline{b}^{\prime\prime})\right)e^{-it\underline{\Gamma}^{\prime}\cdot\underline{b}^{\prime}}.\end{split}$
(12)

##  IV Temperature sensitivity near the critical point

The theoretical model at hand, we will now assess the maximum sensitivity of a
phase probe to the temperature of an Ising lattice sample of $N$ spins at
thermal equilibrium. Different approximate treatments will be compared to an
exact simulation. Our aim is to scrutinize and quantify our expectation that
the temperature sensitivity will reach a pronounced peak, and even diverge for
$N\to\infty$, at the critical temperature of the ferromagnetic-paramagnetic
phase transition.

We start with the extreme case of probing the whole sample magnetization
($n=N$) and introduce the instructive Curie-Weiss model, which serves as a
form of mean-field approximation of the Ising lattice. The model can also
represent a scenario of its own: a fully and uniformly connected spin network,
see Fig. 2.

Subsequently, we will evaluate the exact results for the Ising lattice and
varying cluster size $n$ and compare them to the sensitivity predicted by the
conventional mean-field approximation. Moreover, we will consider a high-
temperature expansion to model the thermometric precision in the paramagnetic
regime.

###  IV.1 Curie–Weiss model

Figure 2: In the Curie–Weiss model, every spin interacts with all the others
uniformly, as can be represented by a complete graph. In this graph, each spin
is represented by a vertex. An edge connecting a pair of vertices $i$ and $j$
represents the coupling contribution
$-(J/N)\hat{\sigma}_{z}^{i}\otimes\hat{\sigma}_{z}^{j}$ to the total energy,
where $J$ represents the interaction strength and $N$ denotes the total number
of spins (here: $12$).

The Curie-Weiss model is based on the idea that, in a ferromagnetic lattice,
all spins align in response to a collective magnetization or to an external
magnetic field [45]. Rather than nearest-neighbour or more general site-
specific couplings, each spin $\hat{\sigma}_{z}^{i}$ in the lattice shall
react only to the collective spin of the whole lattice,
$\sum_{j}\hat{\sigma}_{z}^{j}$. For consistency with the microscopic two-spin
coupling model of strength $J>0$, we divide by $N$ and thus couple each spin
to a lattice-averaged macroscopic spin density, described by the Hamiltonian

$\hat{H}_{s}^{\rm
CW}=-\frac{J}{2N}\sum_{i,j=1}^{N}\hat{\sigma}_{z}^{i}\otimes\hat{\sigma}_{z}^{j}-h\sum_{i=1}^{N}\hat{\sigma}_{z}^{i}.$
(13)

The factor $1/2$ prevents double counting of spin pairs, and the diagonal
terms ($i=j$) in the double sum merely give an irrelevant energy offset. One
can envision this model as a spin network defined by a complete graph with $N$
vertices, wherein every pair of distinct vertices is connected by an edge, as
illustrated in Fig. 2. We will restrict to the zero-field case ($h=0$) later,
but all the following derivations also apply for $h\neq 0$.

At thermal equilibrium, the partition function of the model is

$\displaystyle\mathcal{Z}$ $\displaystyle=$
$\displaystyle\sum_{m_{1}=-1}^{1}\ldots\sum_{m_{N}=-1}^{1}e^{(J/2N)\beta\left(\sum_{i=1}^{N}m_{i}\right)^{2}+h\beta\sum_{i=1}^{N}m_{i}}$
(14) $\displaystyle=$
$\displaystyle\sum_{M=-N}^{N}\binom{N}{\frac{N+M}{2}}e^{(J/2N)\beta
M^{2}+h\beta M},$

where $M=\sum_{i=1}^{N}m_{i}$ is the total dimensionless spin or
magnetization. For a large lattice, $N\to\infty$, we can approximate the
discrete sum over $M$ by an integral over a continuous magnetization density
$m=M/N\in[-1,1]$ per number of spins,

$\mathcal{Z}\approx\int_{-1}^{1}e^{-Nf(m)}dm,$ (15)

where the free energy density is approximately given as

$f(m)\approx\frac{1+m}{2}\ln\frac{1+m}{2}+\frac{1-m}{2}\ln\frac{1-m}{2}-\frac{\beta
Jm^{2}}{2}-\beta hm,$ (16)

by virtue of the Stirling formula. The integral in Eq. (15) is dominated by
the vicinity of the saddle point $m_{0}$ of the function $f$ where
$f^{\prime}(m_{0})=0$ [45], which satisfies

$\tanh{\left(J\beta m_{0}+h\beta\right)}=m_{0}.$ (17)

We can thus perform the saddle point approximation,

$\mathcal{Z}\approx
e^{-Nf(m_{0})}\sqrt{\frac{2\pi}{Nf^{\prime\prime}(m_{0})}},$ (18)

with $f^{\prime\prime}(m_{0})=1/(1-m_{0}^{2})-J\beta$.

The saddle point solution $m_{0}$ of (17) also constitutes the mean
magnetization,
$\langle\sum_{i}\hat{\sigma}_{z}^{i}\rangle=\beta^{-1}\partial_{h}\ln\mathcal{Z}\approx
Nm_{0}$, and thus describes the spontaneous symmetry breaking at the
ferromagnetic-paramagnetic phase transition. For $h=0$, spontaneous finite
magnetization occurs at the critical point $J\beta_{c}=1$ [63, 45].

#### Phase thermometry

For temperature sensing, we consider a qubit probe interacting with all $N$
lattice spins with uniform coupling rate $g$. The decoherence factor (10) can
then be given consistently in the continuum and saddle-point approximation for
$N\to\infty$ as

$\displaystyle r(t,\beta)$ $\displaystyle=$
$\displaystyle\sum_{M=-N}^{N}\binom{N}{\frac{N+M}{2}}e^{(J/2N)\beta
M^{2}+(h\beta-igt)M}$ (19) $\displaystyle\approx$ $\displaystyle
e^{-i\tilde{g}tm_{0}-(\tilde{g}t)^{2}/[2Nf^{\prime\prime}(m_{0})]},$

with $\tilde{g}=gN$. The effect of the lattice is a frequency shift
proportional to the mean magnetization and, to second order in $\tilde{g}t$, a
dephasing due to the thermal spin fluctuations with a characteristic
($e^{-1/2}$) coherence decay time scale
$\tilde{g}\tau(\beta)=\sqrt{Nf^{\prime\prime}(m_{0})}$. The coherence decays
more rapidly the greater the lattice, as $\tau(\beta)\propto 1/\sqrt{N}$.

To illustrate the behaviour, consider the probe initialized in the
X-eigenstate $|+\rangle$ (which is the optimal initial state for phase
thermometry [21]) and an observation of the spin-X component as a function of
time. The expectation value of this coherence observable shows oscillations
and the so-called free induction decay (FID) [64],

$\mathrm{FID}(t)={\rm
tr}\\{\hat{\sigma}_{x}^{p}\rho_{p}(t)\\}=\cos(\tilde{g}tm_{0})e^{-t^{2}/2\tau^{2}(\beta)}.$
(20)

We plot this quantity at zero field and $N=400$ for various temperatures
around the critical point in Fig. 3(a). We observe a Gaussian decay for above-
critical temperatures in the paramagnetic phase ($J\beta<1$) where the mean
magnetization vanishes ($m_{0}=0$). Conversely, in the ferromagnetic phase, we
observe a damped oscillation with a Gaussian envelope. The decay time is
minimal at the critical point (black solid line) and increases further away
from the critical point in both directions, see also Fig. 3(b). In the limit
of zero temperature, we have no decay, $m_{0}\to 1$ and
$\tau(\beta\to\infty)\to\infty$, as expected. In the paramagnetic phase, we
always have $m_{0}=0$, and in the high-temperature limit,
$\tau(0)=1/g\sqrt{N}$. In the ferromagnetic phase close to the critical point,
we can Taylor-expand the $\tanh$ in (17) in small $m_{0}>0$ and
$\epsilon=1-J\beta=1-\beta/\beta_{c}<0$, resulting in
$m_{0}\approx\sqrt{-3\epsilon}$.

Figure 3: (a) Coherence dynamics of the spin probe interacting with the Curie-
Weiss lattice at different temperatures over time, measured in terms of the
free induction decay (20). (b) Coherence decay time versus inverse
temperature, as predicted by the Curie-Weiss model in the thermodynamic limit.
We set the parameters $N=400,h=0$.

#### Temperature precision at criticality

The high sensitivity of the probe coherence to the thermal fluctuations of the
sample around the critical point suggests that temperature could be probed
more precisely. This is quantified by the QFI $\mathcal{F}_{\beta}^{Q}$; the
greater it is at a given temperature, the higher the estimation precision and
the lower the uncertainty around that temperature. Inserting the decoherence
factor (19) into the QFI (11) and several steps of algebra result in

$\displaystyle\mathcal{F}_{\beta}^{Q}$ $\displaystyle=$
$\displaystyle\frac{J^{2}}{4}\frac{[(1-m_{0}^{2})^{2}f^{\prime\prime}(m_{0})-2m_{0}^{2}]^{2}}{(1-m_{0}^{2})^{6}[f^{\prime\prime}(m_{0})]^{4}}\frac{t^{4}/\tau^{4}(\beta)}{e^{t^{2}/\tau^{2}(\beta)}-1}$
(21)
$\displaystyle+NJ^{2}\frac{m_{0}^{2}}{f^{\prime\prime}(m_{0})}\frac{t^{2}}{\tau^{2}(\beta)}e^{-t^{2}/\tau^{2}(\beta)},$

In the paramagnetic phase ($\beta<\beta_{c}=1/J$), where $m_{0}=0$ and
$f^{\prime\prime}(0)=1-\beta/\beta_{c}=\epsilon>0$, the QFI simplifies to

$\mathcal{F}_{\beta<\beta_{c}}^{Q}=\frac{J^{2}}{4\epsilon^{2}}S\left(\frac{\tilde{g}^{2}t^{2}}{N\epsilon}\right),\quad
S(x)=\frac{x^{2}}{e^{x}-1}.$ (22)

Optimal sensitivity is achieved at the time $t$ that maximizes the function
$S$, which yields
$\max_{t}\mathcal{F}_{\beta\leq\beta_{c}}^{Q}=0.162J^{2}/\epsilon^{2}$ when
$\tilde{g}^{2}t^{2}=1.594N\epsilon$.

In the ferromagnetic phase near the critical temperature, the magnetization is
also small, $m_{0}\approx\sqrt{-3\epsilon}\ll 1$, and we have
$f^{\prime\prime}(m_{0})\approx-2\epsilon$. Hence, to leading order in the
temperature deviation $\epsilon<0$,

$\mathcal{F}_{\beta>\beta_{c}}^{Q}\approx\frac{J^{2}}{\epsilon^{2}}S\left(-\frac{\tilde{g}^{2}t^{2}}{2N\epsilon}\right).$
(23)

In total, this leads to the following critical exponents: the optimal QFI
diverges like $\epsilon^{-2}$ and the optimal time before measuring the probe
scales like $t\propto|\epsilon|^{1/2}$ for $\epsilon\to 0$ in the
thermodynamic limit. The values remain finite for finite $N$.

It is instructive to compare the temperature sensitivity of the phase spin
probe to directly probing the local magnetization on any single lattice spin.
This amounts to measuring $\hat{\sigma}_{z}^{i}$ at a random lattice site $i$.
In the thermodynamic limit, the probability to obtain either outcome $\pm 1$
is given by $p_{\pm}=(1\pm m_{0})/2$, and hence the FI of this measurement
with respect to $\beta$ is

$\mathcal{F}_{\beta}^{\rm(lc)}=\frac{(\partial_{\beta}m_{0})^{2}}{1-m_{0}^{2}}\xrightarrow[\epsilon<0]{\epsilon\to
0}-\frac{3J^{2}}{4\epsilon}.$ (24)

Above the critical point, both $m_{0}$ and the local FI vanish, whereas
closely below the critical point, the FI diverges like $|\epsilon|^{-1}$ –
more weakly than the QFI of the phase probe.

Figure 4: Optimal QFI of the phase probe (optimized over the interrogation
time) for the Curie-Weiss model in the thermodynamic limit, plotted as a
function of the relative temperature deviation $\epsilon$ from the critical
point. The right half ($\epsilon>0$) represent the paramagnetic regime of
$\beta<\beta_{c}$ and zero mean magnetization. The left half marks the
ferromagnetic regime of low temperatures, where the direct local measurement
of a single lattice spin yields non-vanishing temperature information. We
assume a probe-lattice coupling strength $\hbar g/J=0.4$ and the lattice
coupling $J=1/4$ in units of an arbitrary reference energy scale.

We illustrate the difference between the local FI and the optimal phase probe
QFI close to the phase transition in Fig. 4. In both cases, only a single spin
degree of freedom is measured, but the phase probe interacts with the whole
lattice and thus picks up correlated thermal fluctuations around the mean
magnetization $m_{0}$ through the dephasing factor (19), leading to a much
higher temperature sensitivity. In the following, we will see a similar,
albeit less pronounced, advantage of the phase probe in the Ising lattice.

###  IV.2 2D Ising model at zero field

We now turn to phase thermometry of the finite ferromagnetic Ising lattice at
zero field, as described by the Hamiltonian (6) with $h=0$ and $J>0$. To
obtain our results and sample the thermally occupied spin configurations, we
performed a Monte Carlo simulation on a $20\times 20$ lattice ($N=400$) with
periodic boundary conditions. The probability $p(\underline{b})$ of each spin
configuration is simulated by the relative frequency over $1.2\times 10^{7}$
instances. The decoherence factor of the phase probe (12) can then be used to
compute the QFI (11) at arbitrary times and its maximum
$\mathcal{F}_{\beta}^{\rm opt}=\max_{t}\mathcal{F}_{\beta}^{Q}(t)$ with
respect to time.

Figure 5: (a) Optimal QFI of the phase probe (optimized over time) interacting
with spin clusters of varying size $n$ in a $20\times 20$ Ising lattice, as a
function of inverse temperature relative to the critical point (vertical
dashed line). (b) FI associated to a direct local measurement of the
configuration of the $n$ cluster spins. (c) Ratio of optimal probe QFI and
direct local FI close to the critical point. We use the probe-cluster coupling
strength $\hbar g/J=0.4$ and set $J=1/4$ (in arbitrary reference units).

Figure 5(a) shows the simulation results for the maximum QFI in relative units
as a function of inverse temperature $\beta$, for various cluster sizes $n$
the probe interacts with. The inverse temperature is plotted relative to the
known critical value $\beta_{c}=\ln{(1+\sqrt{2})}/2J$ for a 2D Ising model in
the thermodynamic limit [65], which matches the QFI peaks in this $20\times
20$ simulation. The visible jump at $\beta\sim 1.2\beta_{c}$ is a numerical
instability detecting the spontaneous symmetry breaking in the ferromagnetic
phase. At this point, the zero-field Ising system may randomly assume one of
two opposite ferromagnetically ordered configurations that are degenerate in
energy [66]. The sudden change from a probabilistic mixture of both options to
one of them directly affects the complex phase of the decoherence factor (19),
resulting in the observed jump.

Although the phase thermometer probes the occupied configurations of an
$n$-spin cluster, the measurement outcome is still a single bit value. Let us
compare this to a direct measurement of the local cluster configurations
$\underline{b}^{\prime}$ in the Ising lattice, thermally occupied with the
marginal probability
$p_{n}(\underline{b}^{\prime})=\sum_{\underline{b}^{\prime\prime}}p(\underline{b}^{\prime},\underline{b}^{\prime\prime})$.
The temperature sensitivity of this $n$-bit measurement is given by the FI

$\displaystyle\mathcal{F}_{\beta}^{\rm(lc)}$ $\displaystyle=$
$\displaystyle\sum_{\underline{b}^{\prime}}\frac{1}{p_{n}(\underline{b}^{\prime})}\left[\frac{\partial
p_{n}(\underline{b}^{\prime})}{\partial\beta}\right]^{2}$ (25)
$\displaystyle=$
$\displaystyle\sum_{\underline{b}^{\prime}}\frac{1}{p_{n}(\underline{b}^{\prime})}\left[\sum_{\underline{b}^{\prime\prime}}E(\underline{b})p(\underline{b})\right]^{2}-\langle\hat{H}_{s}\rangle^{2},$

which we could sample from our Monte-Carlo simulation for moderate $n\lesssim
30$. The results are plotted in Fig. 5(b) for the same temperature range and
cluster sizes as before. We once again observe a peak at the phase transition,
but now a less steep decay into the paramagnetic phase. Naturally, the FI of
this $n$-bit measurement is greater than the QFI of the phase probe, but the
ratio of both quantities close to the critical point grows only sublinearly
with $n$, as can be seen in (c). Deep in the ferromagnetic phase, in
particular, the lattice spins are strongly aligned (i.e., correlated) such
that the phase probe picks up almost all the temperature information in the
cluster. Conversely, in the paramagnetic phase, the correlations are broken
and the disparity between both measurements grows.

Figure 6: Scaling of the optimal probe QFI with the size $n$ of the lattice
spin cluster that it interacts with, using the parameters of Fig. 5. We
compare the scaling at the critical point to the scaling at a $10\%$ higher
temperature ($\beta=0.9\beta_{c}$) and at a $10\%$ lower temperature. The
solid lines are linear fits on the double-logarithmic scale with slopes
$0.73$, $0.71$, and $0.49$, respectively.

The scaling of the optimal QFI with cluster size $n$ is plotted on a double-
logarithmic scale in Fig. 6 at the critical point, as well as $10\,\%$ above
and below. From the average slopes, we obtain a critical scaling with
$n^{0.73}$; weaker in the ferromagnetic phase and about as strong in the
paramagnetic phase. This sublinear scaling illustrates the strong correlations
among the cluster spins.

In the following, we investigate approximate analytical models for the
temperature sensitivity based on the widely used mean-field theory (MFT) and
on a high-temperature expansion (HTE) of the Gibbs state. We will see that
these models yield accurate predictions for temperatures deep in the
ferromagnetic and paramagnetic phase, respectively, but they fail to reproduce
the critical behaviour of our exact numerical simulation.

#### Mean-field approximation

In MFT, one describes a large ferromagnetic Ising lattice by a mean
homogeneous magnetization $m_{0}=\langle\hat{\sigma}_{z}^{i}\rangle$ $\forall
i$ and expands in the microscopic fluctuations around this mean [67].
Inserting the ansatz $\hat{\sigma}_{z}^{i}=m_{0}+\delta\hat{\sigma}_{z}^{i}$
into the Hamiltonian (6) at $h=0$ and neglecting second-order terms in the
fluctuation operators, one obtains

$\displaystyle\hat{H}^{(\mathrm{MFT})}_{s}$ $\displaystyle\approx$
$\displaystyle-J\sum_{\langle
i,j\rangle}m_{0}\left(\hat{\sigma}_{z}^{i}+\hat{\sigma}_{z}^{j}-m_{0}\right)$
(26) $\displaystyle=$
$\displaystyle\frac{NJqm_{0}^{2}}{2}-Jqm_{0}\sum_{i=1}^{N}\hat{\sigma}_{z}^{i},$

with $q$ the number of nearest neighbours per spin, here $q=4$. The
Hamiltonian describes independent spins aligning with respect to an effective
mean field $h_{\rm MFT}=Jqm_{0}$, which itself is determined by the mean spin
alignment. At thermal equilibrium, the partition function and the thermal
occupation probability are products of $N$ identical single-spin terms,

$\displaystyle\mathcal{Z}$ $\displaystyle=$ $\displaystyle e^{-\beta
NqJm_{0}^{2}/2}\left[2\cosh{\left(Jqm_{0}\beta\right)}\right]^{N},$ (27)
$\displaystyle p(\underline{b})$ $\displaystyle=$
$\displaystyle\frac{e^{-2Jqm_{0}\beta\sum_{i}b_{i}}}{(1+e^{-2Jqm_{0}\beta})^{N}}.$
(28)

The corresponding average alignment of each spin leads to the self-consistency
equation $m_{0}=\langle\hat{\sigma}_{z}^{i}\rangle=\tanh[Jqm_{0}\beta]$, which
overestimates the critical temperature here to $J\beta=1/q=0.25$ (as opposed
to $J\beta_{c}\approx 0.44$ for $N\gg 1$). The self-consistency condition has
the same form as in the Curie-Weiss model in Sec. IV.1 if we replace $Jq$ by
an average effective coupling strength, but the key difference is that we have
a collective spin-$N/2$ model there and $N$ individual spins here, reflected
by the distinct partition functions.

The effect on the phase probe over time follows by inserting (28) into the
decoherence factor expression (12),

$r(t,\beta)=\left[1-\frac{1-e^{-igt}}{1+e^{2Jqm_{0}\beta}}\right]^{n}\approx
e^{\left(1-\cos gt-i\sin gt\right)np_{1}},$ (29)

with the single-spin excitation probability $p_{1}=1/(1+e^{2Jqm_{0}\beta})$.
The approximation holds for moderate cluster sizes $n$ in the ferromagnetic
phase. From this, we arrive at an analytic expression for the QFI,

$\displaystyle\mathcal{F}_{\beta}^{(\mathrm{MFT})}$ $\displaystyle=$
$\displaystyle\frac{J^{2}q^{2}n^{2}\left(m_{0}+\beta\partial_{\beta}m_{0}\right)^{2}}{4\cosh^{4}(Jqm_{0}\beta)}e^{-2\gamma(1-p_{1})}$
(30) $\displaystyle\times\left[\sin^{2}(gt)+\frac{(1-\cos
gt)^{2})}{1-e^{-2\gamma(1-p_{1})}}\right].$

As the derivative $\partial_{\beta}m_{0}$ diverges at the (overestimated)
critical temperature, so does the approximated QFI. Above that temperature,
both the mean field and the approximated QFI vanish.

Figure 7: Optimal QFI of the phase probe as a function of inverse temperature
relative to the critical point (vertical dashed line), extending the results
in Fig. 5(a). We separate the four cases corresponding to different cluster
sizes $n$ into different diagrams (a)-(d), comparing the exact Monte-Carlo
simulation data to the approximate predictions of the mean-field theory (solid
lines) and the high-temperature expansion (dashed).

In Fig. 7(a)-(d), we compare the MFT-approximation for the QFI (black solid
line) with the exact Monte-Carlo simulation results for four different cluster
sizes. We find that (30) predicts the temperature sensitivity well for
subcritical temperatures, but fails at and above the phase transition. This is
expected since above the critical point the magnetization vanishes. To study
the QFI above the critical point, we use a high-temperature expansion [63].

#### High-temperature expansion

At high temperatures in the paramagnetic phase, thermal fluctuations are
dominant and the lattice spins tend to be randomly oriented. The properties of
the sample can be studied by performing an expansion of the relevant physical
quantities in terms of the inverse temperature $\beta\ll 1/J$. Specifically,
we must expand up to second order in $\tanh(J\beta)$ to arrive at a nontrivial
result. The partition function expands as

$\displaystyle\mathcal{Z}$ $\displaystyle=$ $\displaystyle{\rm
tr}\left\\{e^{\beta\sum_{\langle
ij\rangle}J\hat{\sigma}_{z}^{i}\otimes\hat{\sigma}_{z}^{j}}\right\\}$
$\displaystyle=$ $\displaystyle{\rm tr}\left\\{\prod_{\langle
ij\rangle}\Big{[}\cosh{(\beta J)}+\sinh{(\beta
J)}\hat{\sigma}_{z}^{i}\otimes\hat{\sigma}_{z}^{j}\Big{]}\right\\}$
$\displaystyle\approx$ $\displaystyle{\rm tr}\left\\{\cosh^{K}(\beta
J)\Big{[}1+\tanh(\beta J)\sum_{\langle
ij\rangle}\hat{\sigma}_{z}^{i}\otimes\hat{\sigma}_{z}^{j}\right.$
$\displaystyle+\left.\tanh^{2}(\beta J)\sum_{\langle ij\rangle<\langle
kl\rangle}(\hat{\sigma}_{z}^{i}\otimes\hat{\sigma}_{z}^{j})(\hat{\sigma}_{z}^{k}\otimes\hat{\sigma}_{z}^{l})\Big{]}\right\\}$
$\displaystyle=$ $\displaystyle 2^{N}\cosh^{K}(\beta J),$

where the $K$ is the total number of pair bonds (i.e., interacting spin pairs)
in the lattice. The double sum over pair bonds in the fourth line counts each
pair only once. To obtain the final form, notice that the first- and second-
order terms in $\tanh(\beta J)$ in the third and fourth line vanish when
taking the trace, since there are always at least two distinct spins
contributing a factor ${\rm tr}\\{\hat{\sigma}_{z}^{i}\\}=0$. Here, the first
non-vanishing contribution would come from the fourth-order summation over
products of four pair bonds forming a closed loop between four adjacent spins
[63]: $\langle ij\rangle,\langle jk\rangle,\langle k\ell\rangle,\langle\ell
i\rangle$.

The phase probe’s decoherence factor (10) expands as

$\displaystyle r(t,\beta)=\frac{1}{\mathcal{Z}}{\rm
tr}\left\\{e^{\beta\sum_{\langle
ij\rangle}J\hat{\sigma}_{z}^{i}\otimes\hat{\sigma}_{z}^{j}}e^{-igt\hat{Z}_{n}}\right\\}$
(32) $\displaystyle\approx\frac{1}{\mathcal{Z}}\cosh^{K}(\beta J)\Bigg{[}{\rm
tr}\left\\{e^{-igt\hat{Z}_{n}}\right\\}$ $\displaystyle\quad+\tanh(\beta
J)\sum_{\langle ij\rangle}{\rm
tr}\left\\{(\hat{\sigma}_{z}^{i}\otimes\hat{\sigma}_{z}^{j})e^{-igt\hat{Z}_{n}}\right\\}$
$\displaystyle\quad+\tanh^{2}(\beta J)\\!\\!\\!\\!\sum_{\langle
ij\rangle<\langle kl\rangle}\\!\\!\\!\\!{\rm
tr}\left\\{(\hat{\sigma}_{z}^{i}\otimes\hat{\sigma}_{z}^{j})(\hat{\sigma}_{z}^{k}\otimes\hat{\sigma}_{z}^{l})e^{-igt\hat{Z}_{n}}\right\\}\Bigg{]}.$

Here, the first- and second-order terms in $\tanh(\beta J)$ no longer vanish
due to the additional unitary generated by the cluster magnetization
$\hat{Z}_{n}$. Writing

$e^{-igt\hat{Z}_{n}}=\cos^{n}(gt)\prod_{m=1}^{n}\left[1-i\tan(gt)\hat{\sigma}_{z}^{m}\right],$
(33)

we see that the first-order term in the third line of (32) adds a term
$\propto\tan^{2}(gt)$ for every pair bond $\langle ij\rangle$ whose both spins
are contained in the cluster. Let $K_{12}$ denote the number of such intra-
cluster pair bonds. For the second-order term in the fourth line of (32), non-
vanishing contributions can come from doublets of pair bonds with at least two
spins inside the cluster, as illustrated in Fig. 8: $K_{22}$ counts the number
of pair bond doublets with two spins inside and one mutual spin outside the
cluster, $K_{23}$ counts doublets of intra-cluster pair bonds with one mutual
spin, and $K_{24}=K_{12}(K_{12}-1)/2-K_{23}$ counts doublets of intra-cluster
pair bonds with no mutual spin. We arrive at

$\displaystyle r(t,\beta)=\cos^{n}(gt)\bigg{\\{}1-\tanh(\beta J)\tan^{2}(gt)$
$\displaystyle\quad\times\left[K_{12}-\tanh(\beta
J)\left(K_{22}+K_{23}-\tan^{2}(gt)K_{24}\right)\right]\bigg{\\}}.$

Figure 8: Construction of pair bond doublets in the Ising lattice. The cluster
$C$ of spins interacting with probe (dark blue) are contained in a circle.
Exemplary red lines symbolize pair bonds due to the Ising interaction between
any two neighbouring lattice spins. From left to right, each marked doublet of
pair bonds represents one possible contribution that is counted by the numbers
$K_{22}$, $K_{23}$, $K_{24}$, respectively. In the depicted configuration
there are $n=13$ cluster spins, $K_{12}=16$ intra-cluster bonds, and one
obtains $(K_{22},K_{23},K_{24})=(8,34,86)$.

The QFI of the phase probe with respect to the inverse temperature can once
again be given analytically. For a compact expression in the high-temperature
regime, we now approximate $\tanh(J\beta)\approx J\beta$ and arrive at

$\displaystyle\mathcal{F}_{\beta}^{(\mathrm{HTE})}$ $\displaystyle=$
$\displaystyle\frac{J^{2}\Lambda^{2}(t)\tan^{4}(gt)\cos^{2n}(gt)}{1-\left[J\beta\Lambda(t)\tan^{2}(gt)-1\right]^{2}\cos^{2n}(gt)},$
(34) $\displaystyle\Lambda(t)$ $\displaystyle=$ $\displaystyle
K_{12}+2J\beta\left[K_{22}+K_{23}-K_{24}\tan^{2}(gt)\right].$

The phase probe thus picks up temperature information through the correlations
between at most four lattice spins. In Fig. 7, the dashed lines depict the HTE
(34) of the QFI, which matches the exact results for $\beta\ll\beta_{c}$ in
the paramagnetic phase.

We conclude that the actual critical behaviour and the peak value of the phase
probe’s temperature sensitivity cannot be accurately described by either
analytical approximation to the QFI. The MFT fails due to the buildup of
strong long-range thermal fluctuations, while the HTE fails due to the buildup
of long-range multi-spin correlations in the critical Ising lattice. Moreover,
the finite peak value of the QFI cannot be captured by continuous models that
describe the phase transition in the thermodynamic limit, necessitating
numerical simulations.

##  V Conclusion

In this work, we considered non-invasive quantum thermometry of finite,
strongly interacting Ising spin systems at temperatures around their critical
phase transition point. We demonstrated that the non-Markovian dephasing of a
quantum probe in the vicinity of such a critical system picks up thermal spin
fluctuations and thus reaches peak temperature sensitivity at the critical
point. Consequently, one can estimate the system temperature by monitoring the
probe’s phase coherence over time, and one achieves the highest estimation
precision (i.e., minimum relative error, as measured by the QFI) at the
critical temperature—applying critical quantum metrology [35, 36, 37, 38] for
temperature estimation.

To exemplify the critical sensitivity enhancement in an instructive analytical
manner, we first performed a case study of phase thermometry on the Curie-
Weiss model of a fully connected spin network in the thermodynamic limit. The
QFI of a single spin probe exhibits a singularity at the critical temperature,
and the probe dynamics changes from an oscillatory to a non-oscillatory decay
of phase coherence as the system transitions from the ferromagnetic to the
paramagnetic phase. The achievable sensitivity clearly surpasses that of a
local single-spin measurement of the lattice.

For our main results, we performed Monte-Carlo simulations on a $20\times 20$
Ising lattice with ferromagnetic nearest-neighbour coupling, showing a finite
critical sensitivity peak of the probe that depends on the number of lattice
spins it interacts with. In the critical range of temperatures, analytic
approximations such as mean-field theory for the ferromagnetic phase and the
high-temperature expansion for the paramagnetic phase do not accurately
capture the thermal spin fluctuations leading to the strongest response of the
probe. Accordingly, phase thermometry constitutes a suitable, non-invasive
metrology scheme to access the theoretically elusive critical regime of strong
fluctuations in interacting many-body systems. Future research could explore
phase probes for temperature estimation in finite samples of Heisenberg-type
spin systems, for example. Bayesian phase thermometry schemes with optimal,
ancilla-augmented probe states could also be considered [68]. Moreover, the
spreading of Fisher information from a target location across the spin lattice
could be considered for remote sensing and parameter estimation tasks.

###### Acknowledgements.

The University of Siegen is kindly acknowledged for enabling our computation
through the OMNI cluster. This work was supported by the Deutsche
Forschungsgemeinschaft (DFG, German Research Foundation, project numbers
447948357 and 440958198), the Sino-German Center for Research Promotion
(Project M-0294), the ERC (Consolidator Grant 683107/TempoQ), and the German
Ministry of Education and Research (Project QuKuK, BMBF Grant No. 16KIS1618K).

## References

  * Degen  _et al._ [2017] C. L. Degen, F. Reinhard, and P. Cappellaro, Quantum sensing, Rev. Mod. Phys. 89, 035002 (2017).
  * Nielsen and Chuang [2002] M. A. Nielsen and I. Chuang, Quantum computation and quantum information, American Journal of Physics 70, 558 (2002), https://doi.org/10.1119/1.1463744 .
  * Ladd  _et al._ [2010] T. D. Ladd, F. Jelezko, R. Laflamme, Y. Nakamura, C. Monroe, and J. L. O’Brien, Quantum computers, Nature 464, 45 (2010).
  * Wu and Widom [1998] J. Wu and A. Widom, Physical limits on the notion of very low temperatures, Phys. Rev. E 57, 5178 (1998).
  * Mehboudi  _et al._ [2015] M. Mehboudi, M. Moreno-Cardoner, G. D. Chiara, and A. Sanpera, Thermometry precision in strongly correlated ultracold lattice gases, New Journal of Physics 17, 055020 (2015).
  * De Pasquale  _et al._ [2016] A. De Pasquale, D. Rossini, R. Fazio, and V. Giovannetti, Local quantum thermal susceptibility, Nature Communications 7, 12782 (2016).
  * Correa  _et al._ [2017] L. A. Correa, M. Perarnau-Llobet, K. V. Hovhannisyan, S. Hernández-Santana, M. Mehboudi, and A. Sanpera, Enhancement of low-temperature thermometry by strong coupling, Phys. Rev. A 96, 062103 (2017).
  * Hovhannisyan and Correa [2018] K. V. Hovhannisyan and L. A. Correa, Measuring the temperature of cold many-body quantum systems, Phys. Rev. B 98, 045101 (2018).
  * Mehboudi  _et al._ [2019] M. Mehboudi, A. Sanpera, and L. A. Correa, Thermometry in the quantum regime: recent theoretical progress, Journal of Physics A: Mathematical and Theoretical 52, 303001 (2019).
  * Stace [2010] T. M. Stace, Quantum limits of thermometry, Phys. Rev. A 82, 011611 (2010).
  * Correa  _et al._ [2015] L. A. Correa, M. Mehboudi, G. Adesso, and A. Sanpera, Individual quantum probes for optimal thermometry, Phys. Rev. Lett. 114, 220405 (2015).
  * Mok  _et al._ [2021] W.-K. Mok, K. Bharti, L.-C. Kwek, and A. Bayat, Optimal probes for global quantum thermometry, Commun. Phys. 4, 1 (2021).
  * Mihailescu  _et al._ [2023a] G. Mihailescu, S. Campbell, and A. K. Mitchell, Thermometry of strongly correlated fermionic quantum systems using impurity probes, Phys. Rev. A 107, 042614 (2023a).
  * Abiuso  _et al._ [2023] P. Abiuso, P. A. Erdman, M. Ronen, F. Noé, G. Haack, and M. Perarnau-Llobet, Optimal thermometers with spin networks (2023), arXiv:2211.01934 [quant-ph] .
  * Hangleiter  _et al._ [2015] D. Hangleiter, M. T. Mitchison, T. H. Johnson, M. Bruderer, M. B. Plenio, and D. Jaksch, Nondestructive selective probing of phononic excitations in a cold bose gas using impurities, Phys. Rev. A 91, 013611 (2015).
  * Jevtic  _et al._ [2015] S. Jevtic, D. Newman, T. Rudolph, and T. M. Stace, Single-qubit thermometry, Phys. Rev. A 91, 012331 (2015).
  * Mancino  _et al._ [2017] L. Mancino, M. Sbroscia, I. Gianani, E. Roccia, and M. Barbieri, Quantum simulation of single-qubit thermometry using linear optics, Phys. Rev. Lett. 118, 130502 (2017).
  * Hofer  _et al._ [2017] P. P. Hofer, J. B. Brask, M. Perarnau-Llobet, and N. Brunner, Quantum thermal machine as a thermometer, Phys. Rev. Lett. 119, 090603 (2017).
  * Cavina  _et al._ [2018] V. Cavina, L. Mancino, A. De Pasquale, I. Gianani, M. Sbroscia, R. I. Booth, E. Roccia, R. Raimondi, V. Giovannetti, and M. Barbieri, Bridging thermodynamics and metrology in nonequilibrium quantum thermometry, Phys. Rev. A 98, 050101 (2018).
  * Bouton  _et al._ [2020] Q. Bouton, J. Nettersheim, D. Adam, F. Schmidt, D. Mayer, T. Lausch, E. Tiemann, and A. Widera, Single-atom quantum probes for ultracold gases boosted by nonequilibrium spin dynamics, Phys. Rev. X 10, 011018 (2020).
  * Mitchison  _et al._ [2020] M. T. Mitchison, T. Fogarty, G. Guarnieri, S. Campbell, T. Busch, and J. Goold, In situ thermometry of a cold fermi gas via dephasing impurities, Phys. Rev. Lett. 125, 080402 (2020).
  * Oghittu and Negretti [2022] L. Oghittu and A. Negretti, Quantum-limited thermometry of a fermi gas with a charged spin particle, Phys. Rev. Res. 4, 023069 (2022).
  * Mitchison  _et al._ [2022] M. T. Mitchison, A. Purkayastha, M. Brenes, A. Silva, and J. Goold, Taking the temperature of a pure quantum state, Phys. Rev. A 105, L030201 (2022).
  * Yuan  _et al._ [2023] J.-B. Yuan, B. Zhang, Y.-J. Song, S.-Q. Tang, X.-W. Wang, and L.-M. Kuang, Quantum sensing of temperature close to absolute zero in a bose-einstein condensate, Phys. Rev. A 107, 063317 (2023).
  * Boeyens  _et al._ [2023] J. Boeyens, B. Annby-Andersson, P. Bakhshinezhad, G. Haack, M. Perarnau-Llobet, S. Nimmrichter, P. P. Potts, and M. Mehboudi, Probe thermometry with continuous measurements (2023), arXiv:2307.13407 .
  * Zhang  _et al._ [2022] N. Zhang, C. Chen, S.-Y. Bai, W. Wu, and J.-H. An, Non-markovian quantum thermometry, Phys. Rev. Appl. 17, 034073 (2022).
  * Zhang and Tong [2022] D.-J. Zhang and D. Tong, Approaching heisenberg-scalable thermometry with built-in robustness against noise, npj Quantum Information 8, 81 (2022).
  * Rubio  _et al._ [2021] J. Rubio, J. Anders, and L. A. Correa, Global quantum thermometry, Phys. Rev. Lett. 127, 190402 (2021).
  * Boeyens  _et al._ [2021] J. Boeyens, S. Seah, and S. Nimmrichter, Uninformed bayesian quantum thermometry, Phys. Rev. A 104, 052214 (2021).
  * Glatthard  _et al._ [2022a] J. Glatthard, J. Rubio, R. Sawant, T. Hewitt, G. Barontini, and L. A. Correa, Optimal cold atom thermometry using adaptive bayesian strategies, PRX Quantum 3, 040330 (2022a).
  * Rubio [2022] J. Rubio, Quantum scale estimation, Quantum Science and Technology 8, 015009 (2022).
  * Fraser [1964] D. A. S. Fraser, On local unbiased estimation, Journal of the Royal Statistical Society: Series B (Methodological) 26, 46 (1964).
  * Goold  _et al._ [2018] J. Goold, F. Plastina, A. Gambassi, and A. Silva,  _Thermodynamics in the Quantum Regime: Fundamental aspects and new directions_ (Springer, 2018).
  * Potts  _et al._ [2019] P. P. Potts, J. B. Brask, and N. Brunner, Fundamental limits on low-temperature quantum thermometry with finite resolution, Quantum 3, 161 (2019).
  * Rams  _et al._ [2018] M. M. Rams, P. Sierant, O. Dutta, P. Horodecki, and J. Zakrzewski, At the limits of criticality-based quantum metrology: Apparent super-heisenberg scaling revisited, Phys. Rev. X 8, 021022 (2018).
  * Garbe  _et al._ [2020] L. Garbe, M. Bina, A. Keller, M. G. A. Paris, and S. Felicetti, Critical quantum metrology with a finite-component quantum phase transition, Phys. Rev. Lett. 124, 120504 (2020).
  * Chu  _et al._ [2021] Y. Chu, S. Zhang, B. Yu, and J. Cai, Dynamic framework for criticality-enhanced quantum sensing, Phys. Rev. Lett. 126, 010502 (2021).
  * Salvia  _et al._ [2023] R. Salvia, M. Mehboudi, and M. Perarnau-Llobet, Critical quantum metrology assisted by real-time feedback control, Phys. Rev. Lett. 130, 240803 (2023).
  * Mihailescu  _et al._ [2023b] G. Mihailescu, A. Bayat, S. Campbell, and A. K. Mitchell, Multiparameter critical quantum metrology with impurity probes (2023b), arXiv:2311.16931 [quant-ph] .
  * De Palma  _et al._ [2017] G. De Palma, A. De Pasquale, and V. Giovannetti, Universal locality of quantum thermal susceptibility, Phys. Rev. A 95, 052115 (2017).
  * Glatthard  _et al._ [2022b] J. Glatthard, J. Rubio, R. Sawant, T. Hewitt, G. Barontini, and L. A. Correa, Optimal cold atom thermometry using adaptive bayesian strategies, PRX Quantum 3, 040330 (2022b).
  * Aybar  _et al._ [2022] E. Aybar, A. Niezgoda, S. S. Mirkhalaf, M. W. Mitchell, D. Benedicto Orenes, and E. Witkowska, Critical quantum thermometry and its feasibility in spin systems, Quantum 6, 808 (2022).
  * Razavian  _et al._ [2019] S. Razavian, C. Benedetti, M. Bina, Y. Akbari-Kourbolagh, and M. G. Paris, Quantum thermometry by single-qubit dephasing, The European Physical Journal Plus 134, 284 (2019).
  * Brattegard and Mitchison [2023] S. Brattegard and M. T. Mitchison, Thermometry by correlated dephasing of impurities in a 1d fermi gas (2023), arXiv:2307.10132 [cond-mat.quant-gas] .
  * Friedli and Velenik [2017] S. Friedli and Y. Velenik,  _Statistical Mechanics of Lattice Systems: A Concrete Mathematical Introduction_ (Cambridge University Press, 2017).
  * Gross and Bloch [2017] C. Gross and I. Bloch, Quantum simulations with ultracold atoms in optical lattices, Science 357, 995 (2017).
  * Bernien  _et al._ [2017] H. Bernien, S. Schwartz, A. Keesling, H. Levine, A. Omran, H. Pichler, S. Choi, A. S. Zibrov, M. Endres, M. Greiner, V. Vuletić, and M. D. Lukin, Probing many-body dynamics on a 51-atom quantum simulator, Nature 551, 579 (2017).
  * Browaeys  _et al._ [2016] A. Browaeys, D. Barredo, and T. Lahaye, Experimental investigations of dipole–dipole interactions between a few Rydberg atoms, J. Phys. B 49, 152001 (2016).
  * Barredo  _et al._ [2018] D. Barredo, V. Lienhard, S. de Léséleuc, T. Lahaye, and A. rowaeys, Synthetic three-dimensional atomic structures assembled atom by atom, Nature 561, 79 (2018).
  * Jelezko  _et al._ [2004a] F. Jelezko, T. Gaebel, I. Popa, A. Gruber, and J. Wrachtrup, Observation of coherent oscillations in a single electron spin, Phys. Rev. Lett. 92, 076401 (2004a).
  * Jelezko  _et al._ [2004b] F. Jelezko, T. Gaebel, I. Popa, M. Domhan, A. Gruber, and J. Wrachtrup, Observation of coherent oscillation of a single nuclear spin and realization of a two-qubit conditional quantum gate, Phys. Rev. Lett. 93, 130501 (2004b).
  * Robledo  _et al._ [2011] L. Robledo, L. Childress, H. Bernien, B. Hensen, P. F. A. Alkemade, and R. Hanson, High-fidelity projective read-out of a solid-state spin quantum register, Nature 477, 574 (2011).
  * Abobeih  _et al._ [2018] M. H. Abobeih, J. Cramer, M. A. Bakker, N. Kalb, M. Markham, D. J. Twitchen, and T. H. Taminiau, One-second coherence for a single electron spin coupled to a multi-qubit nuclear-spin environment, Nat. Commun. 9, 2552 (2018).
  * Degen  _et al._ [2021] M. J. Degen, S. J. H. Loenen, H. P. Bartling, C. E. Bradley, A. L. Meinsma, M. Markham, D. J. Twitchen, and T. H. Taminiau, Entanglement of dark electron-nuclear spin defects in diamond, Nat. Commun. 12, 3470 (2021).
  * Gulka  _et al._ [2021] M. Gulka, D. Wirtitsch, V. Ivády, J. Vodnik, J. Hruby, G. Magchiels, E. Bourgeois, A. Gali, M. Trupke, and M. Nesladek, Room-temperature control and electrical readout of individual nitrogen-vacancy nuclear spins, Nat. Commun. 12, 4421 (2021).
  * Cetina  _et al._ [2015] M. Cetina, M. Jag, R. S. Lous, J. T. M. Walraven, R. Grimm, R. S. Christensen, and G. M. Bruun, Decoherence of impurities in a fermi sea of ultracold atoms, Phys. Rev. Lett. 115, 135302 (2015).
  * Cetina  _et al._ [2016] M. Cetina, M. Jag, R. S. Lous, I. Fritsche, J. T. M. Walraven, R. Grimm, J. Levinsen, M. M. Parish, R. Schmidt, M. Knap, and E. Demler, Ultrafast many-body interferometry of impurities coupled to a fermi sea, Science 354, 96 (2016).
  * Yu  _et al._ [2023] M. Yu, O. Gühne, and S. Nimmrichter, Exact entanglement dynamics of two spins in finite baths, Phys. Rev. A 108, 022415 (2023).
  * Braunstein and Caves [1994] S. L. Braunstein and C. M. Caves, Statistical distance and the geometry of quantum states, Phys. Rev. Lett. 72, 3439 (1994).
  * Cramér [1999] H. Cramér,  _Mathematical Methods of Statistics (PMS-9)_ (Princeton University Press, 1999).
  * Paris [2009] M. G. A. Paris, Quantum estimation for quantum technology, International Journal of Quantum Information 07, 125 (2009).
  * Breuer and Petruccione [2002] H. Breuer and F. Petruccione,  _The Theory of Open Quantum Systems_ (Oxford University Press, 2002).
  * Stanley [1987] H. Stanley,  _Introduction to Phase Transitions and Critical Phenomena_, International series of monographs on physics (Oxford University Press, 1987).
  * Niknam  _et al._ [2021] M. Niknam, L. F. Santos, and D. G. Cory, Experimental detection of the correlation rényi entropy in the central spin model, Phys. Rev. Lett. 127, 080401 (2021).
  * Onsager [1944] L. Onsager, Crystal statistics. i. a two-dimensional model with an order-disorder transition, Phys. Rev. 65, 117 (1944).
  * Beekman  _et al._ [2019] A. J. Beekman, L. Rademaker, and J. van Wezel, An introduction to spontaneous symmetry breaking, SciPost Phys. Lect. Notes , 11 (2019).
  * Kadanoff [2009] L. P. Kadanoff, More is the same; phase transitions and mean field theories, Journal of Statistical Physics 137, 777 (2009).
  * Bavaresco  _et al._ [2023] J. Bavaresco, P. Lipka-Bartosik, P. Sekatski, and M. Mehboudi, Designing optimal protocols in bayesian quantum parameter estimation with higher-order operations (2023), arXiv:2311.01513 [quant-ph] .