# Evidence for an abundant old population of Galactic ultra long period
magnetars and implications for fast radio bursts

P. Beniamini∗1,2 Z. Wadiasingh ∗3,4,5 J. Hare4,5 K. M. Rajwade6 G. Younes4,7
A. J. van der Horst7
1 Department of Natural Sciences, The Open University of Israel, P.O Box 808,
Ra’anana 4353701, Israel
2 Astrophysics Research Center of the Open university (ARCO), The Open
University of Israel, P.O Box 808, Ra’anana 4353701, Israel
3 Department of Astronomy, University of Maryland, College Park, Maryland
20742, USA
4 Astrophysics Science Division, NASA Goddard Space Flight Center, Mail Code
661, Greenbelt, Maryland 20771, USA
5 Center for Research and Exploration in Space Science and Technology,
NASA/GSFC, Greenbelt, Maryland 20771, USA
6 ASTRON, the Netherlands Institute for Radio Astronomy, Oude Hoogeveensedijk
4, 7991 PD Dwingeloo, The Netherlands
7 Department of Physics, George Washington University, 725 21st St NW,
Washington, DC 20052, USA

###### Abstract

Two recent discoveries, namely PSR J0901–4046 and GLEAM-X J162759.5–523504.3
(hereafter GLEAM-X J1627), have corroborated an extant population of radio-
loud periodic sources with long periods (76 s and 1091 s respectively) whose
emission can hardly be explained by rotation losses. We argue that GLEAM-X
J1627 is a highly-magnetized object consistent with a magnetar (an ultra long
period magnetar - ULPM), and demonstrate it is unlikely to be either a
magnetically or a rotationally-powered white dwarf. By studying these sources
together with previously detected objects, we find there are at least a
handful of promising candidates for Galactic ULPMs. The detections of these
objects imply a substantial number, $N\gtrsim 12800^{+19000}_{-10100}$ and
$N\gtrsim 510^{+500}_{-420}$ for PSR J0901–4046 like and GLEAM-X J1627 like
objects, respectively, within our Galaxy. These source densities, as well as
cooling age limits from non-detection of thermal X-rays, Galactic offsets,
timing stability and dipole spindown limits, all imply the ULPM candidates are
substantially older than confirmed galactic magnetars and that their formation
channel is a common one. Their existence implies widespread survival of
magnetar-like fields for several Myr, distinct from the inferred behaviour in
confirmed Galactic magnetars. ULPMs may also constitute a second class of FRB
progenitors which could naturally exhibit very long periodic activity windows.
Finally, we show that existing radio campaigns are biased against detecting
objects like these and discuss strategies for future radio and X-ray surveys
to identify more such objects. We estimate that ${\cal O}(100)$ more such
objects should be detected with SKA-MID and DSA-2000.

###### keywords:

stars: magnetars - fast radio bursts - radio continuum: transients - white
dwarfs

††pagerange: Evidence for an abundant old population of Galactic ultra long
period magnetars and implications for fast radio bursts–Evidence for an
abundant old population of Galactic ultra long period magnetars and
implications for fast radio bursts**footnotetext: These authors contributed
equally to this work

##  1 Introduction

Since the first radio detection of an astronomical object by Karl Jansky,
about 90 years ago, radio astronomy has become a major pillar in astrophysics.
It has revealed classes of objects and types of radiation that were previously
unknown such as pulsars, quasars, the cosmic microwave background, and fast
radio bursts (FRBs: Linscott & Erkes, 1980; Lorimer et al., 2007). Still, much
of the parameter space for radio detections remains as of yet unexplored
(Cordes & McLaughlin, 2003; Pietka et al., 2015), particularly for rare or
intermittent transients that are currently being discovered (e.g., Wang et
al., 2021, 2022). This in large part owes to the many trade-offs that must be
made when designing a radio experiment or survey, e.g. selecting the bands to
be observed, compromising between area of coverage and limiting flux or
compromising between temporal and spectral resolution, sensitivity and
specificity of triggering algorithms, and computational or storage
constraints, just to name a few.

At the faint end, modern day radio experiments routinely measure lower
luminosity emission, i.e. intrinsically faint sources, than probed in any
other region of the electromagnetic spectrum. As an example, we can consider
the electromagnetic observations of the afterglow of the first GW detected NS
merger, GRB 170817A. In that case, the VLA radio limits reached as deep as
$\sim 3\times 10^{-19}\mbox{erg s}^{-1}\mbox{cm}^{-2}$ (Alexander et al.,
2020) while even $\sim 100$ks Chandra $X$-ray limits are more than a thousand
times brighter (Troja et al., 2020). Similarly, while many fast radio burst
(FRB) models predict inherently bright multi-wavelength counterparts (often
expected to be orders of magnitude more luminous than the radio emission), so
far (barring the somewhat special case of the Galactic FRB (Mereghetti et al.,
2020; Tavani et al., 2020; Ridnaia et al., 2020; Li et al., 2020; Margalit et
al., 2020; Lu et al., 2020; Li et al., 2022), only upper limits exist in other
wavelengths, and these are typically far above the expected sources’ output.

In light of this, it is perhaps to be expected that many new and exciting
astrophysical objects could be found in current and future radio surveys. In
the recent months, two particularly interesting such objects were detected,
GLEAM-X J162759.5–523504.3 (henceforth GLEAM-X J1627 for short) detected by
the MWA (Hurley-Walker et al., 2022) and PSR J0901–4046 discovered by the
MeerTRAP program (Caleb et al., 2022). Both objects exhibit periodicity on
long timescales (1091 s and 76 s respectively) and exhibit some magnetar-like
phenomenology***Although neither would presently be classified as AXPs or
SGRs, the original two sources classes which fit the historical high-energy
operational definition of and motivation for “magnetars” (Usov, 1984;
Paczynski, 1992; Duncan & Thompson, 1992; Thompson & Duncan, 1995, 1996;
Kouveliotou et al., 1998; Kaspi et al., 2003) where the present luminosity is
manifestly greater than spin-down power.. Indeed, there have been a variety of
previously detected sources with some similar features as will be discussed in
this work. However, the association of these objects with magnetars is far
from obvious. Galactic magnetars, as historically defined, are observed to
have periods $2\,\rm{s}\\!\lesssim\\!\,P\,\\!\lesssim\\!12\,$s with one
notable exception (§2.3). The dearth of longer period conventional magnetars
is understood to stem from the freezing of the spindown (SD) process as their
(crust-anchored) dipole field $B_{\rm dip}$ decays on a timescale of
$\tau_{d,{\rm dip}}\sim 10^{3}-10^{4}$ yrs (Kouveliotou et al., 1998; Colpi et
al., 2000; Dall’Osso et al., 2012; Viganò et al., 2013; Beniamini et al.,
2019; Jawor & Tauris, 2022). Clearly if the new class of observed objects are
indeed highly magnetized, they require a different dominant spindown mechanism
to drive their evolution, potentially related to a retardation of field decay
in their crusts (with constraints on core field expulsion). In this work we
adopt a broad definition of what a magnetar is – a highly magnetized neutron
star whose present observed luminosity does not necessarily exceed its spin-
down power (but likely did exceed it sometime in its past).

Various lines of evidence (e.g., Spitler et al., 2016; Lyutikov, 2017;
Wadiasingh & Timokhin, 2019; Ravi, 2019), and perhaps most remarkably the
association with a known Galactic magnetar SGR 1935+2154 (Bochenek et al.,
2020; The Chime/Frb Collaboration et al., 2020), connect FRBs with magnetar
sources. Therefore, the discovery that some FRBs have periodic activity
windows (Chime/Frb Collaboration et al., 2020; Pastor-Marazuela et al., 2020;
Rajwade et al., 2020a; Cruces et al., 2021), lasting tens of days and more,
was unexpected. Beniamini et al. (2020) explored the idea that the observed
periods reflect the rotational periods of magnetar-like sources. These authors
have shown that giant flare (GF) kicks, charged particle winds or fallback
accretion can all lead to enhanced spindown and that as a result, a small
fraction of magnetars (that they dubbed ULPMs, for ultra long period
magnetars) may develop extremely long rotation periods. It is still actively
debated what is the source of the observed periodicity seen in FRBs (other
suggestions include binary periodicity and free precession, Lyutikov et al.
2020; Ioka & Zhang 2020; Levin et al. 2020; Zanazzi & Lai 2020). However, FRB
20180916’s chromatic windowing, lack of strong dispersion or rotation measure,
random flux variation with activity phase, and absence of a low frequency
cutoff (Pastor-Marazuela et al., 2020; Pleunis et al., 2021a; Mckinven et al.,
2022; Bethapudi et al., 2022) disfavors simple binary shrouding models. There
are exist deep limits on any underlying shorter periodicity less than 1 ks for
FRB 20121102 (Li et al., 2021; Hewitt et al., 2022) and other prolific
repeaters (e.g., FRB 20180916, Chime/Frb Collaboration et al., 2020).

A possible connection to a Galactic (and nearby) class of objects could
therefore hold an important key towards deciphering the origin of FRB
periodicity. Even if FRB periodicity is unrelated to their rotational periods,
the class of Galactic objects with magnetar-like features and long periods can
be considered as ULPM-candidates and the enhanced spindown mechanisms
discussed in the FRB context to create ULPMs should be revisited in this new
context. There is also tentative evidence that FRB sources may encompass two
or more populations (Pleunis et al., 2021b). ULPMs would naturally fit as a
second, perhaps older population, compared to FRBs from relatively young SGR
1935+2154-like (or younger) magnetars. Evidence of offsets from a star forming
region for the location of FRB 20180916 (Tendulkar et al., 2021) also are
compatible with an older neutron star (NS), disfavoring precession scenarios
for the periodic windowing activity. Lastly, we note the recent discovery of a
repeating FRB 20200120E in a globular cluster of M81 (Bhardwaj et al., 2021;
Kirsten et al., 2022; Majid et al., 2021) that exhibits ‘burst storms’ much
like galactic magnetars (Nimmo et al., 2022) and for which short-term
periodicity signatures can be ruled out in the burst arrival time data (in the
range $1\mbox{ ms}\\!\lesssim\\!P\lesssim\\!25\\!\mbox{ s}$). As pulsars in
globular clusters are at least few Myr old (Michel, 1993; Lyne et al., 1996;
Manchester et al., 2005; Boyles et al., 2011; Lynch et al., 2012) and star
formation is almost nonexistent in such environments, an object much older
than conventional magnetars produced in a low-kick process is favored,
possibly produced by electron capture supernovae, accretion-induced collapse
or white dwarf (WD) collisions rather than a conventional core collapse
supernovae (CCSN) (Lu et al., 2021; Kremer et al., 2021).

In this paper we present the emerging class of Galactic ULPM candidates (§2).
We then discuss empirical evidence supporting episodes of enhanced spindown in
known Galactic magnetars in §3. We then explore the source densities of ULPMs
inferred from observations (§4) use this and other arguments to estimate the
ages of the ULPM candidates, §5, and consider the implications regarding the
temporal evolution of their magnetic field (§6). In §7 we turn to show that
the population of Galactic ULPM candidates can only be related to FRB
periodicities if the latter are extremely strongly beamed and / or only a
small fraction of ULPMs are efficient FRB producers. Finally, in §8 we discuss
prospects and strategies for future detections of ULPMS in radio and X-ray
campaigns. Our results are summarized in §9.

Magnetar | P | $\dot{P}$ | $B_{\rm SD}$ | SD age | Dist. | $F_{X}$ | $kT_{\rm surface}$ | Gal. b | Duty Cycle | Refs.
---|---|---|---|---|---|---|---|---|---|---
candidate | (s) |  (s s${}^{-1})$ | (G) | (yr) | (pc) |  (erg s-1 cm-2) | (keV) |  | (radio) |
GLEAM-X J1627 | 1091 | $<10^{-9}$ | $\lesssim 7\times 10^{16}$ | – | $1300$ | $<5\times 10^{-14}$ | $<0.072$ | $-2.6^{\circ}$ | $2-5\%$ | [1]
PSR J0901–4046 | 76 | $\sim 2.21\times 10^{-13}$ | $2.6\times 10^{14}$ | $5.3\times 10^{6}$ | 327-450 | $<2\times 10^{-13}$ | $<0.064$ | $+3.7^{\circ}$ | $\sim 0.7\%$ | [2]
PSR J0250+5854 | 23.54 | $2.72\times 10^{-14}$ | $2.56\times 10^{13}$ | $1.37\times 10^{7}$ | 1558 | $<2\times 10^{-15}$ | $<0.053$ | $-0.5^{\circ}$ | $0.3-0.4\%$ | [3]
1E 161348–5055 | 24030.42 | $<1.6\times 10^{-9}$ | $<10^{17}$ | – | 3100-4700 |  few$\times 10^{-11}$ | $0.58$ | $-0.4^{\circ}$ | – | [4]
FRB 20180916B (R3) | $1.4\times 10^{6}$ | – | – | – | $4\times 10^{7}$ | – | – | – | $\sim 40\%$ | [5]
IGR J16358–4726 | 5970 | – | – | – | >7000? | – | – | $\sim 0.0^{\circ}$ | – | [6]
4U 2206+54† | 5750 | $6\times 10^{-7}$ | – | – | 3700 |  3$\times 10^{-11}$ | $\sim 1.5$ | $-1.1^{\circ}$ | – | [7]
SGR 0755-2933 | 308 | – | – | – | 3500 | – | – | $-0.6^{\circ}$ | – | [8]
4U 1954+319 | 20500 | – | – | – | 3300 | – | – | $+1.9^{\circ}$ | – | [9]
4U 0114+65 | 9350 | – | – | – | 7000 | – | – | $+2.6^{\circ}$ | – | [10]
GCRT J1745–3009 | 4630 | – | – | – | 8500? | – | – | $-0.54^{\circ}$ | $\sim 13\%$ | [11]
PSR J2251–3711∗ | 12.12 | $1.3\times 10^{-14}$ | $2.5\times 10^{13}$ | $1.47\times 10^{7}$ | 540-1344 | $\lesssim 10^{-14}$ | – | $-62.9^{\circ}$ | $0.6-0.7\%$ | [12]
RX J0420–5022 | 3.45 | $2.76\times 10^{-14}$ | $9.9\times 10^{12}$ | $1.98\times 10^{6}$ | 345 | $4.4\times 10^{-14}$ | 0.046 | $-44.4^{\circ}$ | – | [A]
PSR J0726–2612 | 3.44 | $2.93\times 10^{-13}$ | $3.2\times 10^{13}$ | $1.86\times 10^{5}$ | 100-1000 (<3200) | $1.6\times 10^{-12}$ | 0.074 | $-4.64^{\circ}$ | $2-3\%$ | [B]
RX J0720–3125 | 8.39 | $7\times 10^{-14}$ | $2.45\times 10^{13}$ | $1.9\times 10^{6}$ | 286 | $3.2\times 10^{-11}$ | 0.082 | $-8.17^{\circ}$ | – | [A]
RX J0806–4123 | 11.37 | $5.6\times 10^{-14}$ | $2.54\times 10^{13}$ | $3.24\times 10^{6}$ | 250 | $2.1\times 10^{-12}$ | 0.058 | $-4.98^{\circ}$ | – | [A]
RX J1308+2127 | 10.31 | $1.12\times 10^{-13}$ | $3.44\times 10^{13}$ | $1.46\times 10^{6}$ | 500 | $8.4\times 10^{-12}$ | 0.069 | $+83.1^{\circ}$ | – | [A]
RX J1605+3249 | 6.88 | – | – | – | 390 | $1.1\times 10^{-11}$ | 0.065 | $+48.0^{\circ}$ | – | [A]
RX J1856–3754 | 7.06 | $2.98\times 10^{-14}$ | $1.47\times 10^{13}$ | $3.76\times 10^{6}$ | 156 | $2.7\times 10^{-11}$ | 0.062 | $-17.2^{\circ}$ | – | [A]
RX J2143+0654 | 9.43 | $4.1\times 10^{-14}$ | $2\times 10^{13}$ | $3.65\times 10^{6}$ | 150-430 | $1.3\times 10^{-11}$ | 0.055 | $-33.1^{\circ}$ | – | [A]
Table 1: Candidates and potentially related Objects which demand more
observational scrutiny. Data collected from the literature for ULPM candidates
and (below the double solid line) another class of unusual NSs, XDINS.
† The $\dot{P}$ is likely dominated by accretion torques. ∗ XDIN candidate.
References: [1] Hurley-Walker et al. (2022) [2] Caleb et al. (2022) [3] Tan et
al. (2018) [4] De Luca et al. (2006); De Luca et al. (2008); Esposito et al.
(2011); D’Aì et al. (2016); Braun et al. (2019) [5] Chime/Frb Collaboration et
al. (2020); Pastor-Marazuela et al. (2020); Bethapudi et al. (2022) [6]
Revnivtsev et al. (2003); Kouveliotou et al. (2003); Patel et al. (2004, 2007)
[7] Torrejón et al. (2018) [8] Barthelmy et al. (2016); Harrison et al. (2017)
[9] Enoto et al. (2014); Hinkle et al. (2020); Bozzo et al. (2022) [10] Li &
van den Heuvel (1999); Sanjurjo-Ferrrín et al. (2017) [11] Hyman et al. (2005,
2007); Spreeuw et al. (2009); Kaplan et al. (2008) [12] Morello et al. (2020)
[A] Kaplan & van Kerkwijk (2009); Tetzlaff et al. (2011); Yoneyama et al.
(2019); Potekhin et al. (2020) [B] Burgay et al. (2006); Speagle et al.
(2011); Rigoselli et al. (2019)

##  2 Galactic ULPM Candidates

We list below several Galactic objects which we consider as tenable ULPM
candidates, and briefly summarize the evidence supporting such an association.
Additional, less secure, associations are listed in Appendix A. These are not
part of the analysis presented in the paper and are provided for completeness.
Table 1 summarizes some of the main observations of all the candidates as well
as those of XDINS which share some of the properties with ULPM candidates
(nearby objects, low inferred surface temperatures and magnetar-like fields).

###  2.1 PSR J0901–4046

PSR J0901–4046 is a prominent $P\approx 76$ s pulsar discovered by MeerTRAP
(Rajwade et al., 2020b) with a $\dot{P}\sim 2.21\times 10^{-13}$ s s-1
implying a spindown polar dipolar magnetic field component of $B_{\rm SD}\sim
2.6\times 10^{14}$ G (Caleb et al. 2022; see §2.1.1 below for a more general
limit on $B$ that does not assume dipole SD). A faint radio shell appears to
be spatially coincident with PSR J0901–4046 (Caleb et al., 2022), although, if
associated, it is much too small to be consistent with a supernova remnant.
PSR J0901–4046 displays many of the radio characteristics seen in common radio
pulsars, such as polarization position angle swings and variability within
single pulses, yet with relatively wide duty cycle for its long period (see
Appendix B). It also exhibits unusual features, such as harmonically-spaced
QPOs (see discussion in Appendix C). Moreover, it is remarkably stable in
timing, with timing residuals of $5.7$ ms over 7.4 months – a relative
stability of $\sim 10^{-6}$, and no timing instabilities or glitches in
contrast to what is common within the young magnetar population.

Some single pulses in PSR J0901–4046 also exhibit partial nulling or spikes
within the pulse that occurs on timescales $\ll 1$ ms, less than the temporal
sampling resolution of MeerKAT 1.4 GHz observations. This, along with a
characteristic radio luminosity of sub-components of $\gtrsim 1$ mJy implies a
brightness temperature $T_{b}>10^{16}$ K (assuming a distance 328 pc from the
measured DM), i.e. favoring coherent pulsar-like emission and a NS origin. The
discussion regarding PSR J0901–4046 hosting a thick crust and strong core
field are given in Appendix C and §6.2, respectively.

A potential optical counterpart of PSR J0901–4046 is detected by Gaia, having
a G-band magnitude of 17.5, and is offset from the source by about 3.3′′. The
measured Gaia parallax, places the optical source at a distance of 2.4 kpc
(Bailer-Jones et al., 2021), which is far beyond the inferred dispersion
measure distance for PSR J0901–4046. Additionally, Caleb et al. (2022) suggest
the optical source is a reddened A-type star, based on optical spectroscopic
and photometric observations of the source. For these reasons, it appears that
PSR J0901–4046 does not have any plausible optical counterpart in Gaia.

####  2.1.1 The magnetic field of PSR J0901–4046

The spin down value of $2\times 10^{-13}$ s s-1 for PSR J0901–4046 results in
the inferred magnetar-like dipolar magnetic field estimate, which exceeds
$10^{14}$ G and argues for its magnetar nature. Yet, it may be argued that
this is an overestimate due to multipolar fields or particle winds enhancing
spin down. The former argument can be largely discarded, as the influence of
higher-order multipoles on spindown is negligible given its large light
cylinder at $\sim 10^{4}$ stellar radii (unlike the case for millisecond
pulsars). Particle winds, and enhanced spindown, on the other hand may lower
the magnetization estimate but even this is disfavored as described below. In
order for the wind to dominate the spindown over dipole EM spindown, it must
be powered magnetically. This demands the wind luminosity is limited by
$L_{\rm p}\lesssim E_{\rm B}/\tau\propto B^{2}$. At the same time, the
particle wind luminosity decreases with $B$ according to $L_{\rm p}\approx
6c^{3}I^{2}\dot{P}^{2}/(B^{2}P^{2}R^{6})\propto B^{-2}$ (Harding et al.,
1999). Here $B$ is the dipole field at the pole. Combining these expressions,
we find a lower limit on the required $B$ to sustain a sufficiently strong
particle wind

$B\gtrsim\left(\frac{6c^{3}I^{2}\dot{P}^{2}\tau}{P^{2}R^{9}f^{2}}\right)^{1/4}\approx
10^{13}f^{-1/2}\tau_{\rm Myr}^{1/4}\quad\mbox{G}$ (1)

where $\tau$ is the time over which this luminosity is sustained (normalized
here by the age estimate of the source, $\sim 1$ Myr) and $f\geq 1$ is the
ratio between the internal and dipole field in the star. Eq. (1) is strictly a
lower limit on $B$, as it assumes $100\%$ efficiency in converting the
magnetic energy to a wind luminosity. As we can see a strong field is still
required in this case. This $B$ limit, in turn, corresponds to an upper limit
on the particle wind luminosity, $L_{\rm p}$,

$L_{\rm
p}\lesssim\left(\frac{6c^{3}I^{2}\dot{P}^{2}f^{2}}{P^{2}R^{3}\tau}\right)^{1/2}\approx
6\times 10^{30}f\tau_{\rm Myr}^{-1/2}\quad\mbox{erg s}^{-1}.$ (2)

An additional argument disfavors the existence of any dynamically significant
wind (and in turn increases the lower limit on $B$): PSR J0901–4046 is
extraordinarily stable in its observed timing.

Figure 1: Left: Detectability, with Gaia G band (black) and DECam r-band
(orange), of hot WDs as a function of mass and age. WDs having ages younger
than the respective lines for a given mass would be detectable. Right: Similar
to the left plot, but for effective WD temperature as a function of mass. WDs
hotter than the respective lines for a given mass would be detectable. The WD
evolutionary tracks of Bédard et al. 2020 are adopted here.

###  2.2 GLEAM-X J1627

GLEAM-X J1627 is a bright $20-50$ Jy periodic radio transient discovered at
low ($\sim 70-230$ MHz) frequencies by the MWA (Hurley-Walker et al., 2022).
It exhibits pulsations at $P\approx 1091$ s which are nearly $100\%$ linearly
polarized, as well as unresolved bright bursts that imply a brightness
temperature $\gg 10^{12}$ K for its DM-inferred distance of $1.3\pm 0.5$ kpc,
possibly with $\gg 100$ kJy bursts for the unresolved spikes if they are
intrinsically ${\cal O}(10\rm\,ms)$ in duration (i.e. similar in duration to
extragalatic FRBs at low frequencies after accounting for scattering). At this
distance, its persistent radio luminosity exceeded $10^{31}\mbox{ erg s}^{-1}$
in the limited bandwidth of MWA, and is possibly higher given its relatively
flat radio spectrum. Owing to its long period, high persistent radio
luminosity, and spindown upper limit ($\dot{P}\lesssim 10^{-9}$ s s-1), the
source is unlikely to be a rotation-powered object with moment of inertia
compatible with neutron or quark stars***Of course, these conclusions depend
on the beaming of the radiation as well as the efficiency of conversion
between spindown luminosity and radio emission (see Erkut 2022). . It’s radio
duty cycle is also much too large to be compatible with its physically small
polar cap size or known empirical relations for rotation-powered radio pulsars
(see Appendix B). This suggests it could be of a magnetar origin. If the
observed period is due to rotation, this would imply that this object is an
ULPM. Alternatively, Ekşi & Şaşmaz (2022) have suggested that GLEAM-X J1627 is
a $\sim 0.1$ s period magnetar that is undergoing precession with the $\sim 1$
ks observed period duration. A major difficulty with such a scenario is that
the polarization would be variable and likely depolarized on MWA’s 0.5 second
temporal sampling. There could also be a beat period with the temporal
sampling and the putative short period. In addition, the large toroidal
magnetic field $(\sim 10^{16}$ G) required for sustaining the deformation
needed for precession, and a hot core required for the inhibition of
superfluidity (e.g., Shaham, 1977) generally means the object must be
extremely young which is inconsistent with a large number of observational
constraints (see §5 for details). Other authors (Loeb & Maoz, 2022; Katz,
2022) have suggested that GLEAM-X J1627 may instead be a rotation or magnetic
powered WD. We outline below some challenges that this interpretation
encounters.

####  2.2.1 Optical constraints on a WD or stellar association to GLEAM-X
J1627

We do not find a plausible counterpart to GLEAM-X J1627 in either Gaia or
DECam archival data. Assuming no counterpart to GLEAM-X J1627 is detected by
either Gaia or DECam (see Appendix D), we can place constraints on the
mass/age of a putative WD at 1.3 kpc. Using the evolutionary tracks of Bédard
et al. (2020), we estimate at what age (or temperature) a WD of a given mass
becomes undetectable***The actual spectra for the evolutionary tracks
calculated by Bédard et al. (2020) are not provided, so we use a black-body
model for the WDs.. For Gaia, we use the G band magnitude and a limiting
magnitude of 21, while for DECam, we use the r-band to minimize the effects of
absorption, and a limiting magnitude of 22.8. We account for the extinction in
the direction of GLEAM-X J1627 at a distance of 1.3 kpc using the mwdust
package (Bovy et al., 2016). These estimates are shown in Figure 1. For a
large range of masses ($0.2\\!<\\!M_{\rm WD}\\!<\\!1.2$ $M_{\odot}$) a WD
requires an age $t_{\rm WD}\gtrsim(1\\!-\\!3)\times 10^{8}$ yrs to be
undetected by DECam. Such an old age is a severe constraint to WD models of
GLEAM-X J1627 (see below).

####  2.2.2 Why GLEAM-X J1627 is likely a magnetized neutron or quark star and
not a WD

Magnetically powered WDs

Below we consider some limitations on a magnetically active WD as the source
of GLEAM-X J1627. The energetic radio output of GLEAM-X J1627 in its two-month
active episode was $\sim(10^{37}-10^{39})f_{\Omega}$ erg assuming a typical
flux density of $10-30$ Jy and flat spectral index over bandwidth of $\sim
100$ MHz (at a distance of 1.3 kpc). Hurley-Walker et al. (2022) estimate
$L_{R}\sim$ few $\times 10^{31}$ erg s-1 (isotropic-equivalent) based on the
spectral index and extrapolation to higher frequencies. The $f_{\Omega}<1$
above is the intrinsic (unknown) collimation of emission, which may be
conservatively taken as $f_{\Omega}\sim 10^{-2}$, the observed duty cycle of
pulsations. This disfavors magnetically-powered activity from a highly
magnetized WD similar to Ar Sco (Marsh et al., 2016; Buckley et al., 2017;
Stiller et al., 2018; Pelisoli et al., 2022), as the total field (and
ephemeral free energy, the magnetic field helicity) reservoir corresponds to
$E_{B}\\!\sim\\!B_{8}^{2}R_{8.5}^{3}/6\\!\sim\\!5\times 10^{40}$ erg. Here
$R\\!\sim\\!3\times 10^{8}$ cm is a typical WD radius. The efficiency of
conversion (from magnetic free energy) to coherent radio emission needs to be
$\eta\\!\sim\\!2\\!\times\\!10^{-4}-2\\!\times\\!10^{-2}$ for the energetics
to be justifiable. Yet, the non-thermal radio fluence of the outburst episode
likely severely underestimates the true energetics (e.g. source Poynting flux
and particle kinetic power, transient or persistent), as in known natural
coherent radio emitters. For bright broadband coherent radio emitters such as
rotation-powered pulsars $\eta\ll 1$ (e.g., Figure 4 in Wadiasingh & Timokhin,
2019) except the oldest, most charge-starved***For WDs, as pair production
cannot take place without extreme values of magnetization, and a relativistic
plasma is required for coherent radio emission, significant kinetic power must
go to baryons implying $\eta\ll 1$. ones in the class (e.g., Arzoumanian et
al., 2002). Likewise, the bright radio burst from SGR 1935+2154 was only $\sim
10^{-5}$ the fluence of the associated hard X-ray burst (Margalit et al.,
2020). Then, given the radio luminosity of GLEAM-X J1627, it can maintain this
level of activity up to at most

$\tau_{B}\lesssim\frac{\eta E_{B}}{f_{\Omega}L_{R}}\sim
50B_{8}^{2}R_{8.5}^{3}\eta f_{\Omega}^{-1}L_{R,31.5}^{-1}\quad\mbox{yr,}$ (3)

seven orders of magnitude shorter in age than the lower limits corresponding
to the optical limits from WD cooling (§2.2.1). One may imagine a situation in
which the active lifetime of the sources is shorter than their true ages,
however it would be highly unlikely to detect such an object if its “active
duty cycle" (i.e. active lifetime compared to true age), $\eta_{\rm act}$ were
$\ll 1$. In particular the number of inferred Galactic objects of this type
(see §4.1) would be increased by $\eta_{\rm act}^{-1}$. Independent of
$\eta_{\rm act}$, the low value of $\tau_{\rm B}$ leads to unrealistic minimum
requirements on the formation rate of high $B$ WDs as detailed next. As shown
in §4.1, the distance to GLEAM-X J1627 requires that the number of similar
active objects in the Galaxy be $N_{\rm act}\gtrsim 2500$ (where we have
considered the more conservative limit, that is independent of object type and
we have considered a 3D distribution of the putative objects, as appropriate
for WDs). This corresponds to a formation rate $\dot{r}\gtrsim
2500\tau^{-1}\mbox{ yr}^{-1}$ where $\tau\leq\tau_{\rm B}$ is the active
lifetime of the sources. Combining with Eq. (3), and relating to the Galactic
WD formation rate, $\dot{r}_{\rm WD}\approx 0.25\mbox{ yr}^{-1}$ (Hills,
1978), we find

$\frac{\dot{r}}{\dot{r}_{\rm WD}}\gtrsim
200B_{8}^{-2}R_{8.5}^{-3}f_{\Omega}\eta^{-1}L_{R,31.5}.$ (4)

We see that for $B\lesssim 2\times 10^{9}\mbox{ G}$, the required formation
rate is greater than that of all WDs, which is clearly ruled out. Even for
$B\approx 10^{10}$ G, one needs the formation rate of such extremely magnetic
WDs to be more than 2$\%$ of the total WD formation rate. This is inconsistent
with the observed population (Ferrario et al., 2015, 2020)—much weaker
magnetized WDs with $B\approx 10^{8}\mbox{ G}$ WDs are already limited to at
most a few percent of the total population.

Rotation powered WDs

A typical WD with $B\sim 10^{6}$ G, has a dipole spin power that is much too
low to account for the observed radio luminosity. Thus, even though the
rotational energy is $\sim 10^{45}M_{1M_{\odot}}R_{8.5}^{2}\mbox{ erg}\gg
E_{B}$, spindown cannot power the radio luminosity through dipole spindown,
unless the magnetic field is substantially larger.

An outlandishly high-field magnetized WD (with moment of inertia $I\sim
10^{50}$ g cm2, corresponding to one expected from WDs) rotationally powering
the radio emission is also disfavored. Such a putative WD would require a
global dipolar field $B\gtrsim 10^{11}$ G to power the observed luminosity
(assuming $\eta/f_{\Omega}\sim 1$), a regime in which the total stellar
magnetization is only an order of magnitude less than the gravitational
binding energy (i.e. the Chandrasekhar hydromagnetic limit for stability), and
two to three orders of magnitude larger than Ar Sco or the most highly
magnetized WDs known (Ferrario et al., 2015, 2020). For a typical WD core
conductivity, the resulting ohmic timescale for field decay is $1-3$ Gyr
(Cumming, 2002). This implies a persistent ohmic magnetic dissipation
luminosity of $L_{\rm md}\gtrsim 10^{30.5}B_{11}^{2}R_{8.5}$ erg s-1 and a
minimum effective temperature of $T_{\rm eff}\gtrsim 1.4\times
10^{4}\,B_{11}^{1/2}R_{8.5}^{-1/4}$ K***This estimate is strictly a lower
limit, as this magnetic dissipation is in additional to any residual heat from
formation that may be estimated from appropriate WD cooling models (for
varying age and mass).. At its distance and expected extinction, optical/UV
observations with HST or Webb could readily rule out a WD origin, irrespective
of high nebular emission expected from energetic WDs. Such temperatures are
ruled out by DECam in Figure 1 for masses up to $\sim 0.8M_{\odot}$.

A core difficulty with the WD origin of GLEAM-X J1627, is its small duty
cycle, on the order of a few percent, the likes of which has not been
previously observed in pulsating WDs. For instance, Ar Sco’s pulse duty cycle
exceeds $30\%$ (Buckley et al., 2017). Observing highly collimated radiation
(as needed for a small duty cycle) from a WD, presents also a theoretical
difficulty. To attain a high degree of beaming, one typically needs a
relativistic particle population. However, it is generally difficult to source
and accelerate plasma near the surface of a WD to relativistic velocities as
in NSs, considering that the escape velocity from the surface of the WD is
$\ll c$ so large voltages are not necessary to launch plasma. Furthermore
single-photon magnetic pair production requires high magnetic fields and small
curvature radii. This then sets demands on any voltage and magnetic field an
isolated WD must attain to source its relativistic plasma (in contrast to Ar
Sco, which has a stellar companion).

We may generically quantify the above arguments and constraints at the
measured period of GLEAM-X J1627 for arbitrary stellar compactness
$GM/(c^{2}R)$ versus dipolar surface magnetic field $B$, assuming the moment
of inertia is approximately $2MR^{2}/5$ (a good approximation for degenerate
stars within a factor of two). The following restrictions largely rule out
rotation-powered WDs for GLEAM-X J1627. Broadly, we require (i) a putative
rotation-powered nature of radio emission $\eta\dot{E}/f_{\Omega}>L_{R,\rm
obs}\sim~{}10^{31.5}$ erg s-1 (ii) electromagnetic spindown rate less than the
observed constraint***A more restrictive smaller value, e.g. in Extended
Figure 2 in Hurley-Walker et al. (2022) yields no viable phase space for a
rotation-powered WD scenario in Figure 2. $\dot{P}<P_{\rm limit}\sim 10^{-9}$
s s-1, (iii) hydromagnetic stability $U_{B}\sim B^{2}R^{3}\ll GM^{2}/R\sim
U_{G}$ (iv) compactness above the rotational shedding or breakup limit
$\Omega^{2}R\ll GM/R^{2}$ ($\Omega=2\pi/P$) and (v) curvature photon magnetic
pair cascades for high brightness temperature, highly polarized, pulsar-like
emission in a vacuum gap-like condition. Propagation of such radio emission in
the magnetosphere of the object also sets constraints on the hierarchy of
plasma and cyclotron frequencies demanded (if, for example, observed emission
are ordinary eigemodes of the magnetized plasma, e.g., Arons & Barnard, 1986),
but it turns out that the pair cascade constraint is most limiting for the
scales of interest here. As shown in §E, the gap constraint is essentially a
radio death line (assuming curvature radius $\rho_{c}\sim 10R$), derived using
the condition that the curvature photon gap height is much smaller than the
polar cap rim boundary $\ell_{\rm gap}<r_{\rm pc}\sim\sqrt{2\pi R^{3}/(Pc)}$.
Equivalently the voltage drop in the gap is less than the open field line
potential $\Delta V_{\rm gap}<\Phi_{\rm open}$, or the accelerated primary
luminosity (which cascades into $e^{+}e^{-}$ pairs) $L_{e^{\pm}}\sim 2\pi
A\ell_{\rm gap}^{2}\rho_{\rm gap}^{2}c<\dot{E}_{\rm SD}$. Here $A$ is the
characteristic polar cap area $\pi r_{\rm pc}^{2}$, while $\rho_{\rm
gap}\sim\rho_{\rm GJ}\sim B/(Pc)$ is the Goldreich-Julian charge density.
Finally, for a WD scenario, optical limits and standard WD cooling models
(§2.2.1) imply that the true age and thus also the characteristic spindown age
($\propto P/\dot{P}$) is greater than about $10^{8}$ yr. This is a
conservative estimate as epochs of enhanced spindown only lower the maximum
magnetization allowed.

We plot these constraints for $M=\\{0.4,1.2\\}M_{\odot}$ in Figure 2, along
with a horizontal line indicating a typical compactness value from known WD
mass-radius relations, as well as the mass and magnetization of Ar Sco
(Stiller et al., 2018; du Plessis et al., 2019; Pelisoli et al., 2022). There
are several import features of the constraints in Figure 2 to note. First, the
observed limit $\dot{P}_{\rm limit}\lesssim 10^{-9}$ s s-1 is quite
constraining to models with large magnetic fields, particularly objects with
low compactness or large radii. At the same time, the gap constraints above
and Eq. (11) require a minimum compactness for a given $B$, yielding the
rather restrictive allowed parameter space of a WD with outlandishly high
surface dipole field $B>10^{10}$ G and $\dot{E}_{\rm SD}\sim 10^{33}-10^{34}$
erg s-1. Note that even the narrow $B$ range that is consistent with
conditions (i)-(v) listed above, is ruled out for standard composition (and
mass) cooling limits considering the optical upper limits (§2.2.1). Thus the
arguments presented above largely rule out a rotation-powered WD nature for
GLEAM-X J1627.

Other WD models

We close this discussion with a brief description of other models that have
been suggested to power GLEAM-X J1627. Loeb & Maoz (2022) have suggested that
GLEAM-X J1627 is a hot proto-WD or sub-dwarf. Due to the large ($\sim
0.3R_{\odot}$) radius of the WD in this model (required to explain the
observed period as corresponding to the break-up limit), the spindown
luminosity in this model is extreme ($\sim 10^{38}\mbox{erg s}^{-1}$). The
corresponding effective temperature at the surface would be $\sim
10^{4.5}-10^{5}\mbox{ K}$ (Heber, 2009). The existence of such an object is
strongly ruled out by existing Gaia/DECam observations discussed above (see
Figure 1 and Appendix D). For instance, a hot proto-WD with $T\sim 10^{5}$ K
and other parameters as listed above would have a Gaia G band magnitude of
$\approx 12.5$ at the distance and corresponding absorption of GLEAM-X J1627,
which is about 8.5 magnitudes brighter than Gaia’s limiting magnitude. We
conclude that an ultra-long period magnetar description for GLEAM-X J1627 is a
more natural working hypothesis.

Figure 2: Constraints on a rotation-powered origin for GLEAM-X J1627 ($P=1091$
s) expressed in terms of stellar compactness versus dipolar magnetic field at
the surface, $B$. The red shaded region indicates the allowed region from
spindown energetics $\eta\dot{E}_{\rm SD}>f_{\Omega}L_{R}\sim
f_{\Omega}10^{31.5}$ erg s-1 and the limit $\dot{P}_{\rm SD}<P_{\rm
limit}\approx 10^{-9}$ s s-1 while the blue region is a curvature pair cascade
death line bounded by the hydromagnetic stability limit; both regions are for
mass $1.2M_{\odot}$. Only the region within the blue and red shaded regions is
allowed for rotation powered WDs. Likewise, the second lower panel indicates
similar regions for $0.4M_{\odot}$. The mass shedding or break up limit is
indicated by the gray horizontal line, as well as relevant compactness for
stable WD equation of state at the indicated masses. The thick solid green
segment indicates the allowed region, requiring outlandishly high
magnetization $B~{}>~{}10^{10}-10^{11}$ G for a putative WD. This allowed
parameter region is narrow and largely untenable. Furthermore, even this
narrow region is ruled out for optical limits corresponding to ages from
standard WD cooling models (black dashed line, applicable to the range of
standard WD compactness and composition).

###  2.3 1E 161348–5055

1E 161348–5055 is the pulsating ($P\sim 6.7$ hour) central compact object in
SNR RCW 103, and its magnetar nature is well-established (De Luca et al.,
2006; De Luca et al., 2008; Esposito et al., 2011; Tong et al., 2016; Rea et
al., 2016; D’Aì et al., 2016; Tendulkar et al., 2017; Ho & Andersson, 2017;
Borghese et al., 2018). From detailed modeling of the SNR, Braun et al. (2019)
estimate the age of 1E 161348–5055 to be $880-4400$ yrs.

###  2.4 PSR J0250+5854

Until recently, PSR J0250+5854 was the longest period bona fide pulsar with
$P\approx 23.5$ s and $B\approx 3\times 10^{13}$ G (Tan et al., 2018). It is
however consistent with being rotation-powered and has not yet exhibited any
magnetar-like activity. That being said, it’s spindown age of $14$ Myr is the
largest in Table 1. We note that in this case, the characteristic age may not
be the best estimate of the true age as slow can pulsars show a significant
amount of timing noise leading to spurious measurements of age (see Lower et
al., 2020, for more details). Thus its spindown age can be regarded as an
upper limit.

###  2.5 GCRT J1745–3009

GCRT J1745–3009 is the Galactic center “Burper" discovered serendipitously by
the VLA at 330 MHz (Hyman et al., 2005). It exhibited bright $\sim 1$ Jy
10-minute-wide “pulses" at a cadence of $P\sim 77$ minutes, as well as
sporadic on/off epochs in follow-up observations (Hyman et al., 2007; Spreeuw
et al., 2009). The $\lesssim$ minute-timescale variability implies a
brightness temperature $\gg 10^{12}(d/70\rm pc)^{2}$ K (Hyman et al., 2007;
Spreeuw et al., 2009) consistent with a coherent emission mechanism, if the
source is not associated with a nearby stellar object. Optical observations
(Kaplan et al., 2008) largely rule out any plausible nearby M/brown-dwarf
counterpart within $\sim 0.1-1$ kpc, implying GCRT J1745–3009 is indeed a
distant coherent emitter. The only known bright coherent emitters are NSs, and
if the pulsations are ascribed to spin periodicity***Spreeuw et al. (2009)
disfavors aperiodicity in the Hyman et al. (2005) data-set, suggesting the
source is a relatively stable rotator., the source cannot be rotation-powered,
suggesting a magnetar origin. In this regard, the source is analogous to
GLEAM-X J1627 but whose full window of activity was missed owing to
technological limitations of the time.

##  3 Evidence for enhanced spindown

As shown in §2, there are several ULPM candidates in the Galaxy. The existence
of such systems requires a mechanism for substantially spinning down some
magnetars to spin frequencies well below those that are realized by the
confirmed Galactic magnetar population. Indeed it has been shown that in the
presence of field decay, dipole spindown becomes increasingly inefficient, and
as a result standard magnetars’ spin periods tend to reach a maximum value of
$\sim 12$ s (Dall’Osso et al., 2012; Beniamini et al., 2019). However, the
presence of enhanced spindown mechanisms, operating for limited durations can
be inferred from observed features in the confirmed magnetar population. Below
we list some of this evidence, and discuss how it could result in extremely
prolonged spin periods in a fraction of the magnetar population.

Several known magnetars exhibit enhanced spindown in association with high-
energy activity. For example, SGR 1900+14, has exhibited an enhanced spindown
of $x_{\rm P}\equiv\Delta P/P\approx 10^{-4}$ after its 1998 giant flare (GF).
Since its prolific 2004 GF, SGR 1806$-$20 has increased its $\dot{P}$
considerably relative to its pre-GF value (Younes et al. 2015; although
$x_{\rm P}\lesssim 5\times 10^{-6}$ during the time of the main spike of the
GF, Woods et al. 2007). By 2012, its $P$ increased by $\sim\\!\,2\%$ compared
to the extrapolated value from 1994. Kinematic age constraints for SGR
1806$-$20 and SGR 1900+14 (Tendulkar et al., 2012) that imply effective
braking indices $n\ll 3$ also demand enhanced spindown in the past history of
these magnetars, either as a result of field decay or sporadic stretches with
monopolar ($n\approx 1$) winds (Thompson & Blaes, 1998; Harding et al., 1999).
More generally, it has been observed that immediately during magnetar
outbursts they may suffer ‘glitches’ or ‘anti-glitches’, where their periods
decrease or increase by $|x_{P}|\sim 10^{-9}-10^{-5}$ abruptly (Kaspi &
Beloborodov, 2017). Perhaps more important, is that these phenomena are often
followed by a subsequent episode of enhanced spindown, the effect of which
often dominates the temporal evolution of $P$. An illuminating example is 1E
1048.1-5937, which experienced several glitches over the years, that were
followed by enhanced spindown episodes on years-long timescale (Archibald et
al., 2020). These episodes resulted in $x_{P}\sim 2\times 10^{-4}$, about two
orders of magnitude greater than associated with the glitches themselves (and
with opposite sign). Another interesting recent example is that of Swift
J1818.0–1607 which after an outburst in 2020 has experienced both a candidate
glitch and an anti-glitch, as well as a long-term highly fluctuating (on a
timescale of days) spindown rate (Hu et al., 2020).

Such events of small fractional increase in $P$ can accumulate to produce a
large total increase in $P$, if their number per magnetar lifetime is $N_{\rm
P}\gtrsim\\!x_{\rm P}^{-1}$. Indeed, if the value of $x_{\rm P}$ is
independent of $P$, then $P$ can grow exponentially, $P_{\rm f}=P_{\rm
0}\exp(N_{\rm P}x_{\rm P})$ and $P_{\rm f}\gg P_{\rm 0}$ for $N_{\rm
P}\gtrsim\\!x_{\rm P}^{-1}$. For the inferred GF energy, $E_{\rm GF}\approx
4\times 10^{44}\,$erg and $x_{\rm P}\approx 10^{-4}$ from SGR 1900+14, a
significant increase of $P$ due to the cumulative effect of many GFs requires
an initial magnetic energy reservoir of $E_{B,0}\approx E_{B,{\rm
int},0}>E_{\rm GF}x_{\rm P}^{-1}\sim 4\times 10^{48}\,$erg that can be used to
power GFs. This magnetic energy corresponds to a quite reasonable initial
internal magnetic field strength $B_{{\rm int},0}~{}>~{}5\times 10^{15}$ G
(compare this to the current surface field strength of SGR 1900+14 inferred
from dipole spindown $B_{\rm d}=7\times 10^{14}$ G). Although far from being a
proof, this simple calculation shows the idea that some magnetars (that start
with a sufficiently large energy reservoir) can spin down substantially due to
the accumulated effect of outbursts, is plausible and merits further
investigation.

##  4 Source densities of ULPMs

The local density of sources similar to a given object, can be estimated by
the distance to that object, assuming it to be the closest member of the
group. If due to observational selection effects (for example due to beaming),
there are closer (and missed) similar objects, this analysis will overestimate
the minimum distance, $d_{\rm min}$ and thus underestimate the true number of
sources. The discussion below is divided into two parts. In the first (§4.1)
we provide an approximate, but more straightforward estimate for the minimum
source density, designed to give intuition to the reader. This analysis has
the advantage of being applicable to a general category of underlying sources
(not restricted to magnetars). In the second part (§4.2), we consider a more
elaborate and realistic model that is however, dependent on the ULPM
candidates being NSs. Moderate errors due to these types of analysis are
unavoidable due to small-number Poisson statistics. However, at a 90%
confidence limit, we can still constrain the number of sources to better than
an order of magnitude, which is constraining for the underlying physical
picture.

###  4.1 Simplified estimate - isotropic distributions

We consider the underlying distribution to be approximately isotropic (in
either 3 or 2D, see below for details) on a scale comparable to $d_{\rm min}$
(motivated by the proximity of PSR J0901–4046 and GLEAM-X J1627). The method
applied follows the analysis of Lu et al. (2021). Given a source number
density, $n$, the probability density of the volume enclosing the closest
member, $V_{1}$ is (using Poisson statistics)
$\mbox{dP}/\mbox{d}V_{1}|_{n}\\!=\\!n\exp{(-V_{1}n)}$. Using Bayes theorem
with an uninformative prior, this can be converted to
$\mbox{dP}/\mbox{d}n|_{V_{1}}\\!=\\!V_{1}\exp{(-V_{1}n)}$. Under the
assumptions, we have, by construction, $d_{\rm min}\\!=\\!(3V_{1}/4\pi)^{1/3}$
for a 3D distribution and $d_{\rm min}\\!=\\!(V_{1}/\pi)^{1/2}$ for the 2D
case. We see that $d_{\rm min}$ gives an estimate of the local source density.

Using this approach we turn to estimate the total number of ULPM candidate
sources in the Galaxy. This is attained by first calculating the fraction of
Galactic systems that reside up to a distance $d$ from the Sun. Assuming an
object population that tracks the stellar population of the Milky Way (we
address below the possibility that ULPMs are confined to the disk only), this
fraction is given by:

$f_{d}=\frac{4\pi\int_{0}^{d}dr^{\prime}\int_{0}^{\pi}d\theta^{\prime}\int_{0}^{2\pi}d\phi^{\prime}n_{*}(r,z)r^{\prime
2}}{2\pi\int_{0}^{\infty}dr\int_{0}^{\infty}dzn_{*}(r,z)r}$ (5)

where $r,z,\phi$ are polar coordinates measured from the Galactic center and
$n_{*}(r,z)$ is the (unnormalized) stellar density in the Milky Way. As a
first approximation, the latter can be described by
$n_{*}(r,z)=\exp(-r/h_{r})\exp(-|z|/h_{z})$ with $h_{r}=3.5\mbox{
kpc},h_{z}=0.25\mbox{ kpc}$ (Bahcall & Soneira, 1980) being the scale-lengths
in $r,z$ (more sophisticated distributions are discussed in §4.2). Finally,
$r^{\prime},\theta^{\prime},\phi^{\prime}$ in Eq. 5 are spherical coordinates
as measured from the Sun, i.e

$\displaystyle r^{\prime}=\left[(r\cos\phi-R_{\rm
sun})^{2}+r^{2}\sin^{2}\phi+z^{2}\right]^{1/2}$
$\displaystyle\cos\theta^{\prime}=\frac{z}{r^{\prime}}\quad\&\quad\tan\phi=\frac{r\sin\phi}{r\cos\phi-
R_{\rm sun}}$ (6)

where $R_{\rm sun}=8\mbox{ kpc}$.

With these definitions, the total number of objects in the Galaxy is on
average $\langle N\rangle=\int dn\int dVn\frac{dP}{dn}2\pi
h_{r}^{2}h_{z}n_{*}(r,z)=f_{\rm d}^{-1}$. Taking a distance of $0.4$kpc for
PSR J0901–4046 we estimate a lower limit of
$N\\!\approx\\!3620^{+9000}_{-3000}$ similar objects in the entire Galaxy
(where quoted uncertainties are 90% confidence limits, given Poisson
statistics). Similarly, with a distance estimate of $1.3$ kpc for GLEAM-X
J1627 we estimate a lower limit on the number of similar objects to be
$N\\!\approx\\!225^{+565}_{-190}$. If we take a further conservative
assumption, that all similar objects reside in the Galactic plane (as
magnetars are observed to have small offsets from the plane, Kaspi &
Beloborodov, 2017; Enoto et al., 2017), the numbers are reduced to $N_{\rm
2D}\\!\approx\\!2500^{+5450}_{-2145}$ for PSR J0901–4046 and $N_{\rm
2D}\\!\approx\\!220^{+490}_{-188}$ for GLEAM-X J1627. These lower limits are
likely to be overly conservative as the observation of magnetars being
confined to the Galactic plane could be an observational bias towards
detecting X-ray active and therefore also much younger objects. A similar
analysis can be realized with the other ULPM candidates but becomes less
constraining as the distances get larger. For example, for PSR J0250+5854, we
get $N\\!\approx\\!150^{+450}_{-135}$ ($N_{\rm
2D}\\!\approx\\!140^{+440}_{-130}$). If some or all of these sources belong to
the same category, then the estimate will be dominated by the closer object.

###  4.2 More realistic estimate - population synthesis approach

A more realistic (but model dependent) estimate of the (lower limit on the)
number of Galactic sources needed to explain a given observed $d_{\rm min}$ is
achieved by accounting for the birth rates, locations, velocities and ages of
a hypothetical population. Our analysis here builds upon the work by Faucher-
Giguère & Kaspi (2006) done in the context of the general pulsar population.
We perform a Monte Carlo calculation to derive the birth properties of the
objects. We then propagate the motion of each object according to the Galactic
potential.

The initial location of each system is assumed to follow the spiral arm
structure of the Milky Way (see Faucher-Giguère & Kaspi 2006 for more details)
in either 2 (i.e. confined to the plane) or 3 (i.e. taking into account the
disk scale height) dimensions. As ULPM candidates are most likely NSs (see
discussion in §2.2), the magnitudes of the birth kicks (relative to the
Galactic rotational velocity curve at their birth radius) are taken to follow
the two-component Gaussian favoured by the analysis of pulsar proper motions
(Arzoumanian et al., 2002; Igoshev, 2020). As we show below, the results are
weakly dependent on this assumption. The directions of the kick velocities are
arbitrary relative to the Galactic coordinates. The characteristic lifetime of
ULPM candidates is left as a free parameter, $T_{\rm age}$. The age of a
specific object, $t$, is uniformly distributed in the range
$t\\!\in\\![0,T_{\rm Age}]$. Similarly, the formation rate of systems (systems
formed in the Galaxy per unit time), $r$ is left as a free parameter. At a
given time, the expected number of Galactic systems is $N_{0}\\!=\\!rT_{\rm
age}$, and the actual number is given by a Poisson distribution with a mean of
$N_{0}$.

After initializing the birth properties, we propagate each system’s Galactic
motion for a duration $t$ with galpy***http://github.com/jobovy/galpy (Bovy,
2015) using the MWPotential2014 Galactic potential. Finally, we compare the
systems’ locations to that of the Sun and calculate $d_{\rm min}$, between the
Sun and any of the simulated systems. This process is repeated $10^{3}$ times
for each $\\{r,T_{\rm Age}\\}$. The results are shown in Figure 3. As expected
from §4.1, the results depend mostly on $N_{0}$, and are largely insensitive
to other combinations of $r,T_{\rm Age}$. We see that at a $90\%$ confidence
level, $N\\!=\\!19200_{-13600}^{+40000}$ ($N\\!=\\!590_{-500}^{+1010}$) for
PSR J0901–4046 (GLEAM-X J1627 ) if the initial distribution of systems is
taken in 3 dimensions. Confining the birth locations to the Galactic plane,
only slightly changes these estimates: $N_{\rm
2D}\\!=\\!12800_{-10100}^{+19000}$ ($N_{\rm 2D}\\!=\\!510_{-420}^{+500}$) for
PSR J0901–4046 (GLEAM-X J1627).

Finally, we note that if systems are preferentially born closer to the
Galactic center (e.g. within the Galactic bar), $d_{\rm min}$ increases
dramatically for a given $\\{r,T_{\rm Age}\\}$. The implication is that the
total required number of systems in the Galaxy will increase significantly
compared to the estimates above under such an assumption.

Figure 3: Minimum distance between the Sun and ULPM candidates, assuming the
later are born with rate $r$ and have a characteristic lifetime of $T_{\rm
Age}$. Candidates’ birth locations and velocities are randomized as described
in §4.2. The numbers shown are median estimates of $d_{\rm min}$, averaged
over many realizations with the same underlying properties. The top panel
depicts results for sources initially distributed throughout the Galaxy and
the lower panel for sources being born in the Galactic plane.

##  5 Age estimates

We turn next to show that the Galactic ULPM candidates are likely to be old
systems relative to SGRs or AXPs. This is motivated by several independent
considerations that, combined, allow us to estimate the ages of ULPM
candidates.

###  5.1 Spindown age limits

The period and period derivative of a NS provide us with the “spindown age",
$\tau_{\rm SD}=P/(2\dot{P})$. This is generally an upper limit on the true age
of the NS (and this conclusion is strengthened in the presence of magnetic
field decay, see Dall’Osso et al. 2012; Viganò et al. 2013; Beniamini et al.
2019). For GLEAM-X J1627, the extant upper limit on $\dot{P}$ leads to
$\tau_{\rm SD}>170\mbox{ kyr}$. For PSR J0901–4046, $\dot{P}$ is measured,
yielding $\tau<\tau_{\rm SD}=12\mbox{ Myr}$.

Figure 4:  Left: 3$\sigma$ upper-limit on the thermal X-ray luminosity as a
function of temperature for PSR J0901–4046, GLEAM-X J1627, and PSR J0250+5854.
PSR J0901–4046 is a factor of $\\!\sim\\!3$ and 4 closer than GLEAM-X J1627
and PSR J0250+5854, respectively. Right: Bolometric thermal luminosity versus
age for different types of thermally emitting NSs (figure adopted from Viganò
et al. 2013). The $3\sigma$ upper-limits for PSR J0901–4046 (yellow), GLEAM-X
J1627 (blue), and PSR J0250+5854 (black), assuming the entire NS surface is
thermally emitting, have been placed on the plot for comparison (see Appendix
F for details).

Figure 5: Top: Inferred active lifetime $\tau$ and source density $n$ for FRB
20180916B (aka R3, the closest of the highly active cosmological FRBs), SGR
1935+2154 type FRBs, ULPM candidates discussed in this work and XDINS. $n$ is
inferred from the distance to the nearest object of that class taking into
account Poisson statistics and 90$\%$ confidence limits. For Galactic
populations, we correct for the density of Milky Way-like galaxies taking
$n_{\rm MW}=0.01\mbox{Mpc}^{-3}$. In addition, for SGRs, the density ranges
consider whether it is only SGR 1935+2154 or all Galactic magnetars which are
capable of producing an FRB. Age estimates for the SGR-type population are
informed by Galactic magnetar evolution studies (Dall’Osso et al., 2012;
Beniamini et al., 2019), while for R3 we only take a conservative lower limit
of 10 years, informed by its steady behaviour since it has first been
observed. The figure demonstrates that the two types of FRBs correspond to
vastly diverging source birth rates. The circle sizes corresponds to
$\log_{10}(P[s])$ and shows that the different Galactic sources are unlikely
to be part of a single evolutionary trajectory. Bottom panels depict the
projections on the $n-P$ (left) and $P-\tau$ (right) planes.

###  5.2 Timing stability

PSR J0901–4046 exhibits remarkable timing stability for a magnetar, the
stability is more commensurate with middle-aged pulsars. Over the span of 7
months, PSR J0901–4046 has shown a significantly lower timing noise compared
to other slow radio pulsars (Lower et al., 2020). The timing residuals show no
signature of higher order frequency derivatives (see Caleb et al., 2022, for
more details on the timing properties), suggesting an exceedingly stable
magnetospheric state. This provides valuable information about the age of the
NS. Assuming that increasing stability with age is universal, a relative
timing stability (over 7 months) and glitch size constraint of
$\sim\Delta\nu/\nu\lesssim 10^{-6}$ implies an age $\gtrsim 100$ kyr from the
known glitch activity distribution of young rotation-powered pulsars (Espinoza
et al., 2011; Fuentes et al., 2017; Basu et al., 2022).

###  5.3 Source density age limits

A conservative limit on the local Galactic formation rate is given by the CCSN
rate, $\dot{N}<\dot{N}_{\rm ccSNe}\approx 17\pm 2\mbox{ kyr}^{-1}$. Given this
birth rate and our estimates for $N$, the number of ULPM candidates in the
galaxy (§4), we can estimate the typical lifetime of a given object as
$\tau\approx N/\dot{N}\gtrsim 1.2_{-0.9}^{+0.8}\mbox{ Myr}$ for PSR J0901–4046
($\tau\approx N/\dot{N}\gtrsim 38_{-28}^{+44}\mbox{ kyr}$ for GLEAM-X J1627).
These numbers are slightly modified if the systems are only born in the plane,
in which case we find $\tau\approx N/\dot{N}\gtrsim 0.8_{-0.6}^{+2.7}\mbox{
Myr}$ for PSR J0901–4046 ($\tau\approx N/\dot{N}\gtrsim 35_{-27}^{+33}\mbox{
kyr}$ for GLEAM-X J1627).

###  5.4 Cooling age limits

An independent limit on the age comes from the X-ray flux limits which can be
translated to effective cooling luminosity and surface temperature upper
limits (see Appendix F for details).

Relating the luminosity and temperature to the age of the NS, depends on
various uncertain physical inputs, such as the NS mass, equation of state, the
existence or lack of superfluidity and superconducitivity, the composition of
the core and the envelope and the magnetic field structure throughout the
star. Nonetheless, a data-driven comparison can be made with other NSs with
measured thermal luminosities and age estimates. Such a comparison is shown in
Figure 4. The left panel shows the derived X-ray upper-limits on the
bolometric thermal luminosity of the three long period pulsars. The right
panels shows how these pulsars’ upper-limits on the bolometric thermal
luminosities compare with other pulsars of various classes (e.g., magnetars,
RPPs, XINs, CCOs). Based on the cooling curves (adopted from Viganò et al.
2013), we find that GLEAM-X J1627 is likely to have an age
$\tau\\!\gtrsim\\!10^{5}$ yr. The upper-limits on the thermal luminosities of
PSR J0901–4046 and PSR J0250+5854 are consistent with their inferred
characteristic ages, but they could be as young as $\tau\\!\approx\\!1$ Myr
and still remain undetected in X-rays.

###  5.5 Proper motion

The most nearby supernova remnants to GLEAM-X J1627 are RCW 103, Kes 32, and
G332.0+00.2 (Green, 2019). The first hosts a known NS remnant (§2.3). The
other two have angular offsets of $2.72^{\circ}$ and $2.84^{\circ}$,
respectively. Considering the distance to GLEAM-X J1627, this corresponds to a
physical offset of $\\!\sim\\!60$ pc. The Earth’s proximity to GLEAM-X J1627
implies that a nearby SNR would have remained detectable unless it was
extremely old. Combining this with the fastest space velocities of observed
pulsars ($\sim 2000$ km s-1), we conservatively estimate a minimum kinematic
age of $\gtrsim 3\times 10^{4}$ yr. Indeed $60$ pc is also the approximate
offset between GLEAM-X J1627 and the Galactic plane (see Table 1). Another
possibility is that GLEAM-X J1627 formed in a young (i.e., $<10$ Myr) stellar
cluster. The nearest cluster to GLEAM-X J1627 with a compatible distance and
age is DBSB 154 (Kharchenko et al., 2013), which has an estimated distance of
$d\approx 2.1$ kpc and an age of 5 Myr. GLEAM-X J1627 is offset from this
cluster by a distance of $2.84^{\circ}$, leading to a similar lower-limit on
its kinematic age***This presumes the cluster is at the GLEAM-X J1627
distance, and not at 2.1 kpc. If instead GLEAM-X J1627 is at 2.1 kpc, the
lower-limit on the kinematic age would be $\gtrsim 5\times 10^{4}$ yr..
Therefore, requiring that the NS formed in the plane results in a commensurate
lower limit on the age of the object.

We have also carried out a similar study of the field surrounding PSR
J0901–4046. The closest SNR is G261.9+05.5 (2.14∘ offset) but it is at a
distance of 2.9 kpc (Green, 2019), which is likely too far to be compatible
with PSR J0901–4046’s DM distance. The next closest SNRs (G266.2-01.2 and the
Vela complex) all have known compact objects associated with them. G272.2-03.2
is the next closest SNR, which is 11.4∘ offset from PSR J0901–4046. This
corresponds to a physical offset of about 80 pc and to a limit on the
kinematic age of $\gtrsim 40$ kyr. The closest stellar cluster with a
compatible distance estimate is Trumpler 10, which lies at 420 pc and has an
age of about 24 Myr (Kharchenko et al., 2013). This cluster is $\approx$3.1∘
offset from PSR J0901–4046. This leads to a fairly unconstraining lower-limit
on the kinematic age of $>11$ kyr.

###  5.6 Combined age limits and implications

These different age estimates suggest that systems like PSR J0901–4046 and
GLEAM-X J1627 are rather old, perhaps $\sim\\!10^{5.5}-10^{6}$ yr or so.
Although there are uncertainties involved in each of the limits discussed
above, they are independent from each other, lending overall credence to the
old age of these systems compared to the known magnetar population. We
consider the ULPM candidates’ location on the $n\\!-\\!\tau$ diagram shown in
the top panel of Fig. 5, where we have also supplemented an additional class
of objects, the XDINS, a nearby ($\sim$ 150 – 500 pc) class of cooling NSs
whose origin is still debated but likely related to standard magnetars
(relevant details of these sources are given in Table 1). The tentative
alignment of the confirmed magnetars with the RCW103 source, GLEAM-X J1627,
PSR J0901–4046 and the XDINS along a roughly constant value of the volumetric
rate density, $\dot{n}$, suggests that these objects arise from similar
formation channels, commensurate with common CCSN. At the same time, their
vastly differing periods (changing non-monotonically along this supposed
sequence), imply that these represent different evolutionary sequences, with
similar occurrence rates. It is noteworthy that FRB 20180916B (R3) is a clear
outlier from these objects (see discussion in §7) suggesting it probes a much
rarer population altogether, involves fine-tuning in physical conditions
necessary for FRB emission, or strong observational selection effects against
R3-like objects. Finally, we note that loosely $n\propto P^{-1}$ (even
ignoring R3) and that the SGR source density is possibly an underestimate, as
the sample is likely not as complete as that of XDINS. This will become
important in §7, when we discuss the possible connection with periodic FRBs.

##  6 Magnetic field and spin evolution - Different decay to standard
magnetars

###  6.1 Phenomenological constraints on field evolution

The dipole surface field evolution of confirmed Galactic magnetars has been
characterized through $\dot{B}_{\rm d}\\!\propto\\!B_{\rm d}^{1+\alpha}$
(Colpi et al., 2000). It follows that a magnetar’s age $t$ is related to the
fractional decrease in the dipole magnetic field, $B_{\rm d}/B_{\rm d,0}$ as

$\frac{t}{\tau_{\rm d,0}}=\frac{1}{\alpha}\left(\left(\frac{B_{\rm d}}{B_{\rm
d,0}}\right)^{-\alpha}-1\right).$ (7)

where $\tau_{\rm d,0}$ is a characteristic time on which the dipole field
decays by a significant amount. We see that a lower limit on $B_{\rm d}/B_{\rm
d,0}$ and $t/\tau_{\rm d,0}$ translates to a lower limit on $\alpha$ (see also
Fig. 6). Comparison to Galactic magnetar properties suggests typical initial
decay times of the dipole field of $\tau_{\rm d,0}\sim 10^{3}-10^{4}\mbox{
yr}$ for an initial surface dipole field of $B_{\rm
d,0}~{}\lesssim~{}10^{15}\mbox{ G}$. Motivated by our results in §5, we
consider the age of PSR J0901–4046 to be approximately $1$ Myr. If, in
addition, we take the surface dipole field as estimated from dipole spindown
$B_{\rm d}\sim 2.6\times 10^{14}$ G, Eq. (7) leads to an effective
$\alpha\gtrsim 4.5$ ($\alpha\gtrsim 6.5$) for $t/\tau_{\rm d,0}\gtrsim 100$
($t/\tau_{\rm d,0}\gtrsim 1000$). Such values are inconsistent with the known
Galactic magnetar population that generally imply $-1\lesssim\alpha\lesssim 1$
for $\tau_{\rm d,0}\sim 10^{4}\mbox{ yr}$ (Beniamini et al., 2019) or
$1\lesssim\alpha\lesssim 2$ for $\tau_{\rm d,0}\sim 10^{3}\mbox{ yr}$
(Dall’Osso et al., 2012). They are also inconsistent with certain theoretical
predictions (e.g. $\alpha_{\rm int}=6/5$ for the solenoidal mode of ambipolar
diffusion, Dall’Osso et al. 2012). This, along with the arguments given in §5,
suggests a different field decay evolution for the ULPM candidates than seen
in the confirmed Galactic magnetar population. We note that the ULPM field
decay mechanism may be episodic rather than continuous, in which case the
effective $\alpha$ above describes only an averaged property of this mechanism
on long timescales.

A different possibility is that ULPM candidates have similar $\alpha$ to
confirmed magnetars, but begin their life with much larger values of $B_{\rm
d,0}$. However, for either $\alpha<1,t/\tau_{\rm d,0}\gtrsim 100$ or
$\alpha<2,t/\tau_{\rm d,0}\gtrsim 1000$, one requires $B_{\rm d,0}\gtrsim
10^{16}\mbox{ G}$. Such large dipole fields (and likely even larger internal
fields that accompany them) are difficult to obtain from standard channels.
This problem becomes more severe considering the large implied formation rates
of ULPM candidates (see Figure 5), comparable to those of confirmed magnetars.

Finally, if the current $\dot{P}$ of the ULPM candidates is not dominated by
dipole spindown, then their surface field may be overestimated by
$\\{P,\dot{P}\\}$. As demonstrated in §2.1.1, even if spindown is dominated by
a particle wind, we still demand $B_{\rm d}\gtrsim 10^{13}\mbox{ G}$. Such a
reduced field strength leads to $\alpha>1$ ($\alpha>1.6$) for $t/\tau_{\rm
d,0}\gtrsim 100$ ($t/\tau_{\rm d,0}\gtrsim 1000$) which is marginally
consistent with the confirmed magnetar population, but requires a distinct
spin evolution between ULPM candidates and confirmed magnetars. Furthermore,
the possibility that the current rate of change in $P$ is dominated by a wind
could be tested with future observations by constraining the braking index of
the ULPM candidates, which should approach in this case $n=1$ rather than the
dipole value, $n=3$. We conclude that the ULPM candidates require a distinct
magnetic and / or spin decay mechanisms to standard Galactic magnetars.

###  6.2 Physical origins of long-lived strong fields in ULPMs and powering
transient radio emission

The large phenomenological effective values of $\alpha$ above, and the high
source density of PSR J0901–4036 like objects implies the physics of field
evolution in the crust and core, substantially different from standard
magnetars, ought to be generically realized. Apparently magnetar-like fields,
with dipolar components of $\gtrsim 10^{13}-10^{14}$ G, are commonly produced
and survive for Myr or longer with limited observational consequences beyond
the radio.

For young magnetars, their field evolution is compatible with nonlinear Hall
evolution (in a resistive electron MHD treatment) in their crusts (e.g.,
Goldreich & Reisenegger, 1992; Cumming et al., 2004; Gourgouliatos et al.,
2013; Viganò et al., 2013; Pons & Viganò, 2019; Gourgouliatos et al., 2022) on
a timescale on about $10^{4}$ yr (with significant dependence on crust
thickness), about 2–3 orders of magnitude shorter than Ohmic evolution from
electron scattering. The tentative evidence in favor of a crust in PSR
J0901–4046 (Appendix C) and the high source density disfavors nonstandard
crust scenarios (e.g. quark stars, or very light or massive NSs) for this
class of object. Crust-threading fields (which may possess significant
toroidal components) are thought to power the energetic high-energy activity
of conventional magnetars (e.g., Lander, 2022). The existence of non-dipolar
components is exemplified by the “low-field" magnetar SGR 0418+5729 (with
$P=9.1$ s) whose dipolar field is only $B_{p}\\!\sim\\!\\!1.2\\!\times
10^{13}$ G (Rea et al., 2013); yet evidence points to much stronger local
field strengths up to $10^{14}\\!-\\!10^{15}$ G (Tiengo et al., 2013) as
expected from 3D crustal Hall evolution (e.g., Wood & Hollerbach, 2015;
Gourgouliatos et al., 2016). For electromagnetic spin evolution in quiescence,
especially on long timescales, the dipolar component is most germane, as that
is most influential at the light cylinder.

The X-ray limits on PSR J0901–4046 and GLEAM-X J1627 (see Figure 4) suggest
low magnetic dissipation within the crust, commensurate with the crustal Ohmic
timescale (Cumming et al., 2004) or even lower. A rough guide on the relevant
magnetic dissipation timescale in the current epoch is
$B_{14}^{2}R_{6}^{3}/\tau_{\rm decay}\\!\sim\\!L_{X}\ll 10^{32}$ erg s-1
implying $\tau_{\rm decay}\\!\gg\\!1$ Myr. This is generally beyond the
standard Hall evolution timescale for typical densities of the inner crust,
and suggests a different physical mechanism currently reigning for field
evolution. The crustal nonlinear Hall evolution may be substantially retarded
with “Hall attractor" solutions where the Hall evolution saturates to a level
commensurate with Ohmic evolution (Gourgouliatos & Cumming, 2014a, b; Wood &
Hollerbach, 2015; Gourgouliatos et al., 2016). Such attractor solutions are
predominantly independent of initial conditions. The dipolar component in this
scenario may actually increases with time, implying unconventional tracks in
the $P-\dot{P}$ diagram (see also Popov et al., 2017). However, a major caveat
is that these models assume a crust-dominant field, compatible with Meissner
type-I superconductor for the outer core where the field is entirely expelled.
This core boundary condition is a strong assumption and unrealistic for common
CCSN formation channels (see below). Another assumption in these models is a
time-independent specification for the crust conductivity, rather than coupled
magnetic and thermal evolution. The assumption of a crust dominant field is
also likely inapplicable for old systems which manifest strong global dipolar
components exceeding $10^{13}$ G (see 2.1.1). This perhaps points towards
strong core-bound fields in PSR J0901–4046, GLEAM-X J1627 and similar objects.

Proton superconductivity is expected to set in very soon after the formation
of NSs by CCSN (Baym et al., 1969a; Yakovlev & Pethick, 2004), as the core
temperature drops from neutrino cooling (by either modified URCA or direct
URCA processes, Shapiro & Teukolsky, 1983) below the critical proton
superconducting temperature. In general, two approximate critical fields
exist, $H_{c1}\sim 10^{15}$ G and $H_{c2}\\!\sim\\!10^{16}$ G (which depend on
depth and local properties) where for the magnetic field $H\\!<\\!H_{c1}$ the
protons are in a type-I Meissner state, while for $H\\!>\\!H_{c2}$
superconductivity is destroyed. Numerous works have considered the evolution
of core fields in NS interiors in the context of neutron superfluidity and
proton superconductivity (e.g., Muslimov & Tsygan, 1985; Jones, 2006;
Glampedakis et al., 2011a, b; Henriksson & Wasserman, 2013; Lander, 2013,
2014; Palapanidis et al., 2015; Graber et al., 2015; Sinha & Sedrakian, 2015;
Elfritz et al., 2016; Gusakov & Dommes, 2016; Gusakov, 2016; Passamonti et
al., 2017a, b; Ho et al., 2017; Dommes & Gusakov, 2017; Graber, 2017; Haskell
& Sedrakian, 2018; Sur & Haskell, 2021; Wood & Graber, 2022) following the
pioneering work of Baym et al. (1969a, b). The physics in this multicomponent
and multiscale system is a poorly understood domain of NSs, and still debated.
The initial configuration of the field at formation also influences its
evolution at later times, both in the crust and the core, and there can be
strong differences based on the equation of state, constituents in the core,
and mass of the NS. Among complications are how superfluid vortices interact
with quantized flux tubes (or “fluxoids"), vortex pinning (e.g., Link, 2012a,
b) and entrainment, formation of magnetic domains, and the rate of diffusion
of flux tubes out of the core into the crust. For a recent discussion of
$B\\!\neq\\!H$ in the core, see Rau & Wasserman (2021). Yet, the timescales
for core field evolution are much longer than the crust (e.g., Glampedakis et
al., 2011b; Graber et al., 2015) – strong fields may be retained in the core
for $\gg$Myr in superconducting cores. Diffusion of fields out of the core
into the crust is a slow process. Proton superconductivity in the core is also
expected in all magnetars soon after formation (e.g., Ho et al., 2012, 2017).
Here, the typical diffusive timescale near the outer core is $\sim 0.1-0.3$ km
Myr-1 (e.g., Graber et al., 2015; Elfritz et al., 2016; Ho et al., 2017), or
equivalently, $\tau_{\rm sc}\gtrsim 10\,\ell_{\rm km}^{2}$ Myr where $\ell$ is
a characteristic lengthscale for evolution. As the NS cools, increasingly
deeper regions of the core become superconducting with the outer core type-II
proton superconductor. If the core field is $H\\!<\\!H_{c1}$, a metastable
arrangement (where fluxoids and fields are frozen within interspersed zones of
superconducting and normal matter) is likely as the growth rate of
superconducting zones (from neutrino cooling) is much faster than the flux
expulsion timescale (e.g., Ho et al., 2017). Deeper in the core, where the
Ginsburg-Landau parameter is small, a type-I Meissner state is likely
realized.

The field diffusion out of the outer core into the crust, and its timescale,
may eventually power transient radio emission by stressing by the crust
(possibly beyond the breaking strain for high fields, e.g., Lander et al.,
2015), whose stored energy may be released abruptly (as FRBs) or gradually
with low field external twists (both scenarios requiring currents, particle
accelerate, pair production and concomitant radio emission). Indeed, much of
the dissipation can occur at the core-crust interface (Lander, 2014; Elfritz
et al., 2016) where current sheets may develop. At late times, i.e., $\gtrsim
10^{5}-10^{6}$ yr, weak toroidal components may also develop (e.g., Elfritz et
al., 2016) that may stress the crust and twist external fields.

The metastable superconducting arrangement is stable insofar as the neutrino
cooling timescale is shorter than the flux expulsion timescale. Yet, at around
a Myr, the cooling timescale begins to become commensurate with the flux
expulsion timescale (Ho et al., 2017). The timescale associated with flux
migration from vortex-fluxoid pinning is commensurate with the spin-down
timescale (Chau et al., 1992), i.e. also Myr for PSR 0901–4046, and similar in
value to $\tau_{\rm sc}$. This late-time diffusion and evolution is likely
most germane for radio activity from some Galactic ULPMs or the M81 FRB (FRB
20200120E). Flux expulsion and vortex pinning in NS cores has actually
previously been invoked in deprecated models of GRBs, e.g. (Ruderman, 1991a,
b) that may have renewed qualitative applicability to late-time radio activity
from ULPMs and FRBs.

In all of these physical scenarios, where the fields are in the crust only or
both core and crust, the field evolution generally results in strongly
dominant poloidal fields after $\sim 100$ kyr. For the external magnetosphere,
this is a low-twist state with low persistent charge density – and as such is
highly susceptible to single photon pair cascades even from weak crust
disturbances. Such magnetospheric states are generically suited for generation
and escape of transient coherent radio emission and FRBs (Wadiasingh &
Timokhin, 2019).

Figure 6: Magnetar’s $B$-field strength normalised by its initial value,
$B/B_{0}$, versus its true age $t$ normalized by the $B$-field’s initial decay
time, $t/\tau_{d,0}$ (different curves correspond to $\alpha=5,4,3,2,1,0$ from
top to bottom).

##  7 Are periodic FRBs related to Galactic ULPMs?

Long activity window periods were reported in the two most prolific repeaters,
FRB 20180916B and FRB 20121102A (with periods of $\sim 16.3$ and $\sim 160$
days, respectively) (Chime/Frb Collaboration et al., 2020; Rajwade et al.,
2020a). Beniamini et al. (2020) have suggested that these could be a
manifestation of the extremely slow spin periods achieved by a rare sub-class
of magnetars. These authors have shown that the natural predictions of the two
main competing mechanisms – binary companion or precession – are both
disfavoured by the data. In particular, both mechanisms would require an
underlying shorter period associated with the magnetar spin. However, a
continuous $\sim\\!5\,$hr observation of 20121102A revealed no such
periodicity in burst arrival time, suggesting either a wide beaming cone
(relevant only to the binary model, a precession model necessarily implies
beaming) or else a long spin period of $P\\!\gtrsim\\!5\,$hr. Likewise, FAST
observations of 1652 bursts from FRB 20121102A spanning over a month (Li et
al., 2021) find no periodicity up to at least $10^{3}$ s. Tendulkar et al.
(2021) also have ruled out any stellar companion brighter than O6V for FRB
20180916B, and favor an older ${\cal O}$(Myr) age if its birth site is the
nearest observed clump of star formation. Other evidence disfavouring a binary
companion include the observed lack of an active phase fraction that increases
with period (Pleunis et al. 2021a; this is expected since absorption of radio
waves decreases with increasing frequency) and the lack of change in DM
between bursts in the beginning and the center of the active phase. Moreover,
all observed Galactic magnetars have neither a confirmed binary companion nor
do they show robust signatures of precession on similar timescales.

As FRB 20180916B is much closer to us than FRB 20121102A, it is natural to
expect that it represents a significantly more common type of object, which,
at this stage, dominates the source density estimates of periodic FRBs (Lu et
al., 2021). The source densities estimated for ULPM candidates are at least
one order of magnitude greater than estimated for the Galactic FRB (see §4),
two orders of magnitude greater than estimated for the M81 FRB and seven(!)
orders of magnitude greater than estimated for FRB 20180916B (Lu et al.,
2021). This source density discrepancy may be at least partially resolved if
ULPM sources very rarely produce any FRBs (such as due to peculiar necessary
physical conditions or triggers) and / or the beaming of their FRB emission is
much narrower than that of the radio pulses of the Galactic ULPM candidates.
It is interesting to note that FRB 20200120E (associated with a globular
cluster in M81) represents a class of FRB sources with significantly larger
source densities than either FRB 20180916B or FRB 20121102A. Remarkably,
recent observations of FRB 20200120E reveal a lack of underlying periodicity
with values comparable to those seen in confirmed magnetars or pulsars. At the
same time, while the existence (or lack of) of a much longer active periodic
window (like in FRB 20180916B or FRB 20121102A) cannot yet be established,
there is a hint of periodic activity window with a period of 12.5 days (Nimmo
et al., 2022). If this connection is eventually established, it suggests that
the discrepancy between the ULPM candidate and periodic FRB source densities
are much reduced.

To explore the possibility that FRB beaming accounts for the vast source
density gap, we consider a simplified but generic beaming model. Our basic
assumption is that the periodicity of FRBs is associated with their spin
periods. To allow for dark periods (where no bursts are seen from FRB
20180916B), the direction of burst emission (denoted by $\chi$) must be
removed from the rotation axis by at least $\theta_{\rm b}$ (where
$\theta_{\rm b}$ is the beaming or collimation angle). i.e. $\chi>\theta_{\rm
b}$. Considering that the ‘active phase’ covers $\sim\\!1/4\\!-\\!1/2$ of the
burst period, the region of the sky from which bursts can be produced, can be
very large compared to the emission collimation. The large active fraction
suggests that the size of the region is itself is of the same order as $\chi$.
These considerations are schematically depicted in Figure 7. We denote the
time averaged isotropic equivalent luminosity of bursts from FRB 20180916B by
$L_{\rm iso}$. If bursts are beamed, the true energy of each burst is lower by
$\sim\theta_{\rm b}^{2}/4$. At the same time, a fraction $\pi\chi/\theta_{\rm
b}$ of all bursts produced by the active region are beamed away from us and
not observed. Putting all this together, the energy output of FRB 20180916B is

$E\sim L_{\rm iso}\tau\epsilon_{\rm r}^{-1}\theta_{\rm b}\chi\approx
10^{45}\tau_{10\rm yr}\epsilon_{\rm r,-3}^{-1}\theta_{\rm b}\chi\mbox{ erg}$
(8)

where $\tau$ is the active FRB producing lifetime and $\epsilon_{\rm r}$ is
the intrinsic radio efficiency. We see that a smaller beaming cone relaxes the
energetic requirements on FRB 20180916B. Therefore, energetic requirements
alone do not pose a problem for the $\chi,\theta_{\rm b}\ll 1$ scenario.

Increasing the source density of FRB 20180916B like sources due to beaming,
would require $\chi\\!\ll\\!1$. This would imply that for each FRB
20180916B-like source observed, there are $\sim\chi^{-1}$ sources whose
emitting regions are always beamed away from us and are therefore
missed***Using Eq. (8), we see that the total energy output from FRB
20180916B-like sources is independent of $\chi$.. That is, the true source
density of FRB 20180916B-like sources could be larger by $\chi^{-1}$ compared
to the one mentioned above.

The beaming model discussed above requires in particular
$\chi\\!>\\!\theta_{\rm b}$. In FRB models where the emission takes place far
from the NS surface (outside the light-cylinder), the geometric model proposed
above, works only if both $\theta_{b}\\!\ll\\!1$ and the Lorentz factor of the
emitting material is $\Gamma\\!\geq\\!\theta_{\rm b}^{-1}$. While some beaming
can be accounted for in this way, these models generically find much larger
values of $\theta_{\rm b}$ and smaller values of $\Gamma$ than needed to
explain the Galactic ULPM candidate and periodic FRB source density
discrepancy.

An alternative set of FRB models, consist of radiation from regions close to
the NS surface. In this case, beaming can be the result of radio waves having
to propagate along field lines***This is because in a homogeneous plasma, the
low-altitude waves are in fact the normal modes (e.g., Melrose, 1979; Arons &
Barnard, 1986; Rafat et al., 2019) of the ambient upstream plasma (these are
necessarily superluminal O-modes to avoid Landau damping or nonlinear wave-
particle interactions) and are ducted or refracted along field lines analogous
to an optical fiber, decoupling at higher altitudes transforming into
electromagnetic vacuum modes.. One may invoke the extremely small polar cap
region to explain the small active region size. Indeed, for a period of 16
days, and dipole spin down, the polar cap opening angle is tiny $\theta_{\rm
c}\sim 10^{-5}$. Empirically, the radio in known pulsars appears to decouple
around $10-100$ NS radii (see Appendix B), which suggests an effective beaming
angle that is at least $3-10\theta_{\rm c}$ (due to field line flaring with a
dipolar structure). Furthermore, for such long periods as discussed above, the
open field line region is very easily increased much beyond this value due to
outflows (Beniamini et al., 2020) or due to any small perturbance to the
dipole field line structure. In addition, intrinsically nonneglibile surface
area is expected to participate in FRBs and magnetar oscillations (e.g.,
Wadiasingh et al., 2020; Wadiasingh & Chirenti, 2020), particularly quasi-
polar closed zones as observed in SGR 1935+2154 (Younes et al., 2020, 2021).
The field line flaring in closed zones away from the pole can be more extreme
compared to the flaring in the near polar cap region (small angle
approximation), substantially increasing the area of the sky where the radio
could decouple from the field and escape. Finally, if we were to associate
observed FRB periodicities with $\chi\ll 1$, it would require a physical
mechanism that causes small $\chi$ to be associated with a higher activity
rate, while at the same time preventing $\chi=0$, which would wipe out the
observed periodicity in the FRB signal. This circumstance is perhaps unnatural
and requires fine-tuning, unless the mechanism which removes angular momentum
to produce a ULPM also leads to small but non-vanishing magnetic obliquity.

Figure 7: Schematic of the generic beaming model discussed in §7. The region
from which bursts can be observed is separated by $\chi$ from the principal
rotation axis $\boldsymbol{\Omega}$ and has a typical angular size of the same
order (to account for $\sim 0.25-0.5$ active phase observed in periodic FRBs).
The beaming of each particular burst has an angle $\theta_{\rm b}$, which may
be much smaller than $\chi$.

In the Timokhin-Arons model for pulsar radio emission (Timokhin, 2010;
Timokhin & Arons, 2013), a crude estimate of $\theta_{\rm b}$ resulting from
an infinitesimally narrow active region, is the angle of the pair formation
front. This quantity couples the characteristic gap and field radius of
curvature length scales. This field curvature is essential in the Timokhin-
Arons model for the generated waves to be electromagnetic (rather than
electrostatic) in character (Philippov et al., 2020; Tolman et al., 2022). The
collimation may be much wider than this estimate if large surface areas
(involved in global magnetar oscillations) is tied to flux tubes where pair
cascades occur, or the field lines flare significantly where the waves
decouple. Thus, the pair formation front angle estimate constitutes a
conservative lower limit on collimation. For at-threshold pair cascades, the
angle is approximately $\vartheta\sim h_{\rm gap}/\rho_{c}$ where $h_{\rm
gap}$ is defined in Eq. (30) of Wadiasingh et al. (2020). This yields,

$\vartheta\sim
10^{-3}\,\rho_{c,7}^{-5/7}B_{14}^{-3/7}\nu_{2}^{-3/7}\xi_{2}^{-3/7}\lambda_{5}^{3/7}$
(9)

which may be significantly wider for smaller amplitudes $\xi$ viable for FRBs
in charge-starved ULPMs (Beniamini et al., 2020). Thus, we conclude that
$\theta_{b}\gg\vartheta$, while we generically require $\chi\geq\theta_{\rm
b}$ in our toy model. Thus beaming is likely insufficient to explain the much
lower source density of FRB 20180916B over ULPMs, and could suggest the FRB
mechanism may require special or contrived conditions (e.g., particular
magnetosphere states or excitation of specific global oscillation modes,
Wadiasingh & Timokhin, 2019) to transpire frequently in a rare subset of
ULPMs. Such peculiarity is already suggested by SGR 1935+2154, where only one
of its thousands of X-ray short bursts in April 2020 resulted in a
contemporaneous bright radio burst (Younes et al., 2020, 2021).

Finally, the observed period of GLEAM-X J1627 is three orders of magnitude
smaller than that of FRB 20180916B. It is plausible that the number of sources
decreases rapidly enough with period to accommodate their respective source
densities, see bottom left panel of Figure 5. This scenario will be
constrained with future discoveries of objects with long periods (§8).

##  8 Observing Strategy

Figure 8: Left Panel: Fraction of S/N recovered after boxcar convolution
during a single pulse search as a function of boxcar-trails for the measured
duty cycles and periods of PSR J0901–4046 and GLEAM-X J1627. The pink and
purple shaded regions depict the space covered by the boxcar trials used in
current surveys for PSR J0901–4046 and GLEAM-X J1627, respectively. Right
Panel: Fraction of recovered S/N as a function of spin period for fixed pulse
duty cycles (taken to be equal to those of PSR J0901–4046 and GLEAM-X J1627).

###  8.1 Radio

The recent discovery of ULPM candidates like GLEAM-X J1627 and PSR J0901–4046
has demonstrated that there might be large population of ultra-long period
radio-loud NSs that has remained hidden to date. This dearth of radio-loud
ULPMs can mostly be attributed to a selection bias of all current time-domain
searches for radio transients. Typically, a time-domain survey for pulsars and
other transients spends 10–20 minutes on any given part of the sky (for e.g.
see Manchester et al., 2001). This limited amount of integration time
immediately precludes detection of long-term periodic sources similar to
GLEAM-X J1627. Furthermore, current real-time single pulse search pipelines at
various facilities are strongly biased against single pulses with large widths
($>$few hundred ms). Single radio pulses from ultra-long period NSs can be
fairly wide depending on the duty cycle. For example, a NS with a period of 5
minutes and 1$\%$ duty cycle can have a pulse as wide as 3 seconds, which
would be completely missed by single pulse search pipelines. To demonstrate
this, we estimate the reduction in the recovered signal to noise ratio (S/N)
of a single pulse after the convolution with a series of boxcars with varying
widths. Following the analytical treatment from Morello et al. (2022), for a
Gaussian pulse with a width $\sigma$ (we ignore any sub-pulse structure in the
pulse to simplify our computations) and a boxcar width $w$, both measured in
units of pulse phase (going from 0 to 1), the efficiency is

$\eta=\left(\frac{4\pi\sigma^{2}}{w^{2}}\right)^{1/4}~{}{\rm
erf}\left(\frac{w}{2\sigma\sqrt{2}}\right),$ (10)

where $\rm erf(x)$ is the error function. We compute the efficiency for
various values of boxcar width trials for the duty cycle of PSR J0901–4046 and
GLEAM-X J1627, and show this in the left panel of Figure 8. For a pulse with a
width of 500 ms or more, the maximum boxcar trial will only cover $\lesssim
30\%$ of the total width. Thus, the recovered S/N in the best case scenario
will be only as high as 40–50$\%$ of the true S/N of the pulse for any wide
single pulses from ULPMs. Furthermore, the right panel of Figure 8 shows the
recovered S/N as a function of period for the duty cycles of PSR J0901–4046
and GLEAM-X J1627, clearly indicating that the efficiency of current surveys
reduces very rapidly for longer periods and that they are insensitive to the
majority of the ULPM population in the Galaxy.

Figure 9: Distance horizon for putative ULPMs with the same radio luminosity
as PSR J0901–4046 for SKA-MID and DSA-2000 time domain surveys as a function
of spin period.

Therefore, it is important to modify observing strategies or develop new
techniques in order to search for these sources in the future. One way of
overcoming the bias is the use of Phased Array Feeds (PAF) on single dish
radio telescopes. PAFs can increase the effective field of view of the
telescope manifold which can enable larger dwell times for surveys, thus
becoming more sensitive to long period sources. On the other hand, performing
searches in radio images for these sources is an effective way to combat the
selection effects. Image domain surveys have large dwell times (few hours)
with large integration times (a few seconds). Since the width of pulses from
ULPMs is expected to be large, radio images can be used effectively to search
for radio-loud ULPMs over a large parameter space without hindering the
sensitivity of the telescope due to the large integration times.

This provides a framework for predicting a horizon for radio-loud ULPMs
assuming they have a radio luminosity comparable to PSR J0901–4046. To
illustrate this, we assume a time-domain radio survey with the SKA-MID
(Dewdney et al., 2017) and the DSA-2000 (Hallinan et al., 2019) array. We
assume an integration time of 3600 seconds per pointing and use the parameters
for the telescope (for e.g. gain, system temperature, bandwidth) as provided
in the literature (Keane et al., 2015; Hallinan et al., 2019, Lina Levin,
priv. comm.). Using the radiometer equation, for a signal to noise ratio
threshold of 8, we obtain a sensitivity of 5.7 $\mu$Jy for the duty cycle of
PSR J0901–4046. We use this sensitivity to compute the distance horizon for
detecting ULPMs with properties similar to PSR J0901–4046 as a function of
spin period. Figure 9 shows the result of such a simulation. We can in turn
use these estimates to put a limit on the detection rate of ULPMs similar to
PSR J0901–4046. Using the estimates in §4.2 for the number of ULPMs of this
type in the entire Galaxy, we estimate that $N_{\rm SKA}=21_{-5}^{+4}$ similar
objects will be detectable by all-sky surveys with SKA-MID and $N_{\rm
DSA}=130_{-13}^{+11}$ with DSA-2000. These numbers, are only slightly
decreased in case the initial birth locations of ULPMs are confined to the
Galactic plane: $N_{\rm SKA2}=16_{-4}^{+4},N_{\rm DSA2}=89_{-10}^{+10}$.

###  8.2 X-ray

The detection of a persistent high energy counterpart to any ULPM candidate
would provide a crucial step forward to understanding their origin and
evolutionary path. As discussed in Section 5, older ULPMs such as PSR
J0901–4046 are expected to be $\gtrsim\\!10^{5}$ years old. Under any B-field
decay scenario, their surface temperature is not expected to exceed the $100$
eV (e.g., Potekhin et al., 2020) with thermal bolometric luminosities
comparable to those of the XDIN population, $L_{X}\lesssim 10^{31}$ erg s-1
(Kaplan, 2008; Kaplan & van Kerkwijk, 2009). The faintness of these sources,
coupled with their cool surface temperatures poses a challenge for their
detection. For instance, an absorption column density $N_{\rm H}$ larger than
$10^{22}$ cm-2 will severely attenuate the intrinsic flux below $\sim 2$ keV.
From He et al. (2013), this limit is reached at around 4 kpc, implying that
the current fleet of X-ray telescopes are practically insensitive to such
population beyond this horizon. From the current population of ULPM
candidates, PSR J0901–4046 is the prime source for deep X-ray counterpart
searches due to its DM-implied distance of 330 pc and $N_{\rm H}\lesssim
10^{21}$ cm-2; both are the smallest among ULPMs and comparable to the XDIN
population. The latter are easily detectable with currently operating X-ray
telescopes, e.g., XMM and NICER.

Prospects of independently discovering ULPM candidates in X-rays must rely on
large field-of-view imaging soft X-ray instruments, such as eROSITA. From
extensive simulations of the Galactic NS population and the eROSITA all sky
survey scans, Pires et al. (2017) derived a detection estimate of about 95
XDINs within a distance of 1.7 kpc. Such sources could then be targeted for
deeper exposures with XMM-Newton and/or NICER to search for slow X-ray
pulsations to identify their nature. Given our source density estimate of
ULPMs and XDINs (Figure 5), we surmise that $\\!\sim\\!10$ ULPMs may hide
within the total eROSITA XDIN population. Another potential route to uncover
ULPMs is through relatively deep (reaching, e.g., $\\!<\\!10^{-15}$ erg s cm-2
observed flux in the 0.2-2 keV band) observations of large swatch of the
Galactic plane. This is achievable with sensitive, large FOV ($\\!>\\!1$
degree2) soft X-ray telescopes, e.g., STAR-X (Zhang et al., 2022) with larger
grasp than any current mission. In such cases, search is more sensitive to
farther distances, balancing the smaller FOV of such instruments compared to
eROSITA.

##  9 conclusions

The recent discovery of several intriguing Galactic radio sources with long
periods compared to the pulsar population reveals a potentially large and
unexplored class of objects. We have considered the observations of these
objects and scrutinized the two latest additions to this list: PSR J0901–4046
with $P\approx 76$ s and GLEAM-X J1627 with $P\approx 1091$ s. We have shown
that the periodicity in these objects is highly unlikely to be associated with
a binary companion and that in both cases, the observations are hardly
explained by a WD nature. Instead we find that PSR J0901–4046 and GLEAM-X
J1627 are naturally explained as ultra-long period highly-magnetized NSs
(ULPMs). Together with 1E 161348-5055 ($P\\!\approx\\!2.4\times 10^{4}$ s),
PSR J0250+5854 ($P\\!\approx\\!23.5$ s) and GCRT J1745-3009
($P\\!\approx\\!4620$ s), they form our list of ULPM candidates.

Evidence based on the confirmed magnetar population suggests they undergo
episodes of enhanced spindown associated with their outbursts. The accumulated
effect of outburst associated spindowns is sufficient to drive, at least some
magnetars, to ultra-long periods.

We have shown that the nearby distances of the most recent ULPM candidates
require that there are thousands to tens of thousands such objects in the
Galactic field. This in turn constrains the formation rate and typical ages of
these objects. Indeed the ages of ULPM candidates, are further constrained by
several other considerations, including dipole spindown ages, high timing
stability of the radio pulsations, upper limits on the X-ray fluxes (requiring
the objects to have cooled down sufficiently) and their offsets from SNR and
stellar clusters. The conclusion is that GLEAM-X J1627 and PSR J0901–4046
should have ages of $\sim 10^{5}-10^{6}$ yr, much older than the ages of the
confirmed Galactic magnetar population. These results suggest that ULPMs
represent a distinct evolutionary sequence from confirmed magnetars and have
undergone a distinct evolution of their poloidal fields.

These long magnetic field decay timescales in a common magnetar sub-
population, are a sizeable fraction of the time to merger in the rapid channel
associated with a large fraction of binary neutron star systems (Beniamini &
Piran, 2019). This suggests that magnetar strength-fields might still be
tenable in binary neutron star systems, just before their merger and might
lead to detectable precursor emission from those systems (Lyutikov, 2019).

Ultra-long period magnetars have also been suggested as an explanation for
many days long periodic activity windows seen in some repeating FRBs. We
explored the possible connection between the Galactic ULPM population and
these periodic FRBs. We found that the source density of the latter is between
five and seven orders of magnitude lower than that of the former. It is on the
one hand fortunate that we have more than enough long period objects than
required by FRB observations. On the other hand, the large source density
discrepancy is puzzling and needs to be understood if the connection between
these two phenomena is a real one. Only a fraction of this discrepancy could
be ascribed to beaming. Other factors which could be important are (i) the
$\sim 3$ order of magnitude gap between the period of GLEAM-X J1627 and that
of periodic FRBs (considering the apparent reduction in the density of sources
with increasing periods, Figure 5) and (ii) some selection process causing the
longest period magnetars to be preferential producers of bursts (see e.g.
Beniamini et al. 2020). It therefore remains to be determined whether periodic
FRBs are directly connected to the Galactic ULPM population.

Finally, the recent observations of Galactic ULPM candidates serves as a vivid
demonstration of the large amount of discovery parameter space that remains
unexplored beyond the standard curvature pair cascade pulsar death line. We
have shown that current time-domain searches for radio transients are heavily
biased against finding long period systems, and many more exciting and nearby
long-period objects therefore likely remain undetected. We advocate the need
to explore different observational strategies in order to truly open up the
window to the detection of this potentially large population of objects hidden
in our Galactic backyard.

## acknowledgements

We especially thank Kostas Gourgouliatos, Sam Lander and Wynn Ho for useful
discussions regarding models for field evolution in the crust and core. We
also thank Jonathan Granot, Alice Harding, Oleg Kargaltsev, Demos Kazanas,
Cole Miller, Sk Minhajur Rahaman, and Andrey Timokhin for interesting
discussions. PB’s research was supported by a grant (no. 2020747) from the
United States-Israel Binational Science Foundation (BSF), Jerusalem, Israel.
The material is also based upon work supported by NASA under award number
80GSFC21M0002. This research has made use of NASA’s Astrophysics Data System.
KMR acknowledges support from the Vici research programme “ARGO” with project
number 639.043.815, financed by the Dutch Research Council (NWO). JH
acknowledges support from an appointment to the NASA Postdoctoral Program at
the Goddard Space Flight Center, administered by ORAU. JH thanks Brad Cenko
for helpful discussions related to the access and use of DECam data.

## Data Availability

Data produced in this study will be shared upon reasonable request.

##  Appendix A Magnetar Candidates in Binaries

For completeness, we here detail magnetar candidates in binaries listed in
Table 1 (also see Popov, 2022, for a recent review) excluding pulsating ULXs.
The magnetar nature of these sources is less secure, and chiefly based on
arguments on source luminosity, and spin equilibrium in the propeller regime,
and their position in the Corbet diagram (Corbet, 1984, 1986). The rarity of
such exceptional high mass X-ray binary systems (including pulsating ULXs)
implies they cannot be the dominant channel for production of ULPMs.

###  A.1 IGR 16358–4726

IGR 16358–4726 (Revnivtsev et al., 2003) is an unusual periodic ($P\sim 6000$
s) source originally detected as a nonthermal hard-spectrum INTEGRAL transient
that exhibited a K$\alpha$ iron florescence line characteristic of accreting
compact objects (Kouveliotou et al., 2003; Patel et al., 2004). The hard
spectrum disfavor it being a compact low-mass X-ray binary (LMXB), and thus
orbital modulation origin for the reported periodicity. A subsequent large
spin-up $\Delta\nu/\nu\\!\sim\\!10^{-2}$ over a baseline of $10^{6}$ s
reported by Patel et al. (2007) also disfavor an orbital origin for the
periodicity. Interpreting such a large spin-up as orbital contraction that
would imply a huge transient companion mass loss (that also does not result in
significant spectral or flux changes) of
$\dot{m}\\!\gtrsim\\!10^{-8}m_{c}\\!\sim\\!10^{25}$ g/s
($m_{c}\\!\sim\\!0.1-1M_{\odot}$ is the companion mass in a putative compact
LMXB). Thus a NS origin of the compact object and association of periodicity
with spin is preferred. A disk accretion scenario for the spin-up implies a
large magnetic moment for the NS for a large lever arm (the Alfvén radius). A
standard calculation implies $B_{*}\\!\sim\\!10^{14}d^{-6}_{7}$ G for a
distance of 7 kpc. A lower limit to distance of $\\!\sim\\!7$ kpc is also
inferred if the source flux is not in the propeller regime. These results
suggest for any plausible distance within the galaxy, IGR 16358–4726 possesses
a strong magnetic field $B_{*}\\!\gtrsim\\!10^{13}$ G.

###  A.2 SGR 0755-2933

SGR 0755-2933 is a candidate magnetar discovered by magnetar-like short bursts
by Swift/BAT (Barthelmy et al., 2016). The error region of BAT contains a HMXB
harboring a NS with spin period $P\\!\approx\\!308$ s (Harrison et al., 2017).
Although Doroshenko et al. (2021) suggest the candidate HMXB may be
unassociated with the BAT burst, we will argue against this possibility in a
future dedicated work. The parallax distance to the counterpart is $3.5$ kpc.

###  A.3 4U 1954+319

4U 1954+319 is a peculiar HMXB hosting a 12-50 Myr M supergiant (Hinkle et
al., 2020). It was suggested by Enoto et al. (2014) to host a magnetar with
spin period $P\\!\approx\\!5.7$ hr. Since Enoto et al. (2014), Gaia has
attained distance of $3.3\pm 1$ kpc (Bailer-Jones et al., 2018; Hinkle et al.,
2020) for the counterpart. Stellar evolution and lack of a supernova remnant
constrains the age of the NS to $\\!<\\!43$ Myr. Recently, Bozzo et al. (2022)
argue a highly-magnetized (or magnetar-like) magnetic field is demanded given
the recent Gaia association, in either a wind accretion or spin-equilibrium
propeller scenario.

###  A.4 4U2206+54

First suggested by Ikhsanov & Beskrovnaya (2010, 2013) as a possible wind-fed
magnetar, 4U2206+54 is a HMXB with NS spin period $P\\!\approx\\!5750$ s
(Torrejón et al., 2018). A high spin frequency derivative of
$\dot{\nu}\\!\approx\\!-1.8^{-14}$ Hz s-1 implies a relatively large level arm
and magnetization. Gaia parallax implies a distance of $3.7\pm 0.4$ kpc for
the optical counterpart. Torrejón et al. (2018) suggest that a wind settling
accretion scenario requires a NS with magnetar-like field.

###  A.5 4U0114+65

Suggested by Li & van den Heuvel (1999) as a potential magnetar, 4U0114+65 has
a spin period of $P\approx 9350$ s and positive frequency derivative
$\dot{\nu}\approx 10^{-14}$ (Sanjurjo-Ferrrín et al., 2017). The companion
distance is estimated $7\pm 3$ kpc.

##  Appendix B GLEAM-X J1627 and PSR J0901–4046, and the Radio Duty Cycles of
Rotation-Powered Pulsars

The radio duty cycles of GLEAM-X J1627 and PSR J0901–4046 are generally
inconsistent empirical relations of period dependence of pulsar beam width and
frequency in young rotation-powered pulsars. This is independent of standard
(rotational) theoretical death lines for pair formation tenability, for which
both GLEAM-X J1627 and PSR J0901–4046 are beyond without significant fine-
tuning. This suggests the radio emission is not rotation-powered in both
objects, and that they are entirely (GLEAM-X J1627) or perhaps partially (PSR
J0901–4046) magnetically powered.

More specifically, in a fixed band, rotation-powered pulsar beam full widths
(empirically) scale as $\Delta\phi\sim C/\sqrt{P}$ for a “core" component
aligned with the magnetic axis, and hollow “conal" components. This is the
expected scaling for the corotating magnetospheric open zone and polar caps
where single photon pair production occurs. The absolute value of the beam
width, however, depends on the decoupling height of the radio emission and
other variables such as the pulsar magnetic obliquity. Empirically,
$C\\!\sim\\!6^{\circ}$ s1/2 for “conal" components and
$\\!\sim\\!1-1.5^{\circ}$ s1/2 for “core" Gaussian components. Such core-cone
components (Rankin, 1983, 1990, 1993) are expected as the pair formation front
angle and gap physics varies between the conducting boundary of the polar cap
rim and changing curvature photon angle with respect to the magnetic field.
These empirical relations constitute the maximum observed pulse width expected
from a rotation-powered pulsar, since the observer may have a glancing impact
with the radio beam. The core component angular size scales with the square
root of decoupling height (Kijak & Gil, 1998), $\rho_{\rm
core}\\!\approx\\!1.2^{\circ}\sqrt{h/(h_{\rm kg,0}P)}$, where the height is
normalized to $h_{\rm
kg,0}\\!\sim\\!40R_{6}\dot{P}^{0.07}_{-15}P^{0.3}\nu^{-0.26}_{\rm GHz}$ (Kijak
& Gil, 2003). This empirical relation is an explicit radius-to-frequency
mapping (Cordes, 1978). For a putative ULPM magnetar with $P\\!\sim\\!10^{2}$
s such as PSR J0901–4046, naive extrapolation of these rotation-powered
relations yields $\rho_{\rm core,ulpm}\\!\sim\\!0.4^{\circ}$ or a duty cycle
of $0.1\%$, about $20-50$ times smaller than that observed in the L-band. For
GLEAM-X J1627, the situation is even more at odds. Thus, ULPM emission is
unlikely to be a mere extrapolation of the behaviour of regular pulsars, and
different physical effects are involved. The core component is where most of
the flux density resides in rotation-powered pulsars, and the presence of a
conal components are revealed by a double-peaked pulse profiles. However, note
conal components begin to dominate over core components at larger periods
(Arzoumanian et al., 2002), up to the death line where single photon pair
cascades become untenable. Yet, neither GLEAM-X J1627 or PSR J0901–4046
exhibit features (pulse profile or polarization) typically associated with
conal components.

##  Appendix C Radio QPOs and The presence of a crust in PSR J0901–4046

PSR J0901–4046 is a “pulsar inside a pulsar" exhibiting many QPOs modulating
emission within the envelope of its 76-second period radio pulsations,
including $\\!\sim\\!50$ ms period QPOs consistent with low-order nonmagnetic
$l=2$ or $l=3$ ($n=0$) crustal shear modes (i.e. torsional modes). The QPOs in
PSR J0901–4046 appear to persist for several hundred milliseconds or longer
and appear in about 30% of pulse envelopes (Caleb et al., 2022). That is, the
covering fraction is $\dot{q}{\cal{T}}\\!\sim\\!1/3$ where ${\cal{T}}$ is the
mode damping time and $\dot{q}$ is the excitation rate. If continually but
randomly excited this implies $\dot{q}\\!<\\!1/(3\Delta t)$ ($\Delta
t\\!\sim\\!0.5$ seconds the typical pulse width) for ${\cal{T}}\\!>\\!\Delta
t$. This rate of QPO excitations is much larger than rare observed QPOs in
high-energy magnetar bursts (e.g., Li et al., 2022), which may be detectable
only if the burst is bright enough and the corresponding amplitude of
oscillations is large. Thus the radio observations in PSR J0901–4046 are
likely sampling a population of lower in amplitude (but much more numerous)
crustal vibrations. This suggests a hierarchy in magnitude of vibrations
versus commonality, much like the logN-logS power-law event size distribution
of magnetar bursts or earthquakes (Cheng et al., 1996). Note that if the
effective power (total elastic energy divided by ${\cal{T}}\\!\sim\\!1$ s) in
these magneto-elastic toroidal oscillations (which appear significantly
modulate the radio emission) are comparable to or greater than the radio
luminosity, then their required amplitude is larger than only about
$\xi\\!\gtrsim\\!10^{-4}$ cm. The required screening charge density for such
small amplitude oscillations is much smaller than the persistent one from
corotation ($\propto B/(cP)$), and ought not to lead to dramatic modification
of the pair cascades and radio luminosity as in FRBs, but rather apparently
only the wave/plasma dynamics and inhomogeneity along open field lines (see
also, Timokhin et al., 2000; Morozova et al., 2010; Wadiasingh & Timokhin,
2019).

If these are indeed crustal magneto-elastic modes, they may damp by coupling
to a continuum of Alfvénic modes associated with the magnetized core of the NS
(Levin, 2006). The QPOs appear inconsistent with bouncing magnetospheric
Alfvén modes along field lines, as the pulse profile is commensurate with
polar cap radio emission and has a characteristic “S" shape curved in its
polarization position angle. This would also imply field lines involved in
radio emission attain very high altitudes (or are open, sampling a
characteristic scale of the light cylinder at few $\\!\sim\\!10^{11}$ cm) and
thus the observed harmonic frequencies (and potential damping) timescales are
inconsistent with low-amplitude Alfvénic modes in a bounded cavity as
sometimes invoked for QPOs in magnetar bursts and giant flares (e.g., Castro-
Tirado et al., 2021).

The presence of crustal modes rules out PSR J0901–4046 as a strange (quark)
star (Witten, 1984; Alcock et al., 1986; Weber, 2005; Alford et al., 2008), as
strange quark matter phase transitions and formation scenarios generally
preclude the existence of a thick ionic crust bound by Coulomb interactions
(e.g., Usov, 2004; Niebergal et al., 2006). Although nonbare or hybrid star
models exist, their predicted quasinormal frequencies are much higher into the
kHz regime (e.g., Mannarelli et al., 2015).

It is presently unknown if the QPOs of PSR J0901–4046 are long lived (i.e. if
${\cal{T}}\\!\gg\\!\Delta t$). However, it appears that $P\\!>\\!{\cal{T}}$
since QPOs do not usually persist across consecutive single pulses, nor are
they coherent across pulses. The expected damping time of crustal shear modes
is quite model and magnetic field configuration dependent. Configurations with
(MHD-mode) spectral gaps, tangled magnetic fields, or fields confined largely
to the crust can significantly enhance the longevity of modes (e.g., Sotani &
Kokkotas, 2009; Colaiuda & Kokkotas, 2011, 2012; Bretz et al., 2021). In the
simplest models with a core-bound global dipolar field (Gabler et al., 2012)
that is not entirely expelled by type II proton superconductivity, the damping
(i.e. e-folding) time is $\eta_{\rm damp}\\!\sim\\!0.1$ Alfvén crossing times
and yields an upper bound on a core-tied dipolar magnetic field: $B_{\rm
core}\overset{?}{\lesssim}2\times 10^{13}\,(\eta_{\rm
damp}/0.1)\,{\cal{T}}_{-0.5}^{-1}R_{6}\sqrt{\rho_{14}}$ G. Here
${\cal{T}}\\!\sim\\!0.3-0.5$ s is the minimum damping time,
$\rho\\!\sim\\!10^{14}$ g cm-3 the core density at the crust boundary and
$R\\!\sim\\!10^{6}$ cm the characteristic spatial scale. This together with
the constraints on $B_{\rm d}$ found in §2.1.1 leaves a relatively narrow
parameter space for the simplest picture of the mode structure. The
implication could be that a more complicated scenario is required, such as the
core field structure is significantly modified by proton superconductivity,
and that there are spectral gaps in the MHD continuum of the core field.

##  Appendix D Constraints on optical counterparts to GLEAM

The positional uncertainty of GLEAM-X J1627 was reported as 2′′ (Hurley-Walker
et al., 2022). A search in the Gaia EDR3 catalog returns one optical source
within this positional uncertainty, which is offset from the radio position by
0.94′′, making it a potential optical counterpart to the GLEAM-X J1627. The
Gaia source has an estimated distance of 5.4 kpc, suggesting the source is not
related to GLEAM-X J1627 which has an estimated distance of 1.3 kpc (Hurley-
Walker et al. 2022; Bailer-Jones et al. 2021). However, the parallax
($\bar{\omega}=0.63\pm 0.68$ mas) is not well measured, making this distance
estimate unreliable and warranting additional investigation. The Gaia source
does have well measured proper motion, suggesting that this optical source is
Galactic.

To place a lower-limit on the distance to this Gaia source, the following
procedure was undertaken: First, we extracted all Gaia sources within 10′ of
GLEAM-X J1627’s position. Then we removed all sources with G-band magnitudes
$<20.0$, to keep only those sources with a comparable brightness to this
potential optical counterpart (i.e., $G=20.25$). Next, we cross-matched these
Gaia sources with the parallax distance catalog of Bailer-Jones et al. (2021).
We then plotted the derived distances from Bailer-Jones et al. (2021) (i.e.,
rgeo) versus the measured parallax divided by the parallax error (i.e.,
$\bar{\omega}/\sigma_{\bar{\omega}}$). This allows us to determine out to
roughly what distance sources of a similar brightness to the potential GLEAM-X
J1627 counterpart will have meaningful parallax measurements (i.e., not due to
statistical noise), thus allowing us to place a lower-limit on the distance to
the potential Gaia counterpart. The results of this procedure are shown in
Figure 10. It appears that sources still have well measured parallaxes (i.e.,
$\bar{\omega}/\sigma_{\bar{\omega}}\gtrsim 2$) out to distances of about 3
kpc. Therefore, it is likely that this Gaia source lies far beyond GLEAM-X
J1627 and is an unrelated background source.

Figure 10: Parallax inferred distances from the catalog of Bailer-Jones et al.
2021 plotted versus the measured parallax divided by the parallax error (i.e.,
$\bar{\omega}/\sigma_{\bar{\omega}}$) for all Gaia sources within 10′ of the
potential GLEAM-X J1627 optical counterpart and having a similar G band
magnitude (i.e., $>20.0$). The orange point shows the location of the
potential Gaia counterpart in this plot. The red line shows roughly where
sources begin to have poorly measured parallaxes
($\bar{\omega}/\sigma_{\bar{\omega}}\lesssim 2$), corresponding to a distance
of $\sim 3$ kpc.

There also exists more sensitive optical observations of this field from the
Dark Energy Camera Plane Survey (DECaps; Schlafly et al. 2018). The DECaps
observed in five bands, $grizY$, and goes much deeper than Gaia. The Gaia
G-band limiting magnitude is about 21, while the DECaps images of this field
go to 5$\sigma$ limiting magnitudes of 23.7, 22.8, 22.8, 22.4, and 21.5 in the
five bands, respectively. There are four sources detected in the DECaps images
in the 2′′ error circle of GLEAM-X J1627. Two of these sources have less than
three secure detections and have magnitudes very close to the limiting
magnitude, suggesting that they may be spurious. Additionally, their
magnitudes are not well constrained so we exclude these sources, leaving only
two sources. The brightest of these sources is the source detected by Gaia,
which we already ruled out as a potential counterpart. The second source is
1.3′′ offset from the GLEAM-X J1627 position and has magnitudes $g=<22.8$,
$r=22.8$, $i=21.8$, $z=21.0$, $Y=20.7$.

This leaves three potential options, the first is that the optical source is
the counterpart to GLEAM-X J1627. The second is that GLEAM-X J1627 may be in a
binary with the optical counterpart being the binary companion. Lastly, the
source could simply be a background star, unrelated to the GLEAM-X J1627. To
constrain the source properties in the first two scenarios, we assume the
counterpart lies at a distance of 1.3 kpc, and account for extinction using
the mwdust package (Bovy et al., 2016) with the Fitzpatrick (1999) reddening
law tabulated in the Dark Energy Camera bands by Schlafly & Finkbeiner (2011).
We find that the $E(B-V)=0.21$ at 1.3 kpc in the direction of GLEAM-X J1627.
After accounting for distance and extinction, we find the potential optical
counterpart has an absolute magnitude $M_{r}=11.8$ and color $r-z=1.6$.

To explore the first scenario, we use the WD atmospheric models of Bédard et
al. (2020) to compare the color and absolute magnitude of this potential
counterpart. We find that the absolute magnitude of the potential optical
counterpart is too large for a source with such a red color to be consistent
with a WD, thus ruling out this scenario. If, however, the source is instead a
binary companion to GLEAM-X J1627, it would be more consistent with a late-
type star. Given the absolute magnitude and color, the source would be most
consistent with an M4-6 type M-dwarf (see Figures 15 and 16 in Kiman et al.
2019). This spectral type is also consistent with this sources absolute
$M_{i}\\!=\\!10.9$ and $M_{z}\\!=\\!10.2$ magnitudes and non-detection of the
source by 2MASS as all NIR magnitudes are expected to be $\gtrsim 18$. Further
follow-up observations to obtain a source spectrum could be carried out to
confirm or refute this possibility.

The final possibility is that the source is a background star that is
coincident with GLEAM-X J1627 by chance. The probability of having one or more
DECaps sources in the $\delta r\\!=\\!2^{\prime\prime}$ positional uncertainty
of GLEAM-X J1627 can be calculated as $P_{\rm DECaps}\\!=\\!1-e^{\rho_{\rm
DECaps}\pi\delta r^{2}}$, where $\rho_{\rm DECaps}$ is the density of sources
with three or more good detections in the DECaps survey and $\delta r$ is the
positional uncertainty of GLEAM-X J1627 (see e.g., Hare et al. 2018). This
calculation yields a large chance coincidence probability of
$\\!\approx\\!90\%$, suggesting that this scenario is the most likely. If this
is the case, it is difficult to place any constraints on the spectral type as
the source distance and extinction are unknown.

##  Appendix E Gap Voltage and the Curvature Pair Death Line for Arbitrary
Stellar Radii

The three conditions ($\ell_{\rm gap}<r_{\rm pc}$, $\Delta V_{\rm
gap}<\Phi_{\rm open}$, and $L_{e^{\pm}}<\dot{E}_{\rm SD}$) yield identical
parameter dependence of the “death line" on $B$ and compactness within factors
of unity, verifying the self-consistency of the constraints. Requiring the
pair luminosity is greater than the uncorrected radio luminosity
$L_{e^{\pm}}\\!>\\!L_{R,\rm obs}$ applies only to very compact stars (e.g.
NSs) in an untenable part of the parameter space and is not shown in Figure 2.
For these gross estimates, we adapt the curvature photon gap height $\ell_{\rm
gap}\\!\sim\\!(B_{\rm cr}^{4}\lambda\hbox
to0.0pt{\hss--}^{2}c^{3}/\pi^{3})^{1/7}\rho_{c}^{2/7}P^{3/7}B^{-4/7}$
(Timokhin & Harding, 2015) where $B_{\rm cr}$ is the quantum critical or
Schwinger field, $\lambda\hbox to0.0pt{\hss--}$ the reduced Compton wavelength
and $\rho_{c}$ is the characteristic field curvature radius near the polar
cap. We omit weak order unity factors associated with the relative gap speed
and pair attenuation exponential of the Erber formula (Erber, 1966). The open
field line voltage is $\Phi_{\rm open}\sim\Omega^{2}BR^{3}/c^{2}$ while for a
vacuum gap $\Delta V_{\rm gap}\sim(\Omega/c)B\ell_{\rm gap}^{2}$. These all
yield the death line,

$R\gtrsim 4\times
10^{9}\left(\frac{\varrho_{c}}{10}\right)^{4/17}P_{3}^{13/17}B_{9}^{-8/17}\quad\rm
cm$ (11)

where we adopt dimensionless $\varrho_{c}=\rho_{c}/R=10$ for the
characteristic field curvature radius. This constraint on $R$ then may be
readily expressed in terms of constraint on compactness for a given mass. A
more realistic choice for $\varrho_{c}\gg 10$ commensurate with the expected
small polar cap size for a $P\sim 10^{3}$ s rotator would eliminate any
allowed parameter space in Figure 2.

##  Appendix F X-ray Upper-limits

Source | Observatory | ObsID | Start Date | Exposure
---|---|---|---|---
|  |  | MJD | ks
PSR J0901–4046 | Swift-XRT | 00014019002 | 59245.83 | 3.87
| Swift-XRT | 00014019003 | 59352.04 | 1.35
| Swift-XRT | 00014019004 | 59353.29 | 2.70
GLEAM-X J1627 | Chandra | 26228 | 59601.87 | 20.09
| Chandra | 26282 | 59602.25 | 10.11
PSR J0250+5854 | XMM-Newton | 0844000201 | 58692.15 | 5.20
| XMM-Newton | 0844000301 | 58696.12 | 3.40
| XMM-Newton | 0844000401 | 58694.17 | 2.50
| XMM-Newton | 0844000501 | 58698.12 | 8.00
| XMM-Newton | 0844000601 | 58728.04 | 4.90
| XMM-Newton | 0844000701 | 58730.03 | 4.50
| XMM-Newton | 0844000801 | 58481.59 | 5.30
| XMM-Newton | 0844000901 | 58734.02 | 2.80
| XMM-Newton | 0844001001 | 58741.98 | 5.80
| XMM-Newton | 0844001101 | 58743.98 | 5.10
| XMM-Newton | 0844001201 | 58089.03 | 2.70
Table 2: Observations of 3 long period NSs. Exposure times for XMM-Newton are
calculated after removing instances with high particle background flares

PSR J0901–4046, PSR J0250+5854, and GLEAM-X J1627, were observed by Swift-XRT,
XMM-Newton (PI: Tan), and Chandra (PI: Hurley-Walker), respectively. Table 2
lists a log of all of the observations for each source. The Swift-XRT data
were reduced and analyzed using version 6.29 of the HEAsoft software package,
while the XMM-Newton data were reduced using SAS version 20.0. For Chandra, we
used CIAO version 4.14. We followed the standard procedures for reducing and
cleaning each data set. We note that the exposure times for XMM-Newton listed
in Table 2 are calculated after removing the times with high particle
background flares. Next, we calculated the number of counts at each source’s
radio position using radii of $30^{\prime\prime}$, $15^{\prime\prime}$ and
$2^{\prime\prime}$ for Swift-XRT, XMM-Newton, and Chandra respectively. Then
we calculated the net counts for each source by subtracting the total number
of expected background counts from large source free regions in each image
after correcting for the difference in area of extraction regions. PSR
J0901–4046 and GLEAM-X J1627 both had zero net source counts detected, while
PSR J0250+5854 had 15$\pm 10$ net source counts detected. We converted these
net source counts into $3\sigma$ upper-limits on the count rates in each
observatory using Table 1 from Gehrels (1986) and dividing by the total
exposure, while also correcting for the energy containment fraction and
exposure map. This lead to $3\sigma$ upper-limit count rates of $10^{-3}$,
$9\times 10^{-4}$, and $2\times 10^{-4}$ cts s-1 for PSR J0901–4046, PSR
J0250+5854, and GLEAM-X J1627, respectively. These count rates were calculated
in the 0.5-10 keV bands for Swift and XMM-Newton, and 0.5-8 keV for Chandra.

We used the count rates to estimate the maximum temperature that each NS could
be emitting from its entire surface, yet still remain undetected. We placed
each source at its nominal dispersion measure estimated distance to calculate
the blackbody model normalizations (i.e., bbodyrad in Xspec). The absorbing
column density for each source was estimated from the $N_{\rm H}$-DM relation
of He et al. (2013). Since it is quite common for NSs to only emit from a
small fraction of their surface (e.g., from hot spots), we repeated this
exercise by shrinking the NS emitting area and recalculating the largest
temperature that would remain undetected at the $3\sigma$ level. The curves
shown in the left panel of Figure 4 were calculated adopting this procedure.

## References

  * Alcock et al. (1986) Alcock C., Farhi E., Olinto A., 1986, ApJ, 310, 261
  * Alexander et al. (2020) Alexander K. D., Hajela A., Margutti R., Bright J., Eftekhari T., Kathirgamaraju A., Berger E., 2020, GRB Coordinates Network, 29053, 1
  * Alford et al. (2008) Alford M. G., Schmitt A., Rajagopal K., Schäfer T., 2008, Reviews of Modern Physics, 80, 1455
  * Archibald et al. (2020) Archibald R. F., Scholz P., Kaspi V. M., Tendulkar S. P., Beardmore A. P., 2020, ApJ, 889, 160
  * Arons & Barnard (1986) Arons J., Barnard J. J., 1986, ApJ, 302, 120
  * Arzoumanian et al. (2002) Arzoumanian Z., Chernoff D. F., Cordes J. M., 2002, ApJ, 568, 289
  * Bahcall & Soneira (1980) Bahcall J. N., Soneira R. M., 1980, ApJS, 44, 73
  * Bailer-Jones et al. (2018) Bailer-Jones C. A. L., Rybizki J., Fouesneau M., Mantelet G., Andrae R., 2018, AJ, 156, 58
  * Bailer-Jones et al. (2021) Bailer-Jones C. A. L., Rybizki J., Fouesneau M., Demleitner M., Andrae R., 2021, AJ, 161, 147
  * Barthelmy et al. (2016) Barthelmy S. D., et al., 2016, The Astronomer’s Telegram, 8831, 1
  * Basu et al. (2022) Basu A., et al., 2022, MNRAS, 510, 4049
  * Baym et al. (1969a) Baym G., Pethick C., Pines D., 1969a, Nature, 224, 673
  * Baym et al. (1969b) Baym G., Pethick C., Pikes D., 1969b, Nature, 224, 674
  * Bédard et al. (2020) Bédard A., Bergeron P., Brassard P., Fontaine G., 2020, ApJ, 901, 93
  * Beniamini & Piran (2019) Beniamini P., Piran T., 2019, MNRAS, 487, 4847
  * Beniamini et al. (2019) Beniamini P., Hotokezaka K., van der Horst A., Kouveliotou C., 2019, MNRAS, 487, 1426
  * Beniamini et al. (2020) Beniamini P., Wadiasingh Z., Metzger B. D., 2020, MNRAS, 496, 3390
  * Bethapudi et al. (2022) Bethapudi S., Spitler L. G., Main R. A., Li D. Z., Wharton R. S., 2022, arXiv e-prints, p. arXiv:2207.13669
  * Bhardwaj et al. (2021) Bhardwaj M., et al., 2021, ApJ, 910, L18
  * Bochenek et al. (2020) Bochenek C. D., Ravi V., Belov K. V., Hallinan G., Kocz J., Kulkarni S. R., McKenna D. L., 2020, Nature, 587, 59