Algorithm 407 DIFSUB for Solution of Ordinary Differential Equations

This subroutine integrates a set of up to N ordinary differential equations one step of length H, where H may be specified by the user, but is controlled by the subroutine to control the estimated error within a specified tolerance, if possible. A maltistep predictor corrector method is used whose order is automatically chosen by the subroutine as the integration proceeds. Either an Adams' method or methods suitable for stiff equations can be selected. The starting procedure is automatic and the information retained by the program about previous steps is stored in such a way as to make the interpolation to a nonmesh point straightforward. (See the description of the parameter Y in the subroutine.) The methods used are described from a mathematical point of view in the papers referenced in [1 ]. The programs may call on up to three subroutines. They are