Tight bounds on the mutual coherence of sensing matrices for Wigner D-functions on regular grids

Many practical sampling patterns for function approximation on the rotation group utilizes regular samples on the parameter axes. In this paper, we relate the mutual coherence analysis for sensing matrices that correspond to a class of regular patterns to angular momentum analysis in quantum mechanics and provide simple lower bounds for it. The products of Wigner d-functions, which appear in coherence analysis, arise in angular momentum analysis in quantum mechanics. We first represent the product as a linear combination of a single Wigner d-function and angular momentum coefficients, otherwise known as the Wigner 3j symbols. Using combinatorial identities, we show that under certain conditions on the bandwidth and number of samples, the inner product of the columns of the sensing matrix at zero orders, which is equal to the inner product of two Legendre polynomials, dominates the mutual coherence term and fixes a lower bound for it. In other words, for a class of regular sampling patterns, we provide a lower bound for the inner product of the columns of the sensing matrix that can be analytically computed. We verify numerically our theoretical results and show that the lower bound for the mutual coherence is larger than Welch bound. Besides, we provide algorithms that can achieve the lower bound for spherical harmonics.


Introduction
In many applications, the goal is to recover a function defined on a group, say on the sphere S 2 and the rotation group SO(3), from only a few samples [1][2][3][4][5]. This problem can be seen as a linear inverse problem with structured sensing matrices that contain samples of spherical harmonics and Wigner D-functions.
In the area of compressed sensing and sparse signal recovery, there are recovery guarantee results for random sampling patterns [4,5] based on proving Restricted Isometry Property (RIP) of the corresponding sensing matrix. Regular deterministic sampling patterns are, however, more prevalent in practice due to their easier deployment [1,2]. RIP based results cannot be used for analyzing deterministic sensing matrices, because it has been shown in [6,7] that verifying RIP is computationally hard. For deterministic sampling patterns, the mutual coherence is widely used as performance indicator, which measures the correlation between different columns in the sensing matrix. For the case of sparse recovery on S 2 and SO(3), the product of two Wigner D-functions and the product of two spherical harmonics appear in the coherence analysis.
The product of orthogonal polynomials, i.e., Legendre and Jacobi polynomials, is sought-after in mathematics, and it is related to hypergeometric functions. Several works to obtain a closed-form and simplified version of these products have been presented, for instance in [8][9][10][11][12]. In those articles, the authors attempt to derive a compact formulation of the product and represent it as a linear combination of a single orthogonal polynomial with some coefficients. This representation reveals interesting properties that can also be applied in quantum mechanics [13,14], geophysics [15], machine learning [16], and low-coherence sensing matrices [1][2][3]. In quantum mechanics, these coefficients are used to calculate the addition of angular momenta, i.e., the interaction between two charged particles, and such coefficients are called the Clebsch-Gordan coefficients or Wigner 3j symbols. Specifically, these coefficients appear when we want to determine the product of Wigner D-functions and spherical harmonics. Since Jacobi and Legendre polynomials are specifically used to express the Wigner D-functions and spherical harmonics, the product of those polynomials in terms of Wigner 3j symbols is also commonly used in angular momentum literature [13,14]. The product of spherical harmonics also emerges in the spatiospectral concentration or Slepian's concentration problem on the sphere, as discussed in [15]. In this case, the goal is to maximize the concentration of spectrum frequencies, i.e., spherical harmonic coefficients, given a certain area on the spherical surface. The problem is similar to finding relevant eigenvalues from a matrix that consists of the product of spherical harmonics. As a consequence, Wigner 3j symbols appear as a tool to analyze the problem. In the area of machine learning, Wigner D-functions are used for analyzing group transformations of the input to neural networks and to implement equivariant architectures for rotations. The authors in [16] develop Clebsch-Gordan nets to generalize and to improve the performance of spherical convolutional neural networks for recognizing the rotation of spherical images, 3D shape, as well as predicting energy of the atom. The contribution of this article includes the implementation of Wigner D-functions to perform the transformation of a signal in the Fourier domain and to tailor a representation of the product in terms of Clebsch-Gordan coefficients that will primarily support theoretical analysis of the networks.
In this work, we are interested in coherence analysis and sparse recovery tasks. It has been shown in [1] that a wide class of modular symmetric regular sampling patterns, like equiangular sampling, yield highly coherent sensing matrices and thereby perform poorly for signal recovery tasks. Besides, the coherence was shown to be affected by the choice of elevation sampling pattern independent of azimuth and polarization sampling patterns. It was numerically shown that for regular sampling points on the elevation for Wigner D-functions and spherical harmonics, the mutual coherence is lower bounded by the inner product of columns with zero orders and two largest degrees 1 , which are then equal to Legendre polynomials. This bound is not contrived because one can show that this bound is achievable by optimizing azimuth angle φ ∈ [0, 2π). Consequently, the resulting deterministic sampling points can be implemented into a real-system to carry out measurements on the spherical surface, as discussed in [3,17]. In this article, we confirm mathematically the numerical findings of [1]. Our proof relies on using results for angular momentum analysis in quantum mechanics and properties of Wigner 3j symbols. To the best of our knowledge, this work is the first to provide the coherence analysis of a sensing matrix using the tools from angular momentum in quantum mechanics.

Related Works
The construction of sensing matrices from a set of orthogonal polynomials is widely studied in the area of compressed sensing. For instance, in [18] the authors show that the sensing matrix construction from random samples of Legendre polynomials with respect to Chebyshev measure fulfils the RIP condition and thus performs a robust and stable recovery guarantee to reconstruct sparse functions in terms of Legendre polynomials by using 1 -minimization algorithm. The extensions for random samples spherical harmonics and Wigner D-functions are discussed in [1,4,5,19]. The key idea in those articles bears a strong resemblance to the technique discussed in Legendre polynomials' case, which is based on carefully choosing random samples with respect to different probability measures and preconditioning techniques to keep the polynomials uniformly bounded. Despite the recovery guarantees with regard to the minimum number of samples, it has been discussed in [20] that these theoretical results is seemingly too pessimistic. Practically, one can consider a smaller number of samples and still achieve a very good reconstruction by using 1 -minimization algorithm. Therefore, a gap between theoretical and practical settings exists. Another concern in antenna measurement system is designing a smooth trajectory for robotic arms to acquire electromagnetic fields, which causes a practical obstacle in using random samples, as mentioned in [3,17].
One of the most prevalent application of orthogonal polynomials is in the area of interpolation, where those polynomials are used to approximate a function within a certain interval. Recently, the 1 -minimization-based technique is tailored to interpolate a function which has sparse representation in terms of Legendre and spherical harmonic coefficients, as discussed in [19]. In this case, random samples of Legendre and spherical harmonics are used to construct sensing matrices. Similar to the results in compressed sensing, the RIP plays a pivotal role in showing that a particular number of samples is required to achieve certain error approximation.
Another construction of sensing matrices from some orthogonal polynomials related to the sparse polynomial chaos expansion is investigated in [21]. In order to optimize the sensing matrices, the authors adopt the minimization of coherence sensing matrices from several random samples of Legendre and Hermite polynomials. Using Monte Carlo Markov Chain (MCMC), the authors also derive the coherence of optimal-based random sampling points to employ 1 -minimization algorithm. Additionally, they also derived the coherence bound for a matrix constructed from those polynomials.
In contrast to all aforementioned results which rely on random samples from certain probability measures, in this work, we derive the coherence bound for deterministic sampling points and utilize properties of Wigner d-functions and their products and their representation in terms of Wigner 3j symbols.

Summary of Contributions
In [1], it was conjectured that the lower-bound on the mutual coherence is tight. In other words, the inner product of columns with equal orders is dominated by the inner product of two columns with zero order and highest degree. In this paper, we prove this conjecture and derive a set of related corollaries. The main contributions and some of the interesting conclusions of our paper are as follows: • We show that the product of Wigner D-functions can be written as a linear combination of single Legendre polynomials and Wigner 3j symbols. For equispaced sampling points on the elevation and using symmetry of Legendre polynomials, we show in Section 2.3 that only even degree Legendre polynomials contribute in the analysis, which in turn simplifies the problem formulation.
• In Section 4 and 5, we provide various inequalities and identities for sum of Legendre polynomials and Wigner 3j symbols. We establish monotonic properties of these terms as a function of degree and orders of Wigner D-functions. These results establish a certain ordering between inner products of the columns of the sensing matrix. Particularly we show that, under some conditions, the inner products have a specific order as a function of degrees and orders. As a corollary, we also present that the inner product of columns of equal orders are decreasing with orders and increasing with degree. The result can be used to obtain a lower bound on the mutual coherence. Proofs of main theorem, supporting lemmas, and propositions are given in Section 8 and 9 • We numerically verify our results and show that our bound is larger than the Welch bound. Therefore, the desideratum of regular sampling pattern design should be this lower bound rather than the Welch bound. We also extend the sampling pattern design algorithm of [1] to a gradient-descent based algorithm and show that it can achieve the lower bound for spherical harmonics 2 . f (θ, φ, χ)g(θ, φ, χ)dν(θ, φ, χ).

Notation
Similar to the Fourier basis, which can be seen as the eigenfunctions of Laplace operator on the unit circle, we can also derive eigenfunctions on the rotation group SO(3). These functions are called Wigner D-functions, sometimes also called generalized spherical harmonics. They are orthonormal basis for L 2 (SO(3)). It can be written in terms of Euler angles θ ∈ [0, π], and φ, χ ∈ [0, 2π) as follows D k,n l (θ, φ, χ) = N l e −ikφ d k,n l (cos θ)e −inχ (1) where N l = 2l+1 8π 2 is the normalization factor to guarantee unit L 2 -norm of Wigner D-functions. The function d k,n l (cos θ) is the Wigner d-function of order −l ≤ k, n ≤ 1 and degree l defined by The function P (ξ,λ) α is the Jacobi polynomial. Wigner D-functions are equal to complex spherical harmonics when the order n is equal to zero, namely where Y k l (θ, φ) := N k l P k l (cos θ)e ikφ . The term N k l := 2l+1 4π (l−k)! (l+k)! is a normalization factor to ensure the function Y k l (θ, φ) has unit L 2 -norm and P k l (cos θ) is the associated Legendre polynomials. Wigner D-functions form an orthonormal basis with respect to the uniform measure on the rotation group, i.e., sin θdθdφdχ.
Similarly, Wigner d-functions are also orthogonal for different degree l and fix order k, n: The factor 2 2l+1 can be used for normalization to get orthonormal pairs of functions. As discussed earlier, the Wigner d-functions, or weighted Jacobi polynomials, are generalization of hypergeometric polynomials including associated Legendre polynomials, where the relationship between those polynomials can be expressed as The associated Legendre polynomials of degree l and order −l ≤ k ≤ l is given as where P l (cos θ) is Legendre polynomial, and using Rodrigues formula it can be written as Similar to Wigner d-functions, Legendre polynomials are also orthonormal. The important properties of associated Legendre polynomials that are necessary in this article are also presented in Supplementary Material in Section 10.2.

Problem Formulation
In many signal processing applications, it is desirable to study properties of matrices constructed from samples of Wigner D-functions.

4
Suppose that we take m samples of this function at points (θ p , φ p , χ p ) for p ∈ [m]. The samples are put in the vector y, and the goal is to find the coefficients g. This is a linear inverse problem formulated by y = Ag, with A given as For index column q ∈ [N ], we denotes the degree and orders of the respective Wigner D-function by l(q), k(q) and n(q) 3 . The column dimension of this matrix is given by . Using these functions, the elements of this matrix are given by In compressed sensing, the sensing matrix A with lower mutual coherence are more desirable for signal reconstruction [22]. The mutual coherence, denoted by µ(A) is expressed as where we adopt the following convention For the rest of the paper, we focus mainly on the inner product between the samples. It can be numerically verified that the 2 -norm of columns do not affect the coherence value. We will comment later on how these norms scale.
Although the closed-form derivation of mutual coherence is in general difficult, the authors in [1] observed empirically that the mutual coherence of sensing matrices with equispaced sampling points on the elevation angle is indeed tightly bounded by a single term under certain assumptions. This is because the inner products of Wigner D-functions are ordered in a regular way as a function of their orders and degrees. In this paper, we provide theoretical supports for these observations. In other words, we provide simple analysis of mutual coherence for specific sampling patterns. Central to our analysis is a set of combinatorial identities about the sum and product of Wigner D-functions. We will focus on the equispaced sample on θ p for p ∈ [m], which are chosen such that This means that −1 = cos θ 1 < cos θ 2 < . . . < cos θ m−2 < cos θ m−1 < cos θ m = 1. There are multiple reasons for using this sampling pattern. First of all, this sampling pattern has been shown to be beneficial in spherical near-field antenna measurement [2,3] where the robotic probe can acquire the electromagnetic field samples and move in the same distance. Second of all, this sampling pattern induces orthogonal columns in the sensing matrix between even and odd degree polynomials as discussed in [1,Theorem 5]. Interestingly, fixing the sampling patterns on the elevation imposes a lower bound on the mutual coherence, which is tight in many cases. In this paper, we study the mutual coherence for this elevation sampling and arbitrary sampling patterns on φ and χ.

Product of Wigner D-functions
In the expression for mutual coherence, product of Wigner D-functions appears constantly. This product appears also in the study of angular momentum in quantum mechanics and can be written as linear combination of a single Wigner D-functions with coefficients, called Wigner 3j symbols [13,14,23]. Using this representation, the discrete inner product of Wigner D-functions can be simplified as follows.
The parameters of Wigner 3j symbols are non-zero only under certain conditions known as the selection rules. The selection rules state that Wigner 3j symbols • The summation of all k i should be zero: • Triangle inequality holds for l i 's: • The sum of all l i 's should be an integer l 1 + l 2 + l 3 ∈ N.
Note that, from the selection rules, if l 1 + l 2 +l is an odd integer, then l 1 l 2l 0 0 0 is zero. Further properties and the exact expression of the Wigner 3j symbol will be included in the Supplementary Material in Section 10.4. It is trivial to derive the product of same orders Wigner d-function, i.e., k 1 = k 2 = k and n 1 = n 2 = n, as follows An interesting conclusion of the above identities is that the sampling pattern affects the inner product through the sum of individual functions, for instance Legendre polynomials.

Main Results
The starting point of bounding the coherence is the following trivial inequality, which holds in full generality: This is obtained by choosing the column with equal orders of k and n, and using Definition 1. In other words, for any sampling pattern, regardless of the choice of φ p and χ p , the mutual coherence is lower bounded by merely choosing θ p . This indicates the sensitivity of mutual coherence to the sampling pattern on the elevation. The following theorem shows that the maximum inner product in above expression has a simple solution for the sampling pattern (12) if m is moderately large.
Intuitively, this theorem states that if one considers equispaced samples on the elevation, the maximum inner product occurs at the zero-order Wigner d-functions, which are Legendre polynomials. Additionally, we can obtain similar results for maximum discrete inner product of associated Legendre polynomials, since for k = 0 or n = 0 the Wigner d-functions are the associated Legendre polynomials. The result is given in the following corollary.
where C k l = (l−k)! (l+k)! is the normalization factor.
Corollary 1 implies that the maximum product of two associated Legendre polynomials for different degrees and same orders is also attained at degrees B − 1 and B − 3 and order k = 0.
Remark 1. A byproduct of Theorem 1 is that for a fixed number of measurement numbers m, the inner product of Wigner d-functions with degree less that √ 10m are ordered. This surprising behavior in presented in numerical experiments.
The proof of the main result follows from a sequence of inequalities and identities. In the next sections, we provide some of them that are of independent interest. Interestingly, the well-ordered behaviour of the inner products of Wigner d-functions is linked to orders in the summand of (15). The proof leverages mostly classical inequalities and identities, e.g., Abel partial summation, and the orders between Wigner d-functions. All the proofs for the following sections are given in Section 9.

Finite Sum of Legendre Polynomials
The starting point of the proof is to use (15) for the product of Wigner d-functions. The proof strategy is based on establishing inequalities for each term in the sum (15), and then using them to bound the final sum. In this section, we provide a set of results for sum of Legendre polynomials. The following lemma provides an identity for the sum of equispaced samples of Legendre polynomials. Lemma 1. Suppose we have equispaced samples as in (12), then the sum of sampled Legendre polynomials for even degrees l > 0 is given by is a zeta function. For odd degrees l, the sum is equal to zero.
Lemma 1 shows that we can simplify the summation of equispaced samples Legendre polynomials with respect to number of samples m, degree of polynomials l, and the residual R l (m).
Remark 2. If l = 0, then the sum of equispaced samples Legendre polynomials is equal to m, since P 0 (cos θ) = 1. In addition, for l = 2, the summation is equal to 1 + 1 m−1 and we do not have any residual. If we take a number of samples larger than the degree l, it is obvious that the summation converges to 1.
The residual R l (m) is important in the summation in Lemma 1. We provide upper and lower bounds on the residual R l (m). The next proposition provides a bound on this summation. Proposition 2. Suppose we have m ≥ (l+1) 2 10 + 1 with even degree l ≥ 4. The residual R l (m) is therefore bounded by −0.463 < R l (m) < 0.
The previous proposition shows that the residual, conditioned on m ≥ (l+1) 2 10 +1, is inside the interval, and most importantly is always negative. We present the numerical evaluation of this bound in Section 6. Using this property, not only that the summation m p=1 P l (cos θ p ) is non-negative, but also it is monotonically increasing for an increasing even degree l.

Inequalities for Wigner 3j Symbols
To prove the main theorem, we establish that there is a similar order between Wigner 3j symbols. Some of these properties of Wigner 3j symbols are given in Section 2 and in Supplementary Material in Section 10. 4. In what follows, we will have some combinatorial identities and inequalities related to Wigner 3j symbols. Despite ample investigation of authors, it is not clear whether these results bear interesting implications for other areas particularly angular momentum analysis in quantum physics. For compressed sensing, however, these are quite interesting as they show that the sensing matrix from samples Wigner D-functions possesses a lot of structures and symmetries. We start with the following lemma.
By reducing the indices one at a time, this result shows that the maximum of Wigner 3j symbols for zero order is achieved at zero degree. Using property in (14), one can directly obtain that the maximum is equal to 1. The summation of Wigner 3j symbols for k 1 = k 2 is also essential to prove the main result and can be decomposed as the summation of odd and even degrees, as given below.
Lemma 4. Suppose that l 1 = l 2 ∈ N and |l 1 − l 2 | ≤l ≤ l 1 + l 2 , then for −min(l 1 , l 2 ) ≤ k = n ≤ min(l 1 , l 2 ), we have l ,even Furthermore, for k = n = τ and 1 ≤ |τ | ≤ min(l 1 , l 2 ) we have l ,even Remark 3. In this lemma, we do not include the condition τ = 0 in the second property because it is obvious that we have l 1 l 2l 0 0 0 2 and from the selection rules in Section 2 the sum l 1 + l 2 +l should be an even integer, which means that for evenl and odd l 1 + l 2 , the Wigner 3j symbols value is zero. From the orthogonal property of Wigner 3j symbols as in (14), we have l ,even (2l + 1) l 1 l 2l 0 0 0 2 = 1, for l 1 + l 2 is even. 8 The last result of this section is related to the product of Wigner 3j symbols andl(l + 1), where |l 1 − l 2 | ≤l ≤ l 1 + l 2 and l 2 = l 1 + 2. As discussed in Lemma 1, the sum of equispaced samples Legendre polynomials can be expressed as m p=1 P l (cos θ p ) = 1 + l(l+1) 6(m−1) + R l (m), where R l (m) is the residual. The following lemma gives an expression of the inner product between Wigner 3j symbols and l(l + 1).
Since we consider degree 0 ≤ l 1 < l 2 ≤ B − 1 and order k = n = 0, the previous lemma has a direct implication for l 1 = B − 3 and l 2 = B − 1, as given in the following corollary. This corollary is important to determine the product of Wigner 3j symbols and the sum of equispaced samples of Legendre polynomials. Since the latter can be expressed byl(l + 1), as shown in Lemma 1 , we can directly apply Lemma 5 to estimate the product.

Experimental Results
In what follows, we conduct a series of experiments to verify some of the results in the paper, as well as applications in compressed sensing.

Numerical verification of theoretical results
In Lemma 1, we can express the sum of equispaced samples Legendre polynomials as 1, where it can be seen that by considering m = (l+1) 2 10 + 1, represented with the black line, the residual is restricted within the interval by the blue and red colors, respectively. In other words, the obtained constants are indeed tight in Proposition 2.
Next, we numerically evaluate Theorem 1. In Figure 2, we show for which pairs of (m, B), the identity of Theorem 1 holds using the color red. We have furthermore included a black line indicating the number of samples as in Theorem 1 for different B. It can be seen for small B, the condition is tight. However, it seems that it can be improved for larger values of B. In Figure 2, the numerical experiments are performed by considering the normalization with respect to the 2 -norm. Thereby, we can numerically verify that the normalization does not affect the inequality in Theorem 1. In order to design sensing matrices from spherical harmonics and Wigner D-function, we choose points on azimuth and polarization angles φ, χ ∈ [0, 2π) using the optimization problem below: where f q,r (θ, φ, χ) is given as For spherical harmonics, one can generate the problem from a relation in (3) without constraint on the polarization angle χ and order n.
This optimization problem is a challenging min-max problem with non-smooth objective function and generally non-convex. In [1], we have used a search-based method for optimization, which turns out to be difficult to tune and more time-consuming. We introduce a relaxation of the above problem and use gradient-descent based algorithms for optimizing it. Derivative of spherical harmonics and Wigner D-functions are given in Section 10.1.1.
Using property of p -norm, we can write the objective function as One can choose large enough p and calculate the gradient. Therefore, we use gradient descent algorithms to solve this problem, as given in Algorithm 1. 10 Figure 3 shows the coherence of a sensing matrix from Wigner D-functions using sampling points generated from several stochastic gradient algorithms. Although there are no sampling points that reach the lower bound in Theorem 1, it can be seen that Adam algorithm [26] yields the best sampling points.
Bandwidth of Wigner D-functions is B = 4 or the column dimension N = 84. Algorithm 1 can be tailored for spherical harmonics, by only considering parameter on azimuth φ and using relation between Wigner D-functions and spherical harmonics in (3). It can be seen that most of the gradient descent based algorithms converge to the lower bound for spherical harmonics, as shown in Figure 4. Therefore, we can provide sampling points on the sphere with mutual coherence that can achieve the the lower bound in Theorem 1. In this case, the spherical harmonics are generated with the bandwidth B = 10 or equivalently we have the column dimension N = 100. 10

Conclusions and Future Works
We have established a coherence bound of sensing matrices for Wigner D-functions on regular grids. This result also holds for spherical harmonics, which is a special case for Wigner D-functions. Estimating coherence involves non-trivial and complicated product of two Wigner D-functions for all combination of degrees and orders, yielding an obstacle to derive a simple and compact formulation of the coherence bound. Using the tools in the area of angular momentum in quantum mechanics disentangles this problem and represents the product as a linear combination of single Wigner D-functions and angular momentum coefficients, so called Wigner 3j symbols. In this paper, we derive some interesting properties of these coefficients and finite summation of Legendre polynomials to obtain the coherence bound. We also provided numerical experiments in order to verify the tightness of this bound. For practical application, it is also necessary to provide sampling points on the sphere and the rotation group that can achieve the coherence bound. We have shown that, for spherical harmonics, one can generate points to achieve this bound by using class of gradient descent algorithms. Convergence analysis of these algorithms and construction of deterministic sampling points to achieve this bound will be relegated for future works.
Additionally, the product of Legendre polynomials in left hand side can be written as Pl(cos θ p ).
The strategies are divided into two parts, which are for order k = n = 0 and for other conditions of order k, n.
Let consider the first case, for zero order k = n = 0. From (6), we know that d 0,0 l (cos θ) = P l (cos θ). Thus, it is equivalent to prove the maximum product of two Legendre polynomials is attained at l 1 = B−3 and l 2 = B − 1, i.e., max l1 =l2 m p=1 P l1 (cos θ p )P l2 (cos θ p ) = m p=1 P B−1 (cos θ p )P B−3 (cos θ p ) . We show that for an even l 1 + l 2 , if (l 1 , l 2 ) is increased to either (l 1 + 1, l 2 + 1) or (l 1 + 2, l 2 ), the sum of the product of two Legendre polynomials increases. It is enough to consider these two situations since from any pair (l 1 , l 2 ), one can use a sequence of inequalities to arrive at (B − 3, B − 1). For an odd l 1 + l 2 , it is implied that the inner product is zero.
From Lemma 3 we have al ≥ bl for both cases and additionally l al = 1 and l bl = 1, as pointed out in (14). By using Abel's partial summation formula as stated in (52), we have where Al = l j,even=|l1−l2| a j and Bl are defined accordingly for bl. Since Al ≥ Bl and cl is increasing, it is clear that l alcl ≤ l blcl, which establishes desired result by increasing the degrees until their reach l 1 = B − 3 and l 2 = B − 1.
For k = n and k = n = 0, we want to show that the inequality also holds. Let us first define a variable for Wigner 3j symbols αl = (2l + 1) l 1 l 2l −n n 0 and we can define an independent variable κ = c |l 1 −l 2 | +c l 1 +l 2 2 . The upper bound of the product m p=1 d k,n l1 (cos θ p )d k,n l2 (cos θ p ) is given by The first inequality is derived by using the triangle inequality and the increasing property of the sum of equispaced samples of Legendre polynomials in Lemma 2, i.e., c 2 ≤ c 4 ≤ c 6 ≤ . . . ≤ c 2B−4 . The last inequality also holds by the Cauchy-Schwarz inequality of αl = (2l + 1) l 1 l 2l −n n 0 To be precise, for −min (l 1 , l 2 ) ≤ k = n ≤ min (l 1 , l 2 ) one can write l1+l2 l =|l 1 −l 2 | l,even (2l + 1) l 1 l 2l −n n 0  (23) and (24), it is enough to consider the upper bound in (24). We then need to show m p=1 We have the product of two Legendre polynomials where the equality is derived by using expansion of Legendre polynomials in (46) and substitute x = 1. For even degree l and k = 0, we have P l (1) = 1 and P l (0), respectively. Identity of P l (0) is derived in (47). Hence, the second sum is equal to 14 The final sum can be obtained as where we use the fact that m = m−1 2 . The remainder term of the summation, R l (m), can be expressed as For j ≥ 3, the Bernoulli number is non-zero only for even values of j, as discussed in Table (10.3.1). By changing the summation index, we have R l (m) as l j=4 j,even From (48), (49) and the relation between Bernoulli number and zeta function in (51) we have Now consider j = 2 for the last equation with the value of ζ(2) = π 2 6 , then the equation is equal to which completes the claim. The same approach can be used to derive the result for even m.

Proof of Proposition in Section 4
Proof of Proposition 2. From Lemma 1, the sum of equispaced samples Legendre polynomials can be written as where upper bound is derived from the fact that the zeta function is decreasing for an increasing even j, that is ζ(j) > ζ(j + 2). In order to show the decreasing property, it should be enough by showing that the ratio above is upper bounded by 1, which is accomplished by considering m − 1 ≥ (l+1) 2 10 for 4 ≤ j ,even ≤ l and l ≥ 4. Now, we want to show the lower bound of R l (m). For an even l 2 , we will have R l (m) as − It is trivial to show that (j−1)!(l−j+1)!(2π) j (m−1) j−1 converges to zero for sufficiently large samples m compared to the degree l. Hence, giving the upper bound of the residual.
Proof of Lemma 2. For l = 2, it is proven in Lemma 1 that m p=1 P l (cos θ p ) = 1 + l(l+1) 6(m−1) > 0. Therefore, we have m p=1 P 4 (cos θ p ) − m p=1 P 2 (cos θ p ) = (j−1)!(l−j+1)!(2π) j . Thus, we only need to prove for even l ≥ 4. For this reason, the increasing property of the summations, m p=1 P 2 (cos θ p ) < m p=1 P 4 (cos θ p ) < . . . < m p=1 P l,even (cos θ p ), directly implies a non-negative property of the sum of equispaced samples Legendre polynomials. Since we compare l +2 to l, then the number of sample m−1 ≥ (l+3) 2 10 should be considered, which means that for even l ≥ 4, it is enough to show m p=1 P l+2 (cos θ p ) − m p=1 P l (cos θ p ) ≥ 0. By using the result from Lemma 1, the condition is equal to 2l+3 Let observe the residual (m − 1) R l (m) − R l+2 (m) and write as l j=4 j,even First of all, we show that Second of all, the sequence (S j l+2 −S j l ) (m−1) j−2 is decreasing if we increase even j, or we have . Therefore, the expression is equivalent to showing that the ratio (S j+2 l+2 −S j+2 l ) (S j l+2 −S j l )(m−1) 2 < 1 holds for m − 1 ≥ (l+3) 2 10 . Let expand this ratio by using (31) as From (29), we know the ratio S j+2 l S j l (m−1) 2 . Hence, (32) can be expressed as The upper bound holds from the fact that ζ(j) > ζ(j + 2) and samples m − 1 ≥ (l+3) 2 10 . Thereby, it proves that for increasing even j, Summarizing the results in (31) and (33) we Therefore, for even l 2 , we write l j=4 j,even The upper bound follows because the subtractions in the brackets are positive. It should be noted that the upper bound also holds for an odd l 2 , where instead of having two terms for l − 2 and l, we only have l which is negative. Finally, collecting all these results, we can bound (30) with (m−1) l . We want to show this upper bound is smaller than 2l+3 3 . Let first observe the upper bound of and using the ratio in (31), we obtain < ζ(4)(24)10 2 3!(2π) 4 (l + 2)(l + 1)l(2l + 3) (l + 3) 3 where the constant C 1 = 0.2778 is derived from the fact that ζ(4) = π 4 90 as in Table ( (m−1) l can be determined by .

Proofs of Lemmas in Section 5
Proof of Lemma 3. In (55), we have the exact expression for Wigner 3j symbols . To have the above ratio be greater than one, the following condition should be satisfied l 2 ≥ l 1 + 1 > l 1 . This condition does not change the assumption in this chapter since we want to find the maximum for l 1 = l 2 .
As discussed earlier, if we choose θ = π, then we have l ,even (2l + 1) However, we know that the sum of squared Wigner 3j symbols for alll is 1, due to the orthogonal property in (14), i.e., For l 1 = 1, we have 2 + 2(l 1 + 2)(l 1 + 1) = 14. This summation becomes complicated since we have two different valuesl for the Wigner 3j symbols, From this relation, we can write