This paper provides an efficient recursive approach of the spectral Tau method, to approximate the solution of system of generalized Abel-Volterra integral equations. In this regards, we first investigate the existence, uniqueness as well as smoothness of the solutions under various assumptions on the given data. Next, from a numerical perspective, we express approximated solution as a linear combination of suitable canonical polynomials which are constructed by an easy to use recursive formula. Mostly, the unknown parameters are calculated by solving a low dimensional algebraic systems independent of degree of approximation which prevent from high computational costs. Obviously, due to singular behavior of the exact solution, using classical polynomials to construct canonical polynomials, leads to low accuracy results. In this regards, we develop a new fractional order canonical polynomials using Müntz-Legendre polynomials which have a same asymptotic behavior with the solution of underlying problem. The convergence analysis is discussed, and the familiar spectral accuracy is achieved in L ∞ -norm. Finally, the reliability of the method is evaluated using various problems.Keywords The recursive approach of the Tau method • System of generalized Abel-Volterra integral equations • Müntz-Legendre polynomials • Fractional vector canonical polynomials • convergence analysis.