ARPIST: Provably Accurate and Stable Numerical Integration over Spherical Triangles

Numerical integration on spheres, including the computation of the areas of spherical triangles, is a core computation in geomathematics. The commonly used techniques sometimes suffer from instabilities and significant loss of accuracy. We describe a new algorithm, called ARPIST, for accurate and stable integration of functions on spherical triangles. ARPIST is based on an easy-to-implement transformation to the spherical triangle from its corresponding linear triangle via radial projection to achieve high accuracy and efficiency. More importantly, ARPIST overcomes potential instabilities in computing the Jacobian of the transformation, even for poorly shaped triangles that may occur at poles in regular longitude-latitude meshes, by avoiding potential catastrophic rounding errors. We compare our proposed technique with L'Huilier's Theorem for computing the area of spherical triangles, and also compare it with the recently developed LSQST method (J. Beckmann, H.N. Mhaskar, and J. Prestin, GEM - Int. J. Geomath., 5:143-162, 2014) and a radial-basis-function-based technique (J. A. Reeger and B. Fornberg, Stud. Appl. Math., 137:174-188, 2015) for integration of smooth functions on spherical triangulations. Our results show that ARPIST enables superior accuracy and stability over previous methods while being orders of magnitude faster and significantly easier to implement.