This subroutine computes a forward-difference approximation to the n by n Jacobian matrix associated with a specified problem of n functions in n variables. If the Jacobian has a banded form, then function evaluations are saved by only approximating the nonzero terms.