General orthogonal polynomials are introduced in order to analyse and approximate the solution of a class of delay systems. Using the operational matrix of integration, together with the operational matrix of linear transformation, the dynamical equation of a delay system is reduced to a set of simultaneous linear algebraic equations. The coefficient vectors of the general orthogonal polynomials can be determined recursively by use of the derived algorithm. Examples are given to demonstrate the validity and applicability of the orthogonal-polynomial approximations.