This paper presents a method which applies pseudospectral tau approximation for retarded functional differential equations (RFDEs). The goal is to construct a system of ordinary differential equations (ODEs), which provides a finite dimensional approximation of the original RFDE. The method can be used to determine approximate stability diagrams for RFDEs. Thorough numerical case studies show that the rightmost characteristic roots of the ODE approximation converge to the rightmost characteristic roots of the original RFDE. Application of the method to time-periodic RFDEs is also demonstrated, and the convergence of the stability boundaries is verified numerically. The method is compared with recently developed highly efficient numerical methods: the pseudospectral collocation (also called Chebyshev spectral continuous-time approximation), the spectral Legendre tau method and the spectral element method. The comparison is based on the stability analysis of three linear autonomous RFDEs. The efficiency of the methods is measured by the convergence rate of stability boundaries in the space of system parameters, by the convergence rate of the rightmost characteristic exponent and by the computation time of the stability charts.