The Search for Gravitational Waves using Pulsar Timing

The Search for Gravitational Waves using Pulsar Timing R. N. Manchester CSIRO Astronomy and Space Science Sydney Australia Summary • Brief introduction to pulsars and timing • Detection of GW using pulsar timing • Pulsar Timing Array (PTA) projects • Current status and future prospects

What are pulsars? • Pulsars are rotating neutron stars radiating beams of emission that we see as pulses as they sweep across the Earth • Neutron stars are formed in supernova explosions by the collapse of the core of a massive star • Millisecond pulsars are old neutron stars that have been recycled by accretion of matter from a binary companion • Almost all known pulsars lie within our Galaxy • Neutron stars have a mass about 1. 4 MSun but a radius of only ~15 km • Their surface gravity is enormous, ~1012 g. Earth, hence their break-up rotation speed is very high, ~1 k. Hz

A Pulsar Census • Currently 2405 known (published) pulsars • 2162 rotation-powered Galactic disk pulsars • 243 in binary systems • 353 millisecond pulsars • 142 in globular clusters • 8 X-ray isolated neutron stars • 21 magnetars (AXP/SGR) • 28 extra-galactic pulsars Data from ATNF Pulsar Catalogue, V 1. 52 (www. atnf. csiro. au/research/pulsar/psrcat) (Manchester et al. 2005)

Pulsars as Clocks • Because of their large mass and small radius, NS spin rates - and hence pulsar periods – are extremely stable • For example, in 2001, PSR J 0437 -4715 had a period of : 5. 757451924362137 +/- 0. 00000008 ms • Although pulsar periods are very stable, they are not constant • Pulsars are powered by their rotational kinetic energy • They lose energy to relativistic winds and low-frequency electromagnetic radiation (the observed pulses are insignificant) • Consequently, all pulsars slow down (in their reference frame) • Typical slowdown rates are less than a microsecond per year • For millisecond pulsars, slowdown rates are ~105 smaller

Measurement of pulsar periods • Start observation at a known time and average 103 - 105 pulses to form a mean pulse profile • Cross-correlate this with a standard template to give the arrival time at the telescope of a fiducial point on profile, usually the pulse peak – the pulse time-of-arrival (To. A) • Measure a series of To. As over days – weeks – months – years • Transfer To. As to an inertial frame – the solar-system barycentre • Compare barycentric To. As with predicted values from a model for the pulsar – the differences are called timing residuals. • Fit the observed residuals with functions representing errors in the model parameters (pulsar position, period, binary period etc. ). • Remaining residuals may be noise – or may be science!

Sources of Pulsar Timing “Noise” Ø Intrinsic noise • Period fluctuations, glitches • Pulse shape changes Ø Perturbations of the pulsar’s motion • Gravitational wave background • Globular cluster accelerations • Orbital perturbations – planets, 1 st order Doppler, relativistic effects e d i w a f o • Wind from binary companion s e b o r a • Variations in interstellar dispersion ful p n e r m o n • Scintillation effects powe e h p e r l a a c s Ø Perturbations of the Earth’s motion i r s a y s l h p o P • u. Gravitational wavef background r t s a • Errors in the Solar-system ephemeris o e g n a Ø Clockrerrors Ø Propagation effects • Timescale errors • Errors in time transfer Ø Instrumental errors • Radio-frequency interference and receiver non-linearities • Digitisation artifacts or errors • Calibration errors and signal processing artifacts and errors Ø Receiver noise

. The P – P Diagram P = Pulsar period. P = d. P/dt = slow-down rate . • For most pulsars P ~ 10 -15. • MSPs have P smaller by about 5 orders of magnitude • Most MSPs are binary, but few normal pulsars are. • τc = P/(2 P) is an indicator of pulsar age • Surface. dipole magnetic field ~ (PP)1/2 Great diversity in the pulsar population! (Data from ATNF Pulsar Catalogue V 1. 50) Galactic Disk pulsars

Pulsars and Gravitational Waves • Pulsars have already shown that GW exist • Double-neutron-star systems emit GW at 2 forb • Energy loss causes slow decrease of orbital period • Predict rate of orbit decay from known orbital parameters and NS masses using GR • For PSR B 1913+16, ratio of measured value to predicted value = 0. 997 +/- 0. 002 ØFirst observational evidence for gravitational waves! PSR B 1913+16 Orbit Decay (Weisberg , Nice & Taylor 2010)

The Double Pulsar R: Mass ratio w: periastron advance g: gravitational redshift r & s: Shapiro delay. Pb: orbit decay WSO: Spin-orbit coupling • Six measured relativistic parameters(!) + mass ratio • Five independent tests of GR • Now limits deviations from GR to 0. 02%. • Pb, GW now more precise than for PSR B 1913+16 PSR J 0737 -3039 A/B (Kramer et al. 2015)

Pulsar Timing Arrays (PTAs) • A PTA consists of many pulsars widely distributed on the sky with frequent high-precision timing observations over a long data span • Can in principle detect gravitational waves (GW) Ø GW passing over the pulsars are uncorrelated Ø GW passing over Earth produce a correlated signal in the TOA residuals for all pulsars • Most likely source of GW detectable by PTAs is a stochastic background from super-massive binary black holes in distant galaxies • PTAs are sensitive to long-period GWs, f. GW ~ 1/(data span) ~ few nanohertz – complementary to laser-interferometer systems • Detection requires precision timing of ~20 MSPs over ~10 years • PTAs can also detect instabilities in terrestrial time standards – establish a pulsar timescale Idea first discussed by Hellings & Downs (1983), Romani (1989) and Foster & Backer (1990)

Correlated Timing Signals Ø Clock errors All pulsars have the same TOA variations: Monopole signature Ø Solar-System ephemeris errors Dipole signature Hellings & Downs GW correlation function Ø Gravitational waves Quadrupole signature Can separate these effects provided there is a sufficient number of widely distributed pulsars (Hobbs et al. 2009)

Major Pulsar Timing Array Projects Ø European Pulsar Timing Array (EPTA) • Radio telescopes at Westerbork, Effelsberg, Nancay, Jodrell Bank, Cagliari • Timing ~40 MSPs, 5 - 18 years (Kramer & Champion 2013) Ø North American pulsar timing array (NANOGrav) • Data from Arecibo and Green Bank Telescope • Timing ~30 MSPs, data spans 3 - 28 years (Mc. Laughlin 2013) Ø Parkes Pulsar Timing Array (PPTA) • Data from Parkes 64 m radio telescope in Australia • Timing 22 MSPs, data spans 3 – 19 years (Manchester et al. 2013)

MSPs suitable for PTAs All (published) MSPs not in globular clusters

PPTA Three-Band Timing Residuals • Regular observations in three bands: 50 cm ~700 MHz 20 cm ~1400 MHz 10 cm ~3100 MHz 50 cm 10 cm 20 cm • Frequency-dependent variable delays: due to dispersion variations resulting from ISM density fluctuations • Delay Δt ~ ν-2 (Manchester et al. 2013)

PPTA Timing Residuals (DR 1) • Timing data for 22 pulsars • Data spans to 19 years • Best 1 -year rms timing residuals about 40 ns, most < 1μs • Low-frequency (red) variations significant in about half of sample – some due to uncorrected DM variations in early data, e. g. , J 1045 -4509 • Some of the best-timing pulsars nearly white (e. g. , J 1713+0747, J 1909 -3744) (Hobbs 2013)

Effect of GW on pulsar timing residuals Redshift z of observed pulsar frequency is difference between the effects of GW passing over the pulsar and the GW passing over the Earth: where the GW is propagating in the +z direction and α, β, γ are direction cosines of the pulsar from Earth. is the difference in GW strains at the pulsar and the Earth (A = polarisation). is the “pulsar is the “Earth term” (d is the pulsar distance) and term”. The observed timing residuals are the integral of the net redshift: (Detweiler 1979, Anholm et al. 2009)

Pulsar timing “antenna pattern” • Timing response is maximum for a pulsar (roughly) in the direction of the GW source • Response is zero for pulsars exactly in both source and propagation direction • For different pulsars of a PTA, the pulsar terms are uncorrelated, both because of the retarded time and the different GW environment, but the Earth terms are correlated. • For an isotropic GW background, many different GW sources – have to integrate antenna pattern over different GW propagation directions • Hellings & Downs (1983) showed that the correlated Earth-term response was a function only of the angle θ between pulsars: z GW source (Chamberlin & Siemens 2012) where • Uncorrelated pulsar term produces “noise” on HD curve

The Isotropic GW Background • Strongest source of n. Hz GW waves is likely to be a background from orbiting super-massive black holes in distant galaxies Rate of energy loss to GW: Chirp mass: Frequency in source frame: • Number of mergers per co-moving volume based on model for galaxy evolution (e. g. Millennium simulation) plus model formation and evolution of SMBH in galaxies (Sesana 2013)

• Simple parameterisation of characteristic strain for cosmological distribution of circular BH binaries: A is the characteristic strain at a GW frequency f = 1/(1 yr) • GW modulates the observed pulsar frequency • Modulation spectrum of observed timing residuals: GW power ~ f -13/3 • Spectrum is very steep – Cholesky spectral estimation methods needed to avoid leakage (Coles et al. 2011) • Strength of GWB often expressed as a fraction of closure energy density of the Universe: (Phinney 2001; Jenet et al. 2006)

Distribution of SMBBH contributing to the GWB • Most of background from SMBBH in galaxies at z = 1 – 2 • Biggest contribution from largest BH masses: 108 – 109 Msun (Sesana et al. 2008)

PPTA Limit – Shannon et al. 2013 • Based on six best PPTA pulsars + archival Arecibo data for one pulsar • Data spans 5 – 25 years – effective f. GW ~ 2. 8 n. Hz • 95% confidence limit based on detection statistic + simulations with GWB injection: 95% Confidence Limits Probability of a GW signal in the PPTA data (Millenium sim) Hydro sim. (Kuiler et al. 2013) (Sesana 2013) (Mc. Williams et al. 2012) A < 2. 4 x 10 -15 Ωgw(f. GW) < 1. 3 x 10 -9 (Shannon et al. , 2013, Science)

Single Sources • Likely that many SMBH binary systems are highly eccentric • GW spectrum may be dominated by a strong individual source First GW detection by PTAs could be a single source with period <~ 1 year! (Sesana 2012)

PPTA Limits on Single GW Sources • Searched PPTA data set for arbitrary time series of GW amplitudes A+, Ax that satisfy the quadrupolar constraints for a given source direction • Searched these time series for sinusoidal components that could result from orbital motion of a SMBHB • Repeated for all 1 possible source directions yr 2 yr -1 -1 Spectral strain sensitivity at different sky locations Distance (Mpc) to which a SMBBH with chirp mass 109 MSun emitting GW at 10 -8 Hz can be detected (Zhu et al. 2014)

Localisation of a GW Source • PPTA data set • Simulated source in Fornax cluster, S/N ~ 5 (100 ns amplitude) • Correlation statistic, Earth term only • Max detection statistic (o) close to source position • Localisation poor ~ 400 sq deg • In principle, could make use of pulsar term as well – pulsar-Earth interferometer • Fringe spacing θf = λGW/[2πd(1+γ)] – few arc minute (γ is z direction cosine) • Localisation ~ θf /(S/N) – sub-arc-min possible • However, pulsar distances must be known to better than λGW/(1+γ) • So far, only one pulsar (J 0437 -4715) has sufficiently accurate dist. (d = 156. 79± 0. 25 pc) (Reardon et al. 2015) (Zhu et al. (2015) (Lee et al. 2011)

GW “Memory” Events • Coalescence of two super-massive black holes creates an offset in space-time which propagates away at the speed of light • When such an offset crosses the pulsar or the Earth it causes a step frequency change (glitch) in the observed pulsar period • Events crossing the Earth have the same quadrupolar signature as a normal GW and so can be detected using similar techniques • PPTA data set used to limit amplitude and rate of GWM events • Long way from predictions! • Detection will take 50 MSPs, 50 ns To. A precision, 25 yr data span – needs FAST and SKA (Wang et al. 2014) PPTA data set 20 -yr “PPTA” data set Predictions based on GWB models

A New Limit on the Nanohertz GW Background • Based on PPTA 10 cm data sets for the J 1909 -3744 four best-timed pulsars • 11 -year data span – 3. 5 yrs longer than J 0437 -4715 Shannon et al. (2013) PPTA data set • Inclusion of 20 cm and 50 cm data J 1713+0747 reduces the sensitivity – DM corrections J 1744 -2304 are not sufficiently accurate Also fixed some small instrumental phase jumps • Bayesian analysis method • New limit: A < 1. 0 x 10 -15 (Shannon et al. , Science, submitted) -10 -1 ΩGW < 2. 3 x 10 at f = 0. 2 yr (6 n. Hz) Essentially rules out all current models for GWB!

Comparison with GWB Models • Observational limits (95% conf): N 13: NANOGrav (Demorest et al. ) E 15: European PTA (Lentati et al. ) P 13: Parkes PTA (Shannon et al. ) P 15: PPTA 10 cm (Shannon et al. ) • GWB models: K 15: Kulier et al. ✖ (99. 5%) S 13: Sesana ✖ (91%) M 14: Mc. Williams et al. ✖ (99. 8%) M 14 e: Mc. Williams et al. ✖ (96%) R 14: Ravi et al. (✔) (20%) R 15: Ravi et al. ✖ (95%) • Various different assumptions about galaxy merger rates and SMBBH formation and evolution (Shannon et al. , Science, submitted) • M 14 e and R 14 include strong gas/stellar-interaction-driven evolution at late times

Scenarios for SMBBH Evolution • Fiducial: Standard model • Slow merger: Dynamical friction x 10 -4 • Recoil: 500 km s-1 kick after merger – SMBBH escapes galaxy • Stalled: stellar/gas interaction x 10 -2 (“last parsec” case) • Environment: stellar/gas (Shannon et al. 2015) interaction at late times x 102 • GW: Phase of evolution where GW emission dominates over other loss processes Slow merger, recoil, or stalled merger unlikely to be important Most probable reason for lack of detection is orbit decay due to stellar/gas interactions at late times – less energy into GW Also implies GWB detectability is greater at higher frequencies (f < ~1 yr-1) – change in observational strategy needed

GW Detection Sensitivity • Simplified model, based on correlation analysis of M pulsars in PTA • White timing noise only, all pulsars with same rms noise σ • At low signal levels (S/N < 1): S/N ~ MT 13/3/σ2 • At higher signal levels (S/N > 1): S/N ~ MT 1/2/σ3/13 - effect of uncorrelated pulsar term (Siemens et al. 2013) For GWB detection, increased observation time and higher precision To. As don’t help so much Most effective strategy is to increase M! Real life is more complicated – pulsars have different σ, cadence and data spans, complex GW spectra, etc

The International Pulsar Timing Array • The IPTA is a consortium of consortia, formed from existing PTAs • Currently three members: EPTA, NANOGrav and PPTA • Aims are to facilitate collaboration between participating PTA groups and to promote progress toward PTA scientific goals • The IPTA organises annual Student Workshops and Science Meetings – 2015 meetings will be in Australia, 20 -31 July • The IPTA has organised Data Challenges for verification of GW detection algorithms • Main current focus is making available the combined data sets. Difficult because of different instruments and signal processing methods in different PTAs

IPTA Pulsars • 50 MSPs with good timing data among the three PTAs • 12 observed by all three PTAs, 10 by two PTAs

IPTA Data Sets: 50 MSPs Bands Black: 70 cm Red: 50 cm Green: 35 cm Blue: 20 cm Aqua: 15 cm Red: 10 cm 1984 2013

Future Developments • Improved wideband feeds and receivers for existing telescopes – higher sensitivity and better DM corrections • Improved signal processor systems – better RFI rejection and coherent dedispersion over wider bandwidths ATNF 0. 7– 4. 0 GHz feed design • Algorithmic development – optimal combination of data sets, detection algorithms for GW background and single sources, separation of correlated processes, etc. • Improved predictions of signals affecting pulsar timing data • Bigger telescopes! – more MSP discoveries and improved timing, especially for weaker sources

LOFAR MWA Soon FAST Now SHAO 65 m Now 2017 SKA Mid-Frequency Array Meer. KAT 2016 2020+

PTAs sensitive to nanohertz GW Complement other GW detection methods Current PTAs working at design sensitivity – but no detection yet The time till detection depends on many factors, with the most uncertain being the predicted GW signal levels

The Gravitational Wave Spectrum

Parkes UWB Receiver Design • New feed design with dielectric loading and 4 -port balanced output – excellent polarisation performance and cooled to reduce losses • Preamplifiers based on MMIC technology for high sensitivity – expect 20 K system temperature across most of band • Triplexer – split each polarisation band into three octave sections to improve linearity and avoid spurious second-harmonic responses • Antenna(s) and amplifier chain(s) for RFI reference signals • 5 GSamp/sec 12 -bit digitisers for receiver and RFI reference signals – in focus cabin for stability – needs triple-shielded RFI enclosure • Single-mode fibers to Tower receiver room • FPGA band combiner, RFI mitigation and CPU/GPU distributor 38 The Parkes UWL Receiver - XAO 2014 | Dick Manchester

DM Correction • Observed To. As are sum of frequency-independent “common-mode” terms t. CM (e. g. , clock errors, GW, etc) and interstellar delays t. DM – assume ~ λ 2 • The interstellar term t. DM is noise – want to minimise it Ø Observe at ~zero wavelength, i. e. , X-ray or γ-ray Ø Observe at two or more wavelengths, λ 1 and λ 2 (with λ 1 > λ 2) • Can then solve for t. DM and t. CM:

Effect of CM Term • If CM term not included in fit, power is extracted from freq-independent variations and coupled into DM variations • With CM term included, all freq-independent power (e. g. , GW signal, clock errors) is contained in CM values No CM Pre-fit (Wh+GW+DM) Pre-fit (Wh+GW) Post-fit 1 yr-1 10 yr-1 CM incl. (Keith et al. 2013) Power spectra of timing residuals

CASS UWL Feed • 0. 7 – 4. 0 GHz, dielectrically-loaded quad-ridge horn Design • Aperture efficiency 0. 70 +/- 0. 05 across whole band • Cross-polar coupling less than -20 db across whole band • Feed/OMT to ring slots cooled to 70 K – feed/OMT loss 1 -3 K in Tsys 0. 7 GHz 1. 7 GHz 4. 0 GHz Design by Alex Dunning 41 | The Parkes UWL Receiver - XAO 2014 | Dick Manchester

MSP Polarisation • Grand-average profiles for PPTA pulsars, very high S/N • New profile features, complex polarisation properties (Yan et al. 2010)

TT(PPTA 2011) – Relative to TAI • PPTA extended data set, 19 pulsars • Timing referenced to TT(TAI) • Clock term sampled at ~1 yr intervals • Constrained to have no quadratic or annual terms • Compared with BIPM 2010 – TAI with quadratic removed PPTA BIPM 2010 First realisation of a pulsar timescale with stability comparable to that of current atomic timescales! (Hobbs et al. 2012)