A matrix series method for the integration of rate equations in a reaction network. An alternative to Runge-Kutta methods

Abstract

A matrix method is presented for obtaining numerical solutions to the differential equations pertaining to reaction networks. Neither matrix inversions nor diagonalizations are required. The kinetic equations are east in matrix and vector form and solved as an exponential matrix, which is approximated by a truncated series. The series matrix. S( t) is in terms of the rate constant matrix K that characterizes the reaction network. Accuracy for a given step size is determined by the order of the approximating series. Optimum coefficients for the series are given through the fifth order. Concentrations of all species and sensitivity coefficients are calculated for each time interval. The procedure for first order networks is particularly simple when all times are equally spaced, since S remains fixed throughout the problem and only the concentration vector changes. A second order process is approximated over a small time interval by a first order expression, but both S and the concentration vector change at each step. The procedure is approx. 2–3 times more rapid than conventional fourth order Runge-Kutta integration to the same level of accuracy.