# Super-Resolution on the Two-Dimensional Unit Sphere

Frank Filbir, Kristof Schröder, Anna Veselovska

(*  This is an extended version of our submitted manuscript for publication
titled as ”Recovery of Atomic Measures on the Unit Sphere”. The present
manuscript will not be submitted for publication, but it is rather meant to
provide details missing in the submitted version.  )

###### Abstract

We study the problem of recovering an atomic measure on the unit 2-sphere
$\mathbb{S}^{2}$ given finitely many moments with respect to spherical
harmonics. The analysis relies on the formulation of this problem as an
optimization problem on the space of bounded Borel measures on
$\mathbb{S}^{2}$ as it was considered by Y. de Castro & F. Gamboa and E.
Candés & C. Fernandez-Granda in 2013. We construct a dual certificate using a
kernel given in an explicit form and make a concrete analysis of the
interpolation problem. Numerical examples are provided and analyzed.

###### Contents

  1. 1 Introduction
  2. 2 Unit Sphere and Super-Resolution
    1. 2.1 Analysis on the Sphere
    2. 2.2 Spherical harmonics
    3. 2.3 Super-Resolution on the Unit Sphere
  3. 3 Localized kernels
  4. 4 Dual Certificate on the Sphere
    1. 4.1 Solution of the Interpolation Problem
    2. 4.2 Bounds for the Interpolant
  5. 5 Numerical Solution
    1. 5.1 Semidefinite Formulation of the Optimization Problem
    2. 5.2 Discretization of the Optimization Problem
    3. 5.3 Numerical Experiments
      1. 5.3.1 Semi-definite Program
      2. 5.3.2 Discretization of Semi-Definite Program
  6. A Proof of Lemma 3.0.2
  7. B Proof of Lemma 3.0.3
  8. C  Properties of the Cross Product

##  1 Introduction

Signals on a spherical manifold arise in various applications, staring from
medical imaging [9, 11, 26], computer graphics [14, 18] and sound recording
[10, 13] to astrophysics [22] and topography [20]. As mentioned in [28], in
many of these settings inverse problems arise, where an unknown signal need to
be recovered from linear measurements acquired through a convolution process.
In this work, we consider the case, when a spatially highly resolved signal is
modelled as a weighted sum of Dirac measures $\mu^{\star}$ on the two-
dimensional Euclidean $\mathbb{S}^{2}$, and the information one can access is
only the convolved version of the signal $\mu^{\star}\ast D_{N}$ with the
Dirichlet kernel $D_{N}$ on the sphere, for possibly low $N$. The problem of
recovery of a spatially highly resolved signal from its coarse scale
information is called the super-resolution problem or, in other words, the
de-convolution problem.

In general, this sort of problem has been treated in different geometric
settings and with respect to different systems of functions. Exact measure
recovery in the classical Fourier setting, i.e. when the underlying space is
the torus $\mathbb{T}^{d}$, has a very long history starting with the initial
work by G. R. de Prony in 1795 [1] in the univariate case, and then moving on
to different one-dimensional [7, 6, 21] and multi-dimensional Prony-based
techniques [31, 37, 41, 47, 48] that have stabilized and generalized the
method in various directions. Recently, the Prony’s method has also been
extend to the $d$-dimensional sphere $\mathbb{S}^{d}$.

Less than ten years ago, considering the one-dimensional de-convolution
problem in the light of convex optimization, several authors [25, 32] have
proposed a variational recovery approach that is to minimize the total
variation over the set of all finite complex measures supported on
$\mathbb{T}$, given the convolved version of the measure. It has been shown
that if the support atoms of the measure are well separated then the target
measure $\mu^{\star}$ is the unique solution of the minimization problem,
therewith a sufficient criteria for $\mu^{\star}$ being the unique solution is
the existence of so-called dual certificate, that is, a polynomial of degree
$N$ whose sup-norm is reached at the points of the measure support.

These two fundamental papers have ushered in new ways of treating the super-
resolution problem. Beside the fact that this variational recovery method does
not need the number of unknowns points beforehand, it also has proved to enjoy
stability in the case when the low frequency information of corrupted by noise
[27, 35, 36]. Hereupon, the generalization to higher dimensions on the torus
has been considered in [39, 44, 43]. Another big advantage of considering the
super-resolution as the minimization problem is its adaptivity to different
geometric settings, namely, semi-algebraic domains in higher dimensions [38],
or compact smooth Riemannian manifolds such es the rotation group $SO(3)$ [40]
and two-dimensional sphere [33, 34]. Lately, this procedure has been
generalized to short-time Fourier measurements [42].

In this work, we consider de-convolution problem on the two-dimensional
Euclidean $\mathbb{S}^{2}$ as a total variation minimization problem. Although
to prove uniqueness of an optimal solution, we follow a general idea from [32,
40]. the actual construction of a dual certificate requires localization
estimates for interpolation kernels and its derivatives on $\mathbb{S}^{2}$
with explicit constants, and heavily depends on special behaviour of the
geodesic distance with respect to the boundedness of the derivatives and on
the fact that there is no nowhere not vanishing vector field on the sphere,
due to the Hairy ball theorem. From numerical point of view, we consider two
approaches. First, we use the dual formulation of the minimization problem
solving it via a single semi-definite program (SDP), going along the same line
as [40, 34]. In the second approach we discretize the primal problem
beforehand, then solve the corresponding finite-dimensional optimization
problem. To analyze the convergence of the discretization process we build on
results stated in [29]. We also would like to mention the for the nonnegative
total variation minimization problem an alternative way of construction of a
dual certificate that involves some algebraic techniques has been proposed in
[48], and in pure compressed sensing setting a recovery of sparse signal on
the two-dimensional sphere has been considered in [23]

The outline of this paper is as follows. In Section 2, briefly the necessary
analytical tools on the sphere including spherical harmonics are introduced,
and the problem of super-resolution is stated. In Section 3 provides the
localization results for the chosen interpolating kernel that are essential
for construction of a dual certificate. The actual construction of a dual
certificate as a solution of the Hermite-type interpolation problem is the
content of Section 4. Finally, we finish by presenting the numerical solution
and the discretization of the problem in Section 5.

##  2 Unit Sphere and Super-Resolution

In this section we briefly summarize analytical tools on the two dimensional
sphere $\mathbb{S}^{2}$ and state the super-resolution problem.

###  2.1 Analysis on the Sphere

The unit sphere is an embedded sub-manifold of the three-dimensional Euclidean
space $\mathbb{R}^{3}$ given by

$\mathbb{S}^{2}=\\{x\in\mathbb{R}^{3}\colon\|x\|_{2}=1\\},$

where the $l_{2}$-norm is $\|x\|_{2}=\sqrt{x^{\mathrm{T}}x}.$ Such embedding
provides a very simple definition of the tangent space $T_{x}\mathbb{S}^{2}$
at a point $x\in\mathbb{S}^{2}$, that namely it is given as the orthogonal
complement of the linear subspace $\mathrm{span}\\{x\\}$, i.e.

$T_{x}\mathbb{S}^{2}=\\{y\in\mathbb{R}^{3}\colon\langle{x,y}\rangle=0\\},$

where $\langle{x,y}\rangle=x^{\mathrm{T}}y$ is the standard inner product.
Induced by the Riemannian metric of the ambient space $\mathbb{R}^{3}$, the
Riemannian metric on the sphere $g_{\mathbb{S}^{2}}:T_{x}\mathbb{S}^{2}\times
T_{x}\mathbb{S}^{2}\to\mathbb{R}$ is given for all $x\in\mathbb{S}^{2}$ by

$g_{\mathbb{S}^{2}}(v,w):=\langle{v,w}\rangle,\quad v,w\in
T_{x}\mathbb{S}^{2}.$

In this metric the geodesic distance between two points $x,y\in\mathbb{S}^{2}$
is given by the great-circle distance

$d(x,y)=\mathrm{arccos}(\langle{x,y}\rangle),\quad x,y\in\mathbb{S}^{2}.$

Let us also shortly describe the differential structure on the sphere. We will
use two different local coordinates on $\mathbb{S}^{2}$. First, for each
starting point $x\in\mathbb{S}^{2}$ and a direction $v\in
T_{x}\mathbb{S}^{2}$, there exists the unique geodesic $\gamma_{x,v}$ such
that $\gamma_{x,v}(0)=x$ and $\gamma_{x,v}^{\prime}(0)=v$, and the equation of
such geodesic reads as

$\gamma_{x,v}(t)=\cos{(\|v\|_{2}t)}x+\sin{(\|v\|_{2}t)}\frac{v}{\|v\|_{2}}.$

Then at a point $x\in\mathbb{S}^{2}$ the exponential map ${\exp_{x}{v}:\
T_{x}\mathbb{S}^{2}\to\mathbb{S}^{2}}$ is given by

$\exp_{x}{v}=\gamma_{x,v}(1)$ (2.1)

Now, let us fix an orthonormal basis $\eta_{1}^{x},\eta_{1}^{x}\in
T_{x}\mathbb{S}^{2}$ such that $\eta_{2}^{x}=\eta_{1}^{x}\times x$ and
$\eta_{1}^{x}=x\times\eta_{2}^{x}$. It is always possible to do so, although
we can not choose a local bases in a continuous way, as there is no continuous
nowhere not vanishing vector field on the sphere, due to the Hairy ball
theorem.  This is a special property of $\mathbb{S}^{2}$ is in contrast to
some other manifolds, e.g. the rotation group $SO(3)$, where the tangent space
is basically a translation of the tangent space at the identity. It was a
classical problem to determine which of the spheres
$\mathbb{S}^{n}=\\{x\in\mathbb{R}^{n}\colon\|x\|_{2}=1\\}$ are parallelizable,
then it has been shown that along with $\mathbb{S}^{0}$ and the unit circle
$\mathbb{S}^{1}$ parallelizable are only $\mathbb{S}^{3}$ and $\mathbb{S}^{7}$
[4].

One way to obtain a local bases $\eta_{1}^{x},\eta_{1}^{x}\in
T_{x}\mathbb{S}^{2}$ is to choose a point $z\in\mathbb{S}^{2}$ and an
orthonormal basis $\eta_{1}^{z},\eta_{1}^{z}\in T_{z}\mathbb{S}^{2}$ and then
set

$\eta_{i}^{x}=\left\\{\begin{matrix}\mathrm{e}^{d(x,z)\cdot[\frac{z\times
x}{\sin{(d(x,z))}}]\eta_{i}^{z}},&x\neq-z,\\\
-\eta_{i}^{z},&x=-z,\end{matrix}\right.$ (2.2)

where for a vector $v\in\mathbb{R}^{3}$

$[v]=\begin{pmatrix}0&-v_{3}&v_{2}\\\ v_{3}&0&-v_{1}\\\
-v_{2}&v_{1}&0\end{pmatrix}$

is the corresponding skew-symmetric matrix in the algebra $\mathfrak{so}(3)$.
In other words, we rotate the local coordinate system at $z$, which is
continuous for all points but the antipodal point $-z$.

Combination of the coordinates of a vector $v\in T_{x}\mathbb{S}^{2}$ in the
basis $\eta_{1}^{x},\eta_{1}^{x}\in T_{x}\mathbb{S}^{2}$ with the exponential
map (2.1) yields the normal coordinates centered at $x\in\mathbb{S}^{2}$, i.e.
we parametrize a neighborhood of $x$ by

$\varphi(v_{1},v_{2})=(\cos{\|v\|_{2}t})x+\sin{(\|v\|_{2}t)}\frac{v_{1}\eta_{1}^{x}+v_{2}\eta_{2}^{x}}{\|v\|_{2}}.$
(2.3)

and the inverse parametrization for $y\in\mathbb{S}^{2}$ in a neighborhood of
$x\in\mathbb{S}^{2}$ is given by

$v_{i}(y)=\frac{d(y,x)}{\sin{(d(x,y))}}\langle{y,\eta_{i}^{x}}\rangle,\quad
i=1,2.$ (2.4)

Moreover, the vectors

$\frac{\partial}{\partial
v_{1}}\varphi(v_{1}(y),v_{2}(y)),\quad\frac{\partial}{\partial
v_{2}}\varphi(v_{1}(y),v_{2}(y))$ (2.5)

form a basis of $T_{y}\mathbb{S}^{2}$. One can show that in the center of the
normal coordinates, the derivatives of the basis vectors have the following
properties

$\frac{\partial^{2}}{\partial^{2}v_{i}}\varphi(v_{1}(x),v_{2}(x))=-x,\quad\frac{\partial^{2}}{\partial
v_{j}\partial v_{j}}\varphi(v_{1}(x),v_{2}(x))=0.$ (2.6)

Since in a normal coordinate system centered at $x\in\mathbb{S}^{2}$, the
Christoffel symbols vanish at the point $x\in\mathbb{S}^{2}$, the gradient of
a differentiable function $f:\mathbb{S}^{2}\to\mathbb{C}$ at the point
$x\in\mathbb{S}^{2}$ has the representation

$\nabla f(x)=\begin{pmatrix}X_{1}f(x)\\\ X_{2}f(x)\end{pmatrix},$ (2.7)

where the differential operators $X_{i}$, $i=1,2$, are defined by

$X_{i}f(x)=\frac{\partial}{\partial
v_{1}}(f\circ\varphi)(v_{1}(x),v_{2}(x))=\lim_{t\to
0}t^{-1}(f(\gamma_{x,\eta_{i}^{x}}(t))-f(x)).$ (2.8)

The Hessian matrix of a twice differentiable function
${f:\mathbb{S}^{2}\to\mathbb{C}}$ in the center $x\in\mathbb{S}^{2}$ of the
normal coordinates, has also a special representation, namely

$Hf(x)=\begin{pmatrix}X_{1}X_{1}f(x)&X_{1}X_{2}f(x)\\\
X_{2}X_{1}f(x)&X_{2}X_{2}f(x)\end{pmatrix}.$ (2.9)

However, the representations (2.7) and (2.9) are only true in the center of
the normal coordinates, since as mentioned before the Christoffel symbols
vanish. For different points, we would have to compute the Christoffel symbols
with respect to the normal coordinates, which becomes quite complicated.
Alternatively, we introduce a second set of coordinates, such that the
computation of the Christoffel symbols is much more simpler. For a point
$z\in\mathbb{S}^{2}$, we parametrize the set $B_{\pi}(0)\setminus\\{0\\}$ by

$v(r,\theta)=r(\cos(\theta)\eta_{1}^{z}+\sin(\theta)\eta_{2}^{z})$ (2.10)

for $(r,\theta)\in(0,\pi)\times[0,2\pi)$. Combining (2.10) with the
exponential map (2.1), i.e.

$\varphi^{\mathrm{pol}}(r,\theta)=\exp_{z}(v(r,\theta)),$ (2.11)

yields the polar coordinates centered at $z\in\mathbb{S}^{2}$, which
parametrize ${\mathbb{S}^{2}\setminus\\{z,-z\\}}$. To give an example, for
$z=(0,0,1)^{\mathrm{T}}$, these are the usual spherical coordinates on the
sphere, given by

$\varphi^{\mathrm{pol}}(r,\theta)=\begin{pmatrix}\sin{r}\cos{\theta}\\\
\sin{r}\sin{\theta}\\\ \cos{r}\end{pmatrix}.$

Following the line, for a vector $x\in\mathbb{S}^{2}\setminus\\{z,-z\\}$, the
inverse parametrization is given by

$\begin{split}r(x)&=\mathrm{arccos}(\langle{x,z}\rangle)=d(x,z),\\\
\theta(x)&=\mathrm{arctan}_{2}(\langle{x,\eta_{2}^{z}}\rangle,\langle{x,\eta_{1}^{z}}\rangle),\end{split}$

where $\mathrm{arctan}_{2}(x,y)$ denotes the $\mathrm{arctan}$ of
$\frac{y}{x}$ with respect to the different branches of the tangent function,
which means

$\cos{(\theta(x))}=\frac{\langle{x,\eta_{1}^{z}}\rangle}{\sin{(d(x,z))}},\quad\quad\sin{(\theta(x))}=\frac{\langle{x,\eta_{2}^{z}}\rangle}{\sin{(d(x,z))}}.$

For each vector $x\in\mathbb{S}^{2}\setminus\\{z,-z\\}$, the vectors

$\gamma_{1}^{x}=\frac{\partial\varphi^{\mathrm{pol}}}{\partial
r}(r(x),\theta(x)),\quad\gamma_{2}^{x}=\frac{1}{\sin{(r(x))}}\frac{\partial\varphi^{\mathrm{pol}}}{\partial\theta}(r(x),\theta(x)),$
(2.12)

form an orthonormal basis of the tangent space $T_{x}\mathbb{S}^{2}$. Notably,
we have $\gamma_{2}^{x}=\gamma_{1}^{x}\times x$ and
$\gamma_{1}^{x}=x\times\gamma_{2}^{x}$. We remark, due to the singularities at
the poles $z,-z$ there is no basis in spherical coordinates of the
corresponding tangent spaces.

In the polar coordinates (2.11), the Riemannian metric takes the form

$g(r,\theta)=\begin{pmatrix}1&0\\\ 0&\sin^{2}{(r)}\end{pmatrix},$ (2.13)

and the Christoffel symbols in these coordinates are therefore given by

$\Gamma^{r}(r,\theta)=\begin{pmatrix}0&0\\\
0&-\sin{r}\cos{r}\end{pmatrix},\quad\Gamma^{\theta}(r,\theta)=\begin{pmatrix}0&\cot{r}\\\
\cot{r}&0\end{pmatrix}.$ (2.14)

For a twice differentiable function ${f:\mathbb{S}^{2}\to\mathbb{C}}$ and a
vector ${x\in\mathbb{S}^{2}\setminus\\{z,-z\\}}$, the Hessian matrix with
respect to the polar coordinates centered at $z$, i.e. with respect to the
basis (2.12), is represented by

$Hf=\begin{pmatrix}\frac{\partial^{2}f\circ\varphi^{\mathrm{pol}}}{\partial
r^{2}}&\frac{1}{\sin{r}}\frac{\partial^{2}f\circ\varphi^{\mathrm{pol}}}{\partial
r\partial\theta}\\\
\frac{1}{\sin{r}}\frac{\partial^{2}f\circ\varphi^{\mathrm{pol}}}{\partial\theta\partial
r}&\frac{1}{\sin^{2}{r}}\frac{\partial^{2}f\circ\varphi^{\mathrm{pol}}}{\partial\theta^{2}}\end{pmatrix}-\frac{1}{\sin^{2}{r}}\frac{\partial
f\circ\varphi^{\mathrm{pol}}}{\partial
r}\Gamma^{r}-\frac{1}{\sin{r}}\frac{\partial
f\circ\varphi^{\mathrm{pol}}}{\partial\theta}\Gamma^{\theta}.$ (2.15)

###  2.2 Spherical harmonics

Now, let us describe the involved basis functions, known as spherical
harmonics. For a detailed overview see [24].

Let consider the space $L^{2}(\mathbb{S}^{2})$ of all functions
${f:\mathbb{S}^{2}\to\mathbb{C}}$ such that

$\|f\|_{2}=\left(\int_{\mathbb{S}^{2}}|f(x)|^{2}\mathrm{d}\Omega(x)\right)^{\frac{1}{2}}=\left(\int_{0}^{2\pi}\int_{0}^{\pi}|f(r,\theta)|^{2}sin{(r)}\,\mathrm{d}r\mathrm{d}\theta\right)^{\frac{1}{2}}<\infty,$

where $\Omega$ is the Riemannian volume form on $\mathbb{S}^{2}$ for the
metric $g$, and the second part of the equality represents $\|f\|_{2}$ in the
spherical coordinates. It is well known, that the space
$L^{2}(\mathbb{S}^{2})$ can be decomposed into an orthogonal sum

$L^{2}(\mathbb{S}^{2})=\mathrm{cl}_{\|\cdot\|_{L_{2}}}\bigoplus\limits_{l=0}^{\infty}H_{l}$

where $H_{l}$ is the eigenspace to the eigenvalue $\lambda_{l}=-l(l+1)$ of the
Laplace-Beltrami operator on $\mathbb{S}^{2}$ with $\mathrm{dim}(H_{l})=2l+1$.
Looking for an orthonormal system of eigenfunctions of the Laplace-Beltrami
operator on the unit sphere leads to the spherical harmonics. Namely, the
spherical harmonics $Y_{m}^{l}:\mathbb{S}^{2}\to\mathbb{C}$ of degree
$m\in\mathbb{N}$ and order $l$ are functions

$Y_{m}^{l}(x(r,\theta))=N_{\ell
m}P_{\ell}^{m}(\cos{(r)})\mathrm{e}^{im\theta},$

where $-l\leq m\leq l$, the normalization constant $N_{\ell m}$ is given by

$N_{\ell m}=\frac{1}{\sqrt{2\pi}}\sqrt{\frac{2l+1}{2}\frac{(l-m)!}{(l+m)!}},$

$P_{l}^{m}(t)$ are associated Legendre polynomials, $r\in[0,\pi]$,
$\theta\in[0,2\pi)$ are the inclination and azimuth respectively. With this,
we have that the system

$\\{Y_{m}^{l}\colon l\in\mathbb{N},\;-l\leq m\leq l\\}$

constitutes an orthonormal basis of $L^{2}(\mathbb{S}^{2})$. Moreover, the
following addition theorem is valid

$P_{l}(\langle{x,y}\rangle):=P_{l}^{0}(\langle{x,y}\rangle)=\frac{4\pi}{2l+1}\sum\limits_{m=-l}^{l}Y_{m}^{l}(x)\overline{Y_{m}^{l}(y)}.$
(2.16)

The space of all finite linear combinations of spherical harmonics with degree
less or equal to $N$ will be denoted as

$\Pi_{N}(\mathbb{S}^{2}):=\mathrm{span}\\{Y_{m}^{l},-l\leq m\leq l,l\leq N\\}$

and will be called generalized polynomials of degree $N$. The projection
operator onto the set of generalized polynomials $\Pi_{N}(\mathbb{S}^{2})$, is
given by

$\mathcal{P}_{N}\colon L^{2}(\mathbb{S}^{2})\to
C(\mathbb{S}^{2})\quad\text{with}\quad\mathcal{P}_{N}f(x)=\int_{\mathbb{S}^{2}}f(y)D_{N}(x,y)\mathrm{d}\Omega(y),$
(2.17)

where the Dirichlet kernel on the sphere reads as

$D_{N}(x,y)=\sum\limits_{l=0}^{N}\frac{2l+1}{4\pi}P_{l}(\langle{x,y}\rangle).$

###  2.3 Super-Resolution on the Unit Sphere

In the following, we will introduce the problem of super-resolution, or in
another words, so called deconvolution problem, on the unit sphere sphere
$\mathbb{S}^{2}$, i.e. exact recovery of Dirac measure from its moments with
respect to the spherical harmonics up to a degree $N$.

To this aim, let us consider a weighted superposition of spikes

$\mu^{\star}=\sum\limits_{i=1}^{M}c_{i}\delta_{x_{i}}$

where $M\in\mathbb{N}$, $\delta_{x_{i}}$ is a Dirac measure centered at
pairwise distinct $x_{i}\in\mathbb{S}^{2}$, $c_{i}\in\mathbb{R}$ are real
valued amplitudes. We assume that all parameters $M,c_{i},x_{i}$ are unknown
and we can only access

$\mathcal{P}_{N}^{*}\mu^{\star}(x)=\int_{\mathbb{S}^{2}}D_{N}(x,y)\mathrm{d}\mu^{\star}(y)$
(2.18)

for possibly low degree $N$. The super-resolution problem is to recover the
unknown locations
$\mathcal{X}=\mathrm{supp}{(\mu^{\star})}=\\{x_{i}\\}_{i=1}^{M}\subset\mathbb{S}^{2}$
and the coefficients $c_{i}$ from the low frequency information (2.18). Due to
consideration of the total variation norm as being the continuous analog of
the $\ell_{1}$ norm, in the light of convex optimization, see [25, 27, 32, 33,
40], the super-resolution problem can be formulated as the following
minimization problem

$\min\limits_{\mu\in\mathcal{M}(\mathbb{S}^{2},\mathbb{R})}\|\mu\|_{\mathrm{TV}},\quad\mbox{
subject to }\quad\mathcal{P}_{N}^{*}\mu=\mathcal{P}_{N}^{*}\mu^{\star},$ (RP)

where the minimization is carried out over the set of all finite measures
$\mu$ supported on $\mathbb{S}^{2}$, and the total variation for a signed
Borel measure is defined by

$\|\mu\|_{\mathrm{TV}}=|\mu|(\mathbb{S}^{2})=\sup\sum_{j}|\mu(B_{j})|,$

and the supremum is taken over all partitions $B_{j}$ of $\mathbb{S}^{2}$.

The main ingredient to ensure the fact $\mu^{\star}$ is the unique minimizer
of the convex program (RP) is the existence of a dual interpolating polynomial
${q\in\Pi_{N}(\mathbb{S}^{2})}$, or so-called dual certificate, that for each
sign sequences ${u_{i}\in\\{-1,1\\}}$ satisfies

$\begin{split}q(x_{i})&=u_{i},\quad x_{i}\in\mathcal{X}\\\ |q(x)|&<1,\quad
x\in\mathbb{S}^{2}\setminus\mathcal{X}.\end{split}$ (2.19)

The connection between the uniqueness of optimal solution and the existence of
a dual certificate has been exploited in different settings, see e.g. [25, 30,
32, 33, 40]. In order to fulfill the interpolating conditions (2.19) let us
consider the following Hemite-type interpolation problem

$\begin{split}q(x_{i})&=u_{i},\\\ X_{1}q(x_{i})&=X_{2}q(x_{i})=0\end{split}$
(2.20)

for $x_{i}\in\mathcal{X}$, where $X_{k}$ are the differential operators
defined in (2.8). This means, we ask not only the interpolation, but also for
local extrema at all the interpolation points. In order to approach the
interpolation problem (2.20), we consider a kernel
$J_{N}\colon\mathbb{S}^{2}\times\mathbb{S}^{2}\to\mathbb{C}$, such that the
kernel $J_{N}(\cdot,y)$ itself, and its derivatives $X_{k}^{y}J_{N}(\cdot,y)$,
where the superscript indicates the action of the differential operators on
the second variable, are generalized polynomials of order $N$ in the first
variable, i.e. $J_{N}(\cdot,y)$,
$X_{k}^{y}J_{N}(\cdot,y)\in\Pi_{N}(\mathbb{S}^{2})$ for all vectors
$y\in\mathbb{S}^{2}$. The dual certificate we construct is of the form

$q(x)=\sum\limits_{i=0}^{M}\alpha_{0,i}J_{N}(x,x_{i})+\alpha_{1,i}X_{1}^{y}J_{N}(x,x_{i})+\alpha_{2,i}X_{2}^{y}J_{N}(x,x_{i}).$
(2.21)

Thus, is easily follows that due to the construction
$q\in\Pi_{N}(\mathbb{S}^{2})$. Applying the interpolation conditions (2.20)
leads to the linear system of equations

$K\alpha=\begin{pmatrix}J_{N}&X_{1}^{x}J_{N}&X_{2}^{x}J_{N}\\\
X_{1}^{y}J_{N}&X_{1}^{x}X_{1}^{y}J_{N}&X_{2}^{x}X_{1}^{y}J_{N}\\\
X_{2}^{y}J_{N}&X_{1}^{x}X_{2}^{y}J_{N}&X_{2}^{x}X_{2}^{y}J_{N}\end{pmatrix}\begin{pmatrix}\alpha_{0}\\\
\alpha_{1}\\\ \alpha_{2}\end{pmatrix}=\begin{pmatrix}u\\\ 0\\\
0\end{pmatrix}.$ (2.22)

where the block $J_{N}$ in the matrix $K$ correspond to matrix of the form
$J_{N}=(J_{N}(x_{i},x_{j}))_{i,j=1}^{M}$, and the same we have for the
derivatives. The entries in the vectors are given by
$\alpha_{k}=(\alpha_{k,j})_{j=1}^{M}$, $k=0,1,2$, and $u=(u_{j})_{j=1}^{M}$.
To find the vector of coefficients $\alpha$, we need to show that the matrix
$K$ is invertible. To this aim, we follow the general idea provided [32, 40].
Namely, due to the block structure of $K$ one can proof the invertibility of
$K$ using an iterative block inversion and the fact that a matrix $A$ is
invertible if $\|I-A\|_{\infty}<1,$ where
$\|A\|_{\infty}=\max_{i}\sum_{j}|a_{i,j}|$. In this case the norm of the
inverse matrix is bounded by

$\|A^{-1}\|_{\infty}<\frac{1}{1-\|I-A\|_{\infty}}.$

Thus, to show the invertibility of $K$, we need to make use of localization
properties  of the entries of $K$, i.e. we need to bound the expressions
$|J_{N}(x_{i},x_{j})|$ $|X_{k}^{y}J_{N}(x_{i},x_{j})|$ and
$|X_{n}^{x}X_{k}^{y}J_{N}(x_{i},x_{j})|$. Namely, we are looking for the
estimation for the kernel $J_{N}$ of the form

$|J_{N}(x_{i},x_{j})|<\frac{c}{((N+1)d(x_{i},x_{j}))^{s}},$

with some constants $s$ and $c$, and for equivalent bounds for the
derivatives. Using these estimates we can find explicit bounds on the supremum
norm of the coefficients $\alpha_{k}$. Once the coefficients are bounded, we
need to show that $|q(x)|<1$ for $x\in\mathbb{S}^{2}\setminus\mathcal{X}$. To
this aim, we need to show convexity property of $q$, that can be done by
estimating the entries of the Hessian matrix of $q$. Since by construction
(2.21) $q$ alredy includes first derivatives of $J_{N}$, it will put into
consideration third mixed derivatives of $J_{N}$ in the Hessian. Therefore, we
also need some estimates of the derivatives of third order. All this together
results in the topic of the next sections, where we choose a specific kernel
and show needed locality estimations.

##  3 Localized kernels

The aim of this section is to show localization results for the constituents
of dual certificate (2.21) for a particularly chosen interpolation kernel
$J_{N}$.

Let us star discussing the choice of the interpolation kernel. First, the
interpolation kernel need to have an expansion in terms of generalized
polynomials. Therefore we choose a kernel that owns the representation

$J_{N}(x,y)=\sum_{l=0}^{N}\widetilde{w}_{l}\sum\limits_{m=-l}^{l}Y_{m}^{l}(x)\overline{Y_{m}^{l}(y)},$

hence, due to the constriction we have $J_{N}(\cdot,y),$
$X_{n}^{y}J_{N}(\cdot,y)\in{\Pi}_{N}(\mathbb{S}^{2})$ for all
$y\in\mathbb{S}^{2}$. Moreover, using the addition theorem (2.16) for
spherical harmonics leads to

$J_{N}(x,y)=\widetilde{J}_{N}(d(x,y))=\sum_{l=0}^{N}\frac{2l+1}{4\pi}\widetilde{w}_{l}P_{l}(\langle{x,y}\rangle)$

where $\widetilde{J}_{N}$ is a trigonometric polynomial. This means that
$J_{N}$ is a zonal function, i.e. its value only depends on the distance
between $x$ and $y$. Such a property of $J_{N}$ lends us a hand in deriving
estimates of the interpolation kernel $J_{N}$, since on condition localization
estimates can derived from localization principles for the trigonometric
polynomial $\widetilde{J}_{N}$.

Relying on the discussed above, as an interpolation kernel we choose the
specific kernel, given by

$J_{N}(x,y)=\widetilde{J}_{N}(d(x,y))=\frac{1}{(\lfloor{N/2}\rfloor+1)^{4}}\frac{\sin^{4}\big{(}(\lfloor{N/2}\rfloor+1)d(x,y)/2\big{)}}{\sin^{4}(d(x,y)/2)},$
(3.1)

i.e. the classical Jackson kernel evaluated at the distance between
$x,y\in\mathbb{S}^{2}$.

######  Lemma 3.0.1.

The Jackson kernel $J_{N}(x,y)$ has an expansion of the form

$J_{N}(x,y)=\sum_{l=0}^{N}\frac{2l+1}{4\pi}\widetilde{w}_{l}P_{l}(\langle{x,y}\rangle)$
(3.2)

with positive Legendre coefficients
$\displaystyle{\widetilde{w}_{l}=2\pi\int_{-1}^{1}P_{l}(t)\widetilde{J}_{N}(\mathrm{arccos}(t))\,\mathrm{d}t,}$
and thus $J_{N}(\cdot,y),X_{n}^{y}J_{N}(\cdot,y)\in\Pi_{N}(\mathbb{S}^{2})$
for all $y\in\mathbb{S}^{2}$.

###### Proof.

For $t\in[-1,1]$, let us consider the function $\widetilde{F}_{n}$ defined as

$\widetilde{F}_{n}(\mathrm{arccos}(t))=\frac{1}{n+1}\frac{\sin^{2}((n+1)^{2}\mathrm{arccos}(t)/2)}{\sin^{2}(\mathrm{arccos}(t)/2)}$

i.e. Fejér kernel evaluated at $\mathrm{arccos}{(t)}$. Then for the
trigonometric polynomial $\widetilde{J}_{N}$ it holds

$\widetilde{J}_{N}(\mathrm{arccos}{(t)})=\widetilde{F}_{n}^{2}(\mathrm{arccos}(t)),\quad
n=\lfloor{N/2}\rfloor$

As it was shown in [16], the Fejér kernel can be represented in terms of
Legendre polynomials as

$\widetilde{F}_{n}(\mathrm{arccos}(t))=\sum_{l=0}^{M}\frac{2l+1}{4\pi}\widetilde{v}_{l}P_{l}(t),$
(3.3)

which shows that
$\widetilde{F}_{n}(d(x,\cdot)),\widetilde{F}_{n}(d(\cdot,y))\in\Pi_{n}(\mathbb{S}^{2})$,
and therefore $\widetilde{J}_{N}(d(x,\cdot))$,
$\widetilde{J}_{N}(d(\cdot,y))\in\Pi_{N}(\mathbb{S}^{2})$. The positivity of
$\widetilde{w}_{l}$ follows from the positivity of the linearization
coefficients of a product of two Legendre polynomials

$P_{n}(t)P_{m}(t)=\sum_{l=0}^{\min(m,n)}\frac{2m+2n-4l+1}{2m+2n-2l+1}\frac{A(m-l)A(l)A(n-l)}{A(n+m-l)}P_{m+n-2}(t),$

where $A(m)=\frac{1\cdot 3\cdot 5\cdot\,\cdots\,\cdot(2m-1)}{m!},$ for details
see [2]. ∎

Before we proceed with stating the necessary localization estimates, we would
like to discuss the behaviour of the Jackson kernel. The geodesic distance of
the sphere behaves in a special way with respect to the boundedness of the
derivatives. Namely the derivatives of the geodesic distance $d(x,y)$ have
true poles at $x=y$ and $x=-y$. Whereas the first case can be handled using
the point-wise estimates, the second case can not be covered in the same way,
since the sign of $\cos(kd(x,y))$ and $\sin(kd(x,y))$ alternates with $k$ in a
neighborhood of $x=-y$. Consequently, we need to deal with the singularity at
$x=-y$ induced by the derivatives of the geodesic distance in different way.
Especially, we need to bound the following trigonometric expressions

$\begin{split}G_{1}(\omega)&=\frac{\widetilde{J}^{\prime}_{N}(\omega)}{\sin\omega},\quad\quad\quad\quad\quad
G_{2}(\omega)=\widetilde{J}^{\prime\prime}_{N}(\omega)\frac{\widetilde{J}^{\prime}_{N}(\omega)\cos\omega}{\sin\omega},\\\
G_{3}(\omega)&=\frac{\widetilde{J}^{\prime\prime}_{N}(\omega)\sin\omega-\widetilde{J}^{\prime}_{N}(\omega)\cos\omega}{\sin^{2}\omega},\end{split}$
(3.4)

that appear in the spherical derivatives. As we see, knowing the asymptotic
estimates for $\widetilde{J}^{\prime}_{N}$ and
$\widetilde{J}^{\prime\prime}_{N}$ is not sufficient, so we have to consider
the difference in the closed form, and to achieve this, we use the closed form
expression (3.1) of the Jackson kernel.

######  Lemma 3.0.2.

Let $\widetilde{J}_{N}$ be a trigonometric polynomial defined in (3.1), then
for some $\omega\neq 0$ and for $n=\lfloor{N/2}\rfloor$, the following
estimations holds true

$\begin{split}|\widetilde{J}_{N}(\omega)|&\leq\frac{\pi^{4}}{(n+1)^{4}|\omega|^{4}}\quad|\widetilde{J}^{\prime}_{N}(\omega)|\leq\frac{3\cdot\pi^{4}}{(n+1)^{3}|\omega|^{4}}\quad|\widetilde{J}^{\prime\prime}_{N}(\omega)|\leq\frac{12.5\cdot\pi^{4}}{(n+1)^{2}|\omega|^{4}}\\\
|G_{1}(\omega)|&\leq\frac{2\cdot\pi^{4}}{(n+1)^{2}|\omega|^{4}}\quad|G_{2}(\omega)|\leq\frac{14.5\cdot\pi^{4}}{(n+1)^{2}|\omega|^{4}}\quad|G_{3}(\omega)|\leq\frac{8\cdot\pi^{4}}{(n+1)|\omega|^{4}}\\\
|\widetilde{J}^{\prime\prime\prime}_{N}(\omega)|&\leq\frac{50.5\cdot\pi^{4}}{(n+1)^{4}|\omega|^{4}},\end{split}$

where the functions $G_{1},G_{2},G_{3}$ are defined in (3.4). Furthermore, for
$\displaystyle|\omega|\leq\frac{\pi}{4(n+1)}$ we have

$\begin{split}|\widetilde{J}^{\prime\prime}_{N}(\omega)-\cos(\omega)G_{1}(\omega)|&\leq\frac{\widetilde{J}^{(4)}_{N}(0)}{2}|\omega|^{2},\;|\widetilde{J}^{\prime\prime}_{N}(0)-\widetilde{J}^{\prime\prime}_{N}(\omega)|\leq\frac{\widetilde{J}^{(4)}_{N}(0)}{2}|\omega|^{2}\\\
|G_{2}(\omega)|&\leq\frac{\widetilde{J}^{(4)}_{N}(0)}{2}|\omega|^{2},\;|G_{3}(\omega)|\leq
0.52\cdot\widetilde{J}^{(4)}_{N}(0)|\omega|,\\\
\left|\frac{\widetilde{J}^{\prime}_{N}(\omega)-\cos(\omega)\sin(\omega)\widetilde{J}^{\prime\prime}_{N}(\omega)}{\sin^{2}(\omega)}\right|&\leq
0.52\cdot\widetilde{J}^{(4)}_{N}(0)|\omega|,\quad|\widetilde{J}^{\prime\prime\prime}_{N}(0)|\leq\widetilde{J}^{(4)}_{N}(0)|\omega|,\end{split}$

and for $\omega=0$ we have

$\begin{split}\widetilde{J}_{N}(0)&=1,\quad\quad\widetilde{J}^{\prime}_{N}(0)=\widetilde{J}^{\prime\prime\prime}_{N}(0)=0,\quad\quad\widetilde{J}^{\prime\prime}_{N}(0)=-\frac{n(n+2)}{3},\\\
\widetilde{J}^{(4)}_{N}(0)&=\frac{1}{30}n(n+1)(9n(n+2)-2).\end{split}$

For the proof of Lemma 3.0.2 please see Appendix A.

Now, using these bounds we can formulate the localization property of the
spherical derivatives of the Jackson kernel. More precisely, the bounds on the
derivatives in normal coordinates are stated in the Theorem 3.0.1 and the
those in polar coordinates are derived in Theorem 3.0.2. In the proofs, we use
several identities regarding the cross product, see Appendix C.

######  Theorem 3.0.1.

The Jackson kernel (3.1) fulfills for $x\neq y$ and $n=\lfloor{N/2}\rfloor$

$\displaystyle|J_{N}(x,y)|$
$\displaystyle\leq\frac{\pi^{4}}{(n+1)^{4}d(x,y)^{4}},\quad|X_{n}^{y}J_{N}(x,y)|\leq\frac{3\cdot\pi^{4}}{(n+1)^{3}d(x,y)^{4}},$
$\displaystyle|X_{i}^{x}X_{n}^{y}J_{N}(x,y)|$
$\displaystyle\leq\frac{16.5\cdot\pi^{4}}{(n+1)^{2}d(x,y)^{4}},\quad|X_{i}^{x}X_{n}^{x}J_{N}(x,y)|\leq\frac{16.5\cdot\pi^{4}}{(n+1)^{2}d(x,y)^{4}}$
$\displaystyle|X_{j}^{x}X_{i}^{x}X_{n}^{x}J_{N}(x,y)|$
$\displaystyle\leq\frac{101\cdot\pi^{4}}{(n+1)d(x,y)^{4}}$

and for the case $x=y$ we have

$\displaystyle J_{N}(x,x)$ $\displaystyle=1,\quad
X_{i}^{x}X_{i}^{y}J_{N}(x,x)=-X_{i}^{x}X_{n}^{x}J_{N}(x,x)=-\widetilde{J}^{\prime\prime}_{N}(0)$
$\displaystyle X_{n}^{y}J_{N}(x,x)$
$\displaystyle=X_{i}X_{n}^{y}J_{N}(x,x)=X_{j}^{x}X_{i}^{y}X_{n}^{y}J_{N}(x,x)=0.$

###### Proof.

The first estimate follows directly from Lemma 3.0.2. For the second estimate,
we first calculate the derivative of $J_{N}$ using its representation (3.2) to
get

$X_{n}^{y}J_{N}(x,y)=\sum\limits_{l=0}^{N}\displaystyle\frac{2l+1}{4\pi}\widetilde{w}_{l}P_{l}^{\prime}(\langle{x,y}\rangle)\langle{x,\eta_{n}^{y}}\rangle$
(3.5)

which immediately yields $X_{n}^{y}J_{N}(x,x)=X_{n}^{y}J_{N}(x,-x)=0$. In case
$x\in\mathbb{S}^{2}\setminus\\{y,-y\\}$, we have for $X_{n}^{y}J_{N}(x,y)$ the
following representation

$\displaystyle X_{n}^{y}J_{N}(x,y)$
$\displaystyle=-\widetilde{J}^{\prime}_{N}(d(x,y))\frac{\langle{x,\eta_{n}^{y}}\rangle}{\sin(d(x,y))}=\pm\widetilde{J}^{\prime}_{N}(d(x,y))\frac{\langle{x,\eta_{i}^{y}\times
y}\rangle}{\sin(d(x,y))}$
$\displaystyle=\mp\widetilde{J}^{\prime}_{N}(d(x,y))\frac{\langle{\eta_{i}^{y},x\times
y}\rangle}{\sin(d(x,y))}=\mp\widetilde{J}^{\prime}_{N}(d(x,y))\langle{\eta_{i}^{y},n_{x,y}}\rangle,$

where $n_{x,y}$ denotes the unique unit vector perpendicular to $x$ and $y$.
This together with Lemma 3.0.2 yields the second estimate. For the next
estimates, let us compute the second derivatives using the representation
$X_{n}^{y}J_{N}$ in terms of Legendre polynomials (3.5) to get

$\displaystyle X_{i}^{x}X_{n}^{y}J_{N}(x,y)$
$\displaystyle=\sum\limits_{l=0}^{N}\displaystyle\frac{2l+1}{4\pi}\widetilde{w}_{l}\Big{(}P_{l}^{\prime\prime}(\langle{x,y}\rangle)\langle{x,\eta_{n}^{y}}\rangle\langle{y,\eta_{i}^{x}}\rangle+P_{l}^{\prime}(\langle{x,y}\rangle)\langle{\eta_{i}^{x},\eta_{n}^{y}}\rangle\Big{)}$
$\displaystyle X_{i}^{x}X_{n}^{x}J_{N}(x,y)$
$\displaystyle=\sum\limits_{l=0}^{N}\displaystyle\frac{2l+1}{4\pi}\widetilde{w}_{l}\Big{(}P_{l}^{\prime\prime}(\langle{x,y}\rangle)\langle{\eta_{n}^{x},y}\rangle\langle{y,\eta_{i}^{x}}\rangle-\delta_{in}P_{l}^{\prime}(\langle{x,y}\rangle)\langle{x,y}\rangle\Big{)}.$

For $i\neq n$, it follows
$X_{i}^{x}X_{n}^{y}J_{N}(x,x)=X_{i}^{x}X_{n}^{y}J_{N}(x,-x)=X_{i}^{x}X_{n}^{x}J_{N}(x,x)=X_{i}^{x}X_{n}^{x}J_{N}(x,-x)=0$.
In case $i=n$, we have

$\begin{split}-X_{i}^{x}X_{i}^{x}J_{N}(x,x)&=X_{i}^{x}X_{i}^{y}J_{N}(x,x)=\sum\limits_{l=0}^{N}\displaystyle\frac{2l+1}{4\pi}\widetilde{w}_{l}P_{l}^{\prime}(1)\\\
&=\lim\limits_{t\to
1}-\frac{\widetilde{J}^{\prime}_{N}(\mathrm{arccos}{(t)})}{\sqrt{1-t^{2}}}=\lim\limits_{\omega\to
0}-\frac{\widetilde{J}^{\prime}_{N}(\omega)}{\sin{(\omega)}}=-\widetilde{J}^{\prime\prime}_{N}(0).\end{split}$
(3.6)

Following the same line, we obtain

$\begin{split}-X_{i}^{x}X_{i}^{x}J_{N}(x,-x)=X_{i}^{x}X_{i}^{y}J_{N}(x,-x)&=\sum\limits_{l=0}^{N}\displaystyle\frac{2l+1}{4\pi}\widetilde{w}_{l}P_{l}^{\prime}(-1)\\\
&=\lim\limits_{\omega\to\pi}\frac{\widetilde{J}^{\prime}_{N}(\omega)}{\sin{(\omega)}}=-\widetilde{J}^{\prime\prime}_{N}(\pi).\end{split}$
(3.7)

For $x\in\mathbb{S}^{2}\setminus\\{y,-y\\}$, using the closed form of
$X_{n}^{y}J_{N}$, the second derivatives reads as

$\begin{split}X_{i}^{x}X_{n}^{y}J_{N}(x,y)&=\widetilde{J}^{\prime\prime}_{N}(d(x,y))\frac{\langle{x,\eta_{n}^{y}}\rangle\langle{\eta_{i}^{x},y}\rangle}{\sin^{2}{d(x,y)}}\\\
&-\widetilde{J}^{\prime}_{N}(d(x,y))\left(\frac{\langle{x,\eta_{n}^{y}}\rangle\langle{\eta_{i}^{x},y}\rangle\cos{(d(x,y))}}{\sin^{3}{(d(x,y))}}+\frac{\langle{\eta_{i}^{x},\eta_{n}^{y}}\rangle}{\sin{(d(x,y))}}\right)\\\
&=\langle{\eta_{i}^{y},n_{x,y}}\rangle\langle{\eta_{n}^{x},n_{x,y}}\rangle
G_{2}(d(x,y))-\langle{\eta_{i}^{x},\eta_{n}^{y}}\rangle
G_{1}(d(x,y)),\end{split}$

and analogically we get

$X_{i}^{x}X_{n}^{x}J_{N}(x,y)=\langle{\eta_{i}^{x},n_{x,y}\langle{\eta_{n}^{x},n_{x,y}}\rangle}\rangle
G_{2}(d(x,y))+\delta_{in}\cos{(d(x,y))}G_{1}(d(x,y)),$

where $G_{1},G_{2}$ are defined in (3.4) and again $n_{x,y}$ denotes the
unique unit vector perpendicular to both $x$ and $y$. This results in

$\begin{split}|X_{i}^{x}X_{n}^{x}J_{N}(x,y)|,|X_{i}^{x}X_{n}^{y}J_{N}(x,y)|&\leq|G_{2}(d(x,y))|+|G_{1}(d(x,y))|\\\
&\leq\frac{16.5\cdot\pi^{4}}{(n+1)^{2}d(x,y)^{4}},\end{split}$

and yields the estimates for the second derivatives. The third derivatives can
be represented as

$\begin{split}X_{j}^{x}X_{i}^{x}X_{n}^{y}J_{N}(x,y)&=\sum\limits_{l=0}^{N}\displaystyle\frac{2l+1}{4\pi}\widetilde{w}_{l}\left(P_{l}^{\prime\prime\prime}(\langle{x,y}\rangle)\langle{\eta_{j}^{x},y}\rangle\langle{\eta_{i}^{x},y}\rangle\langle{\eta_{j}^{y},x}\rangle\right.\\\
&+P_{l}^{\prime\prime}(\langle{x,y}\rangle)\langle{\eta_{i}^{x},y}\rangle\langle{\eta_{n}^{y},\eta_{j}^{x}}\rangle+P_{l}^{\prime\prime}(\langle{x,y}\rangle)\langle{\eta_{i}^{x},\eta_{n}^{y}}\rangle\langle{y,\eta_{j}^{x}}\rangle\\\
&\left.-\delta_{ij}\langle{x,\eta_{n}^{y}}\rangle\big{(}P_{l}^{\prime}(\langle{x,y}\rangle)+P_{l}^{\prime\prime}(\langle{x,y}\rangle)\langle{x,y}\rangle\big{)}\right)\end{split}$

which immediately shows that

$X_{j}^{x}X_{i}^{x}X_{n}^{y}J_{N}(x,x)=X_{j}^{x}X_{i}^{x}X_{n}^{y}J_{N}(x,-x)=0.$

Using the closed form of the second derivatives for
$x\in\mathbb{S}^{2}\setminus\\{y,-y\\}$, the third derivatives can also be
represented as

$\begin{split}X_{j}^{x}X_{i}^{x}X_{n}^{y}J_{N}(x,y)&=\frac{\langle{\eta_{j}^{x},\eta_{n}^{y}}\rangle-\delta_{ij}\langle{x,\eta_{n}^{y}}\rangle\langle{x,y}\rangle}{\sin{(d(x,y))}}G_{3}(d(x,y))\\\
&-\frac{\langle{x,\eta_{n}^{y}}\rangle\langle{\eta_{i}^{x},y}\rangle\langle{\eta_{j}^{x},y}\rangle\langle{x,y}\rangle}{\sin^{3}{d(x,y)}}\left(\widetilde{J}_{N}^{\prime\prime\prime}(d(x,y))+\widetilde{J}^{\prime}_{N}(d(x,y))\right.\\\
&-3\cos{(d(x,y)}G_{3}(d(x,y))\Big{)}\\\
&+\frac{\langle{\eta_{i}^{x},\eta_{n}^{y}}\rangle\langle{\eta_{j}^{x},y}\rangle}{\sin{d(x,y)}}G_{3}(d(x,y))+\delta_{ij}\cos{(d(x,y))}G_{1}(d(x,y)).\end{split}$

Hence, using the estimates from Lemma 3.0.2 , we obtain

$|X_{j}^{x}X_{i}^{x}X_{n}^{y}J_{N}(x,y)|\leq\frac{101\cdot\pi^{3}}{(n+1)d(x,y)^{4}},$

which finishes the proof. ∎

######  Theorem 3.0.2.

Let $x,y,z\in\mathbb{S}^{2}$ be pairwise different and $x\neq-z$,then with
${n=\lfloor{N/2}\rfloor}$ the entries of the Hessian matrix $H$ of
$J_{N}(\cdot,y)$, $X_{n}^{y}J_{N}(\cdot,y)$ in polar coordinates centered at
$z$ fulfills

$\begin{split}|(HJ_{N}(x,y))_{ii}|&\leq\frac{16.5\cdot\pi^{4}}{(n+1)^{2}d(x,y)^{4}},\quad|(HJ_{N}(x,y))_{ij}|\leq\frac{14.5\cdot\pi^{4}}{(n+1)^{2}d(x,y)^{4}}\\\
|(HX_{n}^{y}J_{N}(x,y))_{ii}|&\leq\frac{103\cdot\pi^{4}}{(n+1)d(x,y)^{4}},\quad|(HX_{n}^{y}J_{N}(x,y))_{ij}|\leq\frac{93.5\cdot\pi^{4}}{(n+1)d(x,y)^{4}}.\end{split}$

In case $y=z$ and $d(x,y)\leq\frac{\delta}{(n+1)}$ with $\delta\leq\pi/4$, we
have

$\begin{split}|\widetilde{J}_{N}^{\prime\prime}(0)-(HJ_{N}(x,z))_{ii}|&\leq\frac{3}{20}\delta^{2}(n+1)^{2},\quad|(HJ_{N}(x,z))_{ij}|\leq\frac{3}{20}\delta^{2}(n+1)^{2}\\\
|(HX_{n}^{y}J_{N}(x,z))_{11}|&\leq\frac{3}{10}\delta(n+1)^{3},\quad|(HX_{y}^{n}J_{N}(x,z))_{22}|\leq\frac{1}{5}\delta(n+1)^{3}\\\
|(HX_{n}^{y}J_{N}(x,z))_{ij}|&\leq\frac{1}{5}\delta(n+1)^{3}.\end{split}$

###### Proof.

As it has been discussed at the beginning, in polar coordinates centered at
$z\in\mathbb{S}^{2}$, there is the local parametrization

$\phi^{\mathrm{pol}}(r,\theta)=\cos{(r)}z+\sin{(r)}\big{(}\cos{(\theta)}\eta_{1}^{z}+\sin{(\theta)}\eta_{2}^{z}\big{)},$

where $\eta_{1}^{z}$, $\eta_{1}^{z}$ form an orthonormal basis of
$T_{z}\mathbb{S}^{2}$, such that $\eta_{2}^{z}=\eta_{1}^{z}\times z$, and the
implicit inverse parametrization for $x\in\mathbb{S}^{2}\setminus\\{z,-z\\}$
is given by

$\begin{split}r(x)&=d(x,z),\\\
\cos{(\theta(x))}&=\frac{\langle{x,\eta_{1}^{z}}\rangle}{\sin{(d(x,z))}}=-\langle{\eta_{2}^{z},n_{z,x}}\rangle,\\\
\sin{(\theta(x))}&=\frac{\langle{x,\eta_{2}^{z}}\rangle}{\sin{(d(x,z))}}=\langle{\eta_{1}^{z},n_{z,z}}\rangle.\end{split}$

First, let us compute the partial derivatives of the function
$f_{\xi}(r,\theta)=\langle{\phi^{\mathrm{pol}}(r,\theta),\xi}\rangle$, for a
vector $\xi\in\mathbb{S}^{2}$. Carrying out simple computations results in

$\begin{split}\partial_{r}f_{\xi}(r,\theta)&=-\sin{(r)}\langle{z,\xi}\rangle+\cos{(r)}\big{(}\cos{(\theta)}\langle{\eta_{1}^{z},\xi}\rangle+\sin{(\theta)}\langle{\eta_{2}^{z},\xi}\rangle\big{)},\\\
\partial_{\theta}f_{\xi}(r,\theta)&=\sin{(r)}\big{(}\cos{(\theta)}\langle{\eta_{2}^{z},\xi}\rangle-\sin{(\theta)}\langle{\eta_{1}^{z},\xi}\rangle\big{)},\\\
\partial^{2}_{rr}f_{\xi}(r,\theta)&=-\cos{(r)}\langle{z,\xi}\rangle-\sin{(r)}\big{(}\cos(\theta)\langle{\eta_{1}^{z},\xi}\rangle+\sin(\theta)\langle{\eta_{2}^{z},\xi}\rangle\big{)}\\\
\partial^{2}_{\theta\theta}f_{\xi}(r,\theta)&=-\sin{(r)}\big{(}\cos(\theta)\langle{\eta_{1}^{z},\xi}\rangle+(\sin(\theta)\langle{\eta_{2}^{z},\xi}\rangle\big{)}\\\
\partial^{2}_{r\theta}f_{\xi}(r,\theta)&=\partial^{2}_{\theta
r}f_{\xi}(r,\theta)=\cos{(r)}\big{(}\cos{(\theta)}\langle{\eta_{2}^{z},\xi}\rangle-\sin{(\theta)}\langle{\eta_{1}^{z},\xi}\rangle\big{)}.\end{split}$

Thereafter, inserting the inverse parametrization and simplifying the obtained
expression, we get

$\begin{split}\partial_{r}f_{\xi}(r,\theta)&=\sin{(d(x,\xi))}\langle{n_{x,\xi},n_{z,x}}\rangle,\\\
\partial_{\theta}f_{\xi}(r,\theta)&=\sin{(d(x,z))}\sin{(d(x,\xi))}\langle{x,n_{x,\xi}\times
n_{z,x}}\rangle,\\\ \partial^{2}_{rr}f_{\xi}(r,\theta)&=-\cos{(d(x,\xi))},\\\
\partial^{2}_{\theta\theta}f_{\xi}(r,\theta)&=-\sin{(d(x,z))}\cos{(d(x,z))},\sin{(d(x,\xi))}\langle{n_{x,\xi},n_{z,x}}\rangle\\\
&\quad-\sin^{2}{(d(x,z))\cos{(d(x,\xi))}},\\\
\partial^{2}_{r\theta}f_{\xi}(r,\theta)&=\partial^{2}_{\theta
r}f_{\xi}(r,\theta)=\cos{(d(x,z))}\sin{(d(x,\xi))}\langle{x,n_{x,\xi},n_{z,x}}\rangle.\end{split}$

Now, let us proceed with calculation of the full derivatives. For
abbreviation, we use further the following notation

$f_{\xi}(x)=f_{\xi}(r(x),\theta(x))=\langle{\phi^{\mathrm{pol}}(r(x),\theta(x)),\xi}\rangle.$

First, we start with the kernel $J_{N}(\cdot,y)$ and assume that $x\neq-y$.
For the first element of the Hessian matrix we have

$\begin{split}(HJ_{N}(x,y))_{11}&=(\partial_{r}f_{y}(x))^{2}\left(\frac{\widetilde{J}_{N}^{\prime\prime}(\mathrm{arccos}{(f_{y}(x))})}{(1-f_{y}(x)^{2})}\right.\left.-\frac{\widetilde{J}_{N}^{\prime}(\mathrm{arccos}{(f_{y}(x))f_{y}(x)}}{(1-f_{y}(x)^{2})^{3/2}}\right)\\\
&-\partial^{2}_{rr}f_{y}(x)\frac{\widetilde{J}_{N}^{\prime}(\mathrm{arccos}{(f_{y}(x))}}{(1-f_{y}(x)^{2})^{1/2}}\\\
&=(\langle{n_{x,y},n_{z,x}}\rangle)^{2}G_{2}(d(x,y))+\cos(d(x,y))G_{1}(d(x,y)),\end{split}$
(3.8)

where $G_{1}$ and $G_{2}$ are given in (3.4). Thus, using Lemma 3.0.2, we get

$|(HJ_{N}(x,y))_{11}|\leq\frac{16.5\cdot\pi^{4}}{(n+1)^{2}d(x,y)^{4}}.$

If $x=-y$, we again use the polynomial representation of $J_{N}$, given in the
poof of Lemma 3.0.2, to derive

$\begin{split}(HJ_{N}(x,-x))_{11}&=\sum\limits_{l=0}^{N}\frac{2l+1}{4\pi}\widetilde{w}_{l}\Big{(}P_{l}^{\prime\prime}(\langle{x,-x}\rangle)(\partial_{r}f_{-x}(x))^{2}\\\
&+P_{l}^{\prime}(\langle{x,-x}\rangle)\partial_{rr}^{2}f_{-x}(x)\Big{)}\\\
&=\sum\limits_{l=0}^{N}\frac{2l+1}{4\pi}\widetilde{w}_{l}\Big{(}P_{l}^{\prime\prime}(-1)(\langle{x\times(-x),n_{z,x}}\rangle)^{2}\\\
&-P^{\prime}_{l}(-1)\cos(d(x,-x))\Big{)}=\sum\limits_{l=0}^{N}\frac{2l+1}{4\pi}\widetilde{w}_{l}P^{\prime}_{l}(-1).\end{split}$

Following the same argumentation as in (3.7), we get

$(HJ_{N}(x,-x))_{11}=\widetilde{J}_{N}(\pi).$

In case $y=z$ and $d(x,z)\leq\frac{\pi}{4(n+1)}$, from the representation
(3.8), we get

$(HJ_{N}(x,-x))_{11}=\widetilde{J}_{N}(d(x,z)).$

and combining this together with Lemma 3.0.2 results in

$|\widetilde{J}^{\prime\prime}_{N}(0)-(HJ_{N}(x,z))_{11}|\leq\frac{\widetilde{J}^{(4)}_{N}(0)}{2}d(x,z)^{2}\leq\frac{3}{20}d(x,z)^{2}(n+1)^{4}\leq\frac{3}{20}\delta^{2}(n+1)^{2}.$

For the second diagonal entry of the Hessian matrix, we have

$\begin{split}(HJ_{N}(x,y))_{22}&=\langle{x,n_{x,y}\times
n_{z,x}}\rangle^{2}\left(\widetilde{J}^{\prime\prime}_{N}(d(x,y))-\frac{\widetilde{J}^{\prime}_{N}(d(x,y))\cos{(d(x,y))}}{\sin{(d(x,y))})}\right)\\\
&+\cos{(d(x,y))}\frac{\widetilde{J}^{\prime}_{N}(d(x,y))}{\sin{(d(x,y))}}\\\
&=\langle{x,n_{x,y}\times
n_{z,x}}\rangle^{2}G_{2}(d(x,y))+\cos{(d(x,y))}G_{1}(d(x,y))\end{split}$

which provides the following estimate

$|(HJ_{N}(x,y))_{22}|\leq\frac{16.5\cdot\pi^{4}}{(n+1)^{4}d(x,y)^{4}}$

and gives zero entry $(HJ_{N}(x,-x))_{22}=0$ for $y=-x$.

In case $y=z$, we get

$(HJ_{N}(x,z))_{22}=\cos{(d(x,z))}\frac{\widetilde{J}^{\prime}_{N}(d(x,z))}{\sin{(d(x,z))}},$

that together with Lemma 3.0.2 provides the estimate

$|\widetilde{J}^{\prime\prime}_{N}(0)-(HJ_{N}(x,z))_{22}|\leq\frac{\widetilde{J}^{(4)}_{N}(0)}{2}d(x,z)^{2}\leq\frac{3}{20}d(x,z)^{2}(n+1)^{4}\leq\frac{3}{20}\delta^{2}(n+1)^{2}.$

For the off-diagonal entries of the Hessian matrix, we have

$(HJ_{N}(x,y))_{ij}=\langle{n_{x,y},n_{z,x}}\rangle\langle{x,n_{x,y}\times
n_{z,x}}\rangle G_{3}(d(x,y)),$

and analogous procedure to the one carried out above provides

$\begin{split}|(HJ_{N}(x,y))_{ij}|&\leq\frac{14.5\cdot\pi^{4}}{(n+1)^{2}d(x,y)^{4}},\quad|(HJ_{N}(x,z))_{ij}|\leq\frac{3}{20}\delta^{2}(n+1)^{2},\\\
(HJ_{N}(x,-x))_{ij}&=0.\end{split}$

For the derivative $X_{n}^{y}J_{N}(\cdot,y)$ of the kernel the estimates can
be derived folowing the same lines. Namely, the derivative
$X_{n}^{y}J_{N}(\cdot,y)$ is defined as

$X_{n}^{y}J_{N}(x,y)=-\frac{\widetilde{J}^{\prime}_{N}(d(x,y))}{\sin{(d(x,y))}}\langle{x,\eta_{n}^{y}}\rangle,$

where $\eta_{n}^{y}\in T_{z}\mathbb{S}^{2}$. Then, with $\xi=\eta_{n}^{y}$ and
the abbreviations in (3.4), we have

$\begin{split}\partial_{r}X_{n}^{y}J_{N}(x,y)&=\frac{\langle{x,\xi}\rangle\partial_{r}f_{y}(x)}{\sin(d(x,y))}G_{3}(d(x,y))-\partial_{r}f_{\xi}(x)G_{1}(d(x,y))\\\
\partial_{\theta}X_{n}^{y}J_{N}(x,y)&=\frac{\langle{x,\xi}\rangle\partial_{\theta}f_{y}(x)}{\sin(d(x,y))}G_{3}(d(x,y))-\partial_{\theta}f_{\xi}(x)G_{1}(d(x,y)).\end{split}$

In polar coordinates centered at $z\in\mathbb{S}^{2}$, we have for $y\neq-x$
and $\psi\in\\{r,\theta\\}$

$\begin{split}\partial_{\psi}G_{1}(d(x,y))&=-\frac{\partial_{\psi}f_{y}(x)}{\sin(d(x,y))}G_{3}(d(x,y)),\\\
\partial_{\psi}\left(\frac{G_{3}(d(x,y))}{\sin(d(x,y))}\right)&=-\frac{\partial_{\psi}f_{y}(x)}{\sin^{3}(d(x,y))}\Big{(}\widetilde{J}^{\prime\prime\prime}_{N}(d(x,y))+\widetilde{J}^{\prime}_{N}(d(x,y))\\\
&\quad-3\cos{(d(x,y))}G_{3}(d(x,y))\Big{)}.\end{split}$

Therefore, for the first diagonal entry of the Hessian matrix one has

$\begin{split}(HX_{n}^{y}J_{N}(x,y))_{11}&=\partial_{rr}X_{n}^{y}J_{N}(x,y)\\\
&=\frac{\partial_{r}f_{\xi}(x)\partial_{r}f_{y}(x)}{\sin(d(x,y))}G_{2}(d(x,y))+\frac{\langle{x,\xi}\rangle\partial_{rr}f_{y}(x)}{\sin(d(x,y))}G_{3}(d(x,y))\\\
&+\langle{x,\xi}\rangle\partial_{r}f_{y}(x)\partial_{r}\left(\frac{G_{3}(d(x,y))}{\sin{(d(x,y))}}\right)\\\
&-\partial_{rr}f_{\xi}(x)G_{1}(d(x,y))-\partial_{r}f_{\xi}(x)\partial_{r}G_{1}(d(x,y)).\end{split}$

Inserting the pre-computed derivatives results in

$\begin{split}(HX_{n}^{y}J_{N}(x,y))_{11}&=2\langle{n_{x,\xi},n_{z,x}}\rangle\langle{n_{z,x},n_{x,y}}\rangle\sin{(d(x,\xi))}G_{3}(d(x,y))\\\
&-\frac{\langle{x,\xi}\rangle\cos{(d(x,y))}}{\sin{(d(x,y))}}G_{3}(d(x,y))+\langle{x,\xi}\rangle
G_{1}(d(x,y))\\\
&-\frac{\langle{x,\xi}\rangle\langle{n_{z,x},n_{x,y}}\rangle^{2}}{\sin{(d(x,y))}}\left(\widetilde{J}^{\prime\prime\prime}_{N}(d(x,y))+\widetilde{J}^{\prime}_{N}(d(x,y))\right)\\\
&+\frac{3\langle{x,\xi}\rangle\langle{n_{z,x},n_{x,y}}\rangle^{2}\cos{(d(x,y))}}{\sin{(d(x,y))}}G_{3}(d(x,y)).\end{split}$

Since the following equalities hold

$\langle{x,\xi}\rangle=\langle{x,\eta_{n}^{y}}\rangle=\pm\langle{\eta_{i}^{y},x\times
y}\rangle=\pm\langle{\eta_{i}^{y},n_{x,y}}\rangle\sin(d(x,y)),$

we have that

$\left|\frac{\langle{x,\xi}\rangle}{\sin{(d(x,y))}}\right|\leq 1.$

Then using Lemma 3.0.2, the above derived representation of
$(HX_{n}^{y}J_{N}(x,y))_{11}$ we get

$|(HX_{n}^{y}J_{N}(x,y))_{11}|\leq\frac{103.5\cdot\pi^{4}}{(n+1)d(x,y)^{4}}.$

In case $y=-x$, we again use the polynomial representation of $J_{N}$, which
directly provides $(HX_{n}^{y}J_{N}(x,-x))_{11}=0.$

In case $y=z$ and $d(x,z)\leq\frac{\pi}{4(n+1)}$, we have

$(HX_{n}^{y}J_{N}(x,-x))_{11}=-\frac{\langle{x,\xi}\rangle}{\sin{(d(x,y))}}\widetilde{J}^{\prime\prime\prime}_{N}(d(x,z)),$

and again making use of Lemma 3.0.2, we get the estimate

$|(HX_{n}^{y}J_{N}(x,-x))_{11}|\leq\widetilde{J}^{(4)}_{N}(0)d(x,z)\leq\frac{3}{10}\delta(n+1)^{3}.$

Using an analogous procedure as above, the second diagonal entry can be
represented as

$\begin{split}(HX_{n}^{y}J_{N}(x,y))_{22}&=2\langle{x,n_{x,\xi}\times
n_{z,x}}\rangle\langle{x,n_{x,y}\times
n_{z,x}}\rangle\sin(d(x,\xi))G_{3}(d(x,y))\\\
&-\frac{\langle{x,\xi}\rangle\cos{(d(x,y))}}{\sin{(d(x,y))}}G_{3}(d(x,y))+\langle{x,\xi}\rangle
G_{1}(d(x,y))\\\ &-\frac{\langle{x,\xi}\rangle\langle{x,n_{x,y}\times
n_{z,x}}\rangle^{2}}{\sin{(d(x,y))}}\left(\widetilde{J}^{\prime\prime\prime}_{N}(d(x,y))+\widetilde{J}^{\prime}_{N}(d(x,y))\right),\\\
&+\frac{3\langle{x,\xi}\rangle\langle{x,n_{x,y}\times
n_{z,x}}\rangle^{2}\cos(d(x,y))}{\sin{(d(x,y))}}G_{3}(d(x,y)),\end{split}$

which shows that

$|(HX_{n}^{y}J_{N}(x,y))_{22}|\leq\frac{103.5\cdot\pi^{4}}{(n+1)d(x,y)^{4}}.$

Moreover, we also have $(HX_{n}^{y}J_{N}(x,-x))_{22}=0$, and using Lemma 3.0.2
gives

$\begin{split}|(HX_{n}^{y}J_{N}(x,z))_{22}|&=\left|\frac{\langle{x,\xi}\rangle}{\sin^{2}(d(x,y))}\left(\widetilde{J}^{\prime}_{N}(d(x,z))\right.\right.\\\
&\left.\left.-\cos(d(x,z))\sin{(d(x,z))}\widetilde{J}^{\prime\prime}_{N}(d(x,z))\right)\right|\\\
&\leq
0.52\cdot\widetilde{J}^{(4)}_{N}(0)d(x,z)\leq\frac{1}{5}\delta(n+1)^{3}.\end{split}$

Lastly, similar computation for the off-diagonal entries shows

$\begin{split}(HX_{n}^{y}J_{N}(x,y))_{ij}&=\langle{x,\langle{n_{x,y},n_{z,x}}\rangle(n_{x,\xi}\times
n_{z,x})}\rangle\sin{(d(x,\xi))}G_{3}(d(x,y))\\\
&+\langle{x,\langle{n_{x,\xi},n_{z,x}}\rangle(n_{x,y}\times
n_{z,x})}\rangle\sin{(d(x,\xi))}G_{3}(d(x,y))\\\
&-\frac{\langle{x,\xi}\rangle\langle{x,n_{x,y}\times
n_{z,x}}\rangle\langle{n_{x,y},n_{z,x}}\rangle}{\sin{(d(x,y))}}\left(\widetilde{J}^{\prime\prime\prime}_{N}(d(x,y))+\widetilde{J}^{\prime}_{N}(d(x,y))\right)\\\
&+\frac{3\langle{x,\xi}\rangle\langle{x,n_{x,y}\times
n_{z,x}}\rangle\langle{n_{x,y},n_{z,x}}\rangle\cos{(d(x,y))}}{\sin{(d(x,y))}}G_{3}(d(x,y)),\end{split}$

which results in

$\begin{split}|(HX_{n}^{y}J_{N}(x,y))_{ij}|&\leq\frac{93.5\cdot\pi^{4}}{(n+1)d(x,y)^{4}},\quad(HX_{n}^{y}J_{N}(x,-x))_{ij}=0,\\\
|(HX_{n}^{y}J_{N}(x,z))_{ij}|&=\left|\langle{x,n_{x,\xi}\times
n_{z,x}}\rangle\sin{(d(x,y))}G_{3}(d(x,y))\right|\\\
&\leq\frac{1}{5}\delta(n+1)^{3}.\end{split}$

∎

We finish this section bounding sums of point-wise expression from the
previous theorems.

Let us consider a discrete set $\mathcal{X}\subset\mathbb{S}^{2}$ such that
the minimal separation distance between the set elements is bounded from below
in the following way

$\rho(\mathcal{X})=\min\limits_{x_{i},x_{j}\in\mathcal{X},x_{i}\neq
x_{j}}d(x_{i},x_{j})\geq\frac{\nu}{n+1}.$ (3.9)

Involving classical ringing arguments on the sphere, the following result
holds true.

######  Lemma 3.0.3.

Let $x_{m}\in\mathcal{X}$, where $\mathcal{X}\subset\mathbb{S}^{2}$ is a
discrete set that satisfies a separation condition (3.9). Let a vector
$x\in\mathbb{S}^{2}$ be such that $d(x,x_{m})\leq\varepsilon\frac{\nu}{n+1}$
for $0\leq\varepsilon\leq 1/2$. Suppose that a function
$f\colon\mathbb{S}^{2}\to\mathbb{S}^{2}$ fulfills

$|f(x,y)|\leq\frac{c_{f}}{((n+1)d(x,y))^{s}}$ (3.10)

for $x\neq y$ and some $s\geq 3$, then

$\sum\limits_{x_{k}\in\mathcal{X}\setminus\\{x_{m}\\}}|f(x,x_{k})|\leq\frac{c_{f}a_{\varepsilon}}{\nu^{s}},$

where
$a_{\varepsilon}=\zeta(s-1)\cdot\min\\{9\cdot(1-\varepsilon)^{-s}+25,25\cdot(1-\varepsilon)^{-s}\\}$.
Here $\zeta$ denotes the Riemannian Zeta function.

Now, based on Lemma 3.0.3 we obtain the following bound for the point-wise
summation of the Jackson kernel and its derivatives.

######  Lemma 3.0.4.

Let $x_{m}\in\mathcal{X}$, where $\mathcal{X}\subset\mathbb{S}^{2}$ is a
discrete set, which obeys a separation condition (3.9). Let
$x\in\mathbb{S}^{2}$ such that $d(x,x_{m})\leq\varepsilon\frac{\nu}{n+1}$ for
$0\leq\varepsilon\leq 1/2$, then it follows

$\begin{split}&\sum\limits_{x_{k}\in\mathcal{X}\setminus\\{x_{m}\\}}|J_{N}(x,x_{k})|\leq\frac{\pi^{4}\cdot\tilde{a}_{\varepsilon}}{\nu^{4}},\;\sum\limits_{x_{k}\in\mathcal{X}\setminus\\{x_{m}\\}}|X_{n}^{y}J_{N}(x,x_{k})|\leq\frac{3\cdot\pi^{4}\cdot\tilde{a}_{\varepsilon}(n+1)}{\nu^{4}},\\\
&\sum\limits_{x_{k}\in\mathcal{X}\setminus\\{x_{m}\\}}|X_{i}^{x}X_{i}^{x}J_{N}(x,x_{k})|,\sum\limits_{x_{k}\in\mathcal{X}\setminus\\{x_{m}\\}}|(HJ_{N}(x,x_{k}))_{ii}|\leq\frac{16.5\cdot\pi^{4}\cdot\tilde{a}_{\varepsilon}(n+1)^{2}}{\nu^{4}},\\\
&\sum\limits_{x_{k}\in\mathcal{X}\setminus\\{x_{m}\\}}|X_{i}^{x}X_{i}^{x}X_{j}^{y}J_{N}(x,x_{k})|\leq\frac{101\cdot\pi^{4}\cdot\tilde{a}_{\varepsilon}(n+1)^{3}}{\nu^{4}},\\\
&\sum\limits_{x_{k}\in\mathcal{X}\setminus\\{x_{m}\\}}|(HJ_{N}(x,x_{k}))_{ij}|\leq\frac{14.5\cdot\pi^{4}\cdot\tilde{a}_{\varepsilon}(n+1)^{2}}{\nu^{4}},\\\
&\sum\limits_{x_{k}\in\mathcal{X}\setminus\\{x_{m}\\}}|(HX_{i}^{y}J_{N}(x,x_{k}))_{ii}|\leq\frac{103.5\cdot\pi^{4}\cdot\tilde{a}_{\varepsilon}(n+1)^{3}}{\nu^{4}},\\\
&\sum\limits_{x_{k}\in\mathcal{X}\setminus\\{x_{m}\\}}|(HX_{i}^{y}J_{N}(x,x_{k}))_{ij}|\leq\frac{93.5\cdot\pi^{4}\cdot\tilde{a}_{\varepsilon}(n+1)^{3}}{\nu^{4}},\end{split}$

where
$\tilde{a}_{\varepsilon}=\zeta(3)\cdot\min\\{9\cdot(1-\varepsilon)^{-4}+25,25\cdot(1-\varepsilon)^{-4}\\}$
and $\zeta$ denotes the Riemannian Zeta function.

###### Proof.

The above given estimates follow immediately after combining the estimates
from Theorems 3.0.1 and 3.0.2 and the results of Lemma 3.0.3. ∎

Having established the necessary localization estimates for the Jackson kernel
and its derivatives in this section, we are now able to construct and validate
a dual certificate using the Hermite interpolation in the next section.

##  4 Dual Certificate on the Sphere

###  4.1 Solution of the Interpolation Problem

In this section, we construct a dual certificate as a solution of the Hermite
type interpolation problem. Let us remind, we would like to solve the
interpolation problem

$\begin{split}q(x_{i})&=u_{i},\quad x_{i}\in\mathcal{X}\\\ |q(x)|&<1,\quad
x\in\mathbb{S}^{2}\setminus\mathcal{X},\end{split}$ (4.1)

where $q$ needs to be an element of $\Pi_{N}(\mathbb{S}^{2})$. Following
Section 2.3, we consider the interpolant of the form

$q(x)=\sum\limits_{i=0}^{M}\alpha_{0,i}J_{N}(x,x_{i})+\alpha_{1,i}X_{1}^{y}J_{N}(x,x_{i})+\alpha_{2,i}X_{2}^{y}J_{N}(x,x_{i}).$
(4.2)

and would like to show the the linear system of equations

$K\alpha=\begin{pmatrix}J_{N}&X_{1}^{x}J_{N}&X_{2}^{x}J_{N}\\\
X_{1}^{y}J_{N}&X_{1}^{x}X_{1}^{y}J_{N}&X_{2}^{x}X_{1}^{y}J_{N}\\\
X_{2}^{y}J_{N}&X_{1}^{x}X_{2}^{y}J_{N}&X_{2}^{x}X_{2}^{y}J_{N}\end{pmatrix}\begin{pmatrix}\alpha_{0}\\\
\alpha_{1}\\\ \alpha_{2}\end{pmatrix}=\begin{pmatrix}u\\\ 0\\\ 0\end{pmatrix}$
(4.3)

has a solution for the Jackson kernel $J_{N}$ chosen in Section 3.

We assume that the interpolation points $\mathcal{X}=\\{x_{j}\\}_{j=1}^{M}$
obey a minimal separation distance of the form

$\rho(\mathcal{X})=\min\limits_{x_{i},x_{j}\in\mathcal{X},x_{i}\neq
x_{j}}d(x_{i},x_{j})\geq\frac{\nu}{n+1},$ (4.4)

where $n=\lfloor{N/2}\rfloor$. For abbreviation, we use the following notation

$J_{ij}=X_{i}^{x}X_{j}^{y}J_{N}.$

As have been discussed before, to show the existence of the interpolating
polynomial $q$, we need to show the invertibility of the matrix $K$. Even
more, we need to partially compute the inverse of $K$ to derive bounds on the
coefficients $\alpha$. The following Lemma provides bounds the entries of the
interpolation matrix.

######  Lemma 4.1.1.

Suppose the set of points $\mathcal{X}=\\{x_{j}\\}_{j=1}^{M}$ satisfies the
separation condition (4.4). Then the entries of the interpolation matrix $K$
are bounded in the following way

$\begin{split}\|I-J_{00}\|&\leq\frac{C_{0}}{\nu^{4}},\quad\quad\quad\|-J^{\prime\prime}_{N}(0)I-J_{ii}\|_{\infty}\leq\frac{C_{2}(n+1)^{2}}{\nu^{4}},\\\
\|J_{0i}\|_{\infty},\;\|J_{i0}\|_{\infty}&\leq\frac{C_{1}(n+1)}{\nu},\quad\|J_{ij}\|_{\infty}\leq\frac{C_{2}(n+1)^{2}}{\nu^{4}},\quad\text{for
}i\neq j,\;i,j\neq 0,\end{split}$

where the constants are defined as

$\begin{split}C_{0}=25\cdot\zeta(3)\cdot\pi^{4},\quad
C_{1}=75\cdot\zeta(3)\cdot\pi^{4},\quad
C_{2}=412.5\zeta(3)\cdot\pi^{4}.\end{split}$ (4.5)

If $n\geq 9$ and $\nu^{4}>\frac{3}{0.99}\cdot C_{2}$, then for the inverse of
$J_{00}$ and $J_{ii}$ the following bounds are valid

$\|J_{00}^{-1}\|_{\infty}\leq\frac{1}{1-\frac{C_{0}}{\nu^{4}}},\quad\quad\|J_{ii}^{-1}\|_{\infty}\leq\frac{3}{0.99(n+1)^{2}(1-\frac{3C_{2}}{0.99\nu^{4}})}.$

###### Proof.

The bounds follows directly from Theorem 3.0.1 and Lemma 3.0.2, as well as
from the invertibility of a matrix $A$, if

$\|I-A\|_{\infty}<1,$

and the bounds on the norm of the inverse given by

$\|A^{-1}\|_{\infty}\leq\frac{1}{1-\|I-A\|_{\infty}}.$

For the bound on $J_{ii}^{-1}$, one uses additionally that for $n\geq 9$, so
that

${|J^{\prime\prime}_{N}(0)|=\frac{n(n+1)}{3}\geq\frac{0.99\cdot(n+1)^{2}}{3}}$

. ∎

Relying on the above stated results, we move on to the main Theorem of this
section, which provides a condition on the separation of the interpolation
points to guarantee the invertibility of the interpolation matrix $K$ and
gives bounds on the coefficients $\alpha$.

######  Theorem 4.1.1.

Suppose the separation condition (4.4) is satisfied by the set
$\mathcal{X}=\\{x_{j}\\}_{j=1}^{M}$ for some $\nu\geq 0$ and a constant $b\geq
3$ such that

$\nu^{4}\geq\frac{3}{0.99}\cdot b\cdot C_{2},$

where the constant $C_{2}$ is given in (4.5). Then the matrix $K$ is
invertible and the entries of the coefficient vector $\alpha$ in (4.3) satisfy

$\begin{split}\|\alpha_{0}\|_{\infty}&\leq
1+\frac{1}{45(b-2)-1},\quad\quad\quad\|\alpha_{j}\|_{\infty}\leq\frac{(n+1)^{-1}}{4.5(b-2)-0.1},\quad
j=1,2.\end{split}$

Furthermore, for elements of the vector $\alpha_{0}$ we have the lower bound

$|\alpha_{0,i}|\geq 1-\frac{1}{45(b-2)-1}.$

###### Proof.

In the first place, we partition the matrix $K$ into blocks as

$K=\begin{pmatrix}K_{0}&\widetilde{K}_{1}\\\ K_{1}&K_{2}\end{pmatrix}$

with the blocks given by

$\begin{split}K_{0}&=J_{00}=J_{N},\\\
K_{1}&=\left[J_{01}\;J_{02}\right]^{\mathrm{T}}=\left[X_{1}^{y}J_{N}\;X_{2}^{y}J_{N}\right]^{\mathrm{T}},\\\
\widetilde{K}_{1}&=\left[J_{10}\;J_{20}\right]^{\mathrm{T}}=\left[X_{1}^{x}J_{N}\;X_{2}^{x}J_{N}\right]^{\mathrm{T}},\\\
K_{2}&=\begin{bmatrix}J_{11}&J_{21}\\\
J_{12}&J_{22}\end{bmatrix}=\begin{bmatrix}X_{1}^{x}X_{1}^{y}J_{N}&X_{2}^{x}X_{1}^{y}J_{N}\\\
X_{1}^{x}X_{2}^{y}J_{N}&X_{2}^{x}X_{2}^{y}J_{N}\end{bmatrix},\end{split}$

after that we use an iterative block inversion. For abbreviation, we set

$a_{1}:=\frac{3\cdot C_{2}}{0.99\nu^{4}}$

that represents the quotient of the off-diagonal upper bound and the on-
diagonal lower bound. Thus, the assumption of the theorem reads now as

$a_{1}\leq\frac{1}{b},\quad\quad\text{where}\quad b\geq 3.$

Using the bounds from Lemma 4.1.1, one can shown that

$\left\|I-\frac{K_{2}/J_{22}}{\widetilde{J}^{\prime\prime}_{N}(0)}\right\|_{\infty}\leq\frac{1}{b-1}.$

Due to this, the invertibility of the matrix $K_{2}$ follows easily from

$\|(K_{2}/J_{22})^{-1}\|_{\infty}\leq\frac{3(b-1)}{0.99(n+1)^{2}(b-2)}.$

With the abbreviation $T=K_{2}/J_{22}$, the inverse matrix of $K_{2}$ has the
representation

$K_{2}^{-1}=\begin{pmatrix}T^{-1}&-T^{-1}J_{21}(J_{22})^{-1}\\\
-(J_{22})^{-1}J_{12}T^{-1}&(J_{22})^{-1}+(J_{22})^{-1}J_{12}T^{-1}J_{21}(J_{22})^{-1}.\end{pmatrix}$
(4.6)

In the next step, again relying on Lemma 4.1.1 and on the bound for
$(K_{2}/J_{22})^{-1}$, it can be shown that

$\|I-K/K_{2}\|_{\infty}\leq\frac{1}{45(b-2)},$

which in turn results in

$\|(K/K_{2})^{-1}\|\leq\frac{1}{1-\frac{1}{45(b-2)}}$

and in the invertibility of $K$. Now, using the representation (4.6) of the
inverse of $K_{2}$ and the abbreviation $S=K/K_{2}$, one comes to the solution
of the interpolation problem is given by

$\begin{split}\alpha_{0}&=S^{-1}u\\\
\alpha_{1}&=-T^{-1}(J_{01}-J_{21}(J_{22})^{-1}J_{02})\alpha_{0}\\\
\alpha_{2}&=-J_{22}^{-1}(J_{12}\alpha_{1}+J_{02}\alpha_{0}).\end{split}$

This yields the bounds

$\|\alpha_{0}\|_{\infty}\leq
1+\frac{1}{45(b-2)-1},\quad\quad\|\alpha_{1}\|_{\infty},\,\|\alpha_{2}\|_{\infty}\leq\frac{(n+1)^{-1}}{4.5(b-2)-0.1}.$

Moreover, the absolute value of $\alpha_{0,i}$ has the lower bound

$\begin{split}|\alpha_{0,i}|&\geq\left|\left(\left(I-\left(I-(K/K_{2})^{-1}\right)\right)u\right)_{i}\right|\\\
&\geq\left|1-\left|\left(\left(I-(K/K_{2})^{-1}\right)u\right)_{i}\right|\right|\\\
&\geq\left|1-\left|\left(I-K/K_{2}\right)(K/K_{2})^{-1}\right|\right|.\end{split}$

Since $b\geq 3$, the subtrahend in the last inequality can be bounded as

$\left|(I-K/K_{2})(K/K_{2})^{-1}\right|\leq\|I-K/K_{2}\|_{\infty}\|(K/K_{2})^{-1}\|_{\infty},$

which results in the corresponding bound for $\alpha_{0,i}$, namely

$|\alpha_{0,i}|\geq 1-\frac{1}{45(b-2)-1}.$

∎

Looking for explicit bounds for the coefficients $\alpha$, we obtain the
following result.

######  Corollary 4.1.1.

Suppose the set of interpolation points $\mathcal{X}=\\{x_{i}\\}_{m=0}^{M}$
satisfies a separation distance of

$\rho(\mathcal{X})\geq\frac{19.2\pi}{N}$ (4.7)

for $N\geq 20$. Then the interpolation problem (4.3) has a unique solution,
such that

$\|\alpha_{0}\|_{\infty}\leq 1+6.3\cdot
10^{-3},\quad\quad\|\alpha_{1}\|_{\infty},\,\|\alpha_{2}\|_{\infty}\leq\frac{6.3\cdot
10^{-2}}{(n+1)},$

and for elements of the vector $\alpha_{0}$ we have the lower bound

$|\alpha_{0,i}|\geq 1-6.3\cdot 10^{-3}.$

###### Proof.

For $n=\lfloor{N/2}\rfloor$, we have $(n+1)\geq N/2$ and therefore

$\rho(\mathcal{X})\geq\frac{9.6\pi}{(n+1)}.$

It is easy to check that with $\nu=9.6\pi$ and $b=5.6$, we have

$\nu^{4}\geq\frac{3}{0.99}\cdot b\cdot C_{2}.$

Thus, it follows from Theorem 4.1.1 that

$\begin{split}\|\alpha_{0}\|_{\infty}&\leq 1+\frac{1}{45\cdot 5.6-1}\leq
1+6.3\cdot 10^{-3}\\\ \|\alpha_{j}\|_{\infty}&\leq\frac{(n+1)^{-1}}{4.5\cdot
5.6-0.1}\leq\frac{6.3\cdot 10^{-2}}{(n+1)},\quad j=1,2,\end{split}$

and the lower bound for elements of $\alpha_{0}$ fulfills

$|\alpha_{0,i}|\geq 1-\frac{1}{45\cdot 5.6-1}\geq 1-6.3\cdot 10^{-3}.$

. ∎

###  4.2 Bounds for the Interpolant

Using the derived bounds on the coefficients of the interpolant, we proceed by
showing the upper bound in the absolute value of the interpolating function
$q$ with coefficients of Corollary 4.1.1, i.e.

$|q(x)|\leq 1,$

where $x$ is not an interpolation point. We split the proof into two parts. In
the first place, we consider the points that are close to an interpolation
point. This case in covered by Lemma 4.2.1, where we rely on the convexity
argument via the definiteness of Hessian matrix. Thereafter, we set our sights
on the bounds for points that are well separated from any interpolation
points, which forms the content of Lemma 4.2.2.

######  Lemma 4.2.1.

Suppose the set of interpolation points $\mathcal{X}=\\{x_{i}\\}_{m=0}^{M}$
satisfies the separation condition (4.4) and $x\in\mathbb{S}^{2}$ is such that
$0<d(x,x_{i})\leq\frac{\pi}{6(n+1)}$, for an interpolation point
$x_{m}\in\mathcal{X}$. Then the interpolating function $q$ from Corollary
4.1.1 fulfills

$|q(x)|\leq 1.$

###### Proof.

Let us consider the interpolant of the form (4.2), namely

$q(x)=\sum\limits_{k=0}^{M}\alpha_{0,k}J_{N}(x,x_{k})+\alpha_{1,k}X_{1}^{y}J_{N}(x,x_{k})+\alpha_{2,j}X_{2}^{y}J_{N}(x,x_{k}).$

Without loss of generality, let us assume that for a vector
$x_{m}\in\mathcal{X}$ we have ${q(x_{m})=u_{m}=1}$. To show that $u_{m}=1$ is
a local maximum of $q$, one needs to show that the Hessian matrix of $q$ is
negative definite at $x_{m}$. Considering the normal coordinates at $x_{m}$,
the Hessian matrix of $q$ at the interpolation point $x_{m}$ has the form

$\widetilde{H}q(x_{m})=\begin{pmatrix}X_{1}X_{1}q(x_{m})&X_{1}X_{2}q(x_{m})\\\
X_{2}X_{1}q(x_{m})&X_{2}X_{2}q(x_{m})\end{pmatrix}.$

Now, using the bounds of Theorem 3.0.1, Lemma 3.0.4 and Corollary 4.1.1, we
have

$\begin{split}X_{i}^{x}X_{i}^{x}q(x_{m})&=\sum\limits_{x_{k}\in\mathcal{X}}\alpha_{0,k}X_{i}^{x}X_{i}^{x}J_{N}(x_{m},x_{k})+\alpha_{1,k}X_{i}^{x}X_{i}^{x}X_{1}^{y}J_{N}(x_{m},x_{k})\\\
&+X_{i}^{x}X_{i}^{x}\alpha_{2,k}X_{2}^{y}J_{N}(x_{m},x_{k})\\\
&=\alpha_{0,m}X_{i}^{x}X_{i}^{x}J_{N}(x_{m},x_{m})+\sum\limits_{x_{k}\in\mathcal{X}\setminus\\{x_{m}\\}}\Big{(}\alpha_{0,k}X_{i}^{x}X_{i}^{x}J_{N}(x_{m},x_{k})\\\
&+\alpha_{1,k}X_{i}^{x}X_{i}^{x}X_{1}^{y}J_{N}(x_{m},x_{k})+\alpha_{2,k}X_{i}^{x}X_{i}^{x}X_{2}^{y}J_{N}(x_{m},x_{k})\Big{)}\\\
&\leq-|\alpha_{0,m}|\cdot\frac{n(n+2)}{3}+\Big{(}16.5\cdot\|{\alpha_{0}}\|_{\infty}+202\cdot\|{\alpha_{1}}\|_{\infty}\Big{)}\frac{\tilde{a}_{0}(n+1)^{2}}{(9.6)^{4}}.\end{split}$
(4.8)

Since $N\geq 20$, we have $n+1\geq 10$ and therefore

$|\widetilde{J}^{\prime\prime}_{N}(0)|=\frac{n(n+2)}{3}\geq\frac{n(n+1)}{(n+1)^{2}}\cdot\frac{(n+1)^{2}}{3}\geq
0.99\cdot\frac{(n+1)^{2}}{3}.$ (4.9)

And together with (4.8) it provides the upper bound for the on-diagonal
entries

$X_{i}X_{i}q(x_{m})\leq-0.2248\cdot(n+1)^{2}.$

In the same way, one can bound the off-diagonal entries by

$|X_{i}X_{j}q(x_{m})|\leq 0.1032\cdot(n+1)^{2}.$

Combing these two bounds, we obtain

$\mathrm{tr}(\widetilde{H}q(x_{m}))\leq-0.4496(n+1)^{2},\quad\mathrm{det}(\widetilde{H}q(x_{m}))\geq
0.0398(n+1)^{2},$

which means that the Hessian matrix of $q$ at point $x_{m}$ is strictly
negative definite and, thus, $x_{m}$ is an isolated local maximal point of $q$
and $q(x_{m})=1$ is local maximum.

For $x\neq x_{m}$ with $d(x,x_{m})\leq\frac{\pi}{6(n+1)}$ we use the same
arguments, but instead of the bounds from Theorem 3.0.2 we use the bounds
derived in Theorem 3.0.2 and the corresponding estimates from Lemma 3.0.4.
Thus, for the on-diagonal entries we get

$\begin{split}&(Hq(x))_{ii}=\\\
&=\sum\limits_{x_{k}\in\mathcal{X}}\alpha_{0,k}(HJ_{N}(x,x_{k}))_{ii}+\alpha_{1,k}(HX_{1}^{y}J_{N}(x,x_{k}))_{ii}+\alpha_{2,k}(HX_{2}^{y}J_{N}(x,x_{k}))_{ii}\\\
&\leq\alpha_{0,m}\widetilde{J}^{\prime\prime}_{N}(0)+\|{\alpha_{0}}\|_{\infty}\left(|\widetilde{J}^{\prime\prime}_{N}(0)-(HJ_{N}(x,x_{m}))_{ii}|+\sum\limits_{x_{k}\neq
x_{m}}|(HJ_{N}(x,x_{k}))_{ii}|\right)\\\
&+\|{\alpha_{1}}\|_{\infty}\left(|(HX_{1}^{y}J_{N}(x,x_{m}))_{ii}|+\sum\limits_{x_{k}\in\mathcal{X}\setminus\\{x_{m}\\}}|(HX_{1}^{y}J_{N}(x,x_{k}))_{ii}|\right)\\\
&+\|{\alpha_{2}}\|_{\infty}\left(|(HX_{2}^{y}J_{N}(x,x_{m}))_{ii}|+\sum\limits_{x_{k}\in\mathcal{X}\setminus\\{x_{m}\\}}|(HX_{2}^{y}J_{N}(x,x_{k}))_{ii}|\right).\end{split}$

For the first entry $(Hq(x))_{11}$, the above given inequality results in

$\begin{split}(Hq(x))_{11}&\leq(n+1)^{2}\left(-|\alpha_{0,m}|\frac{0.99}{3}+\|{\alpha_{0}}\|_{\infty}\left(\frac{3}{20}\delta^{2}+\frac{16.5\cdot\tilde{a}_{\delta/\nu}}{(9.6)^{4}}\right)\right)\\\
&+2\|{\alpha_{1}}\|_{\infty}(n+1)^{3}\left(\frac{3}{10}\delta+\frac{103.5\cdot\tilde{a}_{\delta/\nu}}{(9.6)^{4}}\right)\\\
&\leq-0.1552(n+1)^{2},\end{split}$

where $\delta=\pi/6$ and $\nu=9.6\pi$. Following the same line, for
$(Hq(x))_{22}$ and $(Hq(x))_{ij}$ for $i\neq j$ we obtain

$(Hq(x))_{22}\leq-0.1617(n+1)^{2},\quad|(Hq(x))_{ij}|\leq 0.1538(n+1)^{2},$

that together provide us the following bounds

$\mathrm{tr}({H}q(x))\leq-0.3169(n+1)^{2},\quad\mathrm{det}({H}q(x))\geq
0.0004(n+1)^{2},$

Hence, the function $q$ is strictly concave on
$B_{\frac{\pi}{6(n+1)}}(x_{m})\setminus\\{x_{m}\\}$, which shows that
$q(x)<1$. Moreover, the Taylor expansion of the cosine and sine function shows

$\begin{split}J_{N}(x,x_{m})&\geq
1+\frac{\widetilde{J}^{\prime\prime}_{N}(0)}{2}d(x,x_{m})^{2}\geq
1-\frac{(n+1)^{2}}{6}d(x,x_{m})^{2},\\\
|X_{n}^{y}J_{N}(x,x_{m})|&\leq|\widetilde{J}^{\prime}_{N}(d(x,x_{m}))|\leq|\widetilde{J}^{\prime\prime}_{N}(0)|d(x,x_{m})\leq\frac{(n+1)^{2}}{3}d(x,x_{m}),\end{split}$

meaning that for $d(x,x_{m})\leq\frac{\pi}{6(n+1)}$ we have

$\begin{split}q(x)&\geq\alpha_{0,m}J_{N}(x,x_{m})-\Big{(}\|\alpha_{1}\|_{\infty}|X_{1}^{y}J_{N}(x,x_{m})|+\|\alpha_{2}\|_{\infty}|X_{2}^{y}J_{N}(x,x_{m})|\\\
&+\sum\limits_{x_{k}\in\mathcal{X}\setminus\\{x_{m}\\}}\|\alpha_{0}\|_{\infty}|J_{N}(x,x_{k})|+\|{\alpha_{1}}\|_{\infty}|X_{1}^{y}J_{N}(x,x_{k})|\\\
&+\|{\alpha_{2}}\|_{\infty}|X_{2}^{y}J_{N}(x,x_{k})|\Big{)}\geq
0.92.\end{split}$

This shows that $0.92\leq q(x)\leq 1$ for $x$ such that
$d(x,x_{m})\leq\frac{\pi}{6(n+1)}$. In case $q(x_{i})=-1$, the analogous
arguments with changing sings show that $x_{i}$ is an isolated local minimal
point. ∎

######  Lemma 4.2.2.

Suppose the set of interpolation points $\mathcal{X}=\\{x_{i}\\}_{m=0}^{M}$
satisfies the separation condition (4.4) (4.4) and $x\in\mathbb{S}^{2}$ is
such that $d(x,x_{m})\geq\frac{\pi}{6(n+1)}$, for all interpolation point
$x_{m}\in\mathcal{X}$. Then the interpolating function $q$ from Corollary
4.1.1 fulfills

$|q(x)|\leq 1.$

###### Proof.

Here we split the proof into three parts, namely the first case corresponds to
those $x\in\mathbb{S}^{2}$ with $\frac{\pi}{6(n+1)}\leq
d(x,x_{m})\leq\frac{11\cdot\pi}{10(n+1)}$, the second to
$\frac{11\cdot\pi}{10(n+1)}\leq d(x,x_{m})\leq\frac{4\cdot\pi}{(n+1)}$ and the
third to those $x\in\mathbb{S}^{2}$ such that $\frac{4\cdot\pi}{(n+1)}\geq
d(x,x_{m})$ for all interpolation point $x_{m}\in\mathcal{X}$. So we need to
bound

$\begin{split}|q(&x)|=\left|\sum\limits_{k=0}^{M}\alpha_{0,k}J_{N}(x,x_{k})+\alpha_{1,k}X_{1}^{y}J_{N}(x,x_{k})+\alpha_{2,j}X_{2}^{y}J_{N}(x,x_{k})\right|\\\
&\leq\|{\alpha_{0}}\|_{\infty}|J_{N}(x,x_{m})|+\|{\alpha_{1}}\|_{\infty}|X_{1}^{y}J_{N}(x,x_{m})|+\|{\alpha_{2}}\|_{\infty}|X_{2}^{y}J_{N}(x,x_{m})|\\\
&+\sum\limits_{x_{k}\neq
x_{m}}\|\alpha_{0}\|_{\infty}|J_{N}(x,x_{k})|+\|{\alpha_{1}}\|_{\infty}|X_{1}^{y}J_{N}(x,x_{k})|+\|{\alpha_{2}}\|_{\infty}|X_{2}^{y}J_{N}(x,x_{k})|\end{split}$
(4.10)

In the first case, due to the Taylor expansion of the cosine and sine function
and the positivity of the Jackson kernel we have

$|J_{N}(x,x_{m})|=\widetilde{J}_{N}(d(x,x_{m}))\leq
1-\frac{|\widetilde{J}^{\prime\prime}_{N}(0)|}{2}d(x,x_{m})^{2}+\frac{|\widetilde{J}^{(4)}_{N}(0)|}{24}d(x,x_{m})^{4}.$
(4.11)

Then using again the fact for $N\geq 20$, $(n+1)\geq 10$ the inequality (4.9)
holds, and additionally

$|\widetilde{J}^{(4)}_{N}(0)|=\frac{1}{30}n(n+2)(9n(n+2)-2)\leq\frac{3}{10}(n+1)^{4},$

we can rewrite the inequality (4.11) for
$t\in[\frac{\pi}{6},\frac{11\pi}{10}]$ in the form

$\widetilde{J}_{N}(t/(n+1))\leq 1-\frac{0.99}{6}t^{2}+\frac{1}{80}t^{4}.$

The polynomial on the right hand side is monotonic decreasing for
$t\in[\frac{\pi}{6},t_{0}]$, where $t_{0}=\sqrt{\frac{20\cdot 0.99}{3}}\approx
2.5691$ and monotonic increasing for $t\in[t_{0},\frac{11\pi}{10}]$.
Similarly, due to $|\sin{(r\omega)}|\leq k\omega$ we have

$|X_{n}^{y}J_{N}(x,x_{m})|\leq|\widetilde{J}^{\prime\prime}_{N}(d(x,x_{m}))|\leq|\widetilde{J}^{\prime\prime}_{N}(0)|d(x,x_{m})\leq\frac{(n+1)^{2}}{3}d(x,x_{m}).$

Accordingly, for $\frac{\pi}{6(n+1)}\leq d(x,x_{m})\leq\frac{11\pi}{10(n+1)}$,
we can estimate (4.10) using the bounds above and the estimates from Lemma
3.0.4 and Corollary 4.1.1 as

$\begin{split}|q(x)|&\leq\|{\alpha_{0}}\|_{\infty}\left(1-\frac{0.99}{6}t^{2}+\frac{1}{80}t^{4}\right)+2\|{\alpha_{1}}\|_{\infty}\cdot\frac{t}{3}(n+1)\\\
&+\|{\alpha_{0}}\|_{\infty}\cdot\frac{\tilde{a}_{t/\nu}}{(9.6)^{4}}+2\|{\alpha_{1}}\|_{\infty}\cdot\frac{3\cdot\tilde{a}_{t/\nu}\cdot(n+1)}{(9.6)^{4}},\end{split}$

where $t=d(x,x_{m})(n+1)$, which results in the bounds

$|q(x)|\leq\left\\{\begin{matrix}0.9902&\frac{\pi}{6(n+1)}\leq
d(x,x_{m})\leq\frac{t_{0}}{(n+1)},\\\ 0.9682&\frac{t_{0}}{(n+1)}\leq
d(x,x_{m})\leq\frac{11\pi}{10(n+1)},\end{matrix}\right.$

and completes the first case. For the second case. i.e.
$\frac{11\pi}{10(n+1)}\leq d(x,x_{m})\leq\frac{4\pi}{(n+1)}$ we use the bounds
from Theorem 3.0.1

$\begin{split}|J_{N}(x,x_{m})|&\leq\frac{\pi^{4}}{(n+1)^{4}d(x,x_{m})^{4}}\leq\frac{\pi^{4}}{t^{4}},\\\
|X_{n}^{y}J_{N}(x,x_{m})|&\leq\frac{3\cdot\pi^{4}}{(n+1)^{3}d(x,x_{m})^{4}}\leq\frac{3\pi^{4}}{t^{4}}(n+1),\end{split}$

for $t=d(x,x_{m})(n+1)$. For this reason, we have the estimate

$|q(x)|\leq\|{\alpha_{0}}\|_{\infty}\frac{\pi^{4}}{t^{4}}+2\|{\alpha_{1}}\|_{\infty}\frac{3\pi^{4}}{t^{4}}(n+1)\\\
+\|{\alpha_{0}}\|_{\infty}\frac{\tilde{a}_{t/\nu}}{(9.6)^{4}}+2\|{\alpha_{1}}\|_{\infty}\frac{3\tilde{a}_{t/\nu}\cdot(n+1)}{(9.6)^{4}},$

which shows that for $\frac{11\pi}{10(n+1)}\leq
d(x,x_{m})\leq\frac{4\pi}{(n+1)}$ the absolute values of $q$ can be bounded as

$|q(x)|\leq 0.9618.$

Lastly, in case $d(x,x_{m})>\frac{4\pi}{(n+1)}$ for all interpolation points
$x_{m}\in\mathcal{X}$ the set $\mathcal{X}\cup\\{x\\}$ fulfills a separation
distance of $\frac{4\pi}{(n+1)}$. Therefore we can again apply the bounds from
Theorem 3.0.1, Lemma 3.0.3 and Corollary 4.1.1 to bound the absolute values of
$q$. Thus, we have

$\begin{split}|q(x)|&\leq\sum\limits_{x_{k}\in\mathcal{X}}\|\alpha_{0}\|_{\infty}|J_{N}(x,x_{k})|+\|{\alpha_{1}}\|_{\infty}|X_{1}^{y}J_{N}(x,x_{k})|+\|{\alpha_{2}}\|_{\infty}|X_{2}^{y}J_{N}(x,x_{k})|\\\
&\leq\|\alpha_{0}\|_{\infty}\frac{\tilde{a}_{0}}{4^{4}}+2\|\alpha_{1}\|_{\infty}\frac{3\tilde{a}_{0}(n+1)}{4^{4}}\leq
0.1618.\end{split}$

∎

Combining Corollary 4.1.1, Lemma 4.2.1 and Lemma 4.2.2 shows the existence of
a dual certificate $q$. We summarize this result in the following theorem.

######  Theorem 4.2.1.

Suppose that the set $\mathcal{X}=\\{x_{i}\\}_{i=0}^{M}$ fulfills a
separations distance of $\rho(\mathcal{X})\geq\frac{19.2\pi}{N}$ for $N\geq
20$. Then for each sign combinations $u_{i}\in\\{-1,1\\}$, there exists a
generalized polynomial $q\in\Pi_{N}(\mathbb{S}^{2})$ such that

$\begin{split}q(x_{i})&=u_{i},\quad x_{i}\in\mathcal{X},\\\ |q(x)|&<1,\quad
x\in\mathbb{S}^{2}\setminus\mathcal{X}.\end{split}$

The existence of a dual certificate immediately provides the recovery of a
sought measure via the minimization of the total variation.

######  Corollary 4.2.1.

Suppose the support of the singed measure $\mu^{\star}$ fulfills the
separation condition

$\min\limits_{x\neq y}d(x,y)\geq\frac{19.2\pi}{N},\quad
x,y\in\mathrm{supp}{(\mu^{\star})}\subset\mathbb{S}^{2},$

for $N\geq 20$, then the measure $\mu^{\star}$ is the unique solution of the
minimization problem

$\min\limits_{\mu\in\mathcal{M}(\mathbb{S}^{2},\mathbb{R})}\|\mu\|_{\mathrm{TV}},\quad\mbox{
subject to }\quad\mathcal{P}_{N}^{*}\mu=\mathcal{P}_{N}^{*}\mu^{\star}.$

###### Proof.

Theorem 4.2.1 guarantees the existence of a dual certificate. Hence, the
operator $\mathcal{P}^{*}$ has the null-space property with respect to
$\mathrm{supp}{(\mu^{\star})}$ and the therefore $\mu^{\star}$ is a unique
real solution of the minimization problem, which finishes the poof and this
section. ∎

##  5 Numerical Solution

###  5.1 Semidefinite Formulation of the Optimization Problem

At first glance, finding the solution of the total variation minimization
problem (RP) might seem rather complicated since it is an infinite dimensional
optimization problem over the whole measure space
$\mathcal{M}(\mathbb{S}^{2},\mathbb{R})$, therefore numerically not feasible.
The aim of this section is to provide a formulation of the convex program
problem (RP) such that it can be handled by the existing convex optimization
engines. To do so, we follow the ideas in [32, 34, 40].

First, let us switch to the convex pre-dual to (RP) that is given by

$\max\limits_{f\in
L^{2}(\mathbb{S}^{2})}\mathrm{Re}\langle{\mathcal{P}_{N}f,\mathcal{P}_{N}^{*}\mu^{\star}}\rangle\quad\mbox{
subject to }\quad\|{\mathcal{P}_{N}f}\|_{\infty}\leq 1.$

Since $\mathcal{P}_{N}$ in the projection operator given in (2.17),
immediately the equivalent formulation follows

$\max\limits_{f\in\Pi_{N}(\mathbb{S}^{2})}\mathrm{Re}\langle{f,\mathcal{P}_{N}^{*}\mu^{\star}}\rangle\quad\mbox{
subject to }\quad\|{f}\|_{\infty}\leq 1,$ (dRP)

where the constraint imposes that the modulus of the generalized polynomial

$\begin{split}f(x)=f(x(r,\theta))&=\sum_{l=0}^{N}\sum_{m=-l}^{l}f_{l,m}Y_{l}^{m}(r,\theta)=\sum_{l=0}^{N}\sum_{m=-l}^{l}f_{l,m}N_{lm}P_{l}^{m}(\cos(r))\mathrm{e}^{\mathrm{i}m\theta}\end{split}$

is uniformly bounded by $1$ over the whole $\mathbb{S}^{2}$, i.e. for
$(r,\theta)\in[0,\pi)\times[0,2\pi]$. On the grounds that there is a deeply
developed theory for multi-variate trigonometric polynomials providing
numerous condition on their boundedness on a frequency domain, see e.g. [15],
it would be clearly better to represent $f$ as a purely trigonometric
expressions. Since each associated Legendre polynomial $P_{l}^{m}(\cos(r))$
can be uniquely represented in as
$P_{l}^{m}(\cos(r))=\sum_{k=-l}^{l}p_{lmk}\mathrm{e}^{-\mathrm{i}kr}$, we
easily obtain

$\begin{split}f(r,\theta)=\sum_{l=0}^{N}\sum_{m=-l}^{l}\sum_{k=-l}^{l}f_{l,m}N_{lm}p_{lmk}\mathrm{e}^{-\mathrm{i}kr}\mathrm{e}^{\mathrm{i}m\theta}=\sum_{m=-N}^{N}\sum_{k=-N}^{N}\tilde{f}_{mk}\mathrm{e}^{-\mathrm{i}kr}\mathrm{e}^{\mathrm{i}m\theta},\end{split}$

where the coefficients ${f}_{mk}$, see [34], are defined as

$\tilde{f}_{mk}=\sum_{l=0}^{N}\tilde{f}_{mkl},\quad\text{with}\quad\tilde{f}_{lmk}=\left\\{\begin{matrix}f_{l,m}N_{lm}{p}_{lmk},&m,k\in[-l,l],\\\
0,&else.\end{matrix}\right.$

Now, to replace the norm constraint in (dRP) by a finite dimensional
conditions we use so-called the Bounded Real Lemma. To do this, let us first
throw light on the notion of the half-space $\mathcal{H}$. S set $\mathcal{H}$
is called a half-space of $\mathbb{Z}^{2}$ if
$\mathcal{H}\cap(-\mathcal{H})=\\{0\\}$, and
$\mathcal{H}\cup(-\mathcal{H})=\mathbb{Z}^{2}$ and
$\mathcal{H}+\mathcal{H}\subset\mathcal{H}$. A standard way to construct half-
space is given iteratively. We start with $\mathcal{H}_{1}=\mathbb{N}$ and we
say that $\mathbf{k}\in\mathcal{H}$ if either $k_{2}>0$ or $k_{2}=0$ and
$k_{1}\in\mathbb{N}$, of course, such iterative representation is very useful
for numerical purposes. With this preparation, a particular version of the
Bounded Real reads as following.

######  Lemma 5.1.1.

[15] Let $h(r,\theta)$ be a positive orthant polynomial defined by

$h(r,\theta)=\sum_{m=0}^{n}\sum_{k=0}^{n}h_{m,k}\mathrm{e}^{\mathrm{i}kr}\mathrm{e}^{\mathrm{i}m\theta},$

where ${(r,\theta)}\in\mathcal{D}=[-\pi,\pi]^{2}$. Then providing that
$|H(r,\theta)|\leq 1$ on the frequency domain $\mathcal{D}$, there exist a
positive semi-definite matrix
$Q_{0}\in\mathbb{C}^{(2N+1)^{2}\times(2N+1)^{2}}$ and s half-space
$\mathcal{H}\subset\mathbb{Z}^{2}$ such that

$\begin{bmatrix}Q&h\\\ h^{\mathrm{H}}&1\end{bmatrix}\succcurlyeq
0,\quad\mathrm{tr}(\Omega_{\mathbf{k}}Q)=\delta_{\mathbf{k}},\quad\mathbf{k}\in\mathcal{H},\quad-n\leq\mathbf{k}\leq
n,$ (5.1)

where ${h\in\mathbb{C}^{(2N+1)^{2}}}$ is a column stacked vector of
coefficients $h_{m,k}$.

Returning to (dRP), we see that the domain
$(r,\theta)\in[0,\pi]\times[0,2\pi)$ of the generalized polynomial $f$ is
slightly shifted version of the frequency domain $\mathcal{D}$ in Lemma 5.1.1,
moreover $f$ is not a positive orthant polynomial. So to reach the desired
frequency domain and to have the means to apply the Bounded Real Lemma 5.1.1
we can consider, for example, a trigonometric polynomial

$h(r,\theta)=\mathrm{e}^{Nr+N\theta}f(r+\pi,\theta+\pi)=\sum_{m=0}^{2N}\sum_{k=0}^{2N}(-1)^{(m+k)}\tilde{f}_{m-N,k-N}\mathrm{e}^{\mathrm{i}kr}\mathrm{e}^{\mathrm{i}m\theta}.$

The polynomial $h(r,\theta)$ is positive orthant, has the same magnitude as
$f$ does, and is defined on
$(r,\theta)\in[-\pi,0]\times[-\pi,\pi)\subset\mathcal{D}$. Hence, the dual
problem (dRP) is equivalent to

$\max\limits_{f\in\Pi_{N}(\mathbb{S}^{2}),\,Q}\mathrm{Re}\langle{f,\mathcal{P}_{N}^{*}\mu^{\star}}\rangle\quad\mbox{
subject to }\quad\eqref{finite_dim_const},$ (fdRP)

with $h_{m,k}=(-1)^{(m+k)}\tilde{f}_{m-N,k-N}$, for $m,k=0,\dots 2N.$

By strong duality, if $\mu^{\star}$ is a solution of the primal problem (RP)
and $f^{\star}$ is any solution of the pre-dual problem, then it follows

$\mathrm{Re}\langle{f^{\star},\mathcal{P}_{N}^{*}\mu^{\star}}\rangle=\mathrm{Re}\langle{\mathcal{P}_{N}f^{\star},\mu^{\star}}\rangle=\|\mu^{\star}\|_{\mathrm{TV}}.$

In the case $\mu^{\star}$ is a discrete measure, this implies that the
generalized polynomial is exactly equal to the sign of $\mu^{\star}$, when
$\mu^{\star}$ is not vanishing, namely

$f^{\star}(x_{i})=\mathrm{sign}_{\mu^{\star}}(x_{i}),\quad
x_{i}\in\mathrm{supp}{(\mu^{\star})}.$

which in turn means that the supporting points form a subset of the zeros of
the polynomial $1-|f^{\star}(x)|\in\Pi_{2N}(\mathbb{S}^{2})$.

Summing up, we obtain the following computational algorithm.

Data: low frequency information $\mathcal{P}_{N}^{*}\mu^{\star}$ of a measure
${\mu^{\star}\in\mathcal{M}(\mathbb{R},\mathbb{S}^{2})}$

begin

Solve for $f^{\star}\in\Pi_{N}(\mathbb{S}^{2})$,
$Q\in\mathbb{C}^{(2N+1)^{2}\times(2N+1)^{2}}$
$\begin{split}\max\limits_{f\in\Pi_{N},Q}\mathrm{Re}\langle{c,\mathcal{P}_{N}^{*}\mu^{\star}}\rangle\quad\mbox{
s.t. }\quad&\begin{bmatrix}Q&h\\\ h^{\mathrm{H}}&1\end{bmatrix}\succcurlyeq
0,\quad\mathrm{tr}(\Omega_{\mathbf{k}}Q)=\delta_{\mathbf{k}},\\\
&\quad\quad\quad\mathbf{k}\in\mathcal{H},\quad-n\leq\mathbf{k}\leq
n\end{split}$ using interior point method;

Randomly choose $P$ initial points
$\\{x_{1},\dots,x_{p}\\}\subset\mathbb{S}^{2}$

for  _$k=1,\dots,P$_ do

Find $x_{k}^{*}$ via Conjugate Gradient with initial point $x_{l}$;

if  _$(1-|q(x_{k}^{*})|\leq tol)$_ then

$\mathcal{X}=\mathcal{X}\cup\\{x_{k}^{*}\\}$

Set $\nu=\sum_{x_{i}\in\mathcal{X}}c_{i}\delta_{x_{i}}$ with $c_{i}$ such that
$\nu=\mathrm{argmin}\|\mathcal{P}_{N}^{*}\nu-y\|_{2}$

Result: $\nu\in\mathcal{M}(\mathbb{R},\mathbb{S}^{2})$

Algorithm 1 Super-resolution on $\mathbb{S}^{2}$ via SDP

###  5.2 Discretization of the Optimization Problem

To avoid the hing complexity of the the semi-definite program for the higher
order of moments we propose to use a discretization of the primal problem. In
this section we discuss the convergence behavior of the solution of the
discretized problem.

Let us choose a sequence of discrete sets
$\mathcal{G}_{n}\subset\mathbb{S}^{2}$, then the filling distance of
$\mathcal{G}_{n}$ is defined as

$\rho(\mathcal{G}_{n})=\sup\limits_{x\in\mathbb{S}^{2}}\inf\limits_{y\in\mathcal{G}_{n}}d(x,y).$

For the chosen sets $\mathcal{G}_{n}$, we consider the discretized version of
the primal problem (RP), namely

$\min\limits_{\mathrm{supp}{(\mu)}\subset\mathcal{G}_{n}}\|\mu\|_{\mathrm{TV}},\quad\mbox{
subject to }\quad\mathcal{P}_{N}^{*}\mu=\mathcal{P}_{N}^{*}\mu^{\star},$ (RPn)

To discuss the convergence behavior of the solution of the discretized problem
(RPn), we follow ideas from [29], where the convergence for continuously
parameterized dictionaries has been discussed.

######  Theorem 5.2.1.

Let assume that the measure $\mu^{\star}=\sum_{x_{i}\in X}c_{i}\delta_{x_{i}}$
is the unique solution of the primal problem (RP) and sequence
$\mathcal{G}_{n}$ of discretizations is chosen in such way that

$\rho(\mathcal{G}_{n})\to 0,$

when $n\to\infty$. Then each sequence of solutions $\mu_{n}$ of (RPn)
converges to $\mu^{\star}$ in the weak∗-topology. Moreover, there exists
$\varepsilon>0$, such that

$\mu_{n}(B_{\varepsilon}(x_{i}))\to
c_{i},\quad|\mu_{n}|(B_{\varepsilon}(x_{i}))\to c_{i},$

and

$|\mu_{n}|\Big{(}(\cup_{i}B_{\varepsilon}(x_{i}))^{c}\Big{)}\to 0.$ (5.2)

###### Proof.

Following the same line as in [29], we show that each sequence of solutions is
bounded and thus, due to the sequentially Banach-Alaoglu Theorem, admits a
weak∗ convergent subsequence that converges to as solution of the continuous
problem (RP). For the boundedness, let us consider the discretized convex pre-
dual problem to (RPn)

$\max\limits_{f\in\Pi_{N}(\mathbb{S}^{2})}\mathrm{Re}\langle{f,\mathcal{P}_{N}^{*}\mu^{\star}}\rangle\quad\mbox{
subject to }\quad\|{f}\|_{\infty}\leq 1,\,x\in\mathcal{G}_{n}$ (dRPn)

We show that the feasible sets of the dual problem, i.e. the sets of function
$f\in\Pi_{N}(\mathbb{S}^{2})$ such that $|f(x)|\leq 1$ for $x\in G_{n}$, are
bounded and therefore compact. By assumption $\rho(\mathcal{G}_{n})\to 0$,
therefore we have the for large enough $n$,

$\rho(\mathcal{G}_{n})\leq\frac{1}{N}.$

Applying the Marcinkiewicz-Zygmund inequality [16, 17], yields for all
$f\in\Pi_{N}(\mathbb{S}^{2})$

$\|{f}\|_{\infty}\leq(1-N\rho(\mathcal{G}_{n}))^{-1}\max\limits_{x\in\mathcal{G}_{n}}|f(x)|,$

meaning that all feasible sets are uniformly bounded and, thus, compact. This
shows that the discretization problem has a solution, and we denote these
minimizers by $f_{n}$.

The rest of the proof is identical to the proof of Theorem 2 in [29]. Let us
briefly sketch it. Due to the uniform boundedness, it can be shown that the
sequence $f_{n}$ of solutions of the discretized dual problem converges to the
solution $f^{\star}$ of the continuous dual problem. Since strong duality
holds for both the discretized problems and the continuous one, we have

$\mathrm{Re}\langle{f_{n},\mathcal{P}_{N}^{*}\mu^{\star}}\rangle=\|\mu_{n}\|_{\mathrm{TV}}\quad\text{and
}\quad\mathrm{Re}\langle{f^{\star},\mathcal{P}_{N}^{*}\mu^{\star}}\rangle=\|\mu^{\star}\|_{\mathrm{TV}},$

which shows that the sequence $\mu_{n}$ is bounded. Then due to the Banach-
Alaoglu Theorem we get the weak∗ convergence to the minimizer $\mu^{\star}$.
The convergence of the measure of the epsilon balls follows analogically to
Corollary 1 in [29]. ∎

Additionally, as it has been mentioned in [29], for fine enough discretization
the property (5.2) suggests that the support of the solutions $\mu_{n}$ of the
discrete minimization problems cluster around the support of $\mu^{\star}$.

Now let us proceed with the exact digitization procedure via considering the
following grid [8] in the spherical coordinates for inclination $[0,\pi)$ and
azimuth $[0,2\pi)$

$r_{k}=\frac{\pi(2j+1)}{4n},\quad\quad\theta_{k}=\frac{2\pi k}{2n},\quad\quad
j,k=0,\dots,2n-1,$

for some $n\in\mathbb{N}$. These points generate a grid $\mathcal{G}_{n}$ of
$4n^{2}$ points.

Let the matrix ${\bf Y}_{N}$ of spherical harmonics from
$\Pi_{N}(\mathbb{S}^{2})$ evaluated at the grid points,

${\bf
Y_{N}}=\big{(}Y_{m}^{\ell}(x(r_{k},\theta_{k})\big{)}_{\ell=0,|m|\leq\ell}^{N}\quad\in\mathbb{C}^{\dim(\Pi_{N})\times
4n^{2}}.$ (5.3)

Using such a matrix notation, the problem (RPn) can be transformed to the
following finite dimensional basis pursuit problem

$\min\limits_{c\in\mathbb{C}^{4n^{2}}}\|c\|_{1},\quad\mbox{ subject to
}\quad{\bf Y_{N}}\,c=y,$ (dicRPn)

where
$y=\big{(}\langle{Y_{m}^{\ell},\mu^{*}}\rangle\big{)}_{\ell=0,|m|\leq\ell}^{N}$
is the given data, and $c\in\mathbb{C}^{4n^{2}}$ represents the vector of
coefficients of the spherical harmonics. Since we measure is supported only on
a few points of the grid $\mathcal{G}_{n}$, we need to impose sparsity of $c$,
which can be done by minimizing the $\ell_{1}$-norm of $c$.

Due to the well-known basis mismatch phenomenon, a completely sparse solution
of (dicRPn) can not be obtained and one needs a certain threshold in absolute
value keeping only entries that are large enough, i.e. one keeps only points
with solution $c^{\star}\in\mathcal{G}_{n}$ such that $|c^{\star}|>thresh$ .
After that one can see that recover grid points cluster around the support of
$\mu^{\star}$. To find the centers of such cluster, and correspondingly the
support of $\mu^{\star}$ we use an algorithm known as the bivariate kernel
density estimator, for more details see [12], with the normal kernel
$\displaystyle K_{N}({\bf x})=(2\pi)^{-1/2}\exp{\left(-\frac{1}{2}\|{\bf
x}/h\|^{2}_{2}\right)}$ and a scaling parameter $h>0$. Combining all above
described steps results in our second algorithm based on the discretized
optimization problem (RPn).

Data: low frequency information
$y=\big{(}\langle{Y_{m}^{\ell},\mu^{*}}\rangle\big{)}_{\ell=0,|m|\leq\ell}^{N}$
of the measure ${\mu^{\star}\in\mathcal{M}(\mathbb{R},\mathbb{S}^{2})}$

Parameters: grid-parameter $n\in\mathbb{N}$, threshold $thresh>0$, scaling
parameter $h$, $tol>0$

begin

Solve for $c\in\mathbb{C}^{4n^{2}}$
$\min\limits_{c\in\mathbb{C}^{4n^{2}}}\|c\|_{1},\quad\mbox{ subject to
}\quad{\bf Y_{N}}\,c=y$ (dicRPn)

Choose $x_{k}\in G_{n}$ such that  $|c_{k}|>thresh.$

Set $\mathcal{X}=\\{x_{k}\\}_{k}$.

Apply the kernel density estimation algorithm with the scaling parameter $h$
and the tolerance $tol$ to the set $\mathcal{X}=\\{x_{k}\\}_{k}$ and generate
$\mathcal{X}_{mean}=\\{x_{i}^{mean}\\}_{i}$

Set $\nu=\sum_{x_{i}\in\mathcal{X}}c_{i}\delta_{x_{i}}$ with
$c_{i}\in\mathbb{R}$ such that
$\nu=\mathrm{argmin}\|\mathcal{P}_{N}^{*}\nu-y\|_{2}$

Result: $\nu\in\mathcal{M}(\mathbb{R},\mathbb{S}^{2})$

Algorithm 2 Super-resolution on $\mathbb{S}^{2}$ via Discretization

###  5.3 Numerical Experiments

####  5.3.1 Semi-definite Program

Experiment 1 [Noise-free recovery]. The first experiment demonstrates the
performance of Algorithm 1 for the signal reconstruction in the noiseless data
case. A discrete measure

$\mu^{\star}=\sum_{i=1}^{6}c_{i}\delta_{x_{i}}$

with randomly generated support points $x_{i}=x(r_{i},\theta_{i})$ and
amplitudes $c_{i}$, given in the Table 1 is considered. The support fulfills
the minimal separation distance condition

$\min\limits_{i\neq j}d(x_{i},x_{j})\geq 0.549335,$

and the given information about $\mu^{\star}$ is its low frequency information
up to degree $N=6$.

$i$ |  $\theta_{i}$ |  $r_{i}$ |  $c_{i}$
---|---|---|---
1 | 1.366427 | 0.412278 | 4.296273
2 | 1.983298 | 2.591331 | $-1.594284$
3 | 2.589166 | 4.898989 | $-1.058452$
4 | 0.630283 | 3.460063 | 2.005496
5 | 1.294025 | 4.299585 | 4.071419
6 | 3.016381 | 3.455708 | 3.196665
Table 1: The support points and amplitudes of the test measure $\mu^{\star}$.

Having this data at hand, we first compute the solution $f^{*}$ of the convex-
optimization problem (dRP) using CVX-Matlab tools. Then the measure support is
captured by looking for zeros of the function $1-|f^{*}|^{2}$. For this
purpose, the Matlab built-in solver fminunc is applied with $P=20000$ initial
randomly generated points on the sphere. We set the tolerance to
${tol=10^{-8}}$, and identify the resulting minima $x_{i}^{rec}$ that fulfill
$1-|f^{*}|^{2}<tol$. To ensure that there are no several points clustered near
a true support point, we finally apply the kernel density estimator procedure
described in the previous section. The recovery error is measured respectively
for the measure support and intensities by the quantities $\epsilon_{x}$ and
$\epsilon_{c}$ defined as

$\epsilon_{x}=\max\limits_{i}\min\limits_{j}d(x_{i},x_{j}^{rec})\quad\quad\text{and}\quad\quad\epsilon_{c}={\max\limits_{i}|c_{i}-c_{i}^{rec}|}.$
(5.4)

In such way, after running Algorithm 1 we get six recovered support points and
coefficients $\\{x_{j}^{rec}\\}_{j=1}^{6}$ and $\\{c_{i}^{rec}\\}_{i=1}^{6}$
such that

$\epsilon_{x}<2.5895457\cdot
10^{-7},\quad\quad\quad\quad\epsilon_{c}<3.4205068\cdot 10^{-7}$

The results of the experiment are illustrated in Figure 1. As can be seen the
test measure $\mu^{\star}$ is successfully recovered in the noise-free case
via Algorithm 1.

(a) Support of $\mu^{\star}$.  (b) Low-resolution information in longitude-
latitude.

(c) Low-resolution information in longitude-latitude.  (d) Solution
$f^{\star}$ of (dRP) on $\mathbb{S}^{2}$.

(e) Solution $f^{\star}$ of (dRP) in longitude-latitude.  (f) Recovered
support of $\mu^{\star}$.

Figure 1:  Reconstruction scenario of test measure $\mu^{\star}$ given in Tab.
1. The true support $x_{i}$ is marked with blue circles and the reconstructed
support $x_{i}^{rec}$ is indicated by small red crosses.

Experiment 2 [Super-resolution constant]. From the theoretical point of view,
see Corollary 4.1.1, a sufficient criteria for measure recovery is that the
support points $\\{x_{i}\\}$ satisfy the separation condition
$\min\limits_{i\neq j}d(x_{i},x_{j})\geq\frac{\nu}{N},$ for $\nu=19.2\pi$ and
$N\geq 20$. Nevertheless, from numerical point of view the minimal separation
distance may be much smaller. To investigate this issue, we proceed as
follows. Defining for $i=1,\dots,20$ the values $\gamma_{i}=i/10$, we generate
twenty two-point sets $M_{ij}$ on the sphere that enjoy the separation
distance with in the range $[\gamma_{i}-0.05,\gamma_{i}+0.05]$ for every
$j=1,\dots,20$. Each of obtained sets $M_{ij}$ is considered as support of
some measure $\mu_{ij}^{\star}$, while the amplitudes of $\mu_{ij}^{\star}$
are chosen randomly in the rage $[-1,1]\setminus\\{0\\}$. Then Algorithm 1
with $P=20000$ initial points is applied individually to each test measure
$\mu_{ij}^{\star}$.

Figure 2: Results of Experiment 2.

The operation of the algorithm is considered as successful, i.e. the support
of the corresponding measure is considered as recovered, provided that the
maximal recovery error fulfills $\epsilon_{x}<1.056\cdot 10^{-4}.$ Thereafter,
we counted the number of successful runs per separation distance $\gamma_{i}$.
A value $\gamma_{i}$ is considered as the numerical super-resolution constant
once all algorithm runs were successful. As experiments show, see in Fig. 2,
the numerical super-resolution constant $\nu_{{num}}$ is located within the
interval $[1.15\pi,1.33\pi]$.

Remark. When one of the measure support points has spherical coordinates
$(0,0)$, some difficulties by locating zeros of $1-|f^{*}|^{2}$ have been
observed. This might be avoided by using some other solver instead of fminunc.

Experiment 3 [Noisy data]. Assume the noise data scenario, i.e. the obtained
low-resolution information of the measure $\mu^{\star}\in M(\mathbb{S}^{2})$
is of the form

$y=\mathcal{P}_{N}^{\star}(\mu^{\star}+\varepsilon)$ (5.5)

with for some random (and independent of $k$ and $n$) noise term
$\varepsilon$. Here we consider the deterministic noise model assuming that
$\|\mathcal{P}_{N}\varepsilon\|_{L_{2}}\leq\delta$ for some
$\delta\in\mathbb{R}$. The optimization problem (RP), in this case, surely
requires a regularization techniques in order to provide meaningful results.
To this purpose, we consider the corresponding Thikonov-type problem

$\min_{\mu\in
M(\mathbb{S}^{2})}\frac{1}{2}\|\mathcal{P}_{N}^{\star}(\mu^{\star}+\varepsilon-\mu)\|_{L^{2}(\mathbb{S}^{2})}+\tau\|\mu\|_{TV}$
(5.6)

and its semi-definite relaxation

$\max\limits_{f\in\Pi_{N}(\mathbb{S}^{2})}\mathrm{Re}\langle{f,\mathcal{P}_{N}^{\star}(\mu^{\star}+\varepsilon)}\rangle-\tau\|f\|_{L_{2}(\mathbb{S}^{2})}\quad\text{
subject to }\;\|f\|_{\infty}\leq 1.$ (5.7)

Following the same idea as in Subsection 5.1, it can be shown that (5.7) can
be represented as the next finite-dimensional optimization program

$\max\limits_{f\in\Pi_{N}(\mathbb{S}^{2}),\,Q}\mathrm{Re}\langle{f,\mathcal{P}_{N}^{\star}\mu^{\star}}\rangle-\tau\|f\|_{L_{2}(\mathbb{S}^{2})}\quad\text{
subject to\ \eqref{finite_dim_const}}$ (dRPτ)

which we use in this experiment instead of the ordinary optimization program
in Algorithm 1.

To showcase the reconstruction process in the noise corruption scenario we
consider the setting of Experiment 1, but in this case, the data are perturbed
by different levels of noise. $\delta\in\\{10^{-j}\\}_{j=1}^{3}$. Figure 3 and
Table 2 illustrate the results of the experiment.

Noise level $\delta$ |  $\epsilon_{x}$ |  $\epsilon_{c}$ | $\quad\tau\quad$
---|---|---|---
$N=6$
$10^{-3}$ | $0.003281$ | $0.002101$ | ${10^{-2}}$
$10^{-2}$ | 0.038025 | 0.029545 | $10^{-1}$
$10^{-1}$ | 0.927668 | 0.909458 | $10^{-1}$
$N=5$
$10^{-3}$ | 0.014631 | 0.006691 | ${10^{-2}}$
$10^{-2}$ | 0.180824 | 0.091935 | $10^{-1}$
$10^{-1}$ | 0.941574 | 0.762983 | $10^{-1}$
$N=4$
$10^{-3}$ | 0.019439 | 0.020083 | ${10^{-2}}$
$10^{-2}$ | 0.448509 | 0.331496 | ${10^{-1}}$
$10^{-1}$ | 1.154241 | 0.517158 | ${10^{-1}}$
Table 2: Recovery error for the test measure $\mu^{\star}$ via the regularized
optimization problem (dRPτ) given the low-frequency information for $N=4,5,6$.

(a) Solution $f^{\star}$ of (dRP) on $\mathbb{S}^{2}$.  (b) Support of
$\mu^{\star}$ recovered via (dRP).

(c) Solution $f^{\star}$ of (dRPτ) on $\mathbb{S}^{2}$.  (d) Support of
$\mu^{\star}$ recovered via (dRPτ).

Figure 3:  Reconstruction of the test measure $\mu^{\star}$ given the low-
resolution information up to $N=6$ corrupted by deremninistic noise of level
$\delta=10^{-2}$. The true support $x_{i}$ is marked with blue circles and the
reconstructed support $x_{i}^{rec}$ is indicated by small red crosses.

####  5.3.2 Discretization of Semi-Definite Program

Experiment 1 [Grid size]. To showcase the performance of Algorithm 2, here we
consider the discrete measure $\mu^{\star}$ defined in Table 1 as a test
measure once more, assuming the low-frequency information of $\mu$ up to $N=6$
is given. Choosing different sizes $n\in\mathbb{N}$ of the grid
$\mathcal{G}_{n}$ and thresholding with $thresh=0.1$, Algorithm 2 is applied
to recover the support of the mass. The results of such experiments are
plotted in Fig. 4 and the corresponding recovery error is presented in Table
3.

As expected, with increasing grid size the point clusters become denser around
the true support of $\mu^{\star}$. This, in turn, leads to a good support
approximation after when applying the kernel density estimator procedure. And
the recovery error decreases with the grid size $n$ grows.

Grid size $n$ |  $\quad\epsilon_{x}$ |  $\epsilon_{c}$
---|---|---
20 | 0.568226 | 1.923598
40 | 0.020896 | 0.018803
80 | 0.008867 | 0.014393
Table 3: Recovery error for different grid size.

(a) For $n=20$, support clusters.  (b) For $n=20$, recovered support.

(c) For $n=40$, support clusters.  (d) For $n=40$, recovered support.

(e) For $n=80$, support clusters.  (f) For $n=80$, recovered support.

Figure 4: Performance of Algorithm 2 for different grid sizes. The true
measure support is marked with blue circles, clustered solutions are denoted
by magenta crosses and the reconstructed support is indicated by red crosses.

##  Appendix A Proof of Lemma 3.0.2

###### Proof.

For the asymptotic estimates, we will use the property

$|\sin(\omega/2)|\geq\frac{|w|}{\pi}\quad\text{for}\;\omega\in[-\pi,\pi],$
(A.1)

and component-wise estimates of the expressions
$|\omega|\leq\frac{\pi}{4(n+1)}$. The estimate for the kernel itself follows
immediately. For the first derivative, we have the following representation

$\widetilde{J}^{\prime}_{N}(\omega)=2\widetilde{F}_{n}(\omega)\widetilde{F}^{\prime}_{n}(\omega),$

where the second multiplier reads as

$\begin{split}\widetilde{F}^{\prime}_{n}(\omega)&=\frac{1}{2(n+1)^{2}\sin^{2}(\omega/2)}\left((n+1)\sin((n+1)\omega)\right.\\\
&\left.-\frac{2\cos(\omega/2)\sin^{2}((n+1)\omega/2)}{\sin(\omega/2)}\right)\\\
&=\frac{1}{2(n+1)^{2}\sin^{2}(\omega/2)}\Big{(}(n+1)\sin((n+1)\omega)\\\
&-2\cos(\omega/2)\sin((n+1)\omega/2)U_{n}(\cos(\omega/2))\Big{)}\end{split}$

and $U_{n}$ denotes the $n$-th order Chebychev polynomial of the second kind.
Due to $\|U_{n}\|_{\infty}=n$, we get

$|\widetilde{F}^{\prime}_{n}(\omega)|\leq\frac{1.5}{(n+1)\sin^{2}(\omega/2)}.$

Since the Fejér kernel can be written as

$\widetilde{F}_{n}(\omega)=\frac{1}{(n+1)}\left(1+2\sum_{k=1}^{n}\left(1-\frac{k}{n+1}\right)\right)\cos(k\omega),$

we easily get $\widetilde{F}^{\prime}_{n}(\pi)=0$ and therefore
$\widetilde{J}^{\prime}_{N}(\pi)=0$. Moreover, we have

$\begin{split}\frac{\widetilde{F}^{\prime}_{n}(\omega)}{\sin\omega}&=\frac{1}{2(n+1)^{2}\sin^{2}(\omega/2)}\left((n+1)\frac{\sin((n+1)\omega)}{\sin(\omega)}\right.\\\
&\left.-\frac{2\cos{(\omega/2)}\sin^{2}((n+1)\omega/2)}{\sin(\omega)\sin(\omega/2)}\right),\\\
&=\frac{1}{2(n+1)^{2}\sin^{2}(\omega/2)}\left((n+1)U_{n}(\cos(\omega))-(n+1)^{2}\widetilde{F}_{n}(\omega)\right).\end{split}$

Then the boundedness $\|\widetilde{F}_{n}\|_{\infty}=1$ provides the estimates

$\left|\frac{\widetilde{F}_{n}(\omega)}{\sin(\omega)}\right|\leq\frac{1}{\sin^{2}(\omega/2)}\quad\quad\text{and}\quad\quad\lim\limits_{\omega\to
0}\frac{\widetilde{F}^{\prime}_{n}(\omega)}{\sin(\omega)}=\widetilde{F}^{\prime\prime}(0).$

For $|\omega|\leq\frac{\pi}{4(n+1)}$, it can be shown that

$k^{2}\cos(k\omega)-\frac{k\sin(k\omega)\cos(\omega)}{\sin(\omega)}\leq 0,$

which leads to the inequality

$\displaystyle\left|k^{2}-\frac{k\sin(k\omega)\cos(\omega)}{\sin(\omega)}\right|$
$\displaystyle=k^{2}-\frac{k\sin(k\omega)\cos(\omega)}{\sin(d(x,z))}$
$\displaystyle\leq k^{2}(1-\cos(k\omega))\leq k^{4}\frac{\omega^{2}}{2}.$

Since, the first derivative of the Jackson kernel can be written as

$\widetilde{J}^{\prime}_{N}(\omega)=\frac{1}{(n+1)^{2}}\left(-2\sum\limits_{k=1}^{2n}c_{k}k\sin(k\omega)\right)$

with positive Fourier coefficients $c_{k}$, component-wise estimation shows
that

$\left|\widetilde{J}^{\prime\prime}_{N}(0)-\frac{\cos(\omega)\widetilde{J}^{\prime}_{N}(\omega)}{\sin(\omega)}\right|\leq\frac{\widetilde{J}^{(4)}_{N}(0)}{2}|\omega|^{2}$

In the same way, one derives

$|\widetilde{J}^{\prime\prime}_{N}(0)-\widetilde{J}^{\prime\prime}_{N}(\omega)|\leq\frac{\widetilde{J}^{(4)}_{N}(0)}{2}|\omega|^{2}.$

For the function $G_{2}$, we have the following representation

$\begin{split}G_{2}(\omega)&=\widetilde{J}^{\prime\prime}_{N}(\omega)-\frac{\widetilde{J}^{\prime}_{N}(\omega)\cos(\omega)}{\sin{(\omega)}}\\\
&=2\left(\left(\widetilde{F}^{\prime}_{n}(\omega)\right)^{2}+\widetilde{F}_{n}(\omega)\left(\widetilde{F}^{\prime\prime}_{n}(\omega)-\frac{\widetilde{F}^{\prime}_{n}(\omega)\cos(\omega)}{\sin(\omega)}\right)\right).\end{split}$

Since the second derivative of $\widetilde{F}_{n}$ reads as

$\begin{split}\widetilde{F}^{\prime\prime}_{n}(\omega)&=\frac{1}{2(n+1)^{2}\sin^{2}(\omega/2)}\Big{(}(n+1)^{2}\Big{(}(2+\cos(\omega))\widetilde{F}_{n}(\omega)+\cos((n+1)\omega)\Big{)}\\\
&-2(n+1)(1+\cos(\omega)U_{n}(\cos(\omega))\Big{)},\end{split}$

then we get

$\displaystyle\widetilde{F}^{\prime\prime}_{n}(\omega)-\frac{\widetilde{F}^{\prime}_{n}(\omega)\cos(\omega)}{\sin(\omega)}$
$\displaystyle=\frac{1}{2(n+1)^{2}\sin^{2}(\omega/2)}\left((n+1)^{2}\Big{(}2(1+\cos(\omega))\widetilde{F}_{n}(\omega)\right.$
(A.2)
$\displaystyle\left.+\cos((n+1)\omega)\Big{)}-(n+1)\Big{(}2+3\cos(\omega)U_{n}(\cos(\omega))\Big{)}\right)$
$\displaystyle=\frac{1}{2(n+1)^{2}\sin^{2}(\omega/2)}\left((n+1)^{2}\left(2(1+\cos(\omega)\right)\cdot\right.$
$\displaystyle\left(\widetilde{F}_{n}(\omega)-\frac{U_{n}(\cos(\omega))}{n+1}\right)+(n+1)^{2}\cos((n+1)\omega)$
$\displaystyle\left.-(n+1)\cos(\omega)U_{n}(\cos(\omega))\right).$

Using the above derived representations yields the following estimate

$\left|G_{2}(\omega)\right|=\left|\widetilde{J}^{\prime\prime}_{N}(\omega)-\frac{\widetilde{J}^{\prime}_{N}(\omega)\cos(\omega)}{\sin(\omega)}\right|\leq\frac{14.5\cdot\pi^{4}}{(n+1)^{2}|\omega|^{2}}$

and similar result for the second derivative

$\left|\widetilde{J}^{\prime\prime}_{N}(\omega)\right|\leq\frac{12.5\cdot\pi^{4}}{(n+1)^{2}|\omega|^{2}}.$

For $|\omega|\leq\frac{\pi}{4(n+1)}$, the estimate follows by estimating
component-wise in the trigonometric representation, namely

$\left|G_{2}(\omega)\right|=\left|\widetilde{J}^{\prime\prime}_{N}(\omega)-\frac{\widetilde{J}^{\prime}_{N}(\omega)\cos(\omega)}{\sin(\omega)}\right|\leq\frac{\widetilde{J}^{(4)}(0)}{2}w^{2}.$

Furthermore, at $w=0$ we have

$\widetilde{J}^{\prime\prime}_{N}(0)=2\widetilde{F}^{\prime\prime}_{N}(0)=-\frac{4}{(n+1)}\sum\limits_{k=1}^{n}\left(k^{2}-\frac{k^{3}}{(n+1)}\right)=\frac{n(n+2)}{3}.$

Next, we consider the expression

$G_{3}(\omega)=\frac{\widetilde{J}^{\prime\prime}_{N}(\omega)\sin(\omega)-\widetilde{J}^{\prime}_{N}(\omega)\cos(\omega)}{\sin^{2}(\omega)}.$

Observe, that for $\omega=\pi$ this expression vanishes. However, for
$\omega\neq\pi$, we have

$\begin{split}\frac{\widetilde{J}^{\prime\prime}_{N}(\omega)\sin(\omega)-\widetilde{J}^{\prime}_{N}(\omega)\cos(\omega)}{\sin^{2}(\omega)}&=2\left(\left(\frac{\widetilde{F}^{\prime}_{n}(\omega)}{\sin(\omega)}\right)\widetilde{F}^{\prime}_{n}(\omega)\right.\\\
&+\left.\widetilde{F}_{n}(\omega)\left(\frac{\widetilde{F}^{\prime\prime}_{n}(\omega)\sin(\omega)-\widetilde{F}^{\prime}_{n}(\omega)\cos(\omega)}{\sin^{2}(\omega)}\right)\right),\end{split}$

with

$\begin{split}\widetilde{F}_{n}(\omega)&\left(\frac{(\widetilde{F}^{\prime\prime}_{n}(\omega)\sin(\omega)-\widetilde{F}^{\prime}_{n}(\omega)\cos(\omega))}{\sin^{2}(\omega)}\right)=\frac{\sin^{2}((n+1)\omega/2)}{2(n+1)^{4}\sin^{4}(\omega/2)}\times\\\
&\left(2(n+1)^{2}\frac{(1+\cos{(\omega)})}{\sin(\omega)}\left(\widetilde{F}_{n}(\omega)-\frac{U_{n}(\cos(w))}{n+1}\right)\right.\\\
&\quad\quad\left.+(n+1)\frac{(n+1)\cos((n+1)\omega)-\cos(\omega)U_{n}(\cos(\omega))}{\sin(\omega)}\right).\end{split}$

Here the first the right-hand side summand can be represented as

$\begin{split}&\left|\frac{\sin^{2}((n+1)\omega/2)}{(n+1)^{2}\sin^{4}(\omega/2)}\cdot\frac{(1+\cos{(\omega)})}{\sin(\omega)}\left(\widetilde{F}_{n}(\omega)-\frac{U_{n}(\cos(w))}{n+1}\right)\right|\\\
&=\frac{1}{(n+1)\sin^{4}(\omega/2)}\left|\frac{\sin^{2}((n+1)\omega/2)}{(n+1)}\frac{\cos(\omega/2)}{\sin(\omega/2)}\left(\widetilde{F}_{n}(\omega)-\frac{U_{n}(\cos(w))}{n+1}\right)\right|\\\
&=\frac{1}{(n+1)\sin^{4}(\omega/2)}\left|\frac{\sin^{2}((n+1)\omega/2)}{(n+1)\sin(\omega/2)}\right|\times\\\
&\left|\sin((n+1)\omega/2)\cos(\omega/2)\left(\widetilde{F}_{n}(\omega)-\frac{U_{n}(\cos(w))}{n+1}\right)\right|\\\
&=\frac{\sqrt{\widetilde{F}_{n}(\omega)}}{(n+1)\sin^{4}(\omega/2)}\left|\sin((n+1)\omega/2)\cos(\omega/2)\left(\widetilde{F}_{n}(\omega)-\frac{U_{n}(\cos(w))}{n+1}\right)\right|\\\
&\leq\frac{2}{(n+1)\sin^{4}(\omega/2)}.\end{split}$

For the second summand, due to the derivative representation of Chebychev
polynomials, we have

$\frac{(n+1)\cos((n+1)\omega)-\cos(\omega)U_{n}(\cos{(\omega)})}{\sin(\omega)}=-\sin(\omega)U^{\prime}_{n}(\cos(\omega)).$

Using the Bernstein inequality for algebraic polynomials, i.e.

$|P_{n}^{\prime}(t)|\leq\frac{n}{\sqrt{1-x^{2}}}\|P_{n}\|_{\infty},\quad-1<x<1,$

for a polynomial of degree $n$, see e.g. [3], we derive

$\left|\frac{(n+1)\cos((n+1)\omega)-\cos(\omega)U_{n}(\cos{(\omega)})}{\sin(\omega)}\right|\leq
n\|U_{n}\|_{\infty}\leq(n+1)^{2}.$

Hence, the second summand can be estimated as

$\displaystyle\left|\frac{\sin^{2}((n+1)\omega/2)}{(n+1)^{2}\sin^{4}(\omega/2)}\cdot(n+1)\cdot\frac{(n+1)\cos((n+1)\omega)-\cos(\omega)U_{n}(\cos(\omega))}{\sin(\omega)}\right|$
$\displaystyle\leq\frac{0.5}{(n+1)\sin^{4}(\omega/2)}.$

Together, this yields

$\left|\widetilde{F}_{n}\left(\frac{(\omega)(\widetilde{F}^{\prime\prime}_{n}(\omega)\sin(\omega)-\widetilde{F}^{\prime}_{n}(\omega)\cos(\omega))}{\sin^{2}(\omega)}\right)\right|\leq\frac{2.5}{(n+1)\sin^{4}(\omega/2)}$

and consequently we get an estimate for $G_{3}$ of the form

$|G_{3}(\omega)|=\left|\frac{\widetilde{J}^{\prime\prime}_{N}(\omega)\sin(\omega)-\widetilde{J}^{\prime}_{N}(\omega)\cos(\omega)}{\sin^{2}(\omega)}\right|\leq\frac{8}{(n+1)\sin^{4}(\omega/2)}.$

Again, for $|\omega|\leq\frac{\pi}{4(n+1)}$, we estimate the expression
component-wise. Observe, that

$\left|\frac{k^{2}\cos(k\omega)\sin(\omega)-k\sin(k\omega)\cos(\omega)}{/sin^{2}(\omega)}\right|\frac{1}{1+\cos(\omega)}k^{3}\sin(k\omega)\leq\frac{k^{4}|w|}{1+\cos(\omega)}.$

For this reason, we have for $|\omega|\leq\frac{\pi}{4(n+1)}$

$\displaystyle\left|\frac{\widetilde{J}^{\prime\prime}_{N}(\omega)\sin(\omega)-\widetilde{J}^{\prime}_{N}(\omega)\cos(\omega)}{\sin^{2}(\omega)}\right|$
$\displaystyle\displaystyle\leq\frac{\widetilde{J}^{(4)}_{N}(0)|\omega|}{1+\cos(\omega)}\leq\frac{\widetilde{J}^{(4)}_{N}(0)|\omega|}{1+\cos({\pi}/{(4(n+1))})}$
$\displaystyle\displaystyle\leq\frac{\widetilde{J}^{(4)}_{N}(0)|\omega|}{1+\cos(\pi/8)}\leq
0.52\cdot\widetilde{J}^{(4)}_{N}(0)|\omega|$

It can be also shown that the following inequality holds

$\left|\frac{k^{2}\cos(k\omega)\cos(\omega)\sin(\omega)-k\sin(k\omega)}{\sin^{2}(\omega)}\right|\leq\frac{1}{1+\cos(\omega)}k^{3}\sin(k\omega)\leq\frac{k^{4}|\omega|}{1+\cos(\omega)},$

and thus

$\left|\frac{\widetilde{J}^{\prime}_{N}(\omega)-\cos(\omega)\sin(\omega)\widetilde{J}^{\prime\prime}_{N}(\omega)}{\sin^{2}(\omega)}\right|\leq
0.52\cdot\widetilde{J}^{(4)}_{N}(0)|\omega|$

Similarly, we can compute the third derivative of $\widetilde{J}_{N}$, that is

$\widetilde{J}^{\prime\prime\prime}_{N}(\omega)=2\Big{(}3\widetilde{F}^{\prime}_{n}(\omega)\widetilde{F}^{\prime\prime}_{n}(\omega)+\widetilde{F}_{n}(\omega)\widetilde{F}^{\prime\prime\prime}_{n}(\omega)\Big{)}.$

Writing explicitly the first and second summand, we get

$\begin{split}\widetilde{F}_{n}(\omega)\widetilde{F}^{\prime\prime\prime}_{n}(\omega)=\frac{1}{(n+1)\sin^{4}(\omega/2)}\left(-\frac{1}{2(n+1)}\cos(\omega/2)\sin((n+1)\omega/2)\cdot\right.\\\
U_{n}(\cos(\omega/2))\big{(}(5+\cos{(\omega)})\widetilde{F}_{n}(\omega)+3\cos((n+1)\omega)\big{)}+\frac{1}{4}\sin((n+1)\omega)\cdot\\\
\left.\big{(}\widetilde{F}_{n}(\omega)(6+3\cos(\omega))+\cos((n+1)\omega)-1\big{)}\right)\end{split}$

and

$\begin{split}\widetilde{F}^{\prime}_{n}(\omega)\widetilde{F}^{\prime\prime}_{n}(\omega)=\frac{1}{3(n+1)\sin^{4}(\omega/2)}\left(-\frac{3}{2(n+1)}\cos(\omega/2)\sin((n+1)\omega/2)\cdot\right.\\\
U_{n}(\cos(\omega/2))\big{(}(2+\cos(\omega))\widetilde{F}_{n}(\omega)+3\cos((n+1)\omega)+2\big{)}+\frac{3}{4}\sin((n+1)\omega)\cdot\\\
\left.\big{(}\widetilde{F}_{n}(\omega)(4+3\cos(x))+\cos((n+1)\omega)\big{)}\right).\end{split}$

Then using the triangle inequality, we can obtain the estimate for
$\widetilde{J}^{\prime\prime\prime}_{N}$, that is

$|\widetilde{J}^{\prime\prime\prime}_{N}(\omega)|\leq\frac{50.5\cdot\pi^{4}}{(n+1)|w|^{4}}.$

With $\widetilde{J}^{\prime\prime\prime}_{N}(\pi)=0$, we get for
$|\omega|\leq\frac{\pi}{4(n+1)}$,

$|\widetilde{J}^{\prime\prime\prime}_{N}(\omega)|\leq\widetilde{J}^{(4)}_{N}(0)|\omega|.$

The fourth derivative can be written as

$\widetilde{J}^{(4)}_{N}(\omega)=2\Big{(}3(\widetilde{F}_{n}(\omega))^{2}+4\widetilde{F}^{\prime}_{n}(\omega)\widetilde{F}^{\prime\prime\prime}_{n}(\omega)+\widetilde{F}_{n}(\omega)\widetilde{F}^{(4)}_{n}(\omega)\Big{)}$

and therefore we obtain

$\widetilde{J}^{(4)}_{N}(0)=2\Big{(}3(\widetilde{F}_{n}(0))^{2}+\widetilde{F}^{(4)}_{n}(0)\Big{)}.$

Now, since for $\omega=0$ one has

$\widetilde{F}^{(4)}_{n}(0)=\frac{2}{n+1}\sum\limits_{k=1}^{n}\left(k^{4}-\frac{k^{5}}{n+1}\right)=\frac{(n(n+2))}{30}(2n(n+2)-1),$

we easily get the value

$\widetilde{F}^{(4)}_{n}(0)=\frac{1}{30}n(n+2)(9n(n+2)-2).$

∎

##  Appendix B Proof of Lemma 3.0.3

###### Proof.

For $x\in\mathbb{S}^{2}$, with $d(x,x_{j})\leq\varepsilon\frac{\nu}{n+1}$ for
some $x_{j}\in\mathcal{X}$, we define the ring around $x\in\mathbb{S}^{2}$ by

$S_{m}:=\left\\{y\in\mathbb{S}^{2}\colon\frac{\nu m}{n+1}\leq
d(x,y)\frac{\nu(m+1)}{n+1}\right\\}$

for $m\in\mathbb{N}$. As it have been shows in [16], we can estimate the
number of elements in the intersection of $S_{m}$ with the set
$\mathcal{X}\setminus\\{x_{j}\\}$ for $m\geq 1$ by

$\mathrm{card}(\mathcal{X}\setminus\\{x_{j}\\}{\cap}S_{m})\leq 25m.$

Thus, it remains to estimate the number of elements in
$\mathcal{X}\setminus\\{x_{j}\\}{\displaystyle\cap}S_{0}$. We are going to use
the same argument as in [16], namely, we see that for
$x_{i},x_{n}\in\mathcal{X}\in\\{x_{j}\\}{\displaystyle\cap}S_{0}$ one has
${B_{\frac{\nu}{2(n+1)}}(x_{i}){\displaystyle\cap}B_{\frac{3\nu}{2(n+1)}}(x_{n})=\emptyset}$
and

$\displaystyle\cup_{x_{i}\in\\{x_{j}\\}{\cap}S_{0}}B_{\frac{\nu}{2(n+1)}}(x_{i})\subseteq
B_{\frac{3\nu}{2(n+1)}}(x_{j}).$

Since $\varepsilon\leq 1/2$ and the Riemannian volume form is rotation
invariant, we can bound the number of elements by

$\mathrm{card}(\mathcal{X}\setminus\\{x_{j}\\}{\displaystyle\cap}S_{0})\leq
2\frac{\Omega(B_{\frac{3\nu}{2(n+1)}}(e_{3}))}{\Omega(B_{\frac{\nu}{2(n+1)}}(e_{3}))}$

where $e_{3}=(0,01)^{\mathrm{T}}$ is the north pole on the sphere. In polar
coordinates around $e_{3}$, we consequently have the bound

$\begin{split}\mathrm{card}(\mathcal{X}\setminus\\{x_{j}\\}{\displaystyle\cap}S_{0})&\leq\frac{\Omega(B_{\frac{3\nu}{2(n+1)}}(e_{3}))}{\Omega(B_{\frac{\nu}{2(n+1)}}(e_{3}))}=\frac{\displaystyle\int_{0}^{\frac{3\nu}{2(n+1)}}\sin(r)\mathrm{d}r}{\displaystyle\int_{0}^{\frac{\nu}{2(n+1)}}\sin(r)\mathrm{d}r}\\\
&=\frac{1-\cos{\left(\frac{3\nu}{2(n+1)}\right)}}{1-\cos{\left(\frac{\nu}{2(n+1)}\right)}}=\left(1+2\cos{\left(\frac{\nu}{2(n+1)}\right)}\right)^{2}\leq
9.\end{split}$

Due to $d(x,x_{j})\leq\varepsilon\frac{\nu}{n+1}$, we have
$d(x,x_{i})\geq\frac{(1-\varepsilon)\nu}{n+1}$ for
$x_{i}\in\mathcal{X}\setminus\\{x_{j}\\}{\displaystyle\cap}S_{0}$. Using this
and the locality assumption (3.10), we can therefore estimate

$\begin{split}\sum\limits_{x_{i}\in\mathcal{X}\setminus\\{x_{j}\\}}|f(x,x_{i})|&\leq\sum\limits_{x_{i}\in\mathcal{X}\setminus\\{x_{j}\\}{\displaystyle\cap}S_{0}}\frac{c_{f}}{((n+1)d(x,x_{i}))^{s}}\\\
&+\sum_{m=1}^{\infty}\sum\limits_{x_{i}\in\mathcal{X}\setminus\\{x_{j}\\}{\displaystyle\cap}S_{m}}\frac{c_{f}}{((n+1)d(x,x_{i}))^{s}}\\\
&\leq\frac{9c_{f}(1-\varepsilon)^{-s}}{\nu^{s}}+25c_{f}\sum_{m=1}^{\infty}\frac{m}{(m\nu)^{s}}\\\
&\leq\frac{9c_{f}(1-\varepsilon)^{-s}}{\nu^{s}}+\frac{25c_{f}}{\nu^{s}}\sum_{m=1}^{\infty}\frac{1}{(m)^{s-1}}\\\
&\leq\frac{(9(1-\varepsilon)^{-s}+25)c_{f}\zeta(s-1)}{\nu^{s}},\end{split}$

where the last inequality follows by the definition of the Zeta function. On
the other hand, we can define the ring around $x_{j}$ again by

$\widetilde{S}_{m}=\left\\{y\in\mathbb{S}^{2}:\frac{(1-\varepsilon)\nu
m}{n+1}\leq d(x_{j},y)\leq\frac{(1-\varepsilon)\nu(m+1)}{n+1}\right\\}.$

Since $d(x,x_{j})\leq\varepsilon\frac{\nu}{n+1}$, we have
$d(x,x_{j})\leq\varepsilon d(x_{i},x_{j})$ for
$x_{i}\in\mathcal{X}\setminus\\{x_{j}\\}{\displaystyle\cap}\widetilde{S}_{m}$
and therefore $d(x,x_{i})\geq
d(x_{i},x_{j})-d(x,x_{j})\leq\frac{(1-\varepsilon)\nu m}{n+1}$. Using this and
the locality assumption (3.10), for $s\geq 3$ we can estimate the sum as

$\begin{split}\sum\limits_{x_{i}\in\mathcal{X}\setminus\\{x_{j}\\}}|f(x,x_{i})|&\leq\sum_{m=1}^{\infty}\sum\limits_{x_{i}\in\mathcal{X}\setminus\\{x_{j}\\}{\displaystyle\cap}\widetilde{S}_{m}}\frac{c_{f}}{((n+1)d(x,x_{i}))^{s}}\\\
&\leq 25c_{f}\sum_{m=1}^{\infty}\frac{m}{(1-\varepsilon)^{s}(m\nu)^{s}}\\\
&\leq\frac{25c_{f}}{(1-\varepsilon)^{s}\nu^{s}}\sum_{m=1}^{\infty}\frac{1}{m^{s-1}}=\frac{25c_{f}\zeta(s-1)}{(1-\varepsilon)^{s}\nu^{s}}.\end{split}$

Afterwards, we choose
$a_{\varepsilon}=\zeta(s-1)\cdot\min\\{9\cdot(1-\varepsilon)^{-s}+25,25\cdot(1-\varepsilon)^{-s}\\}$
that finishes the proof. ∎

##  Appendix C  Properties of the Cross Product

Let us consider two vectors $a,b\in\mathbb{R}^{3}$, then the cross product is
given by the vector

$a\times b=\begin{pmatrix}a_{2}b_{3}-a_{3}b_{2}\\\ a_{3}b_{1}-a_{1}b_{3}\\\
a_{1}b_{2}-a_{2}b_{1}\end{pmatrix}.$

Then for $a,b,c,d\in\mathbb{R}^{3}$ we have the following relations between
the dot product ”$\cdot$” and the cross product ”$\times$”:

  * (i)

$a\times b=-(b\times a)$,

  * (ii)

$a\cdot(b\times c)=b\cdot(c\times a)=c\cdot(a\times b)$,

  * (iii)

$a\times(b\times c)=b(a\cdot c)-c(a\cdot b)$,

  * (iv)

$(a\times b)\times(a\times c)=(a\cdot(b\times c))a$,

  * (v)

$(a\times b)\cdot(c\times d)=(a\cdot c)(b\times d)-(a\cdot d)(b\times c)$.

## References

  * [1] de Prony B. G. R. (1795). _Essai expérimental et analytique sur les lois de la dilatabilité des fluides élastiques et sur celles de la force expansive de la vapeur de l’eau et de la vapeur de l’alcool á différentes températures._ J. Éc. Polytech, 1 (22): 24–76.
  * [2] Adams, J. (1878). _On the expression of the product of any two Legendre’s coefficients by means of a series of Legendre’s coefficients._ Proc. R. Soc. Lond., 27:63–71.
  * [3] Bernstein, S. (1912). _On the best approximation of continuous functions by polynomials of a given degree._ Comm. Soc. Math. Kharkow, Ser, 2(13):49–194.
  * [4] Milnor, J. (1959). _Differentiable manifolds which are homotopy spheres._ Collected papers of John Milnor, III. Differential topology, 65-88.
  * [5] Cohen-Tannoudji C., Diu B., Laloe F., _Quantum Mechanics_ , 2nd ed. New York: Wiley, 1977, vol. 1.
  * [6] Schmidt, R. (1986). _Multiple emitter location and signal parameter estimation._ IEEE transactions on antennas and propagation, 34(3):276-280.
  * [7] Roy, R., Kailath, T. (1989). _ESPRIT-estimation of signal parameters via rotational invariance techniques._ IEEE Transactions on acoustics, speech, and signal processing, 37(7), 984-995.
  * [8] Healey, D. M., Rockmore, D. N., Moore, S. B. (1996, May). _An FFT for the 2-sphere and applications._ In 1996 IEEE International Conference on Acoustics, Speech, and Signal Processing Conference Proceedings (Vol. 3, pp. 1323-1326). IEEE.
  * [9] Arridge, S. R. (1999). Optical tomography in medical imaging. Inverse problems, 15(2), R41.
  * [10] Meyer, J. (2001). _Beamforming for a circular microphone array mounted on spherically shaped objects._ The Journal of the Acoustical Society of America, 109(1), 185-193.
  * [11] Taguchi, K., Zeng, G. L., Gullberg, G. T. (2001). _Cone-beam image reconstruction using spherical harmonics._ Physics in Medicine & Biology, 46(6), N127.
  * [12] Comaniciu, D., Meer, P. (2002). _Mean shift: A robust approach toward feature space analysis._ IEEE Transactions on pattern analysis and machine intelligence, 24(5), 603-619.
  * [13] Meyer, J., Agnello, T. (2003, October). _Spherical microphone array for spatial sound recording._ In Audio Engineering Society Convention 115. Audio Engineering Society.
  * [14] Ramamoorthi, R., Hanrahan, P. (2004)._A signal-processing framework for reflection._ ACM Transactions on Graphics (TOG), 23(4), 1004-1042.
  * [15] Dumitrescu, B. (2007). _Positive Trigonometric Polynomials and Signal Processing Applications._ Springer.
  * [16] Keiner, J., Kunis, S., and Potts, D. (2007). _Efficient reconstruction of functions on the sphere from scattered data._ Journal of Fourier Analysis and Applications, 13(4):435-458.
  * [17] Filbir, F., Themistoclakis, W. (2008). _Polynomial approximation on the sphere using scattered data._ Mathematische Nachrichten, 281(5), 650-668.
  * [18] Sloan, P. P. (2008, February). _Stupid spherical harmonics (sh) tricks._ In Game developers conference (Vol. 9, p. 42).
  * [19] Jarrett, D. P., Habets, E. A., Naylor, P. A. (2010, August). _3D source localization in the spherical harmonic domain using a pseudointensity vector._ In 2010 18th European Signal Processing Conference (pp. 442-446). IEEE.
  * [20] Audet, P. (2011). _Directional wavelet analysis on the sphere: Application to gravity and topography of the terrestrial planets._ Journal of Geophysical Research: Planets, 116(E1).
  * [21] Filbir, F., Mhaskar, H. N., Prestin, J. (2012). _On the problem of parameter estimation in exponential sums._ Constructive Approximation, 35(3), 323-343.
  * [22] Jarosik, N., Bennett, C. L., Dunkley, J., Gold, B., Greason, M. R., Halpern, M., et al. Larson, D. (2011). Seven-year wilkinson microwave anisotropy probe (WMAP*) observations: sky maps, systematic errors, and basic results. The Astrophysical Journal Supplement Series, 192(2), 14.
  * [23] Rauhut, H., Ward, R. (2011). Sparse recovery for spherical harmonic expansions. arXiv preprint arXiv:1102.4097.
  * [24] Atkinson, K. and Han, W. (2012). _Spherical Harmonics and Approximations on the Unit Sphere: An Introduction. Number 2044 in Lecture Notes in Mathematics._ Springer.
  * [25] De Castro, Y., Gamboa, F. (2012). _Exact reconstruction using Beurling minimal extrapolation._ Journal of Mathematical Analysis and applications, 395(1), 336-354.
  * [26] Johansen-Berg, H., Behrens, T. E. (Eds.). (2013). _Diffusion MRI: from quantitative measurement to in vivo neuroanatomy._ Academic Press.
  * [27] Candéz E.J., Fernandez-Granda C. (2013) _Super-Resolution from Noisy Data._ J. Fourier Anal. Appl., 19:1229-1254
  * [28] McEwen, J. D., Puy, G., Thiran, J. P., Vandergheynst, P., Van De Ville, D., Wiaux, Y. (2013). _Sparse image reconstruction on the sphere: implications of a new sampling theorem._ IEEE Transactions on image processing, 22(6), 2275-2285.
  * [29] Tang, G., Bhaskar, B., and Recht, B. (2013). _Sparse recovery over continuous dictionaries-just discretize._ In Asilomar Conference on Signals, Systems and Computers, 2013.
  * [30] Bredies, K. and Pikkarainen, H. (2013)._Inverse problems in spaces of measures_. ESAIM: Control, Optimisation and Calculus of Variations, 19(1):190–218.
  * [31] Potts, D., Tasche, M. (2013). Parameter estimation for multivariate exponential sums. Electron. Trans. Numer. Anal, 40(204-224), 94.
  * [32] Candéz E.J., Fernandez-Granda C. (2014). _Towards a Mathematical Theory of Super-resolution._ Comm. on Pure and Appl. Math., 67(6): 906–956,
  * [33] Bendory T., Dekel S., and Feuer A. ( 2015). _Exact Recovery of Dirac Ensembles from Projection onto Spaces of Spherical Harmonics._ Constr. Approx., 42(2):183 –207.
  * [34] Bendory T., Dekel S., and Feuer. A. (2015). _Super-Resolution on the Sphere Using Convex Optimization._ IEEE Transactions on Signal Processing, 63(9): 2253–2262.
  * [35] Duval, V., Peyré, G. (2015). _Exact support recovery for sparse spikes deconvolution._ Found. Comp. Math., 15:1315-1355.
  * [36] Duval, V. and Peyré, G. (2015). _Sparse spikes deconvolution on thin grids._ arXiv preprint arXiv:1503.08577.
  * [37] Kunis, S., Peter, T., Römer, T., von der Ohe, U. (2016). _A multivariate generalization of Prony’s method._ Linear Algebra and its Applications, 490, 31-47.
  * [38] De Castro, Y., Gamboa, F., Henrion, D., Lasserre, J. B. (2016). _Exact solutions to super resolution on semi-algebraic domains in higher dimensions._ IEEE Transactions on Information Theory, 63(1), 621-630.
  * [39] Fernandez-Granda, C. (2016)._Super-resolution of point sources via convex programming._ Information and Inference: A Journal of the IMA, 5(3):251-303.
  * [40] Filbir, F., Schröder, K. (2016). _Exact recovery of discrete measures from wigner d-moments._ arXiv preprint arXiv:1606.05306.
  * [41] Cuyt, A., Tsai, M. N., Verhoye, M., Lee, W. S. (2018)._Faint and clustered components in exponential analysis._ Applied Mathematics and Computation, 327, 93-103.
  * [42] Aubel, C., Stotz, D., Bölcskei, H. (2018). _A theory of super-resolution from short-time Fourier transform measurements._ Journal of Fourier Analysis and Applications, 24(1), 45-107.
  * [43] Josz, C., Lasserre, J. B., Mourrain, B. (2019). _Sparse polynomial interpolation: sparse recovery, super-resolution, or Prony?._ Advances in Computational Mathematics, 45(3), 1401-1437.
  * [44] Poon, C., Peyré, G. (2019). _Multidimensional sparse super-resolution._ SIAM Journal on Mathematical Analysis, 51(1), 1-44.
  * [45] Kunis, S., Möller, H. M., von der Ohe, U. (2019). _Prony’s method on the sphere._ The SMAI journal of computational mathematics, 5, 87-97.
  * [46] Prestin, J., Veselovska, H. (2020)._Prony-type polynomials and their common zeros._ Frontiers in Applied Mathematics and Statistics, 6, 16.
  * [47] Cuyt, A., Hou, Y., Knaepkens, F., Lee, W. S. (2020). _Sparse multidimensional exponential analysis with an application to radar imaging._ SIAM Journal on Scientific Computing, 42(3), B675-B695.
  * [48] Kunis, S., Römer, T., von der Ohe, U. (2020)._Learning algebraic decompositions using Prony structures._ Advances in Applied Mathematics, 118, 102044