The well-known central finite difference approximation was combined with the trapezoid quadrature method in this study to provide a numerical solution of the linear system of Volterra integro-fractional differential equations (LSVI-FDEs) of arbitrary orders, where the fractional derivative is described in the Caputo sense and the orders are between zero and one. The method works by first using the central finite difference approximation to approximate the Caputo derivative at any fixed point and then using the trapezoidal rule to obtain a finite difference expression for our fractional equation, while recalling the linear spline approximation for the first steps. This new, more efficient method involves converting sets of equations and conditions into matrix relationships, from which symmetry matrices can be created in some cases. We also present a new approach for error analysis of the discrete numerical scheme and the explicit numerical technique for LSVI-FDEs. The multi-level explicit finite difference approximation’s stability and convergence were explored, and a MatLab application was created to explain the results. Finally, several numerical examples are offered to demonstrate the technique's application.