Trapezoidal rule or similar for integration over sphere (spherical triangle).

Question

I would like to calculate numerically the integral of the function defined on the sphere. Moreover, the sphere is completely covered by non-overlapping spherical triangles, I need the integral to be calculated as sum of integrals over triangles. For simplicity assume, that sphere has radius $1$.

So, if $S$ is a sphere, $T$ is the set of non-overlapping spherical triangles covering this sphere and $t\in T$ is triangle from this set then:

$\int\limits_{\mathbf{x}\in S}f(\mathbf{x})ds=\sum\limits_{t\in T}\int\limits_{\mathbf{x}\in t}f(\mathbf{x})ds$

What I need is the means to numerically calculate integral over spherical triangle $t$.

The simplest way is to use analogues of 'rectangular rule' here, that is

$\int\limits_{\mathbf{x}\in t}f(\mathbf{x})ds\approx f(\mathbf{x}_c)A(t)$

where $\mathbf{x}_c$ is the center of the spherical triangle (whatever point you think is the center of spherical triangle, as long as it is inside this triangle) and $A(t)$ is the spherical triangle's area.

But is it possible to use some more accurate formula which would take, for example, values of $f$ at 3 vertices of the triangle?

Assume, that triangles are of decent shape.

Answer 1

Here is a link to a paper which I believe might interest you

http://homepage.math.uiowa.edu/~atkinson/papers/SphereQuad1982.pdf

It is not certainly not new, but it is free and there is good basic information on centroid rule as well as an isoparametric method which has a higher order of convergence. There are some error bounds as well. I expect that they are pessimistic.

Regardless of how you proceed you will also need a reliable way of estimating the error. Here I recommend that you investigate the existence of an asymptotic error expansion of the form \begin{equation} I - A_h = \alpha h^p + \beta h^q + O(h^r). \end{equation} Here $I$ is your integral, $A_h$ is your approximation of the integral, $h$ is a measure of the size of your spherical triangles, say, the largest arc length used, $\alpha$, $\beta$ are constants independent of $h$, and $0<p<q<r$ are real numbers. Proving the existence of such an expansion is typically hard, but it is easy to verify numerically that it exists, in which case you can use it compute accurate error estimates. If the expansion exists, then you can show that \begin{equation} F_h = \frac{A_{2h} - A_{4h}}{A_h - A_{2h}} \rightarrow 2^p, \quad h \rightarrow 0 \end{equation} and the convergence will be monotone. Here you should view $A_{4h}$ as your initial estimate based on a coarse triangulation, $A_{2h}$ and $A_h$ are improved estimates obtained by recursively cutting each triangle into 4 parts by halving each arc. Observing how $F_h$ behaves as $h$ tends to zero, will allow you to identify $p$, which is not necessarily an integer. Moreover, you will be able to show that the dominance error term $\alpha h^p$ satisfies \begin{equation} \alpha h^p = \frac{A_h - A_{2h}}{2^p - 1} + O(h^q). \end{equation} This will allow you to estimate the error. The error estimates are good, when the fractions are close to $2^p$. There are some numerical issues involving catastrophic cancellation when $h$ becomes small, but this is a minor matter. As long as the computed values of the above mentioned fractions are converging monotonically, then you are in good shape.

The entire technique is known as Richardson extrapolation.

If you functions are nice and smooth, then I expect that you will discover that $p$ is an integer.