Numerical Integration of Systems of Stiff Ordinary Differential Equations with Special Structure

Algorithms for the solution of general systems of stiff differentialequations commonly use implicit integration formulae. The associatednon-linear equations at each step of the integration are efficientlysolved by an iteration such as the parallel chord method, wherethe matrix is an approximation for the Jacobian at a calculatedpoint of the solution. This iteration frequently gives sufficientlyrapid convergence over a number of integration steps beforeupdating and re-inversion of the matrix is required. When thedifferential equations have a special structure, satisfactoryconvergence may be maintained by updating a partition of theJacobian less frequently than the remainder and an efficientcomputational procedure consists in calculating the correspondingupdate of the inverse. Sufficient conditions for local convergencemay be expressed in terms of the difference between the iterationmatrix and the derivative at the solution or in terms of thedifference of the corresponding inverses. Similarly the asymptoticrate of convergence is estimated in terms of the norms of theseperturbations. To assess the effectiveness of updating a partitionof the Jacobian or its inverse we set the corresponding perturbationto zero and evaluate the estimate of the rate of convergence.Variable transformation and "weighting" of equations may beused to obtain more accurate computable estimates of convergencerates and tighter conditions for convergence. This is particularlyrelevant in non-linear stiff systems arising in applicationsfrom physics, chemistry and engineering and associated withfast and slow motions. Such systems exhibit special structurein the Jacobian and higher derivatives related to the sensitivityof the system to the components of the fast motion which makesthem particularly amenable to matrix updating techniques. Anumber of illustrative problems from the literature are cited.A worked example has been solved numerically using an inversewhose partitions are updated irregularly as required for convergenceand a comparison is made of iteration counts and inversion statisticswith updating of the full inverse. Computational savings maysometimes belarge.