Detection of frequency-dependent dispersion measure toward the millisecond pulsar PSR J2241-5236 from contemporaneous wide-band observations

Making precise measurements of pulsar dispersion measures (DMs) and applying suitable corrections for them is amongst the major challenges in high-precision timing programmes such as pulsar timing arrays (PTAs). While the advent of wide-band pulsar instrumentation can enable more precise DM measurements and thence improved timing precision, it also necessitates doing careful assessments of frequency-dependent (chromatic) DMs that was theorised by Cordes et al. (2016). Here we report the detection of such an effect in broadband observations of the millisecond pulsar PSR J2241-5236, a high-priority target for current and future PTAs. The observations were made contemporaneously using the wide-band receivers and capabilities now available at the Murchison Widefield Array (MWA), the upgraded Giant Metrewave Radio Telescope (uGMRT) and the Parkes telescopes, thus providing an unprecedentedly large frequency coverage from 80 MHz to 4 GHz. Our analysis shows the measurable changes in DM that scale with the observing frequency ($\nu$) as $\rm \delta DM \sim \nu^{2.5 \pm 0.1}$. We discuss the potential implications of such a frequency dependence in the measured DMs, and the likely impact on the timing noise budget, and comment on the usefulness of low-frequency observations in advancing PTA efforts.


INTRODUCTION
Pulsar timing array (PTAs) rely on high-precision measurements of pulse arrival times for a large number of pulsars over long time spans for the eventual detection of low-frequency (nanohertz) gravitational waves (GWs; e.g., Foster & Backer 1990;Manchester et al. 2013;Demorest et al. 2013;Bailes et al. 2018;van Haasteren et al. 2011;Janssen et al. 2015). To reach the sensitivity required to detect nanohertz GWs (Cordes & Shan-(uGMRT;120-1450 MHz;Gupta et al. 2017), is highly promising for improved timing precision in PTA observations. However, these large bandwidths also necessitate additional considerations; for example, correcting for the interstellar medium effects on pulsar signals such as dispersion, scintillation and multi-path scattering, all of which scale strongly with the observing frequency.
The dispersive delay is proportional to the integrated column density of free electrons along the line of sight of the pulsar. The time delay for a pulse observed at radio frequency ν, compared to a (hypothetical) pulse at infinite frequency can be approximated as t DM = D DM/ν 2 , where D = e 2 /(2πm e c) is the dispersion constant, e and m e are the charge and mass of the electron, c is the speed of light, and DM is the dispersion measure, traditionally defined as the integral of the electron density (n e ) along the line of sight (LOS). For PTAs, it is a time-varying quantity, primarily due to the pulsar's large space velocities (typically ∼ 50-100 km s −1 ; Hobbs et al. 2005), as a result of which different parts of the interstellar medium (ISM) are sampled by the observations. Since spatial variations in electron density are present on a wide range of scales (∼ 10 6 to 10 13 m or even larger; Lam et al. 2015;Cordes & Rickett 1998), PTA observations require measuring and correcting for DM at every observing epoch. Measuring DMs accurately and applying suitable corrections for their temporal variability has been an integral part of pulsar timing array experiments (Jones et al. 2017;Lam et al. 2016;Keith et al. 2013;Lee et al. 2014;You et al. 2007).
Aside from the time variability in DM, there may also be a frequency dependence; sometimes referred to as 'chromaticity' in DM. The importance of such a subtle effect was first acknowledged by Ramachandran et al. (2006) in their analysis of DM variations of the firstdiscovered millisecond pulsar (MSP) B1937+21. More recently, Cordes et al. (2016) provided a detailed theoretical account of the underlying physics, and discussed possible implications in the wider context of PTAs. The effect arises as a result of multi-path scattering by the electron density fluctuations present in the ISM. Since the radiation received at a given frequency depends on the size of the scattering disk, which scales as the square of the observing wavelength (θ s ∝ λ 2 ), the volume of the ISM (or path lengths) sampled at lower frequencies can be significantly larger than those at higher frequencies, especially when observations are made over very large bandwidths (hundreds of MHz to a few GHz). Such a frequency dependence in DM can give rise to subtle differences in DMs as measured in observations, depending on the observing frequency and the bandwidths over which such measurements are made. For example, as per the formalism presented in Cordes et al. (2016), we may expect an rms DM variation ∼ 10 −5 pc cm −3 for a pulsar at a distance of ∼1 kpc and DM 30 pc cm −3 , which may result in a few hundred nanoseconds of timing noise at the standard timing frequency bands (∼ 1-2 GHz). This can be of the order of up to a few microseconds for higher DM pulsars and over wider observing bandwidths . The investigation of such subtle effects are particularly important in the era of wideband precision pulsar timing, when the observations are routinely made over large bandwidths of the order of hundreds of MHz to a few GHz.
To date, there have been limited observational investigations to study this effect. In their analysis of 20 yr timing data on the MSP B1937+21, Ramachandran et al. (2006) interpreted parts of their analysis (i.e. observations at 327 and 610 MHz) in terms of frequency dependent variations in DM. More recently, Donner et al. (2019) studied observations taken with three German LOng-Wavelength (GLOW) stations, which are part of the LOw-Frequency ARray (LOFAR). They presented 3.5 years of weekly observations of PSR J2219+4754, a long-period pulsar with a DM of 43.5 pc cm −3 . The frequency-dependent DM trends seen towards this pulsar (at 100-200 MHz) were attributed to extreme scattering events (ESEs) caused by a discrete cloud of size ∼20 AU and n e ∼ 10 cm −3 . These results were revisited by Lam et al. (2020), who concluded that the DM variations are due to the turbulent ISM and no ESE occurred along the LOS. They also comment on the suitability of long-period pulsars for investigating frequencydependent DM, given their distinct emission characteristics compared to those of millisecond pulsars that are used for timing-array experiments.
In this paper, we present the measurements of frequency-dependent DMs toward the millisecond pulsar PSR J2241−5236. Its low DM (11.41151 pc cm −3 ) and short pulse period (2.18 ms), along with its brightness and low levels of timing noise make it an important target for PTAs (Keith et al. 2011;Dai et al. 2015;Kaur et al. 2019;Parthasarathy et al. 2021). The pulsar is in a 3.5 hour, almost-circular orbit with a low-mass (∼ 0.012 M ), black-widow type companion, but with no eclipses seen (Keith et al. 2011). For this work, observations were carried out contemporaneously using the Murchison Widefield Array (MWA; Tingay et al. 2013;Wayth et al. 2018), the upgraded Giant Metrewave Radio Telescope (uGMRT), and the new ultra-wide band low-frequency receiver (UWL) at the Parkes radio telescope, thus providing frequency coverage from 80 MHz to 4 GHz. The remainder of the paper is organised as follows. In § 2 we describe the details of our observa-tions. In § 3 and § 4 we summarize data processing and analysis, and in § 5 we present our results. Our conclusions are summarized in § 6.

OBSERVATIONS
Observations of PSR J2241−5236 were made contemporaneously (i.e. within a ∼24-hour duration) using the MWA, uGMRT and Parkes. This pulsar is ideal for our analysis since it has a narrow pulse profile and very little profile evolution across the observed frequency range, as discussed below. The pulsar was observed at multiple (three) epochs with all three telescopes. At two of the epochs, it was observed within a 24 hour block using all three telescopes, while at one of the epochs, only two telescopes (uGMRT and Parkes) were available. Observational details are summarised in Table 1, and are further elaborated below.

The MWA
Observations were made over a total of four epochs, two epochs to sample the full orbit (described in Section 2.1.1) and another two contemporaneously (separated by three weeks) with other telescopes, using the MWA's voltage capture system (VCS; Tremblay et al. 2015). The VCS records 24 (1.28 MHz wide) coarse channels, each of which has been finely channelized to 10 kHz. These data were recorded in both polarizations, from all 128 tiles. The 30.72 MHz observing bandwidth of the MWA can be flexibly placed (e.g., 12 × 2.56 MHz coarse channels) anywhere across the nominal operating frequency range (80-300 MHz). For the observations presented in this paper, a frequency set-up similar to that described in Kaur et al. (2019) was employed, i.e., a distributed 24 × 1.28 MHz channel mode, covering a frequency range from 80 MHz to 220 MHz.

Orbital coverage
Since the pulsar is in a binary system, it was important to examine the data for any orbital DM variations. The ability to sample a full orbit of the pulsar in a single day can help to minimise any temporal variations in DM caused by the ISM. With the VCS, the maximum recording time is however limited to ∼ 90 minutes, which results in a data volume of ∼ 50 TB (data rate = 7.78 GB s −1 ). The data were recorded in multiple 10 minute recordings at ∼15 minute intervals between the recordings. This allowed us to optimally use the available 90 minutes of recording time to sample the full binary orbital period of 3.5 hours. A "picketfence" mode of observation was adopted in order to achieve a large frequency coverage, but with a modified frequency-channel selection where the 30.72 MHz bandwidth was divided into five uneven sub-bands of 2×7.68 MHz and 3×5.12 MHz centred at 83.84 MHz,127.36 MHz,155.52 MHz,187.52 MHz,and 218.24 MHz respectively. This strategy was chosen to increase the sub-bandwidth in order to compensate for the reduced signal-to-noise ratio (S/N) associated with breaking up the recording time into 10 minute observations. The main purpose was to ensure a sufficiently high S/N to detect the pulsar in each of the sub-bands in 10 minute recordings (as compared to the typical observing time of ∼ 1 hour for this pulsar that was used for the work presented in Kaur et al. 2019). In summary, by making optimal use of available resources in bandwidth and time, the full 3.5 hour orbit was sampled in a total of eight 10 minute recording sessions. Note-∆t is the time resolution, T obs is the total observing duration, Nscan is the number of independent observations (scans) each of T obs /Nscan duration, N chan is the number of channels (sub-bands) used to measure the TOAs, and f template is radio frequency of the observations from which an analytic template was derived.

The uGMRT
We conducted observations using the upgraded GMRT (uGMRT; Gupta et al. 2017). Two observations were made in Band 3, which covers from 300 to 500 MHz, and one observation in Band 4, which covers 550 to 750 MHz. Observations were made using the phased array total intensity mode. Depending upon the observing epoch, the number of antennas used varied from 24 to 28. The data were then coherently combined to generate spectral (channelised) voltages time series data. These data were processed in real-time using the coherent dedispersion processing pipeline of the uGMRT wideband backend (GWB; Reddy et al. 2017). For all observations, data were recorded in 512×0.390 MHz coherently-dedispersed filter-bank format with 10.24 µs time resolution. In order to avoid de-phasing due to ionospheric effects or antenna gain changes, we performed re-phasing of the array periodically, at an average interval of ∼ 30 minutes, using a nearby phase calibrator 2225−049, which has a flux density of 15 Jy at P-band (400 MHz). The incoming data were then processed by the GWB, and data were recorded with 256 phase bins across the pulse period, for both the uGMRT Band 3 and Band 4 observations.

Orbital coverage
Given the far southern declination of the pulsar, and the geographic location of the uGMRT, PSR J2241−5236 is visible for a maximum of 2.1 hours on any particular day. This covers approximately 60% of the orbit. In order to span the full 3.5 hour orbit, two successive observations were carried out, five days apart. The first observation on MJD 58794 covered the orbital phase range from ∼ 0.6P b to 0.8 P b , where P b is the orbital period, while the second observation on MJD 58799 covered from ∼ 0.8 to 1.0 and 0.0 to ∼ 0.3 phase ranges. The two observations thus partially covered the orbit, with no coverage for ∼ 0.3 P b to ∼ 0.6 P b in the orbital phase range. Data recording was interrupted for a short duration (∼ 5 minutes) to accommodate the phase calibrators and re-phasing of the array.

Parkes
Pulsar J2241−5236 is regularly monitored at Parkes as part of the ongoing Parkes pulsar timing array (PPTA; Manchester et al. 2013) project. The use of the ultrawideband low-frequency receiver (UWL; Hobbs et al. 2020) provides simultaneous frequency coverage from 704 to 4032 MHz. It records signals in three adjacent radio frequency bands; 704-1344 MHz, 1344-2368 MHz and 2368-4032 MHz. The entire band is sampled by digitisers for both polarisations. Digitised data from the telescope's focus cabin are transported to the signal pre-processing units where the entire band is channelised into 26 contiguous sub-bands, each 128 MHz wide, and for each polarisation.
Parkes observations used the pulsar fold mode where data are synchronously folded modulo the pulse period, with 1024 phase bins in each of the 1 MHz channels (and hence a time resolution ≈2.128 µs), and integrated for 30 s before written out to disk. Data were coherently dedispersed at a DM of 11.41151 pc cm −3 , and the Stokes parameters were recorded.

The MWA
The MWA makes use of the VCS and the postprocessing chain (Bhat et al. 2016;McSweeney et al. 2017). Data calibration and processing were performed on the Galaxy cluster at the Pawsey Supercomputing Centre. The signals from each tile can be either incoherently summed, or combined coherently, to generate a phased-array beam. This process incorporates a model of the polarimetric response and complex gains of each of the tiles, including both cable and geometric delays (Ord et al. 2019). The phase model for the array was determined (the data are calibrated using this model as described earlier) using one of the standard calibrators (e.g., 3C444), recorded in pointed observations prior to pulsar observations. Amplitude and phase calibration solutions were generated for each frequency sub-band using the real time system (RTS; Mitchell et al. 2008). A coherent tied-array beam is produced by phasing up of the signals from all tiles with good calibration solutions, essentially following the procedures as described in Bhat et al. (2016) and Ord et al. (2019). The high time resolution (≈ 1 µs) data were recovered using an enhanced capability of the beamformer software, which performs an inversion of the polyphase filterbank operation, the implementation of which is detailed in McSweeney et al. (2020). Coherently de-dispersed data were then generated using the DSPSR pulsar package (van Straten & Bailes 2011) and average profiles were written in the PSRFITS (or Timer) file format with 0.78 µs time resolution (see also Kaur et al. 2019). All the processing, including calibration and beamforming, was carried out at the Pawsey supercomputing facility.
On MJD 58671 and 58680, standard calibrator data were not available; therefore, for observations made on these dates, we adopted a somewhat non-standard approach. Calibration solutions were constructed by selectively choosing parts of multiple different calibration observations that were recorded using frequency configurations that are not identical to the configuration of our pulsar observations. Even though the related proce-dure was tedious, the data were successfully calibrated, resulting in pulsar detections with signal-to-noise ratios that are comparable to those obtained with a standard calibration process, thereby demonstrating the efficacy of our approach.

The uGMRT
The uGMRT is located in a radio-frequency interference (RFI) hostile environment, because of which data are typically RFI contaminated. Although online RFI mitigation is possible with the uGMRT system, that was not used for our observations. All RFI excision analysis was performed offline using the pazi routine of PSRCHIVE by identifying and flagging the problematic channels. uGMRT Band 3 data are typically more prone to RFI, resulting in ∼ 20% excised data, whereas only ∼ 10% data were excised in Band 4.

Parkes
Parkes data reduction closely followed the procedures described in Hobbs et al. (2020) and Dai et al. (2019). In short, we removed 5 MHz of the bandpass at each edge of the 26 sub-bands to mitigate aliasing. The "Coastguard" package was used to automatically excise RFI by examining data in both time and frequency. Pulsed noise signal recorded before each observation, together with observations of the radio galaxy 3C 218 and PSR J0437−4715, were used to calibrate the differential and absolute gains and polarimetric response of the receiver (cf. van Straten 2004). After the calibration, data were further RFI cleaned using the pazi routine of PSRCHIVE by carefully (manually) flagging the problematic channels. The calibrated and RFI cleaned data were then integrated in time for further timing and DM analysis.

Precursor emission at low frequencies
Some millisecond pulsars (MSPs) are known for their remarkable pulse profile evolution (e.g., Dai et al. 2015;Bhat et al. 2018). Unlike other bright southern MSPs (such as PSRs J0437−4715 and J2145−0750), PSR J2241−5236 has a narrow pulse profile and shows very little profile evolution from 100 MHz to 4 GHz, as described earlier. As seen from Figure 1, overall, the main profile shows very minimal evolution above 500 MHz. At frequencies below 500 MHz, there is clear evidence of the presence of additional precursor components, which were first seen at frequencies below 300 MHz with the MWA (Kaur et al. 2019). Our observations in the uGMRT Band 3 have now unambiguously confirmed the presence of this evolving precursor emission, as evident in Figure 1.  Based on our previous results and a non-detection of precursor emission at Parkes frequencies, we had earlier suggested that the precursor emission may have a steeper spectrum (spectral index α< −3.7, where α is defined as flux density S ∝ ν α , and ν is frequency) com-pared to that of the main pulse (Kaur et al. 2019). Even though our uGMRT data are not flux calibrated, assuming the nominal system parameters, i.e. the system temperature, T sys = 106 K, and a gain G = 0.32 K/Jy (Gupta et al. 2017), we estimate ∼0.1 mJy for the rms noise in uGMRT Band 3, accounting for the number of antennas (26) used in that observation. The precursor amplitude seen in uGMRT Band 3 is ∼ 5% of the main pulse (see Figure 2). This would translate to a precursor peak flux density of ∼3 mJy, which is consistent with < 12 mJy that we would expect based on the extrapolation from the MWA measurements (and assuming α < −3.7; Kaur et al. 2019).
The DM values have been measured using two different techniques, i.e. using the analytic template (narrow-band timing) and wideband timing (the latter is shown in shaded grey). The frequency range covered from each telescope is indicated in the top parentheses. The uncertainties in DMs are given in parentheses and correspond to the least significant digit. For MWA data, the DM precision is of the order of 10 −6 pc cm −3 , whereas it is ∼ 10 −5 pc cm −3 for the uGMRT bands (Bands 3 and 4) and is in the range ∼ 10 −3 to 10 −5 pc cm −3 for the Parkes UWL band. Figure 2) is generally small; in terms of measured pulse width (e.g., quantified using W 50 , at 50% of the peak), we measure a change of ∼ 2% across the two uGMRT bands, ∼5% across the MWA band, whereas it is ∼10% across the large frequency range of Parkes UWL (0.7 to 4 GHz). In light of this, we adopted two different methods for our timing analysis and DM measurements: the first method involved the use of a single template for each observing band, which is akin to "narrow-band" timing, as it does not account for any in-band pulse profile evolution ( § 4.2.1). The second method is essentially wide-band timing based on the pulse portraiture that was originally developed by Pennucci et al. (2014) and further discussed in Pennucci (2019).

Narrow-band timing using analytic templates
For this analysis we used different templates for the MWA, uGMRT Band 3, uGMRT Band 4, and Parkes-UWL, in an effort to incorporate the small degree of observed profile evolution. For the UWL data, we chose the Parkes 20 cm (1.4 GHz) template that we obtained from the PPTA project (Kerr et al. 2020). UWL data were integrated in frequency to 26 channels and then cross-correlated with the template to obtain the TOAs. For the uGMRT data (i.e., Band 3 and Band 4), we chose one of the brightest scans to generate noise free analytic templates by modeling a sum of von Mises functions using the paas utility of PSRCHIVE. The uGMRT data were then integrated in frequency to 32 channels and cross-correlated with respective analytic templates to obtain the TOAs. For the MWA data, we used an analytic template from one of our bright observations recorded in 2017 (also published in Kaur et al. 2019) The analytic template is then cross-correlated with the observed profiles to obtain the TOAs using the PSRCHIVE package.
The TOAs were then analysed using the pulsar timing package TEMPO2 (Hobbs et al. 2006) to determine the DM. We adopted the latest available timing solution of the pulsar (from the PPTA project; Kerr et al. 2020) and fit for only the DM in this analysis. 1 As described in Kaur et al. (2019), for most of our observations with the MWA, we typically achieve a timing precision of the order of ∼ 1 µs. Timing precision of the order of ∼2 µs was achieved with uGMRT Band 3; however, a better timing precision of ∼0.5 µs was obtained with Band 4. For our Parkes analysis, timing precision ranged from 0.9 µs (at the low-frequency end of the UWL band) to 3.5 µs (at the high-frequency end of the UWL band), as expected owing to the typically steep spectrum of radio pulsars (cf. Figure 1). This timing precision is reasonable considering the use of multi-frequency arrival times to determine the DM.

Wideband timing using the pulse portraiture
This analysis made use of the PulsePortraiture package, details of which are described in Pennucci (2019). We obtained a two-dimensional (2D) template from the PPTA project, which was generated by combining 10 very high-quality observations when the pulsar was exceedingly bright and RFI was minimal. For the DM analysis, UWL data were sub-divided into three near-octave sub-bands, using the pazi routine within PSRCHIVE package. DMs and TOAs were calculated for each of these three sub-bands.
For uGMRT Band 3 and Band 4 data, the 2D templates were generated from our own observations. For the MWA, we co-added three of our high-quality VCS observations to generate the 2D template. This was then used to measure the DM and TOAs for each of our MWA observations. The DM measurements from our analysis are summarised in Table 2.
The achievable DM precision strongly depends on the observing frequency band. The low frequencies of the MWA enable the most precise DM measurements, typically of the order of ∼(2-6)×10 −6 pc cm −3 . These are comparable to our previously published results (Kaur et al. 2019), and are still amongst the highest precision DM measurements to date. However, at uGMRT Band 3, DM precision achieved is ∼ 10 −5 pc cm −3 , an order of magnitude lower than that achieved using the MWA. Owing to much lower levels of RFI in Band 4, and due to a fortuitous scintillation brightening, we were able to achieve a comparable precision of the order of ∼ 2 × 10 −5 pc cm −3 . In the Parkes UWL band, the achieved DM precision varies from 10 −5 to 10 −3 pc cm −3 across the low to high segments of the band.

Orbital DM variations
The high DM precision achievable with the MWA offers the prospects of investigating any DM variations as a function of the orbital phase. Figure 3 shows DM measurements across the full 3.5-hour orbital period of the pulsar obtained from our MWA observations (cf. Table  2). The top panel shows the DM measurements over the full orbit from observations on MJD = 58671 and the bottom panel shows similar measurements made nine days later (MJD = 58680). No significant DM variation is seen as a function of orbital phase, except for a noticeable excess near orbital phase ∼ 0.4 -0.5, where a DM change of the order of (1.4 ± 0.6) × 10 −5 pc cm −3 is measured. Figure 4 shows similar measurements for the uGMRT data, which are from observations at two different epochs, separated by five days (MJD = 58794 and 58799). Even though the DM measurements with the uGMRT data are an order of magnitude less precise compared to the MWA measurements, and do not sample the exact same orbital phase ranges, there is a marginal trend for a slightly higher DM around ∼ 0.3 P b to 0.6 P b in the orbital phase. The uGMRT measurements also show larger DM variations on the order of ∼ 10 −4 pc cm −3 . However, unlike the case with MWA observations, these observations were not made on the same day, and therefore these orbital DM variations can be difficult to disentangle from any temporal DM variations given the time span of our observations. On the other hand, the MWA observations were taken on the same day, and therefore give us confidence that the variations are loosely correlated with the binary orbital phase, and are not due to other (temporal) DM variations.

Frequency-dependent DMs
Figure 5 presents DM variations over the large frequency range from 80 MHz to 4 GHz. The three different panels shown here are for different data sets (i.e. from different epochs). For example, the MWA, uGMRT Band 3 and UWL were available for observations made on MJD 58794-58796; the uGMRT Band 3 and UWL were used for observations on MJD 58799-58800; and the MWA, uGMRT Band 4 and UWL were used for observations on MJD 58816-58817.
Our data sets clearly indicate that DM changes measurably across the large frequency range. This frequency dependent behavior is seen consistently at all three epochs. We attempted to model this empirically using a power-law of the functional form δDM ∝ ν x , where δDM = DM meas − DM 0 ; DM meas and DM 0 denote the measured and reference DMs, respectively, and x is the scaling index (i.e. the change in DM relative to a reference value DM 0 , which denotes the asymptotic value near the low end of the frequency range of measurements). We measured scaling indices of 3.2 ± 0.6, 4.1±0.4, 2.9±0.4 from narrow-band timing and 2.2±0.1, 1.5±0.1, 2.4±0.5 from wide-band timing at three epochs, i.e. MJD = 58794, 58799 and 58816, respectively. We note that the fits are dominated by the low-frequency DM measurements which have uncertainties orders of magnitude smaller than the high-frequency DM measurements. An average value of 11.411318(6) pc cm −3 and 11.412391(5) pc cm −3 was obtained for DM 0 from two different timing analysis, compared to the catalog value of 11.41151(2) pc cm −3 (from Kaur et al. 2019).
Evidently, the DM measured at low frequencies ( 1 GHz) is significantly different to (in this case, higher than) those measured at higher frequencies ( 1 GHz). However, it may be possible to measure a predictive scaling, which can be applied to make meaningful DM corrections in PTA observations. Whether such a scaling is pulsar dependent, or epoch dependent, requires further investigation.

DISCUSSION
Our analysis suggests a clear, and quite compelling, evidence that the DM toward the millisecond pulsar J2241−5236, as measured in observations, depends on the observing frequency range spanned. It is seen consistently using two different timing analyses, and across three independent observing epochs spanning about a month. The observations were made contemporaneously using three different telescopes, covering a large Table 3. Summary of DM fit parameters.
frequency range from 80 MHz to 4 GHz. Specifically, in the case of narrow-band timing, we measure a δDM relative to a reference DM, DM 0 , the magnitude of which depends on the frequency range sampled; specifically, we measure δDM ∼ (1.2 ± 0.2) × 10 −4 pc cm −3 across 0.1-0.5 GHz, and δDM ∼ (1.5 ± 0.3) × 10 −2 pc cm −3 across 0.1-3 GHz, whereas the corresponding values from wide-band timing are ∼ (1.1 ± 0.1) × 10 −4 pc cm −3 and ∼ (1.4 ± 0.3) × 10 −2 pc cm −3 , respectively. This change in measured DM scales with the observing frequency (ν) as δDM ∼ ν 3.4±0.3 for narrow-band timing, and δDM ∼ ν 2.2±0.1 for wideband timing, where the quoted scaling index is the mean value of the indices estimated for three different observations. The scaling indices derived from wideband timing measurements are shallower than those using analytic templates. They are also generally consistent within ∼1-2 σ except for those on MJD 58799-58800, when there were no MWA observations, suggesting that lowfrequency measurements are crucial in constraining the scaling indices.
This presents the first reported evidence that the measured DM varies with the observing frequency, and that the change in DM scales with frequency but in a quantifiable way; this is inferred from contemporaneous observations over a large frequency range from 80 MHz to 4 GHz. Since PSR J2241−5236 is amongst the most promising targets for current and future PTAs, this finding has important implications. We consider various plausible sources of DM variations applicable to such high-precision timing, such as temporal variations in DM, or those caused by propagation effects such as scattering, as well those that may arise from intrinsic or environmental properties, and argue that none of them is relevant or can account for the observed changes in DM with frequency.
Temporal variations in DM caused by large space velocities of pulsars (∼80 km s −1 for PSR J2241−5236, Figure 5. DM measurements of PSR J2241−5236 at multiple frequency bands spanning the frequency range from 80 MHz to 4.0 GHz. Observations were made within 24-48 hour observing windows. The dashed line represents an empirical fit to a power law, and the solid gray line represents the reference DM value, DM0, from the fit. The inset plot represents the frequency range ∼100 to ∼500 MHz, to highlight the higher DM precision (of the order of ∼ 10 −6 to 10 −5 pc cm −3 ) obtained using the MWA and uGMRT. Left: DM values measured from narrow-band timing; Right: DM measurements from wide-band timing analysis. Top and middle panels show measurements from the MWA, uGMRT Band 3 and the Parkes UWL receiver, The bottom panel shows measurements from the MWA, uGMRT Band 4, and Parkes UWL. The measured scaling indices and DM0 values are presented in Table 3. based on recent proper motion measurement and assuming a distance of 0.96 kpc; Reardon et al. 2021) is a dominant source of noise in PTA measurements. While there are no published records of longer-term DM variations for this pulsar, our MWA data suggest δDM ∼ 10 −4 pc cm −3 , on timescales of ∼1-2 yr (Kaur et al. 2019). The fact that our observations were made contemporaneously, through careful coordination of the three telescopes within ∼24-48 hour time windows, naturally mitigates this effect.
Since PSR J2241−5236 is in a 3.5-hr binary orbit with a black-widow type companion, another possibility to consider is variation of DM with orbital phase. Even though the pulsar does not show any evidence of eclipse (Keith et al. 2011), in light of the orbital modulations seen in high-energy gamma-ray observations and its interpretation in terms of intra-binary shock emission , it was imperative to consider (and investigate) this effect. However, as discussed in § 4.2.1, the measured δDM from our specially designed MWA observations is rather small, i.e. δDM of ∼ 1.4×10 −5 pc cm −3 , which is ∼1-2 orders of magnitude smaller than measured DM changes across the frequency range. This therefore readily precludes any orbital variations. Even so, this small change in DM will result in a timing delay of ∼ 3 µs at the MWA's observing band (120-220 MHz), which will translate to ∼ 110 ns in the Parkes UWL band (700-4000 MHz).
Another important effect to consider is multi-path propagation that can give rise to a scatter broadening of the observed pulse shape, which varies strongly with frequency. For PSR J2241−5236, our observations do not reveal any visible scatter broadening even at the lowest frequencies of our observations (at 80 MHz). Closer examination of the data presented in Fig. 2 suggests a scintillation bandwidth ∼ 5 MHz at 400 MHz (i.e., uGMRT Band 3), which translates to a pulse broadening time (τ d ) of ∼ 30 ns at this frequency, and ∼ 20 µs at the MWA's 80 MHz (assuming a scaling of τ d ∝ ν −3.9 as per Bhat et al. 2004). Notably, the expected pulse broadening time is 1 µs at frequencies 200 MHz, which is significantly smaller than ∼5 µs dispersive time delay expected for a measured δDM ∼ 1.2 × 10 −4 pc cm −3 across 0.1 -0.5 GHz, and a much larger ∼2.5 ms for δDM ∼ 1.4 × 10 −2 pc cm −3 across 0.1-3 GHz. Even at the lowest frequency (80 MHz), the estimate of τ d is still smaller than the dispersive delay ∼ 75 µs across the MWA band (80-220 MHz) due to the measured change in DM. Therefore, frequency-dependent DMs seen in our data are unlikely due to variable scatter broadening across the frequency range.
The qualitative agreement of results from two distinct timing analyses, as summarised in § 4.2.1 and 4.2.2, clearly suggests that our analysis is not significantly biased by any un-modeled pulse profile evolution. As noted in § 4.2, the observed profile evolution across the large frequency range is comparatively small for this pulsar, with ∼10% change in the measured pulse width (quantified in terms of W 50 ) across 80 MHz -4 GHz spanned by observations. An attempt to fit for FD parameters yielded a marginally significant (∼1.5 σ) value for FD1, with ∼10% improvement in the reduced chi-squared. Furthermore, the measured changes in DM values (i.e. δDM relative to DM 0 ) are significantly higher than the uncertainties in our estimated DM values, as evident in Fig. 5. The wide-band timing analysis method using the PulsePortraiture is designed to model and account for any pulse profile evolution. Given this, the remarkable consistency seen in terms of an implied frequency dependence between the two timing methods, across all three data sets, precludes (or at least significantly diminishes) the possibility that profile evolution is significantly influencing our analysis and results.
While the profile evolution is unlikely to be the cause of an observed frequency dependence in DM measurements, any residual unmodelled profile evolution may give rise to an 'offset' contribution to the measured DM that may still be frequency-dependent. Such an effect will most likely depend on the pulsar, and to a certain extent, also on the instrumentation employed, and even the method used for the DM estimation. In any case, we expect this to be a rather subtle effect for this pulsar.
Another possibility is that there may also be very subtle aberration and retardation (AR) effects which can impact DM measurements. However, AR effects, which depend on the emission height and the radiusto-frequency mapping, are not expected to change over time, as these quantities are functions of the viewing geometry and the magnetic field configuration. Therefore, AR effects cannot be responsible for the difference between our measurements of the DM (for instance), and historical DM measurements at the same frequency. Thus, although AR effects could explain the observed δDM scaling at the reported epoch (at least qualitatively), we reject this explanation on the grounds that the observed effect (e.g., the measured scaling index) appears to be time-dependent.
Finally, there are recent reports of this pulsar exhibiting astonishingly low levels of "jitter noise," which has also been shown to contribute to the variability to DM estimates (∼ 4 ns in an hour of observation; Parthasarathy et al. 2021); this result further strength-ens the confidence in our DM estimates and an argument in support of a frequency dependence in DM. Cordes et al. (2016) present a detailed account of frequency-dependent DMs including their expected magnitudes and scaling in frequency as a function of the observing frequency. According to their Equation (12), we may expect a predicted difference in DM between 150 MHz and 3.2 GHz of the order of ∼ 3 × 10 −6 pc cm −3 , assuming a scattering screen placed half-way between the pulsar and us, and a pulsar distance of ≈1 kpc (Reardon et al. 2021), and scattering strength consistent with our measurement of the diffractive scintillation bandwidth (Fig 2; panel 2). We note that this is however much smaller than our measured difference in DM of ∼ 0.01 pc cm −3 . Furthermore, the difference in DM increases with an increase in frequency, which is in contrast with the theoretical predictions. The analysis and results presented in this paper motivate further detailed studies of frequency-dependence in DM toward other PTA pulsars.

SUMMARY AND CONCLUSIONS
In this paper we have presented the results from contemporaneous observations of PSR J2241−5236 made with the MWA, uGMRT and Parkes at three epochs, spanning about a month. Our observing and analysis strategies have been devised so that effects such as variations of DM with time, or orbital phase, can be either measured or minimised before frequency-dependence in DM can be investigated. This is the first time that this pulsar has been observed over such a wide-frequency range. Our measurements demonstrate that the MWA enables estimation of DMs with precision of the order of 10 −6 pc cm −3 , which is much better than that typically obtained from timing-array observations.
We have reported the first convincing evidence that the DM as measured in wideband observations, depends on the frequency range spanned. This is demonstrated using broadband data obtained for PSR J2241−5236, and is seen consistently at three different observing epochs, and using two different timing techniques. Effects of this kind, if confirmed for other PTA pulsars, will have important implications for timing-array experiments.
Our results cannot be attributed to temporal or orbital variations, or other effects such as profile evolution with frequency, and suggest a bona-fide frequencydependence in measured DMs. The fact that DMs measured at low frequencies differ from those measured at high frequencies suggests that low-frequency observations cannot be straightforwardly used to apply DM corrections in PTA observations. The larger fluctua-tions in DMs at higher frequencies might be attributed to averaging effects caused by the larger sizes of scatter broadened images at lower frequencies Shannon & Cordes 2017).
Our analysis shows that it is possible to empirically estimate the scaling index for frequency-dependent DM, in the same way that it can be done for other frequencydependent effects (e.g., scintillation or pulse broadening). Whether the scaling is constant for a given pulsar, or is epoch dependent, needs further investigation. Regardless, for a given pulsar, if such a scaling can be verified, e.g., via suitably designed monitoring campaigns, low-frequency pulsar observations may still have an important role to play in applying effective DM corrections for pulsar timing arrays.

ACKNOWLEDGMENTS
We thank the referee for their constructive comments that helped improve the quality of the results presented in this paper. This scientific work makes use of the Murchison Radio-astronomy Observatory, operated by CSIRO. We acknowledge the Wajarri Yamatji people as the traditional owners of the Observatory site. This work was supported by resources provided by the Pawsey Supercomputing Centre with funding from the Australian Government and the Government of Western Australia. The uGMRT is operated by the National Centre for Radio Astrophysics of the Tata Institute of Fundamental Research, India. The Parkes radio telescope (Murriyang) is part of the Australia Telescope National Facility which is funded by the Australian Government for operation as a National