\section{An appendix section}

\section*{Acknowledgements}

This work is carried out by InPTA, which is part of the International Pulsar Timing Array consortium. We thank the staff of the GMRT who made our observations possible. GMRT is run by the National Centre for Radio Astrophysics of the Tata Institute of Fundamental Research.  BCJ, PR, AS, SD, LD, and YG acknowledge the support of the Department of Atomic Energy, Government  of India, under project identification \# RTI4002. BCJ and YG acknowledge support from the Department of Atomic Energy, Government of India, under project \# 12-R\&D-TFR-5.02-0700.
AS is supported in part by the National Natural Science Foundation of China grant No. 11988101.

\vspace{-1em}

\section{Conclusion}
\label{conc}

The observing strategy and analysis methods used in the InPTA are useful in informing the design of a PTA experiment with the SKA telescope. A precision better than 2 $\times$ 10$^{-5}$ pc\,cm$^{-3}$ for pulsars such as PSR J1909$-$3744, already demonstrated in our analysis with the InPTA in this paper, is possible with much larger sensitivity and frequency coverage provided by the SKA. We recommend using the SKA-low and SKA-mid as one to four sub-arrays respectively to cover frequencies from 50 to 1800 MHz concurrently based on the InPTA experience. We show that such observing can be used for a 62 pulsar sample three times a month requiring 936 hours similar to observing time required by currently on-going PTA experiments. The development of general relativistic constructs by the InPTA astronomers promise better tuned methods for search for isolated SMBHB GW sources in the PTA data-sets in the SKA  era as well as constrain the astrophysics of these sources with multi-messenger observations. Last, but not the least the InPTA experiment is developing a vibrant Indo-Japanese user community for GW physics with the SKA telescope when this telescope comes online later in the decade.

As the SKA-construction plan gets underway, it is possible to devise an optimised observing plan complimentary to the four major PTA experiments (the EPTA, the InPTA, the NANOgrav and the PPTA) for a higher cadence coverage of an extended sample of PTA pulsars. This involves an optimisation in sensitivity and frequency coverage using different combination of available antennas as they become available in a phased manner from 2024 on-wards. A design of such an optimal strategy is in progress and will be reported in future.

A potential discovery of the nano-Hertz GW  is likely in near future. This can further inform the strategy best suited for post-discovery GW science in the SKA-era. The discussion presented here will hopefully motivate a deeper discussion amongst the major PTA experiments and the IPTA for future directions of exciting GW astrophysics with the SKA, which is likely to revolutionise the radio astronomy with its Square kilo-meter collecting area in the next two decades.

\section{Indian Pulsar Timing Array experiment}
\label{inpta}

\par
Indian Pulsar Timing Array experiment (InPTA) was started in 2015  using the Ooty Radio Telescope \citep[ORT:][]{Swarup1971} and the legacy Giant Metre-wave Radio Telescope \citep[GMRT:][]{Swarup1991} as a pilot experiment. The InPTA now utilises the unique  capabilities of the upgraded Giant Metre-wave Radio Telescope \citep[uGMRT: ][]{gak+17}, which is also a pathfinder telescope for the SKA. The salient features of this experiment are described in this section.

\begin{figure}[H]
\centering
\includegraphics[width=0.5\textwidth]{ToAresiduals_Legacy_plot.pdf}
\caption{Timing residuals for two PTA pulsars with the ORT during and after the legacy phase of the proposal. The wander in the residuals before mid-2016 was traced to instability in the observatory time and frequency standards and associated electronics and was subsequently rectified.}
\label{fig:ORTres}
\end{figure}

\par
The experiment started as a pilot experiment in 2014, first using the ORT for observing two pulsars, PSRs J1713+0747 and J1909-3744, to test the precision timing capabilities of the Indian instruments. A newly commissioned pulsar receiver called PONDER \citep{njmk15} was employed to coherently dedisperse the data over 16 MHz band-pass at 326.5 MHz and generate time-stamped integrated pulse profiles (IPPs) using topocentric pulse frequency calculated from polynomial coefficients with the best known available timing solutions. The folded time-series were analysed using \texttt{TEMPO2}  \citep{Hobbs2006_tempo2,Edwards2006_tempo2} pulsar timing software package to obtain timing residuals. Early observations showed residuals varying over a few milliseconds and with clock jitter (Figure \ref{fig:ORTres}). These were traced to instability in the observatory time and frequency standard and associated electronics. Based on these initial observations, the observatory instrumentation was improved. Further observations showed stable timing residuals demonstrating  the  feasibility of a PTA experiment and also helped to obtain refined timing solutions for the next phase.

\par
Encouraged by these results, a pilot program with a sample of 9 pulsars was initiated using both the legacy GMRT and the ORT in 2015. The multi-frequency data, acquired in this pilot phase, was used to test the feasibility of a PTA with both  instruments for the first time. The data were used (a) to plan improvements in the observatory instrumentation, (b) to assess the timing precision achievable with the legacy GMRT and the ORT, (c) to determine the fixed pipeline delay between data acquisition systems \citep{smj+21}, (d) to build timing solutions for a larger sample of PTA pulsars   and (d) to train new members of the collaboration in observing and data reduction.  Different data reduction techniques using a variety of pulsar off-line processing packages, such as \texttt{PRESTO}\footnote{\url{https://www.cv.nrao.edu/\~sransom/presto/}} \citep{Ransom2011_PRESTO},  \texttt{SIGPROC}\footnote{\url{http://sigproc.sourceforge.net/}} \citep{Lorimer2011_sigproc}, \texttt{PSRCHIVE}\footnote{\url{http://psrchive.sourceforge.net/}} \citep{Hotan2004_psrchive}, and \texttt{DSPSR}\footnote{\url{http://dspsr.sourceforge.net/}} \citep{vanStraten2011_dspsr}  were tried on these data to decide the correct strategy for a bigger PTA experiment in the future.

\begin{table*}
\centering
\begin{tabular}{c|c|c|c|c|c|c}
\hline
 & Name & Period & DM & Median S/N & Median S/N & Obs time \\
& &  (s)   & (pc\,cm$^{-3}$)  & 1460    & 500 & (min) \\ \hline
&J0645+5158   & 0.008853& 18.25& 7&52 & 30 \\
&J1024-0719   & 0.005162&  6.49& 13& 17& 30 \\
&J1455$-$3330 & 0.007987& 13.57& 9&19 & 30 \\
&J1614-2230   & 0.003151& 34.49&12 &18 & 30 \\
Exploratory InPTA &J1640+2224   & 0.003163& 18.43& 20& 34& 30 \\
& J1730$-$2304&0.008123& 9.62 & 61     & 56  & 30 \\
& J1738+0333  & 0.005850&33.77&   7   & 10  & 30 \\
& J2317+1439  & 0.003445&21.91& 8 & 21 & 30\\
& J2302+4442  &0.005192& 13.79& 9  & 24     & 30 \\
\hline
& J1643$-$1224& 0.004622& 62.41& 129  & 423 & 50\\
& J1713+0747  & 0.004570& 15.91& 254 & 174 & 50\\
& J1857+0943  & 0.005362& 13.31& 74  & 134 & 30\\
Classic InPTA& J1909$-$3744& 0.002947& 10.39& 77 & 289 & 50\\
& J1939+2134  & 0.001558& 71.02& 159 & 323 & 20 \\
& J2145$-$0750& 0.016052&  9.00&103 & 946 & 50 \\
& J2124$-$3358& 0.004931& 4.60& 39  & 208 & 50\\
\hline
& J0437$-$4715&0.005757 & 2.64& 958    & 2343   & 15 \\
& J0613$-$0200& 0.003061&38.77&  61      & 119      & 50 \\
& J0751+1807  & 0.003479&30.25&  32      & 59      & 30 \\
Expanded InPTA & J1012$+$5307& 0.005255& 9.02&  259      & 359     & 30 \\
& J1022$+$1001& 0.016452&10.25&  38     & 417     & 20 \\
& J1600$-$3053& 0.003597&52.32&  71      & 67     & 50 \\
&J1744$-$1134 &0.004075 & 3.14&    72   & 343     & 50 \\
\hline \hline
\end{tabular}
\caption{The pulsar sample of the InPTA. The sample labelled "Classic InPTA" was observed from 2018$-$2022, whereas the sample labelled "Exploratory InPTA" was observed only in 2018$-$2019. Likewise, the sample labelled "Expanded InPTA" was mostly observed only after 2021 with high sensitivity. The "Exploratory InPTA" sample also included all pulsars listed in the "Classic InPTA" and "Expanded InPTA" sample. The "Classic InPTA" and the "Expanded InPTA" pulsars are being observed using Band 3 and 5 with 200 MHz bandwidth, whereas "Exploratory InPTA" sample was observed with Band 3, 4 and 5 with 100 MHz bandwidth.}
\label{tab:sample}
\end{table*}

\par
The GMRT was upgraded between 2010 and 2017 \citep{gak+17} providing a seamless frequency coverage from 50 to 1500 MHz with five new wide-band feeds (Band 1 : 50$-$80 MHz; Band 2 : 120$-$250 MHz; Band 3: 250$-$500 MHz; Band 4 : 550$-$850 MHz and Band 5 : 1050$-$1450 MHz). The wide-band data are acquired with new wide-band receivers and back-end, capable of recording data with 400 MHz band-pass, as beam-former time-series from 4 simultaneous beams. After an initial exploratory experiment to assess the timing precision of the uGMRT with 19 pulsars, a strategy to observe each pulsar in more than one frequency band simultaneously with the uGMRT was adopted for the ongoing InPTA experiment since 2018.

\par
The uGMRT is a multi-element interferometer much like the SKA. The large sensitivity of the uGMRT is achieved by compensating the geometric and instrumental phase delays (including ionospheric delays) across the interferometric array, consisting of 45-m diameter antennas, to form a phased array. While a phased array of 20 antennas provides a collecting area of about 180-m diameter single dish, an eight antenna phased array synthesises a 125-m single dish. In addition, the relatively smaller beam of the phased array helps in both reducing the sky background as well as in filtering uncorrelated noise as compared to a single dish telescope.  The second advantage of an interferometer like the uGMRT is a possibility of forming multiple phased arrays using different groups of antennas, called sub-arrays. If different sub-arrays are set up to observe at different frequencies, this allows simultaneous multi-frequency observations covering a large bandwidth. On the other hand, multiple sources can be covered at the same frequencies by multiple sub-arrays. Thus, the uGMRT not only provided a sensitivity similar to the large telescopes used in other PTA experiments, but an opportunity to try out new observations strategies with an interferometer. The ongoing InPTA experiment has explored this flexibility of an interferometer to optimise the precision of the experiment.

\par
In the beginning, the InPTA used three uGMRT bands (400$-$500, 650$-$750 and 1360$-$1460 MHz in Band 3, 4 and 5 of the uGMRT respectively) by grouping five to 15 antennas in three sub-arrays using three out of four beams of the uGMRT. As the uGMRT back-end \citep[GMRT Wideband Backend - GWB : ][]{rkg+17} provides a coherent dedispersion capability for a maximum of 200 MHz bandwidth \citep{dg16}, two of the Bands (Band 3 and 5) were coherently dedispersed to remove dispersive effects of the IISM with 100 MHz band-pass respectively, while Band 4 was observed with a 1024 channel digital filter-bank. While this strategy sampled the 300$-$1500 MHz with three non-contiguous bands, the smaller bandwidth as well as lesser number of antennas in each band limited the achievable sensitivity of InPTA observations.

\par
We addressed this shortcoming since 2019 by using only two sub-arrays with 10 and 15 antennas in Band 3 and 5. As most of PTA pulsars are relatively low DM pulsars (2.6 to 30.0 pc\,cm$^{-3}$), the dispersion smear for these pulsars can be kept smaller than the sampling time employed (40.96 $\mu$s) by using a filter-bank with 1024 sub-bands or more at Band 5. This permitted using the entire 200 MHz bandwidth for coherent dedispersion at Band 3 to eliminate  the significantly larger dispersive smear at frequencies near 300 MHz. Thus, the data were observed in both bands with 200 MHz bandwidth with coherent dedispersion used only in Band 3. The total integration time for all the pulsars was also increased to one hour compared to 30 minute observations in 2018. This strategy provided 4 times larger sensitivity than earlier observations significantly reducing ToA uncertainties as well as improving DM precision as is shown later.

\par
In this new strategy, the monitoring observations were started with 6  bright pulsars and the sample is being gradually increased to 14  pulsars(Table \ref{tab:sample}).  New pulsars were selected from the pulsars which formed a part of the IPTA Data Release I \citep{IPTA_DR1}. Apart from choosing only the pulsars visible with the GMRT (Dec $>$ -50$^\circ$), the  signal-to-noise ratio (S/N) achieved in the exploratory phase of the experiment (labelled "Exploratory InPTA" in Table \ref{tab:sample}) with the selected configuration was used to make a shortlist. As a good sampling of Hellings-Down overlap function \citep{HellingsDowns1983} is required for detection of SGWB, new pulsars were added to optimise a good coverage of this overlap function  with pulsar pairs (Figure \ref{fig:distinpta}).

\begin{figure}
\centering
\includegraphics[scale=0.45]{HD_curve.pdf}
\caption{Histograms for the distribution of pairwise angular separations between the InPTA pulsars. The sample was gradually increased to achieve a more uniform sampling of spatial separations (Hellings-Downs overlap function \citep{HellingsDowns1983}) over the course of the experiment. The red and green histograms represent the sample over two 6 month observing semesters for the uGMRT.   }
\label{fig:distinpta}
\end{figure}

\par
The current data set for InPTA, including the pilot, legacy and exploratory uGMRT phases spans about 8 years. However, the span and cadence of individual pulsars in our sample have varied over the years due to the considerations mentioned above. While most of the recent dataset has been taken using  the upgraded GWB \citep{rkg+17}, the frequency bands and receivers employed in observations prior to 2018 have varied. All the legacy GMRT and the uGMRT data were acquired with two hands of polarisation, but the ORT data were single polarisation data due to the nature of the instrument. Lastly, the pipeline delays as well as other data acquisition features are well characterised for the uGMRT receivers. As a consequence, the collaboration has decided to limit its upcoming first public data release to observations with the uGMRT alone to maintain the homogeneity in data properties.

\begin{figure}
\centering
\includegraphics[scale=0.5]{resfreq2psr.pdf}
\caption{Timing residuals obtained with data observed simultaneously over 300 to 1460 MHz for PSRs J1713+0747 and J1909$-$3744 are shown in this figure after a DM fit over this range without using any frequency dependent FD parameter or DMJUMP. Epoch by epoch fits provide very precise DMs in the InPTA experiment.  }
\label{fig:freqres}
\end{figure}

\par
The pulsar data generated by uGMRT beamformer observations are recorded in a channelised time series binary format by the GWB, together with a timestamp representing the start of the observation. A pipeline named \texttt{pinta}\footnote{\url{https://github.com/abhisrkckl/pinta}} \citep{smj+21}  was developed to process such datasets obtained using the uGMRT for the InPTA experiment. This pipeline performs RFI mitigation using two different softwares, \texttt{RFIClean}\footnote{\url{https://github.com/ymaan4/RFIClean}} \citep{Maan2020_rfiClean} and \texttt{gptool}\footnote{\url{https://github.com/chowdhuryaditya/gptool}}, and produces partially folded PSRFITS archives which can be used for further analysis using standard pulsar software such as \texttt{PSRCHIVE}. After the data reduction offline with this pipeline, the multi-band data were further analysed with post-processing pipelines, \texttt{DMCalc} \footnote{\url{https://github.com/kkma89/dmcalc.git}} and \texttt{PulsePotratiture} \citep{Pennucci+2014} to obtain precision DM measurements. These pipelines were developed and/or tuned for the uGMRT datasets and yield ToAs apart from DM measurements. These were analysed with pulsar timing software, \texttt{TEMPO2}, to obtain the final timing solutions as well as timing residuals, which form the input for GW analysis.

\begin{figure}
\centering
\includegraphics[scale=0.5]{dm2psr.pdf}
\caption{DM measurements using the data from 300 to 1460 MHz obtained for PSRs J1713+0747 and J1909$-$3744 in the InPTA experiment.  }
\label{fig:dmtimeseries}
\end{figure}

\par
The experiment has resulted in the highest precision measurements of DM reported so far by any PTA, primarily due to a combination of simultaneous Band 3 and Band 5 wide-band data covering 300$-$1460 MHz,  which is unique to our experiment. Such a combination is affected by misalignment due to profile evolution and scatter broadening over this frequency range.  As explained later (Section \ref{method}), analysis of ToAs generated using frequency resolved templates showed normally distributed residuals  without using FD parameters or DMJUMP to align Band 3 and Band 5 data   (Figure \ref{fig:freqres}). Consequently, we were able to achieve high precision with the median uncertainty as small as  2 $\times$ 10$^{-5}$ pc\,cm$^{-3}$ for PSR J1909$-$3744. The variation of DM for PSRs J1713+0747 and J1909$-$3744 are shown in Figure \ref{fig:dmtimeseries} as an illustration. The SKA has a similar frequency coverage and is likely to provide these measurements with similar or better precision using a strategy similar to the InPTA. With these measurements, we have also achieved $\mu$s post-fit residuals, comparable to those obtained with the other PTA experiments and the preliminary results from the upcoming public data release of the InPTA are shown in Figure \ref{fig:timres}.

\begin{figure}
\centering
\includegraphics[scale=0.5]{resbat2psr.pdf}
\caption{Timing residuals for PSRs J1713+0747 and J1909$-$3744 with the data obtained in the InPTA experiment are shown as a function of observation epochs. The RMS post-fit residuals achieved are indicated in the legend  }
\label{fig:timres}
\end{figure}

\par
While we did not attempt observations of multiple sources using multiple beams-sub-arrays, the experience from our strategy can be extrapolated to evaluate such an observing strategy. The InPTA observations provide useful experience to understand the phasing constraints from an ionospheric turbulence point of view for different observations. This is important to decide the configurations of interferometer antennas for the sub-arrays to be used, particularly for the most distant antennas in the sub-array.

\section{Lessons from InPTA for GW science with SKA}
\label{inptanska}

\par
The multi-element interferometer nature of SKA makes it an instrument similar to the uGMRT, albeit with a much larger sensitivity. Hence, our experience with the uGMRT, which is also a recognised path-finder for the SKA, can inform the observations and analysis strategy for a future PTA with the SKA. In this section, some possibilities based on the InPTA experience are discussed.

\par
The SKA consists of two telescopes, both of which will be constructed in two phases. Here, we provide a top level description of the capabilities of the phase-1 of the telescope, which will represent about 10 percent of the final SKA collecting area. The SKA-mid will be located in Karoo, South Africa. It will consist of 197 dishes spread over three spiral arms on an area of a diameter of about 150-km. These include 133, 15-m dishes and 64, 13.5-m MeerKat dishes. It is planned to have three frequency bands operational in the first phase. While Band 2 covers 0.95 to 1.76 GHz and is likely to become available first, Band 5 (4.6 $-$ 15.3 GHz) and Band 1 (0.35 $-$ 1.05 GHz) will be commissioned subsequently. The planned back-end consists of a beam-former providing 16 simultaneous beams with coherent dedispersion carried over 300 to 800 MHz band-pass. About 120 antennas, consisting of 46 MeerKat antennas, are located in a cluster of dimensions 2 km by 2 km, whereas the rest of the  77 antennas are spread beyond these baselines. Further, about 20 antennas each from the spiral arms (15 SKA dishes and 5 MeerKat dishes) can  be used as two 3 km by 3 km clusters. While further arm antennas can also be used as phased arrays, the baselines between these antennas are too long for a stable phase solution at a frequency below 2 GHz. These baselines for the outer spiral arms are probably useful as phased arrays at Band 5, where ionospheric propagation effects are much smaller. However, the sensitivity of a Band 5 tied array is smaller than Band 2 for the same number of antennas as the SKA-mid antennas are less sensitive at Band 5. Moreover, pulsars are steep spectrum  sources and fainter at Band 5. Hence, only 120 to 140 antennas are considered here  for high sensitivity PTA observations.

\begin{table}[ht]
\begin{tabular}{|c|c|c|}
\hline
{\bf GMRT } & {\bf Single antenna} & {\bf Subarray }   \\
\hline
Band 2  & & 8 antennas \\
 Aeff/Tsys ($m^2/K$) &	 3.8	& 30 \\
SEFD (Jy) &	727	& 91	    \\
\hline
Band 3 	& &	8 antennas	\\
Aeff / Tsys ($m^2/K$)	& 10.5  & 	84 \\
SEFD (Jy)&	263	&33 	\\
\hline
Band 4 &  & 5 antennas\\
Aeff / Tsys ($m^2/K$)&	9.7	& 48 \\
SEFD (Jy) &	286 & 57		 \\
\hline
Band 5 &  &	15 antennas\\
Aeff / Tsys ($m^2/K$)&	9.2 & 	138   \\
SEFD (Jy) &	300	& 20   \\
\hline
\end{tabular}
\caption{ The system equivalent flux density estimated for the uGMRT sub-array configurations employed in the InPTA experiment}
\label{tab:sefd_gmrt}
\end{table}

\par
The SKA-low telescope consists of  131,072 dual polarised log periodic antenna-based aperture arrays, located  in Western Australia in Murchinson Radio-astronomy Observatory. These dipoles are arranged in 512 stations in the central core and three spiral arms located approximately 65 km apart with each station consisting of 512 dipoles. The core consists of 224 stations with the largest baseline of 1 km with 13 more clusters, each with 6 stations, located within 2 km of the core. The central core can be used as two phased arrays, each proving a physical aperture of 127020 m$^2$ or as a single phased array of two times that area.

\par
SKA-mid uses a real-time digital back-end for pulsar timing, called Pulsar Timing Engine \citep[PST : ][]{k2018}. The PST helps to record high time resolution full-Stokes pulsar data. It receives 16 tied-array beams, with each beam pointing on a pulsar to be timed concurrently. Observation durations can be from 180 to 1800 seconds. Each beam can be coherently dedispersed over 300 to 800 MHz in real-time to provide high time resolution data, which can be further processed to provide frequency-resolved sub-integrated time-stamped IPP in the PSRFITS format \citep{Hotan2004_psrchive}. The multiple beams in both SKA-mid and SKA-low allow multiple targets/bands  to be observed concurrently making both these instruments very flexible for different observing strategies.

\par
The uGMRT, which is a pathfinder telescope for the SKA, is similar to the SKA due to its multi-band, multi-beam and multi-sub-array nature. The 30 antennas of the GMRT can be grouped into 4 sub-arrays to operate over 5 frequency bands using 4 beams provided by the GMRT back-end (See Section \ref{inpta}). Using our experience for such a strategy, a PTA program with the SKA can utilise at least 4 sub-arrays in the SKA-mid and  two sub-arrays in the SKA-low for simultaneous multi-frequency observations as well as concurrent observations of multiple targets. Apart from optimal utilisation of telescope time for the pulsar Key Science Project (KSP), such a strategy will also provide a wide-frequency coverage for separating ISM noise from the precision timing in a way similar to that in the InPTA. Some possibilities are explored below for a suitable strategy for a PTA program with the SKA.

\begin{table*}[ht]
\centering
\begin{tabular}{|c|c|c|c|}
\hline
{\bf } & {\bf 	Single} & {\bf Single}  & {\bf All }   \\
{\bf SKA Low} & {\bf 	station} & {\bf subarray }  & {\bf subarrays  }   \\
{\bf } & {\bf 	 } & {\bf   }  & {\bf  within 1 km }   \\
\hline
&   & & \\
Low Band  & - &	100	& 200\\
(50-350 MHz) & - &	stations	& stations\\
Aeff / Tsys ($m^2/K$)&		1.10  &	110.00 &	220.00	  \\
SEFD (Jy) &	2,510.27 &	25.10 &	12.55	   \\
\hline
\hline
{\bf } & {\bf 	Single} & {\bf Single}  & {\bf All }   \\
{\bf SKA MID} & {\bf 	antenna} & {\bf subarray }  & {\bf subarrays  }   \\
{\bf } & {\bf 	 } & {\bf   }  & {\bf  within 2 km }   \\
\hline
Mid Band 1a  &  & & \\
(0.35-0.65 GHz) & - & 30 antennas & 120\\
 Aeff/Tsys ($m^2/K$) &		2.10 &	63.00 &	252.00 \\
SEFD (Jy) &		1,314.90 &	43.83 &	10.96  \\
\hline
Mid Band 1b &  & & \\
(0.65-1.05 GHz) & - & 30 antennas & 120\\
Aeff / Tsys ($m^2/K$)	& 4.20	& 126.00 &	504.00\\
SEFD (Jy)&	657.45 &	21.92 &	5.48 	\\
\hline
Mid Band 2 &   & & \\
(0.95-1.76 GHz)& - &  30 antennas & 120\\
Aeff / Tsys ($m^2/K$)&		11.00 &	330.0 &	1,320.00  \\
SEFD (Jy) &	251.03 &	8.37 &	2.09   \\
\hline
\end{tabular}
\caption{The system equivalent flux density estimates for the proposed sub-array configurations for a PTA experiment using the SKA Low and SKA Mid telescopes}
\label{tab:sefd_ska}
\end{table*}

\par
A comparison of the system equivalent flux density (SEFD) for the uGMRT arrays employed in the InPTA experiment and those for the SKA-mid and SKA-low sub-arrays,
shown in Tables \ref{tab:sefd_gmrt} and \ref{tab:sefd_ska}, suggests that a 30 antenna SKA sub-array is equivalent to 6, 12 and 35  antenna uGMRT sub-array at band 3, 4 and 5 of the uGMRT respectively. While such  SKA sub-arrays provide marginally smaller sensitivity than the uGMRT Band 3 sub-array used in the InPTA experiment, a higher sensitivity is indicated  for sub-arrays at other bands. With about 120 antennas in the central part of the SKA-mid, this suggests the  use of a maximum 4 sub-arrays using 4 beams with sensitivity similar to or greater than the uGMRT. For weaker pulsars, one could instead use only two or even a single SKA-mid sub-array to achieve almost four times  the uGMRT sensitivity. We recommend a possible configuration of 4 sub-arrays with the SKA-mid telescope as shown in Figure \ref{fig:skamid_sub}. Given that the available bandwidth is 300/800 MHz compared to 200 MHz used in the InPTA program, one would achieve much larger sensitivity than the uGMRT sub-arrays in practice. While 16 sub-arrays are in principle possible in the SKA, this comes at the cost of the sensitivity of each sub-array and is not recommended based on the InPTA experience, where  only 2 sub-arrays are used instead of possible 4 in the current ongoing observations. A similar comparison with the SKA-low configuration suggests that two sub-arrays using each half of SKA-low core stations provide significantly larger sensitivity than Band 2 of the uGMRT (See Tables  \ref{tab:sefd_ska} and \ref{tab:sefd_gmrt}).  We recommend a possible configuration of 2 sub-arrays with the SKA-low telescope as shown in Figure \ref{fig:skalow_sub}. Since SKA-mid and SKA-low are two different telescopes with independent correlator and beam-former, coordinated observations between the telescopes can provide the low-frequency support for higher frequency observations with a SKA-PTA covering an unprecedented frequency range of 50 $-$ 1800 MHz to provide the most sensitive PTA dataset ever.

 \begin{figure*}
 \centering
 \includegraphics[width=0.6\textwidth]{ska_mid_composite_v2_12mar2022.pdf}
 \caption{ The proposed sub-array configuration for a pulsar timing experiment with the SKA-mid Phase-1 Telescope }
 \label{fig:skamid_sub}
 \end{figure*}

\par
The two to four sub-array configurations, discussed above, can be employed in a variety of ways.  One option would be to use the same observing strategy as InPTA, where the different sub-arrays are used for simultaneous multi-frequency observations of the same target pulsar (hereby referred to as OPTION I). This provides a wide frequency coverage, useful for characterising the ISM noise, the jitter noise and bias due to profile evolution across the frequency. The other option would be to use each sub-array to observe a different target pulsar concurrently (OPTION II). This option optimises the use of telescope time to cover a much larger sample of PTA pulsars. A third option is to use a hybrid strategy combining both multi-frequency and multi-target strategies (OPTION III). Finally, fainter pulsars can be covered with the full 120/200 antenna array, one source at a time, with three
observations employing a different frequency band (OPTION IV). It may be noted that the brightness of each pulsar (or equivalently the achievable S/N) determines the precision of ToAs. Consequently, the number of sub-arrays possible depends on the distribution of flux density of the PTA pulsars, which also is a critical factor in deciding the most useful observing strategy. Thus, given the pulsar sample, all four of the observing strategies may be required.

 \begin{figure*}
 \centering
 \includegraphics[width=0.6\textwidth]{ska_low_composite_12mar2022.pdf}
 \caption{ The proposed sub-array configuration for a pulsar timing experiment with the SKA-low Phase-1 Telescope }
 \label{fig:skalow_sub}
 \end{figure*}

\begin{table*}[ht]
\begin{tabular}{|c|c|c|}
\hline
List of pulsars & Number of Pulsars & Observing time required \\
& & Option I/II/III/IV  \\
& & (hrs) ) \\
\hline
\hline
J0437-4715, J0621+1002, J0636+5128, J0900-3144, J1022+1001 & & \\
J1045-4509, J1600-3053, J1603-7202, J1640+2224, J1643-1224 & 24 & 4/2/2 \\
J1713+0747,J1730-2304, J1732-5049, J1747-4036, J1744-1134& &  \\
J1832-0836, B1855+09, J1903+0327, J1909-3744, B1937+21 & & \\
J1944+0907, J2124-3358, J2145-0750, J2234+0944& &  \\
\hline
\hline
J0023+0923, J0030+0451, J0034-0534, J0218+4232, J0340+4130& & \\
J0610-2100, J0613-0200, J0711-6830, J0751+1807, J0931-1902& & \\
J1024-0719, J1453+1902, J1455-3330, J1614-2230, J1721-2457& & \\
J1738+0333, J1741+1351, J1802-2124, J1804-2717, B1821-24A & 38 & 22 \\
J1853+1303, J1910+1256, J1911-1114, J1911+1347, J1918-0642 &&\\
J1923+2515, J1946+3417, J1949+3106, B1953+29, J2010-1323 &&\\
J2017+0603, J2033+1734, J2043+1711, J2129-5721, J2234+0611 && \\
J2214+3000, J2302+4442 J2317+1439& & \\
\hline
\hline
\end{tabular}
\caption{Recommended observation strategy for a PTA program with the SKA telescope. The pulsar sample was divided into two categories. The bright pulsars could be observed for 10 minutes each using Option I (Single Target simultaneous observations with three multi-band  sub-arrays), Option II (Multiple targets with 4 sub-arrays in three non-simultaneous multi-band observations) or Option III(two targets with two band sub-arrays each). For the fainter pulsars in the second category, option IV is recommended with the full array utilised in three multi-band observations)}
\label{tab:obsstrat}
\end{table*}

\par
With the above considerations in mind, we classified the PTA sample into two different categories. The first category consists of the brightest pulsars, where we adopt OPTION I, II or III using the maximum number of sub-arrays available for each telescope optimising the telescope time as well as frequency coverage. The second  category of  pulsars consists of fainter pulsars, where option IV is appropriate. Given the advertised parameters of SKA telescopes, the required observation time  for each pulsar in the respective categories using these options was estimated from the radiometer equation using the flux density estimates from the ATNF pulsar catalog \citep{mhth05} \footnote{\url{https://www.atnf.csiro.au/research/pulsar/psrcat/}}. These calculations provide an estimate for the required telescope time. While details of some parameters of SKA-telescopes are still being worked out in the construction phase and a detailed optimisation is planned in future, some recommendations are made below based on InPTA experience and available information.

\par
A list of target pulsars to be observed using the suggested four options for the SKA-mid telescope is provided in Table \ref{tab:obsstrat}.  The total observing time required for each epoch of observations is also indicated in this table. We find that a sample of 62 pulsars, visible with the SKA telescope, can be observed with about 26 hours with the SKA-mid telescope. Likewise, seven hours are needed for each epoch with the  SKA-low telescope.  Assuming a 10-day cadence, the SKA key science program will require 936 hours for a PTA program with the SKA. This estimate can vary depending on the number of available antennas in the early construction phase of the SKA telescope (which increases progressively from the current 64 MeerKat antennas during AA1 to AA* phase) as well as on achievable SEFD and backend configuration and evolving radio Frequency Interference environment. Finally, the observing proposed here can be further optimised by reserving the SKA observations to weaker pulsars or pulsars which are affected by significant inter-stellar propagation effects. While a more detailed study is planned in future incorporating this information, the recommendation in this section already provides a blue-print for an SKA PTA, probably for the first time.

\par
In conclusion, the experience gained with different configurations of the uGMRT interferometer for frequencies between 300 to 1500 MHz in the InPTA experiment is useful to propose observation configurations that can be used once the SKA phase 1 telescope becomes available. Suggestions for the required observing strategies are described in some detail in this section and can form the basis of future SKA-PTA program of the SKA pulsar key science project.

\section{Introduction}

\par
Gravitational waves (GWs) are propagating ripples in the space-time curvature  as predicted by Einstein's General Theory of Relativity \citep{ein18}.  The routine detection of high-frequency GWs (30$-$800 Hz) by the LIGO-Virgo-KAGRA collaboration\footnote{\url{https://www.ligo.org/science/Publication-O3aCatalog/index.php} } from around 100 merger events that involve black holes (BHs) and neutron stars (NSs) in the last 6 years has opened a new window to the universe which is complementary to the traditional electromagnetic window \citep{Abbot2019_GWTC1,Abbot2020_GWTC-2,Abbot2021_GWTC3}. These observations are providing key insights into a number of astrophysical puzzles that are difficult to tackle by traditional electromagnetic astronomy \citep{Arimoto21}.
This includes the formation channels for stellar-mass BHs, Equation of State for NSs, measuring the Hubble constant by employing standard and dark sirens, and testing general relativity in ultra-strong regimes,  as well as constraining alternative theories which  dispense with dark energy and dark matter    \citep{Arimoto21,Mapelli20,BB17,StdSiren17,Boran18,DrkSiren19,TGR21}.
Therefore, it is reasonable to expect that the other GW frequency windows should provide opportunities to explore many aspects of physics, astrophysics and cosmology \citep{Bailes21}.
Specifically, nanohertz stochastic gravitational wave background (SGWB) formed by the incoherent superposition of GWs from supermassive black hole binaries (SMBHB) is expected to be detected in the coming years \citep{IPTA_GWB_2022}. The discovery of such nanohertz (nHz) GWs will provide unique insights into the formation and evolution of galaxies and their constituent SMBHs, a precise description of solar system, cosmological standard sirens and  fundamental tests of  gravitation \citep{Burke-Spolaor2019}.

\par
These nHz GWs can be detected by employing the precision timing of an ensemble of radio millisecond pulsars \citep[MSPs : ][]{Detweiler1979,S78}.
A pulsar timing array (PTA) is a dedicated  experiment  to time pulsars with the aim of detecting low-frequency GWs \citep[PTAs:][]{FosterBacker1990}. At present, there exist four such established efforts that  employ the world's best radio telescopes along with three more emerging experiments\footnote{\url{http://ipta4gw.org/}}. PTAs require long-term high-cadence monitoring of a large ensemble of MSPs using large highly sensitive radio telescopes. The high sensitivity of the Square Kilometer Array (SKA) telescope  with its multiple-beam design will not only help to discover nHz GW sources, but also will strengthen the  post-discovery nHz GW science with sensitive high-cadence  observations. This is one of the main goals of the pulsar key science project of the SKA \citep{ks15,jhm+15}.
\par
While pulsar astronomy has traditionally been carried out with large single-dish telescopes, SKA is designed to be an interferometer where the signal collected over its large collecting area is recorded using a beam-former. Out of the currently operating four PTAs, only the Indian Pulsar Timing Array (InPTA) primarily uses  an interferometer, the upgraded Giant Metre-wave Radio Telescope (uGMRT), for its observations\footnote{It may be noted that Large European Array for pulsars \citep[LEAP:][]{Bassa2016} also uses a phased array of multi-element telescopes, but these form a subset of European Pulsar Timing Array experiment, which is largely based on single dish observations. Another interferometer, which contributes data to European Pulsar Timing Array experiment is the Westerbok Synthesis Radio Telescope. Recently, the MeerTime experiment \citep{Bailes2018_MeerTime} has started collecting data with MeerKat, which is also an interferometer. It may be noted that these telescopes have not been used as sub-arrays unlike the uGMRT}. The uGMRT is also a pathfinder telescope for the SKA with features similar to the SKA. Hence, the observation strategies uniquely adopted in the InPTA program are very much relevant for initiating the SKA-era PTA program.

\par
One of the unique features of InPTA is characterisation of propagation effects with unprecedented precision, improving the overall precision of the pulsar timing required for eventual nHz GW detection. The techniques developed for the analysis of the low radio frequency data by the InPTA astronomers would allow the Indian community to meaningfully contribute to  PTA work with the SKA. Additionally, InPTA researchers are developing approaches to model nHz GW sources  and associated detection techniques with a particular emphasis on orbital eccentricities. These efforts should provide another pathway for us to contribute to the SKA science efforts, particularly from the multi-messenger perspectives.

\par
We now present the relevance of the InPTA program for the discovery and characterisation of nHz GW sources with the SKA telescope, motivated by the above-mentioned considerations. The structure of the paper is as  follows. A brief review of the current state-of-art in the detection of nHz GWs and the existing pulsar timing array experiments is provided in Section \ref{pta} followed by a description of the  InPTA experiment in Section \ref{inpta}. A comparison of the uGMRT and the SKA from this perspective is discussed in Section \ref{inptanska} outlining features of the InPTA experiment which are relevant to a future SKA program. A discussion of analysis techniques and methods developed for the InPTA program,  which are likely to be carried forward to the SKA program, is presented in Section \ref{method}. Possibilities for  multi-messenger astronomy for individual sources in the context of the SKA are outlined in Section \ref{mmess}.  After a discussion on  contributions of the InPTA on development of the Indian and Japanese neutron star community as future SKA ready users in Section \ref{ucd}, we conclude by listing future directions in Section \ref{conc}.

\section{Analysis methodology from InPTA relevant to SKA}
\label{method}

\par
The reduction and the analysis of the InPTA data, particularly that obtained at the lower frequency band, required the development of new techniques or tuning  the techniques available in the literature. As mentioned in Section \ref{inpta}, the uGMRT raw data were reduced by using an off-line pipeline, called \texttt{pinta}, to obtain  partially folded PSRFITS archives.  While the SKA will employ a real-time pipeline to process pulsar observations, the RFI mitigation and excision techniques (\texttt{RFIClean} and \texttt{gptool}) will be relevant for the SKA pipeline as well, particularly for Band 1 data from SKA-mid and data from SKA-low. The final product of SKA real-time pipeline  is also expected to be in the PSRFITS format, which will allow the InPTA techniques for determination and application of precision DMs to be incorporated in the SKA Science Data processor easily.  The experience relevant to the analysis of the data from the SKA beam-former is described in this section.

\subsection{Online RFI-excision before partially folding the data}
\label{rfi_excision}

\par
Radio frequency interference (RFI) remains a growing concern to obtain quality data and achieve the designed sensitivity. For most of the PTA experiments, the data are processed in real-time to output data that are partially folded over the pulsar spin period. While this real-time processing results in highly reduced data rates as well as requirements of disk-space and computing resources for off-line processing, it is much harder to identify and excise RFI from the resultant partially integrated data. Moreover, a variety of RFI, e.g., the periodic RFI, faint and short time-scale RFI, and temporally varying RFI that mostly causes baseline variations, is not possible to mitigate from the partially integrated data as it requires analysing contiguous temporal sections of data to identify such RFI. The above understanding is reinforced by our InPTA data processing experience, wherein mitigation of RFI in the filterbank data using the excision tools \texttt{RFIClean}\xspace \citep{Maan2020_rfiClean} and \texttt{gptool} greatly aid in the detection significance. In fact, in some particular cases, the pulsars are detected only after RFI is excised using these tools. The SKA will also employ a real-time pipeline to process pulsar observations. To fully utilise the telescope's sensitivity, it is desirable that the above or similar RFI excision tools are employed in the SKA's real-time data reduction pipeline to mitigate RFI before the data are partially folded. The computation time estimates based on the InPTA experiment \citep{smj+21} suggest that these can run in real-time, particularly on GPGPU based hardware as the source code is easily portable to  GPGPU kernels.

\subsection{Narrow-band technique for times-of-arrival and dispersion measure estimation}
\label{narrowband}

\par
Precise estimation of times-of-arrival (ToA) and DM at lower frequencies with the help of traditional narrow-band technique is impeded by the following factors: (1) the measured DM and frequency-dependent profile evolution are intertwined, (2) uncorrected DM variations can introduce smear while frequency-collapsing the data to degrade the S/N, and (3) systematic biases can also get introduced in the estimated DM due to interstellar scatter broadening. A frequency resolved narrow-band analysis method was developed for the post-processing of the InPTA data. This method is called DMcalc and utilises python interface of PSRCHIVE package \citep{kmj+21}.

\par
The ToAs are obtained in a timing experiment by cross-correlating the data profiles with a noise-free template. Traditionally, a noise-free template is obtained by averaging  frequency-collapsed IPP from several epochs. Epoch to epoch DM variations can introduce significant smear when multiple widely separated epochs are averaged to obtain a noise-free template. Additionally, an incorrect alignment of the pulse across the frequencies for pulsars showing significant evolution of profiles with frequency leads to DM offsets. It may be noted that both these effects are significantly larger at lower frequencies.  Indeed, our work shows that DM estimates using a single frequency-collapsed template or an offset fiducial DM used to dedisperse the template itself, can have  constant offsets from the actual DMs \citep{kmj+21}. Traditionally, this is handled by using FD parameters or DMJUMPS in the timing analysis.

\par
In the InPTA narrow-band analysis, we use a frequency-resolved template and an iterative procedure to estimate DM of each epoch. The number of sub-bands used in the frequency resolved template is a trade-off between the profile evolution and the achievable S/N in each sub-band and all DM measurements use frequency resolved ToAs obtained by cross-correlating the data with frequency resolved noise-free templates. The sub-band profiles in the templates for each epoch are aligned by determining the epoch dependent DM in an iterative procedure. Here, two or more sub-integrations Band 3 data of each epoch are analysed using ToAs obtained by cross-correlating time-collapsed  profiles of the epoch aligned first with an assumed DM obtained from \texttt{pdmp} analysis \citep{kmj+21}. This analysis is repeated in subsequent iterations progressively improving the  DM estimate. A final frequency resolved noise-free template is then formed by averaging the epochs aligned using the DM obtained in the final iteration.

\par
The main conclusion from the InPTA analysis with the procedure described here is that FD parameters or DMJUMPS are not needed for any of the pulsars in the InPTA sample. Considering that the sample includes pulsars which show significant profile evolution from 300 to 1460 MHz as well as significant scatter broadening at frequencies around 400 MHz (e.g. PSRs J1643$-$1224 and J2145$-$0750), this suggests that using such frequency resolved analysis appears to be more robust. Please note that our simultaneous multi-frequency observations avoid any possibility of small variations in DM between near-simultaneous high and low frequency observations adopted in other PTA experiments. Additionally, we use  higher precision DM measurements from 300 to 500 MHz. Thus, we recommend following an observing strategy and analysis similar to the InPTA for a PTA experiment with the SKA. Thus, OPTION I and III (See Section \ref{inptanska}) are preferred options provided the required S/N is achievable for the target pulsars.

\par
The post processing pipeline \texttt{DMcalc}{} was developed for application on the combined multi-band data with the noise-free template updated with the fiducial DM obtained in the previous steps. It is straight-forward to incorporate this in the SKA science data processor with a possibility of providing real-time measurements. In turn, these will be useful for a quick follow-up of DM events, such as those reported in PSR J1713+0747 \citep{Lam2018}, as well as updating DM for coherent dedispersion for the subsequent observations.

\subsection{Wideband technique for simultaneous measurement of times-of-arrival and dispersion measure}
\label{wideband}

\par
Pulsar observations are increasingly using wideband receivers (defined as receivers with a fractional band-pass of 100 percent or greater) to acquire more precise TOAs \citep{gak+17,Johnston+2021}. Pulse profiles can evolve significantly over such wide bandwidths. This introduces a bias in arrival times as well as DM measurements when such data are analysed with a single band-collapsed 1-dimensional template.  \citet{Pennucci+2014} devised an algorithm for generating a 2-dimensional template and for using these for simultaneous measurement of DMs and ToAs from wideband pulsar data. They proposed pulsar timing with such wide-band ToAs\citep{Pennucci2019ApJ}. This "wideband timing technique" has been applied on various datasets, such as the NANOGrav 12.5-year data \citep{Alam+2021ApJS} and the CHIME and GBT-L data in the  400$-$800 MHz frequency range \citep{Fonseca+2021ApJ}.  Recently, we applied this technique to the uGMRT data in the 300$-$500 MHz frequency range \citep{Nobleson+2021}. As the profile evolution is much more dominant for many pulsars below 400 MHz apart from a significant smearing due to pulse broadening, this study demonstrated the effectiveness of this technique below 400 MHz for the first time.

\par
The wideband timing technique uses frequency-resolved high S/N IPP to decompose this 2-dimensional data into principal component eigen-profiles which capture the variation of the profile with the observing frequency. A noise-free 2-dimensional \texttt{PulsePortraiture} is generated by fitting spline functions to the coefficients of eigen-profiles  required for modeling the 2-dimensional data \footnote{\url{https://github.com/pennucci/PulsePortraiture}} \citep{2016ascl.soft06013P}.  The phase offset of the wide-band data with the pulse portrait is obtained by least-square minimization of the phase offset of each sub-band pulse profile incorporating the cold-plasma dispersion law. This yields both a high precision ToA and DM estimate. The timing solution and residuals for a PTA  pulsar are then obtained by maximising a likelihood function, implemented in the software package TEMPO \citep{Nice+2015ascl}, with these ToAs and DM estimates as inputs.

\par
Pulsars, such as PSRs J1643$-$1224 and J1939+2134, show significant profile evolution due to scatter-broadening in the InPTA sample, whereas PSR J2145$-$0750 shows an evolving profile. Analysis of Band 3 InPTA data needed three eigen-profiles to model the pulse-portrait of the former two pulsars, while two eigen-profiles were required for the latter \citep{Nobleson+2021}. We also noticed that RFI contributes dominant eigenvalues and corresponding noisy eigen-profiles, so it is important to remove RFI  before the principal component analysis for this technique to work well. The wavelet and  spline fits needed to be appropriately tuned for the technique to work  at our frequency. This experience can easily be adopted for SKA data.

\par
We are currently developing the technique to include the likelihood function used in the TEMPO2 \citep{Hobbs2006_tempo2} software package which will use the wideband DM measurements from the ToAs as priors to estimate the DM model parameters. We also plan to develop a technique to combine multi-band data in the regime of wideband analysis. The standardisation of techniques as part of the InPTA analysis is  directly relevant  to the SKA due to the wide multi-band frequency coverage of the SKA in view of the observations strategy proposed in this paper.

\section{Persistent Multi-Messenger GW Astronomy with the SKA }
\label{mmess}

\par
The LIGO-Virgo-KAGRA collaboration has inaugurated the era of GW astronomy due to their routine detection of hecto-hertz GWs from around 100 stellar-mass BH binaries since 2015 \citep{GWTC2_2020, GWTC3_2021}.
Further, this collaboration heralded the era of multi-messenger GW astronomy by the observations of transient hectohertz GWs from a binary neutron star merger (GW170817) and its counterparts (EM170817) across the electromagnetic spectrum \citep{GW170817_2017, GW170817_EM}.
Strikingly, these multi-messenger observations of
such a $100$ second-long transient GW event
are already providing fundamental contributions to the areas of  physics, astrophysics, and cosmology
\citep{GW170817_2017}.
Recall that multi-messenger astronomy usually involves combining data from  various messengers that arise from spectacular astrophysical events.
\par
During the SKA era, we should be able to practice what we term  the persistent multi-messenger GW astronomy \citep{VDG2021}.
This is influenced by the fact that the rapidly maturing PTAs should be detecting a stochastic nHz GW background created by a population of merging SMBHBs in the coming years \citep{IPTA_GWB_2022}.
Thereafter, continuous nHz GWs from bright individual SMBHBs that stand above this background should be routinely detected with the help of the SKA-era PTA \citep{Burke-Spolaor2019}.
The eventual detection of continuous nHz GWs from individual SMBHBs should allow us to pursue multi-wavelength electromagnetic observations of such GW sources as they are expected to reside in the active galactic nuclei \citep{Burke-Spolaor2019,Xin2021}.
In other words, persistent multi-wavelength electromagnetic observations of such continuous GW sources should allow us to pursue multi-messenger GW astronomy during the SKA era \citep{VDG2021}.
It turns out that the resulting persistent multi-messenger nHz GW astronomy can make profound contributions to the emerging field of multi-messenger GW astronomy \citep{Burke-Spolaor2019}.
\par
The InPTA researchers are well-placed to contribute to the various aspects of this emerging area of astronomy during the SKA era.
This is partly due to our involvement in the general relativity-based predictions of certain  Bremsstrahlung flares from the unique bright blazar OJ~287 and the successful observational campaigns of these predicted flares \citep{Dey2018, Dey2019, Laine2020}.
These investigations and the associated campaigns provided firm observational evidence
for the presence of a nHz GW-emitting SMBH binary in the blazar OJ~287.
The underlying model, developed by Mauri Valtonen and his collaborators, prescribes the observed double-peaked Bremsstrahlung flares from OJ~287 to the repeated impacts of a secondary BH with the accretion disk of the primary BH twice during its general relativistic eccentric orbit.
Notably,  the 2019 {Eddington flare} from OJ~287, predicted in \cite{Dey2018}, was successfully observed with the {Spitzer space telescope} and this campaign even allowed us to provide a parametric constraint on the celebrated BH no-hair theorem \citep{Laine2020}.
\par
These successful observational campaigns of a number of predicted impact flares in optical wavelengths prompted us to explore the possibility of using other electromagnetic wavelengths to probe the SMBH binary central engine description for OJ~287 \citep{VDG2021}.
Fortunately, several high-resolution Very Long Baseline Interferometry (VLBI) observational campaigns have tried to image the radio jet of OJ~287 during the past three decades \citep{Hodgson2017, Cohen2017}.
These observational campaigns reveal that the position angle (PA) of the projected jet on the sky plane shows systematic variations at both millimeter and centimeter wavelengths \citep{Dey2021}.
Therefore, some of us pursued detailed investigations that allowed us to explain the observed PA variations of OJ~287’s radio jet while employing the now established SMBH binary central engine description \citep{Dey2021}.
Further, we have provided specific predictions that should be testable by the GMVA and EHT campaigns on OJ~287, planned for the coming years.
These efforts are ensuring that multi-band EM observations of OJ~287 will be pursued in the coming years to probe various aspects of its SMBH binary central engine description \citep{VDG2021}.
It is expected that
such multi-wavelength and multidisciplinary efforts will be required to pursue multi-messenger nHz GW astronomy with OJ 287 during the SKA era.

\par
We are also developing general relativistic approaches that will be required for extracting nHz GWs from isolated SMBH binaries from the SKA era PTA data sets.
Specifically, efforts are ongoing to develop ready-to-use computationally efficient routines that model the expected pulsar timing residuals induced by GWs from isolated SMBH binaries inspiraling along general relativistic eccentric orbits \citep{Susobhanan2020}.
These approaches employ post-Newtonian(PN) approximation for describing the general relativistic SMBH binary dynamics and this approximation provides general relativity based  corrections to Newtonian  dynamics in terms of a dimensionless parameter that involves SMBH binary orbital velocity \citep{Blanchet2014}.
In \citet{Susobhanan2020}, we employed 3PN-accurate Keplerian parametric solution to model general relativistic eccentric orbits of SMBH binaries while an improved version of GW phasing approach was employed to incorporate the effects of GW emission \citep{MGS2004,DGI2004}.
It may be recalled that these two elements, namely PN-accurate Keplerian type parametric solution and the GW phasing approach,  are employed in the widely used Damour-Deruelle timing formula to time MSPs in general relativistic eccentric orbits \citep{DD1986}.
Further, we are incorporating the spin effects to model pulsar timing residuals induced by inspiral GWs from spinning SMBH binaries in PN-accurate eccentric orbits with the help of \citet{KG2005}

\par
In Figure~\ref{fig:sky_maps}, we show the amplitude of the expected GW-induced pulsar residuals for pulsars across the sky for two SKA era nHz GW emitting
SMBH binaries, namely blazars OJ~287 and PKS~2131-021
while focusing only on the earth term contributions
\citep{ONeill2022}.
The inference of a nHz GW emitting SMBH binary in PKS~2131-021 is attributed to a very recent effort that demonstrated two epochs of strong sinusoidal variation with essentially the same period and phase in the 45.1-year radio light curve \citep{ONeill2022}.
A natural explanation is that this blazar is hosting an SMBH binary with an orbital period of $\sim 4.75$ years (redshifted) in contrast to $\sim 12$ years for OJ~287.
The heat maps in Figure~\ref{fig:sky_maps} reveal that SKA-era PTA efforts should be able to pursue persistent multi-messenger nHz GW astronomy with these two sources that may involve other facilities like the  Thirty Meter Telescope and The Next Generation Event Horizon Telescope \citep{VDG2021}.
\par
 With the help of such general relativistic constructs,

we are working to constrain the presence of isolated eccentric SMBH binaries in PTA data sets, influenced by \citet{N11_isolated}.

 Additionally, we are modelling general relativistic burst with memory PTA signals due to SMBH binaries in PN-accurate hyperbolic orbits with the help of \citet{Cho2018}.
These efforts should allow us to contribute to the ongoing efforts to constrain the presence of such isolated SMBH binaries in the existing IPTA data sets and their updates.
Clearly, such efforts will be crucial to detect routinely  isolated GW sources in the SKA era PTA data sets and therefore pursue persistent multi-messenger GW astronomy during that era.

\begin{figure*}
    \centering
    \includegraphics[width=0.75\textwidth]{OJ287_sky_sensitivity.pdf}\\
    \includegraphics[width=0.75\textwidth]{PKS2131-021_sky_sensitivity.pdf}
    \caption{Distribution of the GW-induced Earth term amplitudes, namely $max - min$ of the Earth term contributions to the timing residuals, for pulsars around two different sky locations associated with two SMBHB candidates: OJ~287 and PKS2131-021. The binary parameters of OJ~287 are taken from \citet{Dey2018}
    while we employ \citet{ONeill2022} for

    PKS2131-021.
    Influenced by \citet{ONeill2022}, we let the total mass $M = 2 \times 10^9$ and mass ratio $q=1$ for the second candidate. The GW-induced residuals are calculated by extending the method of \citet{Susobhanan2020} while incorporating PN-accurate GW emission effects for  tracking the evolution of the SMBHBs. }
    \label{fig:sky_maps}
\end{figure*}
\section{Pulsar Timing Arrays for Nanohertz GW Astronomy}
\label{pta}

\par
SMBHBs with total masses in the $10^{8}-10^{10}\, M_{\odot}$ range  are the primary sources of nHz GWs \citep{Detweiler1979}. Such binaries emit nHz GWs during their GW emission-dominated orbital evolution phases and are expected to have orbital periods of years. This is because the natural (orbital) frequency of such self-gravitating BH binary systems  roughly follow $\omega^2\,r^3 = G\,M$, where $\omega, r$ and $M$ are the orbital angular frequency, the associated separation and the total mass of the system respectively  \citep{SS2009}. For a typical circular SMBHB system,  the frequency of the emitted GWs is given by  \citep{Detweiler1979}

\begin{equation}
    f_{\rm GW} \sim 20\,\text{nHz}\,\left(\frac{200\,M}{r}\right)^{3/2}\,\left(\frac{10^{10}\,M_{\odot}}{M}\right) \,,
\end{equation}
where we have used the fact that the frequency of GWs emitted from a circular binary is twice its orbital frequency.

\par
As noted earlier, the existing PTAs monitor an ensemble of MSPs to detect  nHz GWs emitted by SMBHBs \citep{FosterBacker1990, HobbsDai2017_PTA}. This is because the times of arrival (TOAs) of MSP pulses can be modelled and predicted with very high accuracy due to the extremely stable rotational periods of MSPs \citep{Taylor1992}. It turns out that such predictions require the use of an accurate  \textit{timing model} that incorporates various astrophysical and instrumental delays \citep{Edwards2006_tempo2}. Typically, PTAs compare the observed TOAs from MSPs  with predictions for their arrival times based on a pulsar timing model. The differences between the predictions and the actual measurements are usually referred to as the pulsar `timing residuals'. When a GW passes between the Earth and an MSP, it affects the geodesics along which the pulses arrive at the Earth and therefore modulates the pulsar TOAs \citep{EstabrookWahlquist1975, Anholm2009}. The resulting GW-induced timing residuals at a given epoch depends on the two GW polarization states $h_{+,\times}$ such that \citep{Susobhanan2020}
{\small
\begin{align}
R(t)=&\begin{bmatrix}F_{+} & F_{\times}\end{bmatrix}\begin{bmatrix}\cos2\psi & -\sin2\psi\\
\sin2\psi & \cos2\psi
\end{bmatrix}\nonumber\\
&\qquad\int_{0}^{t}dt'\,\begin{bmatrix}h_{+}(t')-h_{+}(t'+\Delta_{p})\\
h_{\times}(t')-h_{\times}(t'+\Delta_{p})
\end{bmatrix}\,,
\end{align}
}
where $\Delta_p$ is a geometric delay, $F_{+,\times}$ are known as the antenna pattern functions, $\psi$ is the GW polarization angle, and the integration is over the coordinate time measured at the solar system barycentre. The expressions for $\Delta_p$ and $F_{+,\times}$ may be found in, e.g., \citet{Anholm2009}.

\par
The existing PTAs are expected to detect first a stochastic GW background (SGWB) created by the incoherent superposition of  GWs from the whole cosmic population of in-spiralling MBHBs \citep{ptk+21}. Moreover, it is customary to treat the resulting SGWB in the $10^{-9}-10^{-7}$ Hz (GW) frequency window  to be an isotropic Gaussian  stationary signal which is described by a characteristic amplitude $h_c$ that follows a power-law spectrum, $h_c \propto f^{-2/3}$ \citep{Phinney2001}. Further, due to the quadrupolar nature of GWs, the above SGWB is expected to induce the following correlation between the timing residuals of each pulsar pair of PTAs \citep{HellingsDowns1983}:
\begin{align}
c(\theta) &= \frac{3}{2} x \, \ln{x} - \frac{x}{4} + \frac{1}{2} + \frac{1}{2} \delta(x)\,,
\end{align}
where $x = (1 - \cos{\theta} )/2$ for an angle $\theta$ on the sky between two pulsars, and $\delta(x)$ is the  regular Dirac delta function (the resulting curve is usually referred to as the Hellings-and-Downs curve). Therefore, the actual detection of the nHz SGWB, characterised by certain $h_c$, involves detailed investigations that determine how close the observed correlations between the timing residuals of each pulsar pair of a given PTA follow the Hellings-and-Downs curve. An incontrovertible detection of the SGWB due to inspiralling SMBHBs will be made if the timing residuals of all pulsar-pairs in a PTA are shown to follow the Hellings-and-Downs curve with a common red spectrum. In order to achieve these goals,  PTA experiments usually monitor MSPs with a cadence of few weeks and the resulting timing data set now spans a few decades for the three established PTAs \citep{IPTA_DR2}. This ensures that the existing PTA data sets should be sensitive to GWs with periods from weeks to years  corresponding to a nHz GW frequency window typically in the $10^{-8} - 10^{-9}$ Hz range.

\par
Currently, four such experiments are operational with a time-baseline of more than a few years. While Parkes Pulsar Timing Array \citep[PPTA :][]{mhb+13}, North American Nano-hertz Observatory for Gravitational Waves \citep[NANOGrav :][]{abb+18} and European Pulsar Timing Array \citep[EPTA :][]{dcl+16} have been gathering data for more than a decade, an Indian effort that employs the uGMRT, called the Indian pulsar timing array (InPTA),  has been in operation since 2015 \citep{jab+18}. Additionally, PTA efforts are being initiated by employing  two more facilities, namely the Five hundred meter radio telescope (FAST) and the MeerKat \citep{CPTA2016,SAPTA2016}. The PTAs periodically pool their data together to create the International Pulsar Timing Array (IPTA) data releases. Two such data releases have happened during the last decade \citep{IPTA_DR1, IPTA_DR2}. IPTA is  a consortium of consortia that is governed by a steering committee with members from all the above mentioned PTAs.

\par
All the well-established PTAs regularly make public their data sets along with their GW astronomy-related analysis results \citep{Alam2020_NANOGrav12.5_narrowband,Kerr2020_PPTA_DR2,Desvignes2016_EPTA, Chen2021_EPTA_6psr}.
Within the last two years, detailed and independent investigations of the most up-to-date NANOGrav, PPTA and EPTA data sets reveal the presence of  a common red noise process in their pulsars with spectral signature consistent with an SGWB due to in-spiralling SMBHBs \citep{Arzoumanian2020_SGWB,Goncharov2021,Chen2021_EPTA_6psr}.
However, these data sets do not strongly support the presence of the  Hellings-and-Downs spatial correlation, crucial to establish the presence of nHz SGWB that arises from merging SMBHBs. This is being attributed to the presence of various statistical uncertainties in the TOA measurements.  Similar conclusions were drawn from a detailed analysis of more sensitive IPTA DR2 data set formed from a combination of older public data releases from the individual PTAs  \citep{IPTA_DR2,IPTA_GWB_2022}.
At present, the constituent members of IPTA are working in tandem by updating their data sets and employing  multiple analysis algorithms to increase the confidence in a possible detection of SGWB in the coming years.

\par
Unfortunately, the stochastic nature of our astrophysical GWB signal can be mimicked by other random processes that occur in the PTA pulsars. Specifically, the epoch to epoch variations in dispersion measure (DM)\footnote{Dispersion measure is defined as the integral of the column density of electrons over the line-of-sight to the pulsar} and scattering variations need to be accurately incorporated while obtaining TOA residuals from various PTA data sets. These propagation effects are due to the relative motion of the pulsar, the earth and the inhomogenous ionised inter-stellar and inter-planetary medium (IISM and IPM). These  are strong functions of frequency thereby dominating low radio frequency observations of the PTA pulsars \citep{dvt+20,tsb+21}. The wide-band low-frequency coverage in the InPTA experiment provides very high precision measurements of these effects \citep{kmj+21,Nobleson+2021}. A comparable precision is expected from  a similar future  experiment with the SKA. It was argued that one needs to model these variations appropriately to ensure the correctness of the astrophysical interpretation of the expected SGWB signals, extracted from   the timing residuals \citep{Lentati2016}. Work is in progress to build such models for incorporating in a future IPTA data release.

\section{InPTA and the development of SKA user community in Japan and India}
\label{ucd}

The InPTA consortium is playing an active role in building a scientific community that will have expertise to use the SKA products to pursue its key science projects.
Taking a cue from the LIGO collaboration, which recruited many high energy physicists into  gravitational wave astronomy (after the cancellation of the SSC), we have tried to build a community of researchers from diverse areas such as  Cosmology, General Relativity, Nuclear and High Energy Physics, Astronomical Instrumentation, Neutron Star Physics, Solar Physics to list a few examples. These researchers have nicely complemented the original core group, which consisted of
radio astronomers
working on various aspects of pulsar astronomy.

At present, InPTA hosts seven undergrad students, nine graduate students, and ten postdoctoral researchers and we expect this collaboration to grow further as we head into the SKA era. All these young researchers are now well  experienced in performing pulsar observations simultaneously at different uGMRT bands in the sub-array mode, which will be a necessity in the SKA era. There were some other undergraduate students, who were members of InPTA earlier and then moved to different research fields for their Ph.D,   as well as some former postdocs who have moved on to non-academic jobs.
We gather that the training they received as members of InPTA remains  invaluable to them.  Furthermore,  several InPTA members who started with InPTA projects while being based at Indian institutes have continued to remain active after received their under-graduate and graduate degrees, in spite of taking on new positions abroad.

To train new members, InPTA organises training schools and busy weeks (where the whole collaboration members gather for in person interactions).
Further, young InPTA members are encouraged to present research results in national and international levels and they participate in outreach activities.
InPTA members have organised a conference series titled ``Asia SKA Initiative On NS'' twice, first time in 2016 in India and the second time in 2019 in Thailand. Both of the meetings were attended by scientists from India, China, Japan, Thailand, etc. The participants were from various stages of their careers.
InPTA members also played a key role to form a consortium called `Gravitational Radiation and Science with Pulsars' (GRASP) that involve researchers from South Africa and China in addition to InPTA members. Monthly meetings of GRASP is playing a significant role in exchange of knowledge and will help in future international collaborations in the SKA era.
Finally,
InPTA hosted the annual meeting of the IPTA during June 2019 in Pune, India. This meeting included a student training week followed by a conference week. In addition to the 2019 meeting, InPTA members have acted as as members of the organising committees of other IPTA meetings too, as well as delivered lectures in student weeks and presented research works.

In what follows, we provide a brief timeline of our consortium.
InPTA started in 2015 as a pilot program consisting of only five members using the legacy GMRT. It was accompanied by another pilot program using uGMRT during 2016$-$2017. Thereafter, the two programs merged and InPTA program using uGMRT emerged with 19 members during 2018. We were an associate member of the IPTA consortium  until 2021, when we were given the full membership. In 2021, four researchers from Japan (based at Kumamoto University) joined the InPTA consortium. Hence the word ``InPTA'' is now a moniker for an Indo-Japanese  pulsar timing array program using an Indian telescope, but is not confined to only Indian researchers.