# Points on manifolds with asymptotically optimal covering radius

Anna Breger University of Vienna, Department of Mathematics, Oskar-Morgenstern-Platz 1 A-1090 Vienna Martin Ehler University of Vienna, Department of Mathematics, Oskar-Morgenstern-Platz 1 A-1090 Vienna  and  Manuel Gräf University of Vienna, Department of Mathematics, Oskar-Morgenstern-Platz 1 A-1090 Vienna
###### Abstract.

Given a finite set of points on the Euclidean sphere, the worst case quadrature error in Sobolev spaces has recently been shown to provide upper bounds on the covering radius of the point set. Moreover, quasi-Monte Carlo integration points on the sphere achieve the asymptotically optimal covering radius. Here, we extend these results to points on compact smooth Riemannian manifolds and provide numerical experiments illustrating our findings for the Grassmannian manifold.

###### Key words and phrases:
covering radius, worst case integration error, Riemannian manifold, cubature points, quasi-Monte Carlo integration

Corresponding author]martin.ehler@univie.ac.at

## 1. Introduction

Many discretization schemes in numerical analysis are based on finite samples that cover the underlying space sufficiently well, i.e., the sampling points have small covering radius. One way of measuring the covering’s efficiency is by its cardinality in comparison to its covering radius.

Quasi-Monte Carlo integration points have been investigated in [7] for compact smooth Riemannian manifolds. In the special case of the Euclidean sphere, it is shown recently in [9] that the worst case error of integration bounds the covering radius and that thereby quasi-Monte Carlo integration points provide asymptotically optimal covering radii.

Here, we extend these results from the sphere to compact smooth Riemannian manifolds. These theoretical results are derived by combining the ideas in [7] with those in [9].

In the second part of the present note, we numerically construct a sequence of quasi-Monte Carlo integration points for the Grassmannian manifold and illustrate numerically that their covering radii indeed behave in accordance to the theoretical findings, hence, asymptotically optimal.

Our quasi-Monte Carlo integration points are cubatures (in fact designs) in Grassmannians that have been studied in [2, 3, 4, 5] from a theoretical point of view, see [18] for the construction through numerical minimization. For related results on cubatures in more classical settings, see [17, 23, 24, 26, 28, 31, 32] and, for further related results, we refer to [10, 13, 14, 19, 20, 22, 27] and references therein.

The outline is as follows. In Section 2, we briefly discuss the concept of asymptotically optimal covering radii. Low-cardinality cubature points are introduced in Section 3, where we also state our main theoretical result. In Section 4 we state the bound on the covering radius by the worst case error of integration. Section 5 is dedicated to the proof of this bound. In Section 6, we illustrate our theoretical findings by numerical experiments for the Grassmannian manifold .

## 2. Optimal asymptotics of the covering radius

Let be a compact smooth -dimensional Riemannian manifold. We denote its normalized Riemannian measure by and its Riemannian distance by . Given any finite collection of points , the covering radius is

 ρ:=ρ({xj}nj=1):=supx∈Mmin1≤j≤ndist(x,xj).

For denoting the ball of radius centered at , the union covers completely. By compactness of we have We use the notation , meaning the left-hand side is less or equal to the right-hand side up to a positive constant factor. The symbol is used analogously, and means both hold, and . If not explicitly stated, the dependence or independence of the constants shall be clear from the context.

 (1) μ(Br(x))≍rd,x∈M,0

where the constants do not depend on or . Hence, the line of inequalities leads to the lower bound

 (2) n−1d≲ρ.
###### Definition 2.1.

Given a sequence of points , , with , we say that the corresponding sequence of covering radii is asymptotically optimal if the lower bound (2) is matched, i.e., if

According to [29], the expectation of the covering radius of random points , independently identically distributed according to , satisfies

 (3) Eρ≍n−1dlog(n)1d.

Hence, there is an additional logarithmic factor, so that random points do not provide optimal covering radii.

In this brief note, we shall verify that the recently introduced concept of quasi-Monte Carlo systems, cf. [7, 8], lead to point sets with asymptotically optimal covering radii. This generalizes results for the sphere in [9] to compact smooth Riemannian manifolds.

## 3. Optimal coverings from low-cardinality cubatures

Let be the collection of orthonormal eigenfunctions of the Laplace-Beltrami operator on with eigenvalues arranged by . We denote by , , the Banach space of complex-valued -measurable functions on , whose -th power of the absolute value is integrable (with the standard modifications when ).

The space of diffusion polynomials of bandwidth is

 (4) Πt:=span{φℓ:λℓ≤t2}.

For and weights , we say that is a cubature for if

 (5) ∫Mf(x)dμ(x)=n∑j=1ωjf(xj),for % all f∈Πt.

The number refers to the strength of the cubature.

Weyl’s estimates on the spectrum of an elliptic operator yield

 dim(Πt)≍td,

cf. [25, Theorem 17.5.3]. Therefore, any sequence of cubatures of strength must obey .

###### Definition 3.1.

We call a sequence of cubatures for satisfying

 (6) ni≍tdi

with a low-cardinality cubature sequence.

The above definition makes sense since there do exist sequences of cubatures of strength with positive weights and satisfying (6), cf. [16]. We now state our main theoretical result.

###### Theorem 3.2.

The covering radius of any low-cardinality cubature sequence with positive weights is asymptotically optimal.

The remaining part of the present paper is dedicated to prove Theorem 3.2 and to numerically illustrate our findings on the Grassmannian manifold. We conclude this section with a remark concerning the weights.

###### Remark 3.3.

Cubatures of strength , whose weights are all the same, , for , are also called -designs. If is the unit Euclidean sphere, then there are -designs satisfying (6), cf. [6]. By identifying with , the analogous statement holds for the projective space. For general , however, we only know that -designs exist, cf. [30], but it is still an open problem whether or not the asymptotics (6) can be achieved.

## 4. Worst case error of integration

To prove Theorem 3.2, we follow the approach for the sphere in [9]. We shall first introduce the worst case error of integration and shall check that it provides an upper bound on the covering radius. Next, we shall consider the concept of quasi-Monte Carlo points, i.e., points whose worst case error of integration decays sufficiently fast, so that the covering radius is asymptotically optimal. Finally, we shall recapitulate from [7] that low-cardinality cubature points with positive weights are indeed quasi-Monte Carlo points.

The worst case error of integration of points and weights with respect to some Banach space of complex-valued functions on is

 (7) wce({(xj,ωj)}nj=1,H):=supf∈H∥f∥H≤1∣∣∫Mf(x)dμ(x)−n∑j=1ωjf(xj)∣∣.

Although suppressed by our notation, (7) depends on the particular norm , which shall always be clear from the context in the present manuscript. For most parts, we take to be a Sobolev space, which we define next. The Fourier transform of with is

 ^f(ℓ):=∫Mf(x)¯¯¯¯¯¯¯¯¯¯¯¯¯φℓ(x)dμ(x),ℓ=0,1,2,…,

and extends to distributions on . The Sobolev space , for and , is the set of all distributions on with , i.e., with

 (8) ∥f∥Wsp:=∥(I+Δ)s/2f∥Lp=∥∞∑ℓ=0(1+λℓ)s/2^f(ℓ)φℓ∥Lp<∞.

Note that is contained in the space of continuous functions on provided that , cf. [7]. We shall stick to this range throughout the present paper.

It turns out that the covering radius is bounded by the worst case error. The following result has been derived in [9] for being the Euclidean sphere and constant weights .

###### Proposition 4.1.

Let and be given. Then for any points and weights with covering radius it holds

 (9) ρ≲[wce({(xj,ωj)}nj=1,Wsp(M))]1/(s+d/q),

where . The constants may only depend on , , and .

We shall postpone the proof of Proposition 4.1 to Section 5. At this point we turn to sequences of points whose worst case error of integration satisfies decay conditions, which connects to the covering radius via the bound (9). The following definition is due to [7, 8].

###### Definition 4.2.

Given and , a sequence , , of points and weights with is called a quasi-Monte Carlo (qMC) system for if

 (10) wce({(xij,ωij)}nij=1,Wsp(M))≲n−sdi.

According to [7], low-cardinality cubature sequences with positive weights are qMC systems:

###### Proposition 4.3 ([7]).

For and , any low-cardinality cubature sequence with positive weights is a qMC system for .

For , due to (8), the space is continuously embedded into with , for . Thus, if and is a qMC system for , then it is also a sequence of qMC systems for . If is the sphere, then the latter becomes [9, Theorem 4.2]. Therefore, is the strongest requirement among the qMC properties.

###### Remark 4.4.

If is a qMC system for some , then Proposition 4.1 yields that its covering radii are bounded by

 ρi≲n−1d(ss+d(1−1/p))i.

Thus, qMC systems for provide

 ρi≲[wce({(xij,ωij)}nij=1,Ws1(M))]1/s≲n−1di,

so that Theorem 3.2 is a consequence of Propositions 4.1 and 4.3. Hence, to complete the proof of Theorem 3.2, it only remains to verify Proposition 4.1, which is the topic of the subsequent section.

Before we proceed to the proof of Proposition 4.1, we shall discuss a method to compute the worst case error of integration in , the latter being a Hilbert space with inner product

 (11) ⟨f,g⟩Ws2=∞∑ℓ=0(1+λℓ)s^f(ℓ)¯¯¯¯¯¯¯¯¯¯^g(ℓ),f,g∈Ws2(M).

The Bessel kernel given by

 (12) KsB(x,y)=∞∑ℓ=0(1+λℓ)−sφℓ(x)¯¯¯¯¯¯¯¯¯¯¯¯φℓ(y)

is the reproducing kernel for with respect to the inner product (11) provided that . For later reference we consider a slightly more abstract setting. If is a reproducing kernel for some reproducing kernel Hilbert space of continuous functions on , then the worst case error of integration is

 (13) wce({(xj,ωj)}nj=1,HK)2 =n∑j,j′=1ωjωj′K(xj,xj′)−2n∑j=1ωj∫MK(xj,x)dμ(x) +∫M∫MK(x,y)dμ(x)dμ(y),

cf. [21, Theorem 2.7], see also [28]. Note, the norm in is uniquely determined by its reproducing kernel . If has the Fourier expansion

 (14) K(x,y)=∞∑ℓ=0^K(ℓ)φℓ(x)¯¯¯¯¯¯¯¯¯¯¯¯φℓ(y)

with , for , then (13) becomes

 (15) wce({(xj,ωj)}nj=1,HK)2=n∑j,j′=1ωjωj′K(xj,xj′)+^K(0)(1−2n∑j=1ωj),

where we assume without loss of generality . For the Bessel kernel , we observe with

We conclude this section by the following result on the worst case error of uniformly distributed random points, for which constant weights are the natural choice.

###### Proposition 4.5.

If is a reproducing kernel on and are random points on , independently identically distributed according to , then it holds

 √E[wce({(xj,1n)}nj=1,HK)2]=cn−12,

where

 (16) c2=∫MK(x,x)dμ(x)−∫M∫MK(x,y)dμ(x)dμ(y).

By applying (14), the constant (16) is . For the Bessel kernel , the condition implies that , so that on average qMC systems indeed perform better than random points for smooth functions. The proof of Proposition 4.5 follows from the lines in [8] when replacing the sphere by . In fact, Proposition 4.5 is already contained in [21, Corollary 2.8], see also [28].

## 5. Proof of Proposition 4.1

The proof of Proposition 4.1 relies on findings in [7]. We recapitulate the following localization result, which is one of the key ingredients for the proof of Proposition 4.1.

###### Lemma 5.1.

For and being a smooth function, supported in the interval , the kernel

 Ksr(x,y):=∞∑ℓ=0h(λℓr)(1+λℓ)−s/2φℓ(x)¯¯¯¯¯¯¯¯¯¯¯¯φℓ(y)

is bounded by

 (17) ∥∥∫M|Ksr(⋅,y)|dμ(y)∥∥L∞≲r−s,r>0,

where the constant does not depend on .

###### Proof.

According to [7, Lemma 2.8], the estimate

 |Ksr(x,y)|≲rd−s(1+rdist(x,y))−(d+1),r>0,

holds and leads to the requested assertion

 supx∈M∣∣∫M|Ksr(x,y)|dμ(y)∣∣ ≲rd−ssupx∈M∫M(1+rdist(x,y))−(d+1)dμ(y) ≲r−ssupx∈M∫Mrd(1+rdist(x,y))−(d+1)dμ(y) ≲r−s.\qed

We shall make use of Lemma 5.1 to verify the following result.

###### Lemma 5.2.

Let , and be fixed. For any and , there is a function with support in , such that

 (18) ∥f∥Wsp≲r−s+d/p,rd≲∫Mf(x)dμ(x),

where the constants do not depend on or .

###### Proof.

As in [7, Proof of Theorem 2.16], for any radius and , there is a function with support in , such that, for ,

 (19) ∥Δℓf∥L∞≲r−2ℓ,rd≲∫Mf(x)dμ(x),

where by compactness of the constants do not depend on or .

Similarly as in [7, Proof of Theorem 2.16], we bound the norm by using a dyadic decomposition of unity of the Fourier domain, see, for instance, [32]. That is, we set for some smooth function satisfying

 g(x)={1,x≤1,0,x≥2,

and obtain with

 ∞∑m=−∞h(2−mx)={1,x>0,0,else.

Using the Fourier expansion we arrive at

 (I+Δ)s/2f=^f(0)+∞∑ℓ=0∞∑m=−∞h(2−mλℓ)(1+λℓ)s/2^f(ℓ)φℓ.

Fixing an integer and applying the notation of Lemma 5.1 with lead to

 (I+Δ)s/2f =^f(0)+∑2m≤1/r∫MK−s2m(⋅,y)f(y)dμ(y). +∑2m>1/r∫MK2L−s2m(⋅,y)(I+Δ)Lf(y)dμ(y).

We shall now bound the three terms of the right hand side separately. First, the Hölder inequality with (19) for yields

 |^f(0)|=∣∣∫Mf(y)dμ(y)∣∣≤∥f∥L∞μ(M)≲1≲r−s,

where the last estimate is due to being bounded from above and . To bound the second term, we apply Lemma 5.1 and derive

 ∥∑2m≤1/r∫MK−s2m(⋅,y)f(y)dμ(y)∥L∞ ≤∑2m≤1/r∥∫M|K−s2m(⋅,y)|dμ(y)∥L∞∥f∥L∞ ≲∑2m≤1/r∥∫M|K−s2m(⋅,y)|dμ(y)∥L∞ ≲∑2m≤1/r2ms≲r−s.

Lemma 5.1 also leads to a bound on the third term by

 ∥∑2m>1/r∫MK2L−s2m(⋅,y)(I+Δ)Lf(y)dμ(y)∥L∞ ≲∑2m>1/r∥∫MK2L−s2m(⋅,y)(I+Δ)Lf(y)dμ(y)∥L∞ ≲∑2m>1/r∥∫MK2L−s2m(⋅,y)dμ(y)∥L∞∥(I+Δ)Lf∥L∞ ≲∑2m>1/r2−(2L−s)mr−2L ≲r2L−sr−2L=r−s,

where we have also applied (19) and being bounded from above.

Thus, we obtain

 (20) ∥(I+Δ)s/2f∥L∞≲r−s,

where the constants do not depend on or To cover the range , we recall that is supported on , so that the Hölder inequality with (1) yields

 (21) ∥(I+Δ)s/2f∥pLp≲r−psμ(Br(z))≲r−sp+d,

which concludes the proof. ∎

We are now prepared to complete the proof of our main theoretical result.

###### Proof of Proposition 4.1.

For a given point set with covering radius , let be a center of a maximal hole, i.e.,

 (22) ˚Bρ(z)∩{xj}nj=1=∅,

where denotes the interior of . Note that is bounded by the diameter of . Let be as in Lemma 5.2, i.e., and (18) holds with . Since must vanish outside of , (22) implies , for all . Thus, the definition of the worst case error of integration yields

 wce({(xj,ωj)}nj=1,Wsp(M)) ≥∣∣∫Mf(x)dμ(x)−0∣∣∥f∥Wsp≳ρdρ−s+d/p=ρs+d/q,

with , and the constant does not depend on or . ∎

## 6. Numerical experiments for the Grassmannian manifold

This section is dedicated to illustrate the results of the previous sections for the special case of the Grassmannian manifold, i.e., the collection of -dimensional linear subspaces in , which we identify with the set of orthogonal projectors on of rank , denoted by

 (23) Gk,m:={P∈Rm×m:P⊤=P,P2=P,trace(P)=k}.

Hence, the Grassmannian can be considered as a -dimensional submanifold of the Euclidean space . Moreover the Euclidean space induces a Riemannian metric, which in turn yields the canonical probability measure and the canonical geodesic distance on the Grassmannian denoted by and , respectively. In particular, the geodesic distance between is proportional to the -norm of the corresponding principal angles between the subspace associated to and . More precisely, it can be computed by

 (24) distk,m(P,Q)=√2√θ21+…+θ2k,

where and are the -largest eigenvalues of (counted with multiplicities). Note that the factor of accounts for the particular embedding (23) since then

 (25) distk,m(P,Q)=∥P−Q∥F+o(∥P−Q∥2F),P,Q∈Gk,m,

where is the Frobenius-norm of .

### 6.1. Cubature points

Theorem 3.2 tells us that low-cardinality cubature points with positive weights inherit asymptotically optimal covering radii. To illustrate this result, we first construct a sequence of cubature points. It is known that any collection of points satisfies

 (26) 1n2n∑j,j′=1trace(Pj,Pj′)i≥∫Gk,m∫Gk,mtrace(P,Q)idμk,m(P)dμk,m(Q),

for , cf. [5]. According to [11], equality in (26) yields a design of strength , see also [18]. Hence, (26) provides us with a simple approach to numerically compute cubature points by minimization and checking for equality.

Our numerical experiments shall focus on , which has dimension , so that for low-cardinality cubature sequences the number of cubature points must satisfy . Indeed, we have chosen

 ni=⌊13(i+1)2(1+i+12i2)⌋

and computed points , for , by a nonlinear conjugate gradient method on manifolds, cf. [1, 11, 21], such that

 ∣∣∫G2,4f(P)dμ2,4(P)−1nini∑j=1f(Pj)∣∣<10−7,

for all with . Although the worst case error of integration in may not be zero exactly, we shall refer to in the following simply as -designs.

### 6.2. Integration

In view of bounding the covering radius, we may want to provide numerical experiments on the worst case error of integration in Sobolev spaces for . However, determining the worst case error is a tough task in general. For and , on the other hand, we are dealing with reproducing kernel Hilbert spaces, in which we can invoke (15) provided that its reproducing kernel is numerically accessible.

Our first numerical experiments are about integration in , so that the worst case error is indeed given by (15). However, the infinite series of the Bessel kernel in (12) is numerically cumbersome, so that we consider the positive definite kernel

 K1(P,Q)=k1(trace(PQ)),P,Q∈G2,4,

where

 k1(r)=√(2−r)3+2r,0≤r≤2.

Due to a comparison to the Bessel kernel, cf. [7, 11], its reproducing kernel Hilbert space is the Sobolev space equipped with an equivalent norm, i.e., the two norms and are comparable. This implies

 wce({(Pj,1n)}nj=1,W7/22(G2,4)) ≍wce({(Pj,1n)}nj=1,HK1)

where the constants are independent of the point sets. According to Proposition 4.3, the error decays as , for low-cardinality cubature points with positive weights. Hence, we expect a linear behavior with slope in logarithmic error plots.

In a second numerical experiment on integration, we shall consider a second positive definite kernel given by

 K2(P,Q)=k2(trace(PQ)),P,Q∈G2,4,

where

 k2(r)=32exp(−(2−r)),0≤r≤2.

Its reproducing kernel Hilbert space satisfies

 (27) HK2⊂⋂s>2Ws2(G2,4).

According to Theorem 4.3, a low-cardinality cubature sequence with positive weights is a qMC system for any . Due to (27), we expect a super linear behavior of in logarithmic plots.

We now further specify the worst case error in and via (15).

###### Lemma 6.1.

The worst case errors in and are

 (28) wce({(Pj,1n)}nj=1,HK1)2 =1n2n∑j,j′=1K1(Pj,Pj′) −(2+7475√2−25log(1+√2)), (29) wce({(Pj,1n)}nj=1,HK2)2 =1n2n∑j,j′=1K2(Pj,Pj′)−32exp(−1)Shi(1),

respectively, where is the hyperbolic sine integral.

###### Proof.

In view of (15) it remains to compute the -th Fourier coefficients

 ^Ki(0)=∫G2,4Ki(P,Q)dμ2,4(P)dμ2,4(Q),i=1,2.

According to [15, Example 4.3], the orthogonal invariance of and with the variable transformation , where are the principal angles between and , yield

 ∫G2,4Ki(P,Q)dμ2,4(P)dμ2,4(Q)=∫1−1∫1|ξ−|ki(1+ξ−ξ+)dξ+dξ−,i=1,2.

The symmetry of the function leads to

 ^Ki(0)=14∫1−1∫1−1ki(1+ξ−ξ+)dξ+dξ−,i=1,2.

For , we arrive at

 ^K1(0) =14∫1−1∫1−1((1−ξ−ξ+)32+2ξ−ξ++2)dξ+dξ− =2+∫1−1(1+ξ−)52−(1−ξ−)5210ξ−dξ−.

The assertion (28) is then checked by a computer algebra system.

For , we obtain

 ^K2(0) =38∫1−1∫1−1exp(ξ−ξ+−1)dξ+dξ− =38exp(−1)∫1−1exp(ξ−)−exp(−ξ−)ξ−dξ−,

so that (29) follows immediately. ∎

Figure 1 shows logarithmic plots of

 wce({(Pij,1ni)}nij=1,HK1),wce({(Pj,1ni)}nij=1,HK2)

for the cubature points , , from Section 6.1. For comparison, it also depicts the worst case error of integration for random points, independently identically distributed according to , cf. Proposition 4.5. Indeed, we observe the superior integration quality of the -designs over random points. The theoretical results in Propositions 4.3 and 4.5 are in perfect accordance with the numerical experiment. The integration errors of the random points scatter around the expected integration error with rate in both cases, and . Cubature points achieve the optimal rate of for functions in the Sobolev space . Due to , we observe the super linear behavior in the logarithmic plots for .

To conclude this section on numerical integration, we point out that Lemma 6.1 provides analytic expressions for the worst case errors in particular reproducing kernel Hilbert spaces. If we go beyond Hilbert spaces, then numerically computing worst case errors in Sobolev spaces becomes very challenging as is illustrated by the following remark.

###### Remark 6.2.

Let denote the by matrix with two ones in the left upper diagonal entries and zeros elsewhere. We shall consider the integration error with respect to a sequence of low-cardinality cubatures for in . Without loss of generality we assume