1 INTRODUCTION
Pulsars are highly magnetised, rapidly rotating neutron stars that convert their rotational kinetic energy into magnetic dipole radiation. Although their emission can extend to the entire electromagnetic spectrum, they are typically observed at radio wavelengths.
In a simplified picture, pulsar radio emission is generated in proximity of the magnetic poles, and forms radiation beams. If the magnetic and the spin axes are misaligned, then the two beams rotate with the neutron star and sweep through space. An observer whose line of sight crosses one or both of the beams, will observe a pulsed emission, the period of which corresponds to the pulsar’s spin period (a general introduction to pulsar astronomy can be found, e.g., in Lorimer & Kramer Reference Lorimer and Kramer2005). Pulsars are often referred to as ‘cosmic clocks’, because it is possible to predict the arrival time of each pulse at a telescope, sometimes with sub-microsecond precision (Desvignes et al. Reference Desvignes2016; Arzoumanian et al. Reference Arzoumanian2015b; Reardon et al. Reference Reardon2016), through the pulsar timing technique (see Section 3).
Due to their high densities, rapid rotations, strong magnetic fields, and high surface gravity, neutron stars are ideal laboratories for tests of nuclear physics (Lattimer & Prakash Reference Lattimer and Prakash2004), general relativity (GR) (Kramer et al. Reference Kramer2006) and alternative theories of gravity (Shao et al. Reference Shao, Caballero, Kramer, Wex, Champion and Jessner2013) in extreme conditions not feasible in Earth-based laboratories (see also Stairs Reference Stairs2003; Chamel & Haensel Reference Chamel and Haensel2008). In the context of GR tests, Pulsar Timing Array (PTA, Foster & Backer Reference Foster and Backer1990) experiments are among the most exciting projects of the last decade. The primary aim of PTAs is the direct detection of low-frequency gravitational waves (GWs) (Rajagopal & Romani Reference Rajagopal and Romani1995; Wyithe & Loeb Reference Wyithe and Loeb2003; Sesana et al. Reference Sesana, Haardt, Madau and Volonteri2004).
In this article, we describe the basic concepts and approaches of PTA experiments, and we review the recent results from the established PTA experiments. In Section 2, we review the efforts to detect GWs in different parts of the spectrum. In Section 3, we give an overview of the technique of ‘pulsar timing’, used to interpret pulsar data for the purpose of PTAs. In Section 4, we outline the potential sources of GWs at low frequencies, and in Section 5, we describe the basic concepts of PTA experiments. In Section 6, we summarise the latest results from the existing PTA experiments, and in Section 7, we discuss the future prospects of PTAs, also considering the new radio astronomical facilities.
2 THE QUEST FOR GRAVITATIONAL WAVES
2.1. Gravitational waves
GWs were an early prediction of GR (Einstein Reference Einstein1916), and are a consequence of a small (i.e., linearisable) perturbation h μ, ν to an otherwise flat (or Minkowskian) metric ημν of space–time, produced by asymmetric and accelerated mass distributions:
$$\begin{equation}
g_{\mu \nu } = \eta _{\mu \nu } + h_{\mu \nu }.
\end{equation}$$
It is possible to demonstrate that the perturbation h μν propagates in the metric as a transverse wave at the speed of light, and in its propagation, it induces quadrupolar perturbations of space–time. Given a mass distribution (we recall that, in GR, the presence of a mass distribution curves the space–time), it is also possible to demonstrate that such perturbations h μν are generated if the second time derivative of the quadrupole mass moment Q is not zero (Maggiore Reference Maggiore2007).
The amplitude of GWs is typically expressed in terms of the dimensionless strain h, i.e., the fractional change δL induced by GWs over a distance L:
$$\begin{equation}
h = \frac{\delta L}{L}.
\end{equation}$$
GWs can have two polarisations, commonly referred to as ‘plus’ and ‘cross’. If a ‘plus’-polarised GW propagates along the z axis, then it will alternatively stretch and compress space–time along the y and x axes in the orthogonal direction. A ‘cross’-polarised GW will have the same effect, although rotated by 45°.
Pulsar astronomy brought the first indirect confirmation of the existence of GWs, through observations of PSR B1913+16 (Hulse & Taylor Reference Hulse and Taylor1974). This object is a pulsar in a ~7.7-h orbit with another neutron star. By assuming the existence of GWs, GR can predict the rate of orbital decay that can be attributed to GW emission due to the orbital motion of the two neutron stars. The orbital decay of this binary system was found in agreement with the predictions, today to a precision greater than 99.5% (Weisberg & Taylor Reference Weisberg and Taylor1981; Weisberg, Nice, & Taylor Reference Weisberg, Nice and Taylor2010; Weisberg & Huang Reference Weisberg and Huang2016, and see Figure 1).
Figure 1. Figure taken from Weisberg & Huang (Reference Weisberg and Huang2016)Footnote 1 ; orbital phase shift of the double neutron star system that includes PSR B1913+16 versus time. The plot shows the perfect agreement between the observed orbital decay of the (black dots) and the prediction by GR (solid line).
2.2. Searching for GWs in the Cosmic Microwave Background
The cosmological inflation is an epoch in the early history of the Universe, that is conjectured to seed structure formation and primordial GWs. The quasi-exponential expansion of the Universe during this phase is thought to have generated a stochastic background of GWs (Starobinskiǐ Reference Starobinskiǐ1979), that cannot be detected directly with current instrumentation. However, indirect detections may be possible. The inflationary GW background is predicted to have excited both of the polarisation patterns of the Cosmic Microwave Background (CMB); the E-mode pattern (curl-free) and the B-mode pattern (curl). Although the GW-induced E-mode it is not expected to be detectable, the signature in the otherwise quiescent B-mode should be measurable (Polnarev Reference Polnarev1985) as an excess power at large angular scales (the recombination bump at l ~ 100, where l is the multipole moment). A detection of the B-mode would provide crucial information in support of the inflationary model. However, this achievement is challenging because of B-mode contaminations given by the gravitational lensing of the E-mode on small angular scales (l ~ 1000), and the polarised foreground emission (such as from dust and synchrotron radiation) from our Galaxy (Tucci, Martínez-González, Vielva & Delabrouille Reference Tucci, Martínez-González, Vielva and Delabrouille2005) on spatial scales that are searched for the inflationary signature.
Searches for the B-mode of the CMB polarisation are currently ongoing, through experiments such as POLARBEAR (Kermish et al. Reference Kermish and Holland2012), the ongoing observations with the South Pole Telescope (see, e.g., Benson et al. Reference Benson, Holland and Zmuidzinas2014), and the Background Imaging of Cosmic Extragalactic Polarization (see, e.g., Keating et al. Reference Keating, Ade, Bock, Hivon, Holzapfel, Lange, Nguyen, Yoon and Fineschi2003). These experiments, that are focusing on smaller l than the satellites, are expected to detect effects of the inflationary GW background at ultra-low frequencies (below 10−16 Hz, Lasky et al. Reference Lasky2016). No detection of the inflationary-induced B-mode has been made to date, while the E-mode was observed for the first time with the Degree Angular Scale Interferometer (Kovac Reference Kovac, Leitch, Pryke, Carlstrom, Halverson and Holzapfel2002) and followed-up in more detail by the Wilkinson Microwave Anisotropy Probe and the Planck satellites (see, e.g., Komatsu et al. Reference Komatsu2011; Planck Collaboration, et al. 2014).
2.3. Searching for GWs with interferometers
In addition to indirect detections based on pulsars or the CMB, it is possible to make direct detections of GWs. This has been achieved by the Advanced Laser Interferometer Gravitational-wave Observatory (aLIGO), and in the future, direct detections will be possible with the Evolved Laser Interferometer Space Antenna (eLISA), and PTAs.
aLIGO comprises two laser interferometers (in Hanford and Livingston) capable of detecting changes in the length of the interferometer arms induced by GWs. As other ground-based laser interferometers such as Virgo in Italy (Acernese et al. Reference Acernese2015), GEO 600 (Willke et al. Reference Willke2002) in Germany, or KAGRA (Araya et al. Reference Araya2017, not online yet) in Japan, aLIGO explores GW frequencies from approximatively 1 to 103 Hz. The lower limit is set by gravity gradients in the Earth gravitational potential, while the upper limit is given by the shot noise of the laser photons (Aasi et al. Reference Aasi2013). GWs emitted in this frequency range are predicted to be generated by coalescing binary systems of neutron stars or stellar-sized black holes. In 2015 September, aLIGO achieved the first direct detection of GWs (Abbott et al. Reference Abbott2016a) from a coalescing binary of stellar-mass black holes. This detection (followed by other six events and one candidate since 2015, Abbott et al. Reference Abbott2016b, Reference Abbott2017) signed the beginning of the era of GW astronomy.
However, more is required to explore this branch of science—more detections, and a wider range of GW frequencies.
eLISA (eLISA Consortium et al. 2013; Amaro-Seoane et al. Reference Amaro-Seoane2017) is a project of the European Space Agency to deploy a three-body, space-based interferometer with arms 2.5 million km long, that will probe the GW spectrum in a frequency range from 10−1 down to 10−5 Hz. An eLISA pathfinder, a 40-cm one-armed miniature of the future device, was launched in 2015 and reported significantly lower noise levels than expected (Armano et al. Reference Armano2017). This success grew confidence and expectations for the mission, planned for launch in 2034. The frequency range limit of eLISA is given by the measurement accuracies of the free-falling test mass accelerations (Amaro-Seoane et al. Reference Amaro-Seoane2017). The predicted sources of GWs in this frequency range are inspiralling binary systems of white dwarfs and super-massive black holes (SMBH).
The only experiment that can currently provide longer interferometric baselines, on a parsec scale, are PTAs (see Section 5). PTAs explore the frequency range from about 10−6 to 10−9 Hz, where the most likely source of GW emission are coalescing SMBH binaries (SMBHB). Other sources might be cosmic strings and relic GWs from inflation. PTAs are experiments based on the monitoring of an ensemble of selected pulsars, in order to search for spatially correlated deviations in the arrival times of their pulses. A number of phenomena can induce such correlations, included GWs.
The GW frequency bands explored by these three kinds of experiments are complementary, as shown in Figure 2.
Figure 2. Figure taken from Janssen et al. Reference Janssen2015 Footnote 2 ; GW amplitude versus GW frequency, and frequency ranges explored by the interferometric experiments searching for GWs, aLIGO, eLISA, and PTAs. In the ‘PTA band’, the nominal sensitivities for the International Pulsar Timing Array are shown and the Square Kilometer Array, together with a representation of the expected emission from the SMBHB population (solid blue line) in the universe, the emission from ‘GW-loud’ SMBHBs (blue triangles) and from the unresolvable SMBHBs (light purple squares). In the ‘eLISA band’, the nominal eLISA sensitivity curve is shown, together with the expected GW signals from different masses of merging SMBHBs (cyan), a binary with a very high mass ratio (aquamarine), and from the Galactic population of inspiralling white dwarf binaries (yellow). In the ‘aLIGO band’, the sensitivity curve of aLIGO (as of 2015) is shown, together with the expected signals from different inspiralling compact-object binaries (purple). In brown, orange, and red are the GW background expected from inflation and two models of cosmic strings.
3 PULSAR TIMING
As mentioned in Section 1, pulsars are often referred to as ‘cosmic clocks’, as it is possible to predict their times of arrival (ToAs) to high accuracy. Essentially, there are three requirements to enable accurate predictions of the arrival times:
-
1. The pulsar is a stable rotator. As mentioned in Section 1, pulsars lose rotational energy via magnetic dipole radiation, and therefore they spin down. However, the spin-down might suffer from irregularities, in the form of abrupt (‘glitches’, Downs Reference Downs1981) or long-term (‘timing noise’, see Section 3) variations in the spin frequency. For high-precision experiments, pulsars should have predictable spin-evolutions.
-
2. The shape of the integrated pulse profile is stable in time. Pulsars are intrinsically weak sources, with fluxes in the order of a few milliJansky (Lorimer et al. Reference Lorimer, Yates, Lyne and Gould1995; Kramer et al. Reference Kramer, Xilouris, Lorimer, Doroshenko, Jessner, Wielebinski, Wolszczan and Camilo1998). Typically, the individual pulses do not exceed the radiometer noise of the telescope. Thus, many pulsar studies use the integrated pulse profile, i.e., the coherent sum of many thousands of individual pulses. While individual pulses differ often from each other (both in flux distribution and phase), the integrated profile is statistically stable. Known sources of variations in the integrated profile are pulse jitter (Liu et al. Reference Liu2016a), plasma propagation effects (Geyer et al. Reference Geyer2017), or magnetospheric instabilities (Lyne et al. Reference Lyne, Hobbs, Kramer, Stairs and Stappers2010). The long-term temporal stability of the integrated profile is an assumption for high-precision experiments, and efforts are ongoing to mitigate the impact of integrated profile variations in PTA experiments (Lentati et al. Reference Lentati, Kerr, Dai, Shannon, Hobbs and Osłowski2017b).
-
3. The timing model of the pulsar is well known. The timing model of a pulsar (or ephemeris) is a set of parameters that describes the pulsar spin and spin-down, its orbital parameters (if any), its astrometry, and the dispersive influence of the ionised interstellar medium (IISM) along the line of sight to the pulsar. The frequency-dependent dispersive effect of the IISM on radio pulses is quantified by the dispersion measure (DM). The DM is defined as the integrated column density of free electrons along the line of sight:
(3)where d is the distance to the pulsar (pc), and n e is the free electron number density (cm−3).
$$\begin{equation}
\mathrm{DM} = \int _{0}^{{\it d}} n_e \text{d}l,
\end{equation}$$
The first draft of ephemeris for a certain pulsar can be obtained from its discovery, and provides an approximate estimate of the pulsar spin, position, and DM. A precise knowledge of the timing model can be achieved through the technique of pulsar timing (Lorimer & Kramer Reference Lorimer and Kramer2005).
Let us assume that an observing campaign is performed on a pulsar with a given radio telescope. For each observation, we can obtain an integrated pulse profile P, so that the pulse profile is statistically stable and has a suitably high signal-to-noise ratio (S/N). In pulsar timing, average ToAs for each observation are computed via a cross-correlation of each integrated pulse profile with a high-S/N reference template S (Taylor Reference Taylor1992), which is typically a noise-free representation of the pulse profile. This yields a phase shift τ between P and S, if we consider P to be described as
$$\begin{equation}
P(t) = a + bS(t - \tau ) + n(t),
\end{equation}$$
where a, b, and n(t) represent, respectively, an intensity baseline, an intensity normalisation, and the noise level. The topocentric ToA, ToAtopo is the sum of τ to a time stamp associated with the observation.
The topocentric ToAs are then transformed to the (at first-order) inertial reference frame of the Solar System barycentre (SSB). This conversion is based on the parameters included in the timing models, a reference for the time standard, and for the planetary ephemeris (Edwards et al. Reference Edwards, Hobbs and Manchester2006):
$$\begin{equation}
{\rm ToA}_{{\rm SSB}} = {\rm ToA}_{{\rm topo}} + t_{{\rm clk}} - \frac{D}{f^2} + \Delta _{\rm R} + \Delta _{\rm E} + \Delta _{\rm S}.
\end{equation}$$
In this equation, t clk transforms the reference time standard from the (typically) maser-based clock at the observatory to a world-wide recognised time standard such as Terrestrial Time. The third term removes the effects of observing at non-infinite frequency:
$$\begin{equation}
D = \frac{e^2}{2 \pi m_e c} {\rm DM},
\end{equation}$$
where e and m e are the charge and the mass of an electron, and c is the light speed. ΔR is the Roemer delay, that corrects for the difference in travel time between the observatory and the SSB. The Roemer delay is purely based on geometrical considerations, and uses the astrometric parameters of the studied pulsar and the planetary ephemeris. ΔE, the Einstein delay, is based on the planetary ephemeris and corrects for the effects of the gravitational redshift induced by the bodies of the Solar system. ΔS, the Shapiro delay, accounts for the additional time travel required to the light waves for travelling across the gravitational field of the Solar system. Additional corrective parameters are required if the pulsar is part of a binary system.
Once the barycentric ToAs t has been derived, we compute the pulse number N that represents a ‘counter’ for the number of pulsar rotations:
$$\begin{equation}
N (t) = N_0 + \nu _0 (t - t_0) + \frac{1}{2}\dot{\nu }(t - t_0)^2 + ...,
\end{equation}$$
where N
0 is the pulse number at the reference time t
0, and ν0 and
$\dot{\nu }$
are the spin period at t
0 and the spin down rate, respectively. The right-hand side of equation (7) is the Taylor expansion of the pulsar spin.
In the last steps of a timing analysis, the parameters included in the timing model can be varied so that the ToAs are spaced of an integer number of pulsar rotations. The refinement of the timing model can be achieved through dedicated software such as tempo2 (Hobbs et al. Reference Hobbs, Edwards and Manchester2006) and the inspection of the ‘timing residuals’, i.e., the difference between the closest integer number of pulsar rotations and actual number of pulsar rotation among the ToAs. If parameters in the timing model are imprecisely estimated or missing, then we expect to see structures in the timing residuals. For example, we see from equation (7) that an incorrect spin frequency or spin-down rate will show as a linear and a parabolic trend in the timing residuals (see Lorimer & Kramer Reference Lorimer and Kramer2005). If the timing model is sufficiently accurate, then the timing residuals will look ‘white’, i.e., with no correlations (see Figure 3). For reasons that will be explained in Section 5, PTAs are particularly interested in the study of the ‘red noise’. A time series affected by red noise shows long-term correlated structures in the time domain, and an excess in the low-frequency bins of its power spectrum (see Figure 4). Several phenomena can induce red noise, for example, DM variations (You et al. Reference You2007), or intrinsic instabilities in the pulsar spin (better known as ‘spin noise’ or simply ‘timing noise’, Caballero et al. Reference Caballero2016), instrumental imperfections, or GWs.
Figure 3. Timing residuals versus time for PSR J1012+5307. The used observations were obtained at L-band with the Effelsberg radio telescope, and the Effelsberg–Berkeley Pulsar Processor (EBPP) backend. The used ephemeris, with no additional fitting applied, were obtained by Verbiest et al. (Reference Verbiest2016), based on all the available IPTA datasets, including the EBPP one (i.e., the IPTA data release, see Section 6).
Figure 4. Timing residuals versus time for PSR J1939+2134 (upper panel), and corresponding power spectrum versus frequency (lower panel). The used observations were obtained at L-band with the Effelsberg radio telescope, and the Effelsberg–Berkeley Pulsar Processor (EBPP) backend. The used ephemeris, with no additional fitting applied, were obtained by Verbiest et al. (Reference Verbiest2016), based on all the available IPTA datasets, including the EBPP one (i.e., the IPTA data release, see Section 6). The timing residuals of PSR J1939+2134 are clearly affected by red noise, most likely spin noise.
4 SOURCES OF GRAVITATIONAL WAVES AT LOW FREQUENCIES
As mentioned in Section 1, GWs are produced by the second time derivative of the quadrupole moment of the mass distribution that distorts the space–time.
Following a dimensional analysis, the amplitude h of a GW is given by Hughes (Reference Hughes2003):
$$\begin{equation}
h \propto \frac{G}{r c^4} \frac{d^2 Q}{{\rm d}t^2},
\end{equation}$$
where G is the gravitational constant, r is the distance to the GW source, and Q is the quadrupole moment. Due to the factor G/c 4, the peak GW amplitude is expected to be small. GWs are therefore more easily detectable when the second derivative of the quadrupole moment is large (Thorne Reference Thorne, Hawking and Israel1987), as in the case of massive, fast-moving objects. The most likely source of GWs at low frequencies are coalescing SMBHBs, although other potential sources have been identified: GWs from inflation (Grishchuk Reference Grishchuk1974; Starobinskiǐ Reference Starobinskiǐ1979) or cosmic strings (Kibble Reference Kibble1976).
4.1. Super-massive black hole binaries
Observational evidence shows that SMBHs are hosted at the centre of the most or all galaxies (Kormendy & Richstone Reference Kormendy and Richstone1995; Magorrian et al. Reference Magorrian1998). The hierarchical or ‘bottom-up’ scenario (White & Rees Reference White and Rees1978) predicts that larger galaxies are generated via merging of smaller galaxies at high redshifts (z). When two galaxies merge, we then expect that the two SMBHs at their centres form a binary system (Begelman et al. Reference Begelman, Blandford and Rees1980; Volonteri, Haardt, & Madau Reference Volonteri, Haardt and Madau2003). A binary system is characterised by a non-zero value of the second time derivative of its mass quadrupole moment, thus it is a continuous source of GWs. Defining the masses of the individual SMBHs as m
1 and m
2, and assuming a circular orbit for simplicity, with total mass M = m
1 + m
2, reduced mass μ = m
1
m
2/M, and
$\mathcal {M}=(m_1m_2)^{3/5}/M^{1/5}$
(the chirp mass), and using geometricised units such that G = c = 1, the luminosity emitted in GWs (L
gw) by SMBHB is given by (Thorne Reference Thorne, Hawking and Israel1987; Sesana Reference Sesana2013a)
$$\begin{equation}
L_{{\rm gw}} = \frac{32}{5} (\pi \mathcal {M} f)^{10/3},
\end{equation}$$
where f is the observed frequency of the emitted GWs, equal to twice the orbital frequency f B.
The inclination-polarisation averaged amplitude, h, of the radiated GWs is given by (Sesana Reference Sesana2013a)
$$\begin{equation}
h = \sqrt{\frac{32}{5}} \frac{G \mathcal {M}}{r c^4} (\pi f)^{2/3}.
\end{equation}$$
Equations (9) and (10) describe the simple case of a SMBHB in the local universe (i.e., with zero redshift). As pointed out by Vecchio (Reference Vecchio2004), it is possible to ‘move’ the GW source at a different redshift by substituting m x with m x (1 + z), r with r(1 + z) and f with f(1 + z).
Note that both of the expressions for the GW luminosity and strain contain the GW frequency f. During the binary inspiral, f B (and hence f) changes in time. This means that the strain h of the propagated GW might not be the same in two different points in space–time. However, during short time scales over which the change in orbital separation is negligible, the SMBHB can be considered a monochromatic source of GWs (Sesana & Vecchio Reference Sesana and Vecchio2010).
We can identify three stages in the evolution of an SMBHB (Flanagan & Hughes Reference Flanagan and Hughes1998):
-
1. Inspiral, The two SMBHs orbit each other and such orbital separation shrinks due to environmental effects and GW emission.
-
2. Merger, The two SMBHs coalesce, emitting a GW burst that permanently modifies space-time, called a ‘memory’ event.
-
3. Ringdown, the SMBH and the nearby space–time undergo to a relaxation that leads to a spherical configuration.
f gw = 2f B , the low-frequency GW emission is supposed to happen during the inspiral and merger phases, while the ringdown stage is predicted to generate GWs at higher frequencies.
The ‘memory’ phenomenon (Braginskii & Thorne Reference Braginskii and Thorne1987) is a non-oscillatory GW emission that should occur before the actual coalescence of the two BHs. The final stages of the merger generate a net non-zero contribution to the ‘plus’ polarisation mode of the GW emission. Such a DC offset induces a permanent deformation in the metric (Favata Reference Favata2009). Memory events can only be detected during their passage through the detectors, that, in the moments when they operate the metric deformation. The strain of a memory event, h m , is predicted to be (Madison, Cordes, & Chatterjee Reference Madison, Cordes and Chatterjee2014)
$$\begin{equation}
h_m \approx \frac{1-\sqrt{8}/3}{24} \frac{G\mu }{c^2 r} \sin {\mathcal {I}}(17 + \cos ^2{\mathcal {I}})\times [1 + \mathcal {O}(\mu ^2/M^2)],
\end{equation}$$
where
$\mathcal {I}$
is the inclination angle of the binary before the merger. For r = 1 Gpc, and M
1 = M
2 = 109 solar masses, h
m
≈ 10−15 (Madison et al. Reference Madison, Cordes and Chatterjee2014).
The amplitude of a memory event should rapidly increase in the very final stages of the coalescence before the merging, in a time scale τ ≈ 2πR s /c, where R s is the Schwarzschild radius (Cordes & Jenet Reference Cordes and Jenet2012; Madison et al. Reference Madison, Cordes and Chatterjee2014).
5 PULSAR TIMING ARRAYS
Here we outline the expected signatures in pulsar timing data given by the various types of GW emission from SMBHBs.
5.1. Modelling a GW signal from single sources
Two processes are predicted to lead to the orbital shrinking of the SMBHB. First, the shrinking is led by environmental effects until the GW emission reaches approximatively nHz frequencies (see Sesana Reference Sesana2013a for a review). Let us assume a circular, evolving (i.e., shrinking) SMBHB, and a pulsar p, characterised by an angle i between the orbital plane of the SMBHB and the line-of-sight towards the pulsar. Because the GW emission from the SMBHB permeates the entire sky in a quadrupolar fashion, we can expect that ToAs from p will be affected, arriving slightly advanced or slightly delayed than the timing model prediction (Sazhin Reference Sazhin1978). It can be demonstrated that the effect on the timing residuals R(t) of a single pulsar is independent from the travel path of the radiation (Detweiler Reference Detweiler1979). In particular (Babak et al. Reference Babak2016),
$$\begin{equation}
R(t) = R_{\rm E}(t) - R_{\rm p}(t),
\end{equation}$$
R E(t) and R p(t), called the Earth term and the pulsar term, describe the residuals induced by GWs passing over the Earth and the pulsar, respectively. In the quadrupolar approximation (Lommen Reference Lommen2015), we have that (Babak et al. Reference Babak2016)
$$\begin{eqnarray}
R_{\rm E}(t) &=& \frac{h}{\omega } \lbrace (1+{\rm cos}^2{i})F^+[{\rm sin}(\omega t + \Phi ) - {\rm sin}{\Phi }] \nonumber\\[-3pt]
&&+\, 2{\rm cos}{i}F^{\times } [{\rm cos}(\omega t + \Phi ) -{\rm cos}{\Phi }] \rbrace , \nonumber\\[-3pt]
R_{\rm p}(t) &=& \frac{h_{\rm p}}{\omega _{\rm p}} \lbrace (1+{\rm cos}^2{i})F^+[{\rm sin}(\omega t + \Phi + \Phi _{\rm p}) - {\rm sin}(\Phi + \Phi _{\rm p})]\nonumber\\[-3pt]
&&+ \,2{\rm cos}{i}F^{\times } [{\rm cos}(\omega t + \Phi + \Phi _{\rm p}) -{\rm cos}{\Phi + \Phi _{\rm p}}] \rbrace ,
\end{eqnarray}$$
where h and h p are the amplitudes of the GW at the Earth and at the pulsar, ω and ωp are the GW angular frequency at the Earth and at the pulsar, Φ and Φp are the GW phases at the Earth and at the pulsar, and F + and F × are the antenna response functions for the two GW polarisation at the pulsar (i.e., how space–time around the pulsar is affected by the GW). If the SMBHB does not evolve, the power spectrum of the signature described in equation (12) is described by two Dirac delta functions. This signal would be indistinguishable from the signature given by an error in the orbital period and thus not detectable with pulsar timing.
As mentioned in Section 4.1, the last stage of an SMBHB merger lead to a permanent deformation of the metric, similar to a DC offset in space–time. This non-oscillatory phenomenon is propagated, and affects the timing residuals of a pulsar in a way that is equivalent to a change in its rotational spin frequency (vaan Haasteren & Levin Reference van Haasteren and Levin2010). The timing residuals induced by a memory event will be given by (vaan Haasteren & Levin Reference van Haasteren and Levin2010)
$$\begin{equation}
R(t) = h_m B(\theta ,\phi )\times [\Theta (t-t_0) - \Theta (t-t_1)],
\end{equation}$$
where B(θ, ϕ) = 1/2cos (2ϕ)(1 − cos θ), θ is the angular separation between the pulsar and the SMBHB, ϕ is the angular separation between the principal polarisation of the GW signal and the projected line-of-sight to the pulsar onto the plane perpendicular to the GW propagation direction. Θ is the Heaviside function, while t 0 and t 1 are the instants in which the memory event passes the Earth and the pulsar, respectively. In equation (14), it is thus possible to identify an Earth term and a pulsar term, as for the oscillatory contribution to the GW emission shown in equation (12). In particular, the Earth term sensitivity to a memory event is found to increase with the square root of the number of pulsars included in a PTA (Cordes & Jenet Reference Cordes and Jenet2012).
5.2. Modelling a GW background signal
The expected number of SMBHB systems is extremely large, up to 106 depending on the redshift and the mass range of the involved BHs (Sesana, Vecchio, & Colacino Reference Sesana, Vecchio and Colacino2008). The choral GW signal coming from such a population of SMBHBs gives rise to an incoherent superposition of the individual GW signals, that effectively generates a stochastic background of GWs (GWB, Sesana et al. Reference Sesana, Vecchio and Colacino2008; Ravi et al. Reference Ravi, Wyithe, Shannon and Hobbs2015), usually considered isotropic.
The GWB is predicted to induce a red-noise signal in pulsar timing residuals, with a power spectrum P(f) that can be described by a steep power-law (Phinney Reference Phinney2001):
$$\begin{equation}
P(f) = \frac{h^2}{12\pi ^2} \left( \frac{f}{f_{{\rm 1yr}}} \right) ^{2\alpha -3},
\end{equation}$$
where f 1yr normalises the GW frequency at 1/1yr, h is now the amplitude of the GWB, and α is a coefficient whose value is 2/3 in the case of an isotropic and stochastic GWB, thus the spectral index for a GWB is expected to be −13/3. In the case of a GWB, Ωgw(f), the ratio between the energy density ρgw of the GWs (per unit logarithmic frequency) and the critical energy density of the Universe ρ c , is related to the strain h as (Allen & Romano Reference Allen and Romano1999):
$$\begin{equation}
\Omega _{{\rm gw}} (f) = \frac{2 \pi ^2}{3 H_0^2} f^{2} h^2(f),
\end{equation}$$
where H 0 is is the Hubble expansion rate (100 h H km s−1 Mpc−1, with h H being the dimensionless Hubble parameter). In the case of relic GW from inflation and cosmic strings, the spectral index of equation (15) is expected to be −5 (Grishchuk Reference Grishchuk2005) and −16/3 (Damour & Vilenkin Reference Damour and Vilenkin2005), respectively.
5.3. The Hellings and Downs curve
Although SMBHBs are considered the loudest sources of GWs in the universe, the amplitude of such emission is predicted to be extremely tiny (Arzoumanian et al. Reference Arzoumanian2014; Babak et al. Reference Babak2016; Zhu et al. Reference Zhu, Wen, Xiong, Xu, Wang, Mohanty, Hobbs and Manchester2016). The GWB amplitude is expected to exceed the amplitude of the signals from the vast majority of individual SMBHBs. This implies that the detection of a GWB is much more likely than GWs from an individual SMBHB (Rosado, Sesana, & Gair Reference Rosado, Sesana and Gair2015).
As already discussed, the GWB is a stochastic signal that can be described as a red noise process [equation (15)]. Given an individual pulsar, a GWB signal cannot be distinguished unequivocally from other red noise processes such as timing noise, IISM effects, clock noise, or ephemeris errors. Additionally, the pulsar timing procedure absorbs all of the power present in the first two bins of the power spectrum of the timing residuals (corresponding to the spin period and spin period derivative) and in the 1/1yr frequency bin (due to the orbital period of the Earth), effectively decreasing the detectability of a GW signal at this frequency.
However, a GW emission (both from SMBHBs and a background) would affect different pulsars in a correlated fashion depending on their respective sky position. Therefore, detection statistics are based on the correlation among the timing residuals of an array of pulsars (Romani Reference Romani, Ögelman and van den Heuvel1989; Foster & Backer Reference Foster and Backer1990).
The correlation C among the timing residuals of pairs of pulsars perturbed by an isotropic and stochastic GWB was studied by Hellings & Downs (Reference Hellings and Downs1983). Given a pair of pulsars (i, j), separated by an angle θ i, j in the sky, they demonstrated that C takes a specific functional form, known as ‘Hellings and Downs curve’:
$$\begin{equation}
C(\theta _{ij}) = \left[ \frac{3}{2} x\, {\rm log}(x) - \frac{x}{4} + \frac{1}{2} \right](1+\delta _{{\rm {\it i,j}}}),
\end{equation}$$
where x = (1 − cosθi, j)/2. The Hellings and Downs curve is the sky- and polarisation-averaged angular correlation between pairs of pulsars. In the computation, the GWB is assumed to be isotropic (i.e., the power spectrum of the GWB does not have an angular dependence), and the short-wavelength approximation to be valid (i.e., f gw r > >1, that is, the distance between the Earth and the pulsar, and between the pulsars in the array, is large if compared to the wavelength of the GWs). It should be noted that the Hellings and Downs curve is computed using the Earth term only. The pulsar term is estimated to bring a significant contribution only at angular distances close to zero, and only if the pulsar pair is effectively close in space. In this case, the contribution of the pulsar term brings the angular correlation to 1. For pulsar pairs at even smaller angular distances, the contribution of the pulsar term becomes rapidly negligible (Mingarelli & Sidery Reference Mingarelli and Sidery2014). Figure 5 shows the Hellings and Downs curve, without taking into account the additional correlation that would occur at θi, j = 0 when considering the pulsar term.
Figure 5. Angular correlation C given by the Hellings & Downs curve as described by equation (17) (minus the contribution of the pulsar term) versus angular distance θ.
5.4. Aims and characteristics of Pulsar Timing Arrays
PTA experiments aim to detect signals that are angularly correlated across the sky, using the clock-like behaviour of an array of hyper-stable pulsars. The primary goal of PTAs is the detection of low-frequency GWs, and the most likely GW source to emit in the PTA band is coalescing SMBHBs. In this sense, PTAs can be considered as interferometer on Galactic scales, although instead of lasers, PTAs exploit the pulsed radio emission from the pulsars in the array.
To aid our detection prospects, we select pulsars with high rotational stability for PTA analysis (Shannon & Cordes Reference Shannon and Cordes2010). The most rotationally stable pulsars are millisecond pulsars (MSPs; Alpar et al. Reference Alpar, Cheng, Ruderman and Shaham1982). MSPs are pulsars that have been spun up via a transfer of mass and angular momentum by a companion star, which accelerates the neutron star to spin periods in the order of milliseconds (Bhattacharya & van den Heuvel Reference Bhattacharya and van den Heuvel1991). Following the mass transfer, both the magnetic field intensity and the spin-down rate are remarkably lowered. MSPs are much more stable (Verbiest et al. Reference Verbiest2009) than normal pulsars, and are characterised by timing residuals that are typically lower and whiter (Hobbs et al. Reference Hobbs, Lyne, Kramer, Martin and Jordan2004). As such, they are the only class of pulsars that are included in PTA monitoring campaigns.
The sensitivity of PTA experiments lies in the frequency range from approximatively 10−6 to 10−9 Hz. The two boundaries are due to the limits imposed by the Nyquist theorem, and are set by the observing cadence at the higher frequency (assumed to be once per month) and the total timespan at the lower frequency (assumed to be around 20 yrs.
No GW detection has yet been made by PTA experiments. However, the upper limits on the GW amplitudes estimated by PTAs have already given powerful insights in the models for Galaxy formation, aiding to exclude a fraction of them (Sesana Reference Sesana2013b; Shannon et al. Reference Shannon2015).
With the current sensitivities, and amplitude and rate predictions, it is unlikely that PTAs will detect individual SMBHBs in the near future, either in the form of continuous wave or in the form of a memory event (Babak et al. Reference Babak2016; Ravi et al. Reference Ravi, Wyithe, Shannon and Hobbs2015; Wang et al. Reference Wang2015). Concerning a GWB, as mentioned, although its amplitude is thought to be higher than that of individual SMBHB, the peak strain is expected to be very low. Sesana et al. (Reference Sesana, Shankar, Bernardi and Sheth2016), an update of Sesana (Reference Sesana2013b), estimated the spectrum of the GWB amplitude generated by a population of GW-driven, adiabatically inspiralling SMBHBs in quasi-circular orbits, and demonstrated that the majority of the GWB is due to major mergers, where the mass ratio between the two SMBHs is >0.25 within z = 1.5, and for black hole masses larger than 108 M⊙ (Sesana et al. Reference Sesana, Vecchio and Colacino2008). The study also assumes values from different studies available in literature to account for the SMBHB merger rates and masses. Combining the values from observational constraints, the authors generated more than 2 500 realisations of a GWB and computed a distribution for its amplitude. At 3σ, they predict 1.4 × 10−16 < A < 1.1 × 10−15 at 95% confidence. Such uncertainty mainly stems from poorly constrained estimates for the galaxy merger rate and the relation between the mass of the SMBH and the mass of the host galaxy. The influence of the SMBH-host relations is shown in Sesana et al. (Reference Sesana, Shankar, Bernardi and Sheth2016), where the authors compare the GWB predictions obtained by using two of these relations (Kormendy & Ho Reference Kormendy and Ho2013 and Shankar et al. Reference Shankar2016). Kormendy & Ho (Reference Kormendy and Ho2013) is claimed to be biased high, due to a number of overestimated SMBH masses obtained through dynamic measurements, while Shankar et al. (Reference Shankar2016) claims to have corrected the bias. The results, shown in Figure 6, indicate that the GWB predictions based on the two BH-host relations differ by a factor 3. This highlights the importance of the relations in these studies, and the necessity of refining them.
Figure 6. Figure taken from Sesana et al. (Reference Sesana, Shankar, Bernardi and Sheth2016)Footnote 3 ; GWB amplitude versus f gw. The plot shows a comparison between the GWB predictions as based on the BH-host relations from Kormendy & Ho (Reference Kormendy and Ho2013) (left panel) and Shankar et al. (Reference Shankar2016) (right panel). The sensitivity curves for EPTA, PPTA, and NANOGrav are shown in green, blue, and red, respectively, and the differently shaded area represents 99.7, 95, and 65% of probability.
GWs, both in the form of emission from single sources and a background, are not the only signals that can be angularly correlated among pulsars in a PTA. For example, imperfections in the reference time standards and in the planetary ephemeris used to identify the SSB would induce, respectively, a monopolar and a dipolar angular correlation in the timing residuals (Foster & Backer Reference Foster and Backer1990), and the creation of a pulsar-based time reference and the improvement of planetary ephemerides are ongoing projects within the framework of PTAs (Hobbs et al. Reference Hobbs2012; Champion et al. Reference Champion2010). Tiburzi et al. (Reference Tiburzi2016) studied the impact of correlated noise processes other than GWs (such as errors in time standards, planetary ephemeris, and unmodelled effects of instrumentation and the Solar wind) on PTA sensitivity to the GWB. The study demonstrated that, without including mitigation techniques in the detection pipelines, such signals can induce false detections (see also Taylor et al. Reference Taylor, Lentati, Babak, Brem, Gair, Sesana and Vecchio2017). Mitigation is feasible, especially for the monopolar signal, while the dipolar signal is more difficult to subtract without compromising the GWB search.
5.5. Current PTA collaborations
There are currently three well-established collaborations in the world that are leading PTA experiments: the European Pulsar Timing Array (EPTA, Desvignes et al. Reference Desvignes2016) in Europe, the Parkes Pulsar Timing Array (PPTA, Reardon et al. Reference Reardon2016; Manchester et al. Reference Manchester2013) in Australia, and the North American NanoHertz Observatory for Gravitational waves (NANOGrav, Arzoumanian et al. Reference Arzoumanian2015b) in the North America. EPTA, PPTA, and NANOGrav, all based on MSP observations with 100-m class radio telescopes, collaborate as the International Pulsar Timing Array (IPTA, Verbiest et al. Reference Verbiest2016).
5.5.1. EPTA
The EPTA was officially established in 2005, and currently monitors 42 MSPs (Desvignes et al. Reference Desvignes2016) at an approximatively monthly cadence with each of the five largest radio telescopes in Europe: the Effelsberg Radio Telescope (Eff, Germany), the Nançay Radio Telescope (NRT, France), the Westerbork Synthesis Radio Telescope (WSRT, the Netherlands), the Lovell Telescope at Jodrell Bank Observatory (JBO, UK), and the Sardinia Radio Telescope (SRT, Italy). In addition, a special programme within the EPTA, the Large European Array for Pulsars (LEAP), effectively acts as a sixth EPTA telescope (Bassa et al. Reference Bassa2016a).
Besides the dataset collected since its establishment, the EPTA uses archival data of MSPs dating back to the 1990s, and collected under different timing proposals with ‘historical’ backends and receivers (i.e., pulsar instruments now decommissioned).
The current observing setup of the EPTA telescopes is as follows:
-
• EFF. Performs coherently dedispersed observations of MSPs at three different frequencies with the PSRIX backend (Lazarus, Karuppusamy, Graikou, Caballero, Champion, Lee, Verbie Reference Lazarus, Karuppusamy, Graikou, Caballero, Champion, Lee, Verbiest and Kramer2016): 1360 MHz, 2639 MHz, and 4800 MHz.
-
• JBO. Observes MSPs with two backends in parallel, the DFB (incoherent dedispersion, see Manchester et al. Reference Manchester2013) and the ROACH (coherent dedispersion, see Bassa et al. Reference Bassa2016a), at 1532 MHz.
-
• NRT. Performs coherently dedispersed observations of MSPs in two frequency ranges, between 1100 and 1800 MHz, and between 1700 and 3500 MHz with the NUPPI backend (Liu et al. Reference Liu2014).
-
• WSRT. WSRT is currently unavailable for EPTA observations, as a new backend and frontend (ARTS and APERTIF) are being commissioned. The previous setup performed coherently dedispersed observations at 345, 1380, and 2273 MHz with the PuMa II backend (Karuppusamy, Stappers, & van Straten Reference Karuppusamy, Stappers and van Straten2008). The receiver at 345 MHz has been officially decommissioned.
-
• SRT. The first official EPTA observing run commenced in 2016, performing observations between 305 and 410 MHz and between 1300 and 1800 MHz (sometimes simultaneously) with a DFB and a ROACH backend.
-
• LEAP: Performs coherently dedispersed, interferometric observations of MSPs with the five EPTA telescopes at 1396 MHz, using ROACH backends, and collects dual-polarisation baseband data, that are then correlated offline.
5.5.2. PPTA
The PPTA project commenced in 2005, and currently monitors 24 MSPs (Reardon et al. Reference Reardon2016) with the Parkes Radio Telescope (NSW, Australia) every two to three weeks. In addition to observations obtained within the PTA programme, the PPTA uses data sets collected since the 1990s for other timing campaigns. Currently, the observations are carried out at three different frequencies: 3000, 1500, and 600 MHz, using a DFB (incoherent dedispersion) and the CASPSRFootnote 4 (coherent) backends.
5.5.3. NANOGrav
NANOGrav was officially established in 2007, and currently monitors 59 MSPs (a selection that has been expanded after Arzoumanian et al. Reference Arzoumanian2016) with the Arecibo Observatory (AO, Puerto Rico), the Green Bank Telescope (GBT, West Virginia, USA) and the Very Large Array (VLA, New Mexico, USA), every three or four weeks (Arecibo and Green Bank are also carrying out weekly observations of a subset of the monitored MSPs).
The current observing setup of the NANOGrav telescopes is as follows:
-
• AO. Performs coherently dedispersed observations of MSPs at 430, 1410, and 2030 MHz with the PUPPI backend (DuPlain et al. Reference DuPlain, Ransom, Demorest, Brandt, Ford, Shelton, Bridger and Radziwill2008).
-
• GBT. Performs coherently dedispersed observations of MSPs at 820 and 1500 MHz with the GUPPI backend (DuPlain et al. Reference DuPlain, Ransom, Demorest, Brandt, Ford, Shelton, Bridger and Radziwill2008).
-
• VLA. The newest addition to the NANOGrav programme, it observes MSPs between 1000 and 2000 MHz and between 2000 and 4000 MHz since 2017.
5.5.4. New PTA collaborations
Efforts to establish PTA experiment are ongoing in India, China, and South Africa.
The Indian PTA observes MSPs with the Ooty Radio Telescope (ORT) and the Giant Metrewave Radio Telescope (GMRT, both conventional and upgraded). In particular, conventional GMRT (that is timing nine MSPs together with ORT) has 32 MHz of bandwidth available in coherent dedispersion, while updated GMRT (that is timing 18 MSPs) has 200 MHz of bandwidth available in incoherent dedispersion (coherently dedispersion will be available in the near future). ORT has a central frequency of 334 MHz and can observe with coherent dedispersion within 16 MHz of bandwidth (M. Bagchi, private communication).
The Chinese PTA had an inaugural meeting in 2017 May. The Chinese PTA operates several 100-m class telescopes (e.g., NSRT, Kunming, Tianma), but the two most important facilities will be the Five hundred meter Aperture Spherical Telescope (FAST, Peng et al. Reference Peng, Nan, Su, Qiu, Zhu and Zhu2001) and the QiTai Radio Telescope (QTT). Once combined, FAST and QTT will be sensitive to a GWB amplitude of 2 × 10−16 within a few years of observations (Lee Reference Lee, Qain and Li2016).
MeerTIME is an approved proposal dedicated to pulsar timing, that will use the MeerKAT telescope (South Africa, Booth et al. Reference Booth, de Blok, Jonas and Fanaroff2009). MeerKAT is one of the numerous pathfinders for the Square Kilometer Array (SKA) (see Section 7), and is currently under deployment. Among the planned pulsars that will be observed with an increasingly larger number of antennas, are several PTA-relevant sources.
6 RECENT RESULTS FROM PTAS
6.1. EPTA
The current EPTA data set is comprised of 42 MSPs, and is presented in Desvignes et al. (Reference Desvignes2016). It includes updated timing solutions and ToAs spanning more than 15 yrs of data for many of the presented MSPs, besides deepening the astrometric properties of the sources. Caballero et al. (Reference Caballero2016) studied the red noise properties of the EPTA dataset, and found a significant level of red noise in 25 MSPs. Errors in the time standards were estimated to affect for at most 1% of the total noise budget, reducing the sensitivity to the GWB and resolvable SMBHBs.
The six most stable EPTA MSPs were used to derive the upper limits on the GWB amplitude and GWs from individual SMBHB, and to search for anisotropies in the GWB. Lentati et al. (Reference Lentati2015) computed a robust upper limit on the GWB amplitude of A < 3.0 × 10−15, taking into account the presence of other spatially correlated noise (Tiburzi et al. Reference Tiburzi2016). Babak et al. (Reference Babak2016) shows that the highest sensitivity to resolvable sources is reached by EPTA between 5 and 7 × 10−9 Hz, with a strain amplitude limit at 95% between 6 and 14 × 10−15.
Taylor et al. (Reference Taylor2015) assessed that the current EPTA dataset cannot constrain the angular distribution of the anisotropies yet, but their amplitude is 40% of the effect given by the isotropic GWB.
The dataset presented in Desvignes et al. (Reference Desvignes2016) and complemented with historical data was used to carry out individual-pulsar studies of MSPs J1024−0719 (Bassa et al. Reference Bassa2016b), J0613−0200 (McKee et al. Reference McKee2016), and J2051−0827 (Shaifullah et al. Reference Shaifullah2016). The EPTA project LEAP (Bassa et al. Reference Bassa2016a; Smits et al. Reference Smits2017) presented a single pulse analysis of MSP J1713+0747 (Liu et al. Reference Liu2016b), important to assess the impact of pulse jitter on timing precision.
6.2. PPTA
Reardon et al. (Reference Reardon2016) presented an extension to the first PPTA data release (Manchester et al. Reference Manchester2013), which included new timing solutions for 20 MSPs and their red noise analysis based on a new version of the Cholesky method (Coles et al. Reference Coles, Hobbs, Champion, Manchester and Verbiest2011). This study includes the first distance to a pulsar, MSP J0437−4715, measured to sub-parsec precision. A multi-frequency polarisation and spectral analysis of the PPTA MSPs was presented in Dai et al. (Reference Dai2015), finding deviations from the models commonly applied to study pulsar spectra and Faraday rotation in some of the pulsars.
Additional studies of pulse jitter (Shannon et al. Reference Shannon2014), extreme scattering events (Coles et al. Reference Coles2015), differences in measured positions between VLBI, and pulsar timing studies (Wang et al. Reference Wang2017), and variations in the pulse profiles (Shannon et al. Reference Shannon2016) were presented between 2014 and 2016.
The now-established technique of profile-domain pulsar timing (Lentati et al. Reference Lentati, Alexander, Hobson, Feroz, van Haasteren, Lee and Shannon2014) has been expanded to include the frequency evolution of pulse profiles (Lentati et al. Reference Lentati2017a) and the impact on the pulse profile due to the variable scattering effects of the IISM (Lentati et al. Reference Lentati, Kerr, Dai, Shannon, Hobbs and Osłowski2017c).
Shannon et al. (Reference Shannon2015) used the four most stable PPTA sources and placed the most constraining upper limit on the GWB amplitude to date, 1 × 10−15. Madison et al. (Reference Madison2016) developed a new technique to search for individual GW sources, without constrains on the waveform, and the PPTA dataset was searched for individual SMBHBs and memory bursts by, respectively, Zhu et al. (Reference Zhu2014) and Wang et al. (Reference Wang2015). No evidence for GWs was found, and the two studies placed upper limits on the amplitude of the two events. In the case of resolvable SMBHB, an upper limit of A < 1.7 × 10−14 was found at 10−8 Hz, while no burst events with an amplitude lower than 2 × 10−14 could have been detected in the PPTA dataset studied in Wang et al. (Reference Wang2015).
6.3. NANOGrav
The NANOGrav collaboration published its latest data release in 2015 (Arzoumanian et al. Reference Arzoumanian2015b), which included the timing solutions for 37 MSPs obtained from datasets spanning up to 9 yrs. New methods were developed to account for variable DM and profile evolution with frequency, and 10 pulsars were found to be affected by red noise. Fonseca et al. (Reference Fonseca2016) measured the Shapiro delay and masses for 14 MSPs in binary systems in the NANOGrav dataset, while Matthews et al. (Reference Matthews2016) studied the astrometry of the 37 sources, finding the velocity dispersions to be much smaller than for the general pulsar population.
Detailed analyses of the effects induced by the turbulent IISM were also carried out. Levin et al. (Reference Levin2016) analysed the scattering contribution, and concluded that the effect on the ToA errors due to variable multi-path propagation effects is negligible. Jones et al. (Reference Jones2017) studied the DM variations, finding incompatibility with a Kolmogorov spectrum (Armstrong, Rickett, & Spangler Reference Armstrong, Rickett and Spangler1995) in four of the pulsars, but the discrepancies can be explained by the presence of unaccounted trends in the data.
Noise analyses were conducted both on short (Lam et al. Reference Lam2016) and long timescales (Lam et al. Reference Lam2017), and five more pulsars in addition to those identified by Arzoumanian et al. (Reference Arzoumanian2015b) were found to be affected by red noise.
The 9-yr NANOGrav dataset was searched for GWs (Arzoumanian et al. Reference Arzoumanian2016). No evidence of a GWB was found, and an upper limit on the GWB amplitude was set at 1.5 × 10−15. The previous NANOGrav dataset (Demorest et al. Reference Demorest2013) was also searched for GWs from individual sources, in the form of continuous GW emission (Arzoumanian et al. Reference Arzoumanian2014) and memory bursts (Arzoumanian et al. Reference Arzoumanian2015a). No evidence for GWs was found, but upper limits were placed on the amplitude of continuous waves (A < 3.0 × 10−14 at 10−8 Hz) and for the occurrence rate of memory bursts depending on their amplitude (e.g., memory bursts with an amplitude larger than 4 × 10−14 at 6.2 yr−1).
6.4. IPTA
Verbiest et al. (Reference Verbiest2016) and Lentati et al. (Reference Lentati2016) presented the first IPTA data release, based on the combination of the EPTA, PPTA, and NANOGrav datasets for 49 MSPs (see Figure 7). The IPTA dataset consists of the ToAs time series, timing solutions, and noise models for the 49 sources. The noise analysis carried out by Lentati et al. (Reference Lentati2016) showed that the two main sources of red noise are variable DM, and intrinsic timing noise. However, these two sources of noise are often indistinguishable, due to a lack of multifrequency data.
Figure 7. Figure taken from Verbiest et al. (Reference Verbiest2016)Footnote 5 ; Aitoff projection of the IPTA MSPs. The two axes represent Galactic longitude (l) and latitude (b), while the stars represent the position of the IPTA MSPs. The crosses represent the MSPs that have been detected at radio wavelengths which are not part of a globular cluster present in the ATNF Pulsar Catalogue at the time of writing.
A basic search for a GWB was carried out in Verbiest et al. (Reference Verbiest2016), using all of the pulsars in the array. No evidence for GWs was found, and the IPTA placed an upper limit on the GWB amplitude of 1.7 × 10−15. This value, higher than the most stringent upper limit from PTA experiments (1 × 10−15, Shannon et al. Reference Shannon2015), is more constraining than that obtained by the individual PTAs. This indicates that the IPTA as a whole is more sensitive than the individual PTAs by at least a factor of two (Verbiest et al. Reference Verbiest2016).
7 FUTURE PROSPECTS
The main challenge of PTAs progresses lies in identifying and correcting for corrupting effects on the ToAs, many of which were of no importance until recently. As far as current efforts go, the bulk of the ongoing research is dedicated to study several long-period processes affecting the residuals—effects such as inaccuracies in the Solar System ephemerides, the IISM, intrinsic pulsar-timing noise, and instrumental instabilities. In addition to those, there are continuous efforts to increase the number of highly precise MSPs in the arrays, decreasing the levels of white noise, increasing observing baselines, cadence, and frequency coverage, improving analysis methods for multifrequency data, tackling previously intractable issues visible in unprecedented high-S/N data.
The new generation of radio telescopes that are now coming online will greatly increase our sensitivity to low-frequency GWs. The most sensitive instrument will be the SKA (Braun et al. Reference Braun, Bourke, Green, Keane and Wagg2015). The SKA will be built in Western Australia (low-band antennae) and in South Africa (mid- and high-band antennae), and will boost the sensitivity beyond the limits currently set by radio telescopes. The predicted probability of a GW detection after 5 yrs of observations with a SKA-based PTA (even without taking into account the current IPTA dataset, and assuming the original SKA design) is 50% (Janssen et al. Reference Janssen2015).
In preparation for the SKA, several pathfinders have been deployed, such as MeerKAT and the MWA (Western Australia, Tingay et al. Reference Tingay2013), LWA (New Mexico, USA, Ellingson, Clarke, Cohen, Craig, Kassim, Pihlstrom, Rickard Reference Ellingson, Clarke, Cohen, Craig, Kassim, Pihlstrom, Rickard and Taylor2009), and the LOw Frequency ARray (Europe, van Haarlem Reference van Haarlem2013), and they are proving to be vital to tackle some of the mentioned main challenges. For example, the low-frequency facilities such as LOFAR, LWA, and MWA, are fundamental instruments to monitor the turbulent IISM and its effects on pulsar timing, due to the frequency dependence of IISM effects on the propagation of radio waves. In particular, DM variations are one of the main sources of red noise in the ToA time series. IISM studies at low frequencies will be able to provide invaluable insights to improve the red noise models, and to disentangle the IISM contribution from intrinsic timing noise generated from instabilities in the pulsar spin.
ACKNOWLEDGEMENTS
The author is very grateful to Joris Verbiest, Golam Shaifullah, James McKee, Alberto Sesana, Chiara Mingarelli, and Dominik Schwarz, for their availability and helpfulness, for the useful conversations and proofreading.
