An adaptive algorithm for computing an approximation to the integral of each element in a vector function f(x,y) over a two-dimensional region made up of triangles is presented. A FORTRAN implementation of the algorithm is included. The basic cubature rule used over each triangle is a 37-point symmetric rule of degree 13. Based on the same evaluation points the local error for each element in the approximation vector and for each triangle is computed using a sequence of null rule evaluations. A sophisticated error-estimation procedure tries, in a cautious manner, to decide whether we have asymptotic behavior locally for each function. Different actions are taken depending on that decision, and the procedure takes advantage of the basic rule’s polynomial degree when computing the error estimate in the asymptotic case. (Source: http://dl.acm.org/)

This software is also peer reviewed by journal TOMS.