Computing with functions in spherical and polar geometries I. The sphere

TL;DR

This paper introduces efficient, stable algorithms for numerical computation with smooth functions on the sphere, enabling high-precision operations and fast solutions to Poisson's equation with large degrees of freedom.

Contribution

It presents a structure-preserving iterative Gaussian elimination combined with the double Fourier sphere method for stable, efficient computations on the sphere, including an optimal complexity Poisson solver.

Findings

Stable differentiation and integration on the sphere

Efficient solution of Poisson's equation with 100 million degrees of freedom in one minute

Reduction of oversampling near the poles

Abstract

A collection of algorithms is described for numerically computing with smooth functions defined on the unit sphere. Functions are approximated to essentially machine precision by using a structure-preserving iterative variant of Gaussian elimination together with the double Fourier sphere method. We show that this procedure allows for stable differentiation, reduces the oversampling of functions near the poles, and converges for certain analytic functions. Operations such as function evaluation, differentiation, and integration are particularly efficient and can be computed by essentially one-dimensional algorithms. A highlight is an optimal complexity direct solver for Poisson's equation on the sphere using a spectral method. Without parallelization, we solve Poisson's equation with $100$ million degrees of freedom in one minute on a standard laptop. Numerical results are presented throughout. In a companion paper (part II) we extend the ideas presented here to computing with functions on the disk.