The sensitivity of ongoing searches for gravitational wave (GW) sources in the ultra-low frequency regime ($10^{-9}$ Hz to $10^{-7}$ Hz) using Pulsar Timing Arrays (PTAs) will continue to increase in the future as more well-timed pulsars are added to the arrays. It is expected that next-generation radio telescopes, namely, the Five-hundred-meter Aperture Spherical radio Telescope (FAST) and the Square Kilometer Array (SKA), will grow the number of well-timed pulsars to $O(10^3)$. The higher sensitivity will result in greater distance reach for GW sources, uncovering multiple resolvable GW sources in addition to an unresolved population. Data analysis techniques are, therefore, required that can search for and resolve multiple signals present simultaneously in PTA data. The multisource resolution problem in PTA data analysis poses a unique set of challenges such as non-uniformly sampled data, a large number of so-called pulsar phase parameters that arise from the inaccurately measured distances to the pulsars, and poor separation of signals in the Fourier domain due to a small number of cycles in the observed waveforms.

18 pages, 8 figures, appendix A was changed to further clarify AvPhase, submitted to PRD