Access the full text.
Sign up today, get DeepDyve free for 14 days.
Abstract The construction of (near-)minimal cubature formulae on the disk is still a complicated subject on which many results have been published. We restrict ourselves to the case of radial weight functions and make use of a recent connection between cubature and the concept of multivariate spherical orthogonal polynomials to derive a new system of equations defining the nodes and weights of (near-)minimal rules for general degree |$m=2n-1, n \ge 2$|. The approach encompasses all previous derivations. The new system is small and may consist of only |$(n+1)^2/4$| equations when |$n$| is odd and |$n(n+2)/4$| equations when |$n$| is even. It is valid for general |$n$| and has a Prony-like structure. It may admit a unique solution (such as for |$n=3$|) or an infinity of solutions (such as for |$n=7$|). In Section 2, the new approach is described, whereas the new system is derived in Sections 3 and 4. All well-known (near-)minimal cubature rules can be reobtained. Some typical illustrations of how this works are given in Section 5. We expect that this unifying theory will shed new light on the topic of cubature, in particular with respect to the discovery of new bounds on the number of nodes and their connection with the zeros of multivariate orthogonal polynomials. 1. Introduction Let |${\it{\Omega}}$| denote a given region in |$\mathbb{R}^2$| and let |$w(x,y)$| be a non-negative weight function on |${\it{\Omega}}$|. A cubature rule \begin{equation} \iint_{{\it{\Omega}}}f(x,y)w(x,y)\,\,{\rm d}x\,\,{\rm d}y \approx \sum_{i=1}^N A_i\,f(x_i,y_i), \label{numG1} \end{equation} (1.1) is said to be of total degree |$m$|, if it is exact for |$f \in {\mathcal P}_m^2$|, where |${\mathcal P}_m^2$| denotes the space of polynomials of total degree at most |$m$| in the real variables |$x$| and |$y$|. With |$m=2n-1$| or |$m=2n-2$| where |$n \ge 1$|, the number of nodes |$N$| necessarily satisfies (Stroud, 1971, pp. 118–119) \begin{equation} N \ge {n(n+1) \over 2}. \label{trueG} \end{equation} (1.2) A cubature rule of degree |$m$| with |$N$| nodes is called a Gaussian cubature rule when |$N$| attains the lower bound in (1.2). These Gaussian cubature rules rarely exist. Two examples of domains and weight functions that allow Gaussian cubature rules are as follows: |${\it{\Omega}}=\left\{(x,y)\,:\, 1+x+y>0,\, 1-x+y>0,\, x^2>4y\right\}$| (Schmid & Xu, 1994) with $$w_{\alpha,\beta,\pm{1 \over 2}}(x,y)={(1-x+y)}^{\alpha}{(1+x+y)}^{\beta} {(x^2-4y)}^{\pm{1 \over 2}},$$ |${\it{\Omega}}=\left\{(x,y)\,:\, {-3{\left(x^2+y^2+1\right)}^2+8\left(x^3-3xy^2\right)+4}\ge 0\right\}$| (Li et al., 2008) with $$w_{{1 \over 2}}(x,y)={\left(-3{\left(x^2+y^2+1\right)}^2+8\left(x^3-3xy^2\right)+4\right)}^{{1 \over 2}}.$$ If the weight function is centrally symmetric, meaning |$w(x,y) = w(-x,-y)$|, then Gaussian cubature rules of degree |$m=2n-1$| of the form (1.1) do not exist (Xu, 2012). In this case, a stronger lower bound for |$N$| holds (Möller, 1976), namely \begin{equation} N \ge {n(n+1) \over 2} + \left \lfloor {n \over 2} \right\rfloor\!. \label{csymG} \end{equation} (1.3) A cubature rule of degree |$m$| with |$N$| attaining the lower bound in (1.3) is called minimal. When |$N$| is close to this lower bound, the rule is called near-minimal. The latter is the more frequent. From now on, we let |${\it{\Omega}}=B_2$| with $$B_2 = \{(x,y)\in \mathbb{R}^2:\|(x,y)\|_2\le 1\}$$ the unit disk in the Euclidean norm. It is less known that for |$m=2n-1$| and centrally symmetric weight functions, a symbolic (instead of numeric) Gaussian cubature formula exists (Benouahmane & Cuyt, 2000; Cuyt et al., 2004). We denote by |${\mathcal P}_{m, \theta}$| the set of polynomials of degree at most |$m=2n-1$| in the real variable |$z$| with coefficients that are rational functions in |$\cos\theta$| and |$\sin\theta$|. Note that the vector |$(\cos\theta, \sin\theta)$| spans the one-dimensional subspace |$\{(z\cos\theta, z\sin\theta): z\in {\mathbb{R}}\}$| of |${\mathbb{R}}^2$|. Then under some simple and obvious conditions (Benouahmane & Cuyt, 2000; Cuyt et al., 2004, 2011), $$\iint_{B_2} \pi_m(x \cos\theta + y \sin\theta) w\left(\|(x,y)\|_2\right) \; \,{\rm d}x\, \,{\rm d}y = \sum_{i=1}^n A_{i,n}(\theta) \pi_m\left(\xi_{i,n}(\theta)\right)\!,$$ where |$\pi_m {\,\in\,} {\mathcal P}_{m,\theta}$| and the |$\xi_{i,n}(\theta)$| are the zero curves of the polynomial |$P_n {\,\in\,} {\mathcal P}_{m,\theta}$| satisfying the orthogonality conditions \begin{equation} \iint_{B_2} (x\cos\theta + y\sin\theta)^i P_n(x\cos\theta + y\sin\theta) w\left(\|(x,y)\|_2 \right) \;\,{\rm d}x\; \,{\rm d}y = 0, \qquad i=0, \ldots, n-1. \label{MOP} \end{equation} (1.4) So a |$\theta$|-parameterized multivariate polynomial the result also exists in higher dimensions) of degree |$m=2n-1$| is integrated exactly over the unit ball (the result holds in other norms) by a linear combination of |$n$||$\theta$|-parameterized weights and evaluations of the polynomial integrand along zero curves of the spherical orthogonal polynomial of degree |$n$| (Cuyt et al., 2004, 2011). In the next section, we convert this symbolic integration result into a numeric rule of the classical form (1.1), and this for general |$n$|. In this way, we will rediscover all the minimal and near-minimal rules listed in Table 1. So in this article we base the derivation of cubature rules on the new theory of spherical orthogonal polynomials, which has (unfortunately) only been developed after the publication of many rules in separate articles. Table 1. Summary of existing (near-)minimal cubature rules Degree |$2n-1$| Nodes N Configuration References 3 4 A Hammer & Stroud (1958) 5 7 R-Ax-O Radon (1948) 7 12 2B-A Hammer & Stroud (1958) 9 19 3R-Ax-2Ay-O Albrecht (1960) 11 26 4R-2Ax-3Ay Piessens & Haegemans (1975a) 13 35 6R-3Ax-2Ay-O Cools & Haegemans (1988) 13 36 2C-2B-3A Cools & Haegemans (1987) 13 37 2C-2B-3A-O Rabinowitz & Richter (1969) 15 44 2C-3B-4A Rabinowitz & Richter (1969) 17 57 4C-3B-3A-O Cools & Kim (2000) 19 72 6C-2B-4A Kim & Song (1997) Degree |$2n-1$| Nodes N Configuration References 3 4 A Hammer & Stroud (1958) 5 7 R-Ax-O Radon (1948) 7 12 2B-A Hammer & Stroud (1958) 9 19 3R-Ax-2Ay-O Albrecht (1960) 11 26 4R-2Ax-3Ay Piessens & Haegemans (1975a) 13 35 6R-3Ax-2Ay-O Cools & Haegemans (1988) 13 36 2C-2B-3A Cools & Haegemans (1987) 13 37 2C-2B-3A-O Rabinowitz & Richter (1969) 15 44 2C-3B-4A Rabinowitz & Richter (1969) 17 57 4C-3B-3A-O Cools & Kim (2000) 19 72 6C-2B-4A Kim & Song (1997) Table 1. Summary of existing (near-)minimal cubature rules Degree |$2n-1$| Nodes N Configuration References 3 4 A Hammer & Stroud (1958) 5 7 R-Ax-O Radon (1948) 7 12 2B-A Hammer & Stroud (1958) 9 19 3R-Ax-2Ay-O Albrecht (1960) 11 26 4R-2Ax-3Ay Piessens & Haegemans (1975a) 13 35 6R-3Ax-2Ay-O Cools & Haegemans (1988) 13 36 2C-2B-3A Cools & Haegemans (1987) 13 37 2C-2B-3A-O Rabinowitz & Richter (1969) 15 44 2C-3B-4A Rabinowitz & Richter (1969) 17 57 4C-3B-3A-O Cools & Kim (2000) 19 72 6C-2B-4A Kim & Song (1997) Degree |$2n-1$| Nodes N Configuration References 3 4 A Hammer & Stroud (1958) 5 7 R-Ax-O Radon (1948) 7 12 2B-A Hammer & Stroud (1958) 9 19 3R-Ax-2Ay-O Albrecht (1960) 11 26 4R-2Ax-3Ay Piessens & Haegemans (1975a) 13 35 6R-3Ax-2Ay-O Cools & Haegemans (1988) 13 36 2C-2B-3A Cools & Haegemans (1987) 13 37 2C-2B-3A-O Rabinowitz & Richter (1969) 15 44 2C-3B-4A Rabinowitz & Richter (1969) 17 57 4C-3B-3A-O Cools & Kim (2000) 19 72 6C-2B-4A Kim & Song (1997) In the sequel, we assume that the weight function is of the form |$w(\|(x,y)\|_2)$|, such as the popular Gegenbauer weight function $$w_{\lambda} \left(\|(x,y)\|_2 \right)={(1-x^2-y^2)}^{\lambda-1/2},\qquad \lambda\ge 0.$$ When observing the existing numeric near-minimal cubature formulae on the disk, up to degree 31, for weight functions |$w(\|(x,y)\|_2)$|, we notice that the sets of nodes consist of pairs, quadruples or octuples taking one of the configurations in the Figs 1 and 2, namely Ax (nodes on |$x$|-axis), Ay (nodes on |$y$|-axis), B (nodes on bisectors), R (rectangle of nodes) or C (composition of 2 related rectangles), sometimes complemented with the origin O (only in case of odd |$n$|). Fig. 1. View largeDownload slide Configurations (from top to bottom and left to right) O, Ax, Ay, B and C. Fig. 1. View largeDownload slide Configurations (from top to bottom and left to right) O, Ax, Ay, B and C. Fig. 2. View largeDownload slide Configuration R. Fig. 2. View largeDownload slide Configuration R. Moreover, using the system of equations derived here in Section 4, one can see that for general odd |$n=2\nu +1$| a cubature rule with a configuration consisting of |$(\nu-1)\nu/2$| C, |$\nu$| B, |$\nu$| A and O always exists, and so does for general even |$n=2\nu$| a cubature rule with a configuration consisting of |$(\nu-1)\nu/2$| C and |$\nu$| B always exist. So these configurations deserve some special attention. The configurations O, Ax, Ay and B are all special cases of the general configuration R. We call a B configuration rotated over |$\pi/4$| an A configuration (actually a union of Ax and Ay with the same radius). The two rectangles in the configuration C are rotations of one another over |$\pi/2$|. In Table 1, we detail the smallest known configurations satisfying (1.1) up to and including degree 19 (|$n=10$|), for |$w(x,y)=1$|, and with nodes |$(x_i, y_i)$| in the closed unit disk and non-negative weights |$A_i$|. With every entry we only list the oldest known reference and do not mention configurations that are mere rotations. Also, if the cubature rule is not unique, additional constraints may lead to a small reduction of nodes, sometimes with the effect that other nodes move outside of the domain or weights become negative (such as for |$n=5$| (Piessens & Haegemans, 1975b)). Such rules are not listed in the table because we restrict it to rules with non-negative weights and nodes in the closed unit disk. It does, however, by no means imply that they are excluded from our new approach in Sections 3 and 4, as they are then special cases of the general nonunique solution. From degree 21 up to and including degree 31, cubature formulae are presented in Haegemans (1975), but they are not designed for optimality or near-minimality and therefore are not included in our table. The close connection between orthogonal polynomials and quadrature formulae on the one hand, and orthogonal polynomials and Padé approximation theory on the other hand, is well understood in the univariate case (Brezinski, 1980). In the multivariate case, these connections are less established, but remain valid as explained in Benouahmane & Cuyt (2000). In the sequel, we show that they deserve more of our attention, because they can lead to new insights when exploring a difficult topic such as numeric cubature. Making use of this connection, we show in Section 2 how a symbolic cubature rule on |${\it{\Omega}}=B_2$| automatically leads to a numeric cubature rule on the same domain. Once this is established, we describe in Sections 3 and 4 how to obtain numerical rules with a small number of nodes for general degree. In the past, cubature rules were often obtained for a specific degree, following a tailored approach. The current idea is generic and covers all of the rules in Table 1 and more. 2. From symbolic to numeric cubature Let |$c_{jk}=\cos(\theta_{jk})$| and |$s_{jk}=\sin(\theta_{jk})$| with |$\theta_{jk} \in (-\pi/2, \pi/2]$|. A bivariate polynomial of degree |$m$|, \begin{equation} p_m(x,y) = \sum_{j=0}^m \sum_{k=0}^j a_{j-k,k}x^{j-k} y^k \label{pm} \end{equation} (2.1) can be reexpressed as (Cuyt et al., 2012) \begin{equation} p_m(x,y) = \sum_{j=0}^m \sum_{k=0}^j b_{j-k,k} \left( c_{j-k,k}x+s_{j-k,k}y \right)^j, \qquad c_{j-k,k}^2 + s_{j-k,k}^2 = 1, \label{pmr} \end{equation} (2.2) if the |$\theta_{jk}$| satisfy $$\left| \begin{matrix} c_{j,0}^j & c_{j,0}^{j-1}s_{j,0} & \ldots & s_{j,0}^j \\[5pt] c_{j-1,1}^j & c_{j-1,1}^{j-1}s_{j-1,1} & \ldots & s_{j-1,1}^j \\ \vdots & \vdots & & \vdots \\ c_{0,j}^j & c_{0,j}^{j-1}s_{0,j} & \ldots & s_{0,j}^j \end{matrix} \right| \not= 0, \qquad j=0, \ldots, m. $$ Now let |$x_i=z_i\cos(t_i)$| and |$y_i=z_i\sin(t_i)$| with |$z_i \in \mathbb{R}$| and |$-\pi/2 < t_i \le \pi/2$|. Both for |$\theta$| and |$t$| we need not consider the full circle |$]-\pi, \pi]$| because the additional half of this set only brings in a sign change, which can be compensated by |$b_{j-k,k}$| or |$z_i$|. Lemma 2.1 The integral $$\iint_{B_2} (x\cos\theta+y\sin\theta)^j w(\|(x,y)\|_2)\; \,{\rm d}x\, \,{\rm d}y, \qquad j =0, 1, \ldots$$ is zero for |$j$| odd and independent of |$\theta$| for |$j$| even. Proof. Let |$x=z\cos t$| and |$y=z\sin t$| with |$-1 \le z \le 1$| and |$-\pi/2< t \le \pi/2$|. The integral can be written as a product of two one-dimensional integrals, $$\iint_{B_2} (x\cos\theta+y\sin\theta)^j w(\|(x,y)\|_2)\,{\rm d}x\, \,{\rm d}y = \int_{-1}^1 z^j w(|z|) |z|\, \,{\rm d}z \int_{-\pi/2}^{\pi/2} \cos^j(t-\theta) \, \,{\rm d}t.$$ For |$j$| odd, the first factor is zero. For |$j$| even, the second factor equals \begin{equation*} {(j-1)!! \; \pi \over (j/2)! \; 2^{j/2}}. \end{equation*} □ A straightforward connection between the symbolic cubature rules and the numeric ones is obtained from the following. If we can construct a numeric cubature formula of the form (1.1) that is exact for a |$\theta$|-parameterized bivariate polynomial integrand |$q_m(\theta; x, y)$| of degree |$m$| of the form \begin{equation} q_m(\theta; x, y) := q_m(x\cos\theta+y\sin\theta) = \sum_{j=0}^m a_j (x \cos\theta + y \sin\theta)^j, \label{qmparam} \end{equation} (2.3) in other words a cubature formula guaranteeing, with |$A_i$| and |$(x_i, y_i)$| independent of |$\theta$|, that \begin{eqnarray} \iint_{B_2} q_m \left(x \cos\theta + y \sin\theta\right) w(\|(x,y)\|_2) \; \,{\rm d}x \, \,{\rm d}y &= \sum_{i=1}^N A_i q_m \left(x_i \cos\theta + y_i\sin\theta \right)\nonumber\\ &= \sum_{i=1}^N A_i q_m \left( z_i \cos(t_i-\theta) \right)\!, \\ &\qquad \qquad x_i= z_i\cos t_i, y_i = z_i\sin t_i\nonumber \label{numG2} \end{eqnarray} (2.4) then we have at the same time a cubature formula that is exact for the general bivariate polynomial |$p_m(x,y)$| of degree |$m$| given in (2.1). The latter is obtained by combining (2.2), (2.3) and (2.5) into \begin{equation} \begin{split} \iint_{B_2} &p_m(x,y) w(\|(x,y)\|_2)\; \,{\rm d}x \, \,{\rm d}y \\ &= \sum_{j=0}^m \sum_{k=0}^j b_{j-k,k} \iint_{B_2} \left(x \cos\theta_{j-k,k} + y \sin\theta_{j-k,k} \right)^j w(\|(x,y)\|_2)\; \,{\rm d}x \, \,{\rm d}y \\ &= \sum_{j=0}^m \sum_{k=0}^j b_{j-k,k} \sum_{i=1}^N A_i \left( x_i \cos\theta_{j-k,k} + y_i \sin\theta_{j-k,k} \right)^j \\ &= \sum_{i=1}^N A_i p_m(x_i, y_i). \end{split} \label{numG} \end{equation} (2.5) In addition, expression (2.5) is actually independent of |$\theta$| because of Lemma 2.1. In Sections 3 and 4, we now develop some particular formulas (2.5), hich as we know are equivalent to formulas of the form (2.5), with the aim to obtain a unified construction of near-minimal numerical cubature rules. To avoid duplication, we concentrate in these sections on odd |$n$|-values, because the formulas and systems for even |$n$|-values can be obtained in exactly the same way. For completeness, we present the system for even |$n$|-values at the end of the respective Sections 3 and 4 without the derivation. 3. Computing (R, Ax, Ay, O) configured cubature formulae Our aim is to find a combination of nodes on |$H$| rectangles R, |$K_x$| axes Ax and |$K_y$| axes Ay and the origin O, so that \begin{equation} \begin{split} \iint_{B_2} q_m\left( x\cos \theta \right. & \left. + y\sin\theta\right)w \left(\|(x,y)\|_2 \right)\;\,{\rm d}x\,\,{\rm d}y \\ &=\sum_{i=1}^H A_i\left[ q_{m} (z_i\cos(t_i-\theta))+ q_{m}(-z_i\cos(t_i-\theta)) \right. \\ & \qquad \quad \left. +q_{m}(z_i\cos(-t_i-\theta))+ q_{m}(-z_i\cos(-t_i-\theta)) \right] \\ &+\sum_{i=H+1}^{H+K_x} A_i \left[ q_{m}(z_i\cos(0-\theta))+ q_{m}(-z_i\cos(0-\theta)) \right] \\ & +\sum_{i=H+K_x+1}^{H+K_x+K_y} A_i \left[ q_{m}(z_i\cos(\pi/2-\theta))+ q_{m}(-z_i\cos(\pi/2-\theta)) \right] \\ &+A_{H+K_x+K_y+1}\, q_m(0), \label{RAO} \end{split} \end{equation} (3.1) where |$q_m$| is as in (2.3) with |$m=2n-1$|. In the right-hand side, the |$A_i, z_i$| and |$t_i$| are unknown, making a grand total of |$3H+2(K_x+K_y)+1$| values to be determined. For the left-hand side, we denote $$\beta_j = \iint_{B_2} \left(x \cos\theta + y \sin\theta \right)^j w(\|(x,y)\|_2)\;\,{\rm d}x\,\,{\rm d}y,$$ and hence the left-hand side equals \begin{equation} \sum_{j=0}^{n-1} a_{2j} \beta_{2j}. \label{LHS} \end{equation} (3.2) Plugging (2.3) into the right-hand side results in \begin{equation} \begin{split} &\sum_{i=1}^H A_i \left[ \sum_{k=0}^{n-1} 2a_{2k}z_i^{2k} \left(\cos^{2k}(t_i-\theta)+\cos^{2k}(-t_i-\theta)\right)\right] \\ &\quad +\sum_{i=H+1}^{H+K_x} A_i \left[ \sum_{k=0}^{n-1}2a_{2k}z_i^{2k} \cos^{2k}(0-\theta)\right]+ \sum_{i=H+K_x+1}^{H+K_x+K_y} A_i\left[ \sum_{k=0}^{n-1} 2a_{2k}z_i^{2k} \cos^{2k}(\pi/2-\theta)\right] \\ &\quad +a_0A_{H+K_x+K_y+1}. \end{split} \nonumber \end{equation} Furthermore, using the trigonometric identity $$\cos^{2k}(\theta)={1 \over 2^{2k}} \left( {2k \choose k} +2 \sum_{l=0}^{k-1} {2k \choose l} \cos\left((2k-2l)\theta\right) \right)$$ results for the right-hand side in \[\begin{align} & 4{{a}_{0}}\sum\limits_{i=1}^{H}{{{A}_{i}}}+4\sum\limits_{k=1}^{n-1}{\frac{1}{{{2}^{2k}}}}{{a}_{2k}}\sum\limits_{i=1}^{H}{{{A}_{i}}}z_{i}^{2k} \\ & \,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\times \left[ \left( \begin{matrix} 2k \\ k \\ \end{matrix} \right)+2\sum\limits_{l=0}^{k-1}{\left( \begin{matrix} 2k \\ l \\ \end{matrix} \right)}\cos \left( (2k-2l){{t}_{i}} \right)\cos \left( (2k-2l)\theta \right) \right] \\ & +2{{a}_{0}}\sum\limits_{i=H+1}^{H+{{K}_{x}}}{{{A}_{i}}}+2\sum\limits_{k=1}^{n-1}{\frac{1}{{{2}^{2k}}}}{{a}_{2k}}\sum\limits_{i=H+1}^{H+{{K}_{x}}}{{{A}_{i}}}z_{i}^{2k} \\ & \,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\times \left[ \left( \begin{matrix} 2k \\ k \\ \end{matrix} \right)+2\sum\limits_{l=0}^{k-1}{\left( \begin{matrix} 2k \\ l \\ \end{matrix} \right)}\cos \left( (2k-2l)(0-\theta ) \right) \right] \\ & +2{{a}_{0}}\sum\limits_{i=H+{{K}_{x}}+1}^{H+{{K}_{x}}+{{K}_{y}}}{{{A}_{i}}}+2\sum\limits_{k=1}^{n-1}{\frac{1}{{{2}^{2k}}}}{{a}_{2k}}\sum\limits_{i=H+{{K}_{x}}+1}^{H+{{K}_{x}}+{{K}_{y}}}{{{A}_{i}}}z_{i}^{2k} \\ & \,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\times \left[ \left( \begin{matrix} 2k \\ k \\ \end{matrix} \right)+2\sum\limits_{l=0}^{k-1}{\left( \begin{matrix} 2k \\ l \\ \end{matrix} \right)}\cos \left( (2k-2l)(\pi /2-\theta ) \right) \right] \\ & +{{a}_{0}}{{A}_{H+{{K}_{x}}+{{K}_{y}}+1}}. \\ \end{align}\] We further denote $$\gamma_{2j} = {2^{2j} \over {2j \choose j}} \, \beta_{2j}.$$ For |$w(\|(x,y)\|_2)=1,$| for instance, |$\gamma_{2j}=\pi/(j+1)$|. For the Gegenbauer weight function |$w_\lambda\left(\|(x,y)\|_2\right)$|, we find $$\beta_{2j}(\lambda) = {2\pi (2j-1)!! \over (2\lambda+1)(2\lambda+3) \cdots (2\lambda+2j+1)}.$$ Identifying the coefficients of |$a_{2k}, k=0, \ldots, n-1$| in the left- and right-hand side provides |$n$| equations in the unknown values. Doing the same with the coefficients of |$a_{2k}\cos 2\theta, k=1, \ldots, n-1$| leads to an additional |$n-1$| equation and so on with |$a_{2k}\cos 4\theta, \ldots, a_{2n-2}\cos(2n-2)\theta$|. In total, we obtain in this way |$n(n+1)/2$| equations in the |$3H+2(K_x+K_y)+1$| unknowns that we put together below: $$\eqalign{ & \quad \quad \quad \left\{ {\matrix{ {4\sum\limits_{i = 1}^H {{A_i}} + 2\sum\limits_{i = H + 1}^{H + {K_x}} {{A_i}} + 2\sum\limits_{i = H + {K_x} + 1}^{H + {K_x} + {K_y}} {{A_i}} + {A_{H + {K_x} + {K_y} + 1}} = {\gamma _0},} \hfill \cr {4\sum\limits_{i = 1}^H {{A_i}} z_i^2 + 2\sum\limits_{i = H + 1}^{H + {K_x}} {{A_i}} z_i^2 + 2\sum\limits_{i = H + {K_x} + 1}^{H + {K_x} + {K_y}} {{A_i}} z_i^2 = {\gamma _2},} \hfill \cr {\quad \,\,\, \vdots } \hfill \cr {4\sum\limits_{i = 1}^H {{A_i}} z_i^{2n - 2} + 2\sum\limits_{i = H + 1}^{H + {K_x}} {{A_i}} z_i^{2n - 2} + 2\sum\limits_{i = H + {K_x} + 1}^{H + {K_x} + {K_y}} {{A_i}} z_i^{2n - 2} = {\gamma _{2n - 2}}} \hfill \cr } } \right. \cr & \cr & \left\{ {\matrix{ {2\sum\limits_{i = 1}^H {{A_i}} z_i^2\cos (2{t_i}) + \sum\limits_{i = H + 1}^{H + {K_x}} {{A_i}} z_i^2 - \sum\limits_{i = H + {K_x} + 1}^{H + {K_x} + {K_y}} {{A_i}} z_i^2 = 0,} \hfill \cr {\quad \,\,\, \vdots } \hfill \cr {2\sum\limits_{i = 1}^H {{A_i}} z_i^{2n - 2}\cos (2{t_i}) + \sum\limits_{i = H + 1}^{H + {K_x}} {{A_i}} z_i^{2n - 2} - \sum\limits_{i = H + {K_x} + 1}^{H + {K_x} + {K_y}} {{A_i}} z_i^{2n - 2} = 0} \hfill \cr } } \right. \cr & \cr & \quad \,\,\, \vdots \cr & \cr & \left\{ {\matrix{ {2\sum\limits_{i = 1}^H {{A_i}} z_i^{2n - 4}\cos ((2n - 4){t_i}) + \sum\limits_{i = H + 1}^{H + {K_x}} {{A_i}} z_i^{2n - 4} + {{( - 1)}^{n - 2}}\sum\limits_{i = H + {K_x} + 1}^{H + {K_x} + {K_y}} {{A_i}} z_i^{2n - 4} = 0} \hfill \cr {2\sum\limits_{i = 1}^H {{A_i}} z_i^{2n - 2}\cos ((2n - 4){t_i}) + \sum\limits_{i = H + 1}^{H + {K_x}} {{A_i}} z_i^{2n - 2} + {{( - 1)}^{n - 2}}\sum\limits_{i = H + {K_x} + 1}^{H + {K_x} + {K_y}} {{A_i}} z_i^{2n - 2} = 0} \hfill \cr } } \right. \cr & \cr & \left\{ {2\sum\limits_{i = 1}^H {{A_i}} z_i^{2n - 2}\cos ((2n - 2){t_i}) + \sum\limits_{i = H + 1}^{H + {K_x}} {{A_i}} z_i^{2n - 2} + {{( - 1)}^{n - 1}}\sum\limits_{i = H + {K_x} + 1}^{H + {K_x} + {K_y}} {{A_i}} z_i^{2n - 2} = 0.} \right. \cr} $$ (3.3) For |$n$| even, the system looks identical, with the exception that |$A_{H+K_x+K_y+1}$| needs to be removed from the very first equation. 4. Computing (C, B, A, O) configured cubature formulae To reduce the number of unknowns more when |$n$| increases, we apply the same approach to a combination of nodes from |$I$| octuples C, |$J$| bisector sets B, |$K$| axes as in A (Ax and Ay united) complemented with the origin O (remember that |$n$| is odd). To be in line with the previous section |$K=K_x+K_y$|. The gain is in the fact that with |$K_x=K_y$| each couple of Ax and Ay nodes gets assigned the same weight |$A_i$| and the same (signed) distance |$z_i$|. Then the grand total of unknowns |$A_i, z_i$| and |$t_i$| is |$3I+2J+K+1$|. The number of equations in the system is also greatly reduced. We denote |$n-1=2\nu$| and introduce the compact notation $$Q_m(\theta; z_i, t_i) = q_{m}(z_i\cos(t_i-\theta))+ q_{m}(-z_i\cos(t_i-\theta)).$$ The analogue of (3.1) is \begin{equation} \begin{split} \iint_{B_2} & q_m\left( x\cos \theta + y\sin\theta\right)w \left(\|(x,y)\|_2 \right)\;\,{\rm d}x\,\,{\rm d}y \\ &= \sum_{i=1}^{I} A_i\left[ Q_m(\theta; z_i, t_i) + Q_m(\theta; z_i, -t_i) +Q_m(\theta; z_i, \pi/2-t_i) + Q_m(\theta; z_i, \pi/2+t_i) \right] \\ &+\sum_{i=I+1}^{I+J} A_i \left[ Q_m(\theta; z_i, \pi/4) + Q_m(\theta; z_i, -\pi/4) \right] \\ &+\sum_{i=I+J+1}^{I+J+K/2} A_i \left[ Q_m(\theta; z_i, 0) + Q_m (\theta; z_i, \pi/2) \right] +A_{I+J+K/2+1}\,q_m(0). \end{split} \nonumber \end{equation} Note that some of the sums may be empty, as before, and hence do not contribute anything to the expression. The left-hand side remains unchanged and equals (3.2). Now plugging (2.3) into the right-hand side and making use of the trigonometric identity for |$\cos^{2k}(\theta)$| results in \begin{equation} \begin{split} 8a_0 \sum_{i=1}^I A_i &+8 \sum_{k=1}^{n-1} {1 \over 2^{2k}} a_{2k} \sum_{i=1}^I A_i z_i^{2k} \\ &\times\left[ {2k \choose k} + \sum_{l=0}^{k-1} {2k \choose l} \left( 1+(-1)^{k-l} \right) \cos\left((2k-2l)t_i\right) \cos\left((2k-2l)\theta\right) \right] \\ &+4a_0 \sum_{i=I+1}^{I+J} A_i +4 \sum_{k=1}^{n-1} {1 \over 2^{2k}} a_{2k} \sum_{i=I+1}^{I+J} A_i z_i^{2k} \\ &\times \left[ {2k \choose k} + 2 \sum_{l=0}^{k-1} {2k \choose l} \cos\left((2k-2l)\pi/4 \right) \cos\left((2k-2l)\theta \right) \right] \\ &+4a_0 \sum_{i=I+J+1}^{I+J+K/2} A_i +4 \sum_{k=1}^{n-1} {1 \over 2^{2k}} a_{2k} \sum_{i=I+J+1}^{I+J+K/2} A_i z_i^{2k} \\ &\times \left[ {2k \choose k} + 2 \sum_{l=0}^{k-1} {2k \choose l} \left( 1+(-1)^{k-l} \right) \cos\left((2k-2l)\theta \right) \right] \\ &+a_0A_{I+J+K/2+1}. \end{split} \nonumber \end{equation} After identifying the coefficients in the left- and right-hand side again, we arrive at the following system of |$(n+1)^2/4$| equations: \[\begin{align} & \quad \quad \quad \,\,\,\left\{ \begin{array}{*{35}{l}} 8\sum\limits_{i=1}^{I}{{{A}_{i}}}+4\sum\limits_{i=I+1}^{I+J}{{{A}_{i}}}+4\sum\limits_{i=I+J+1}^{I+J+K/2}{{{A}_{i}}}={{\gamma }_{0}}, \\ 8\sum\limits_{i=1}^{I}{{{A}_{i}}}z_{i}^{2}+4\sum\limits_{i=I+1}^{I+J}{{{A}_{i}}}z_{i}^{2}+4\sum\limits_{i=I+J+1}^{I+J+K/2}{{{A}_{i}}}z_{i}^{2}={{\gamma }_{2}}, \\ \,\,\,\,\,\,\,\,\,\vdots \\ 8\sum\limits_{i=1}^{I}{{{A}_{i}}}z_{i}^{2n-2}+4\sum\limits_{i=I+1}^{I+J}{{{A}_{i}}}z_{i}^{2n-2}+4\sum\limits_{i=I+J+1}^{I+J+K/2}{{{A}_{i}}}z_{i}^{2n-2}={{\gamma }_{2n-2}} \\ \end{array} \right. \\\\ & \left\{ \begin{array}{*{35}{l}} 8\sum\limits_{i=1}^{I}{{{A}_{i}}}z_{i}^{4}\cos (4{{t}_{i}})-4\sum\limits_{i=I+1}^{I+J}{{{A}_{i}}}z_{i}^{4}+4\sum\limits_{i=I+J+1}^{I+J+K/2}{{{A}_{i}}}z_{i}^{4}=0, \\ \,\,\,\,\,\,\,\,\,\vdots \\ 8\sum\limits_{i=1}^{I}{{{A}_{i}}}z_{i}^{2n-2}\cos (4{{t}_{i}})-4\sum\limits_{i=I+1}^{I+J}{{{A}_{i}}}z_{i}^{2n-2}+4\sum\limits_{i=I+J+1}^{I+J+K/2}{{{A}_{i}}}z_{i}^{2n-2}=0 \\ \end{array} \right. \\\\ & \,\,\,\,\,\,\,\,\,\,\,\vdots \\\\ & \left\{ \begin{array}{*{35}{l}} 8\sum\limits_{i=1}^{I}{{{A}_{i}}}z_{i}^{4\nu -4}\cos ((4\nu -4){{t}_{i}})+{{(-1)}^{(\nu -1)}}4\sum\limits_{i=I+1}^{I+J}{{{A}_{i}}}z_{i}^{4\nu -4}+4\sum\limits_{i=I+J+1}^{I+J+K/2}{{{A}_{i}}}z_{i}^{4\nu -4}=0 \\ 8\sum\limits_{i=1}^{I}{{{A}_{i}}}z_{i}^{4\nu -2}\cos ((4\nu -4){{t}_{i}})+{{(-1)}^{(\nu -1)}}4\sum\limits_{i=I+1}^{I+J}{{{A}_{i}}}z_{i}^{4\nu -2}+4\sum\limits_{i=I+J+1}^{I+J+K/2}{{{A}_{i}}}z_{i}^{4\nu -2}=0 \\ 8\sum\limits_{i=1}^{I}{{{A}_{i}}}z_{i}^{4\nu }\cos ((4\nu -4){{t}_{i}})+{{(-1)}^{(\nu -1)}}4\sum\limits_{i=I+1}^{I+J}{{{A}_{i}}}z_{i}^{4\nu }+4\sum\limits_{i=I+J+1}^{I+J+K/2}{{{A}_{i}}}z_{i}^{4\nu }=0 \\ \end{array} \right. \\\\ & \left\{ 8\sum\limits_{i=1}^{I}{{{A}_{i}}}z_{i}^{4\nu }\cos (4\nu {{t}_{i}})+{{(-1)}^{\nu }}4\sum\limits_{i=I+1}^{I+J}{{{A}_{i}}}z_{i}^{4\nu }+4\sum\limits_{i=I+J+1}^{I+J+K/2}{{{A}_{i}}}z_{i}^{4\nu }=0. \right. \\ \end{align}\] This system of |$(n+1)^2/4$| equations can be rewritten in the following orm, where the first |$n$| equations are not touched, but the remaining subsystems are: $$\eqalign{ & \quad \quad \left\{ {\matrix{ {8\sum\limits_{i = 1}^I {{A_i}} + 4\sum\limits_{i = I + 1}^{I + J} {{A_i}} + 4\sum\limits_{i = I + J + 1}^{I + J + K/2} {{A_i}} + {A_{I + J + K/2 + 1}} = {\gamma _0},} \hfill \cr {8\sum\limits_{i = 1}^I {{A_i}} z_i^2 + 4\sum\limits_{i = I + 1}^{I + J} {{A_i}} z_i^2 + 4\sum\limits_{i = I + J + 1}^{I + J + K/2} {{A_i}} z_i^2 = {\gamma _2},} \hfill \cr {{\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} \vdots } \hfill \cr {8\sum\limits_{i = 1}^I {{A_i}} z_i^{2n - 2} + 4\sum\limits_{i = I + 1}^{I + J} {{A_i}} z_i^{2n - 2} + 4\sum\limits_{i = I + J + 1}^{I + J + K/2} {{A_i}} z_i^{2n - 2} = {\gamma _{2n - 2}}} \hfill \cr } } \right. \cr & \cr & \left\{ {\matrix{ {8\sum\limits_{i = 1}^I {{A_i}} z_i^{4l}(1 - \cos (4l{t_i})) + 4\sum\limits_{i = I + 1}^{I + J} {{A_i}} z_i^{4l}(1 - \cos l\pi ) = {\gamma _{4l}},} \hfill \cr {8\sum\limits_{i = 1}^I {{A_i}} z_i^{4l + 2}(1 - \cos (4l{t_i})) + 4\sum\limits_{i = I + 1}^{I + J} {{A_i}} z_i^{4l + 2}(1 - \cos l\pi ) = {\gamma _{4l + 2}},} \hfill \cr {{\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} \vdots } \hfill \cr {8\sum\limits_{i = 1}^I {{A_i}} z_i^{2n - 2}(1 - \cos (4l{t_i})) + 4\sum\limits_{i = I + 1}^{I + J} {{A_i}} z_i^{2n - 2}(1 - \cos l\pi ) = {\gamma _{2n - 2}}} \hfill \cr } } \right. \cr & \quad \quad \quad \quad \quad \quad \quad \quad \quad \quad \quad \quad \quad \quad \quad \quad \quad \quad \quad \,l = 1, \ldots ,v. \cr} $$ (4.1) Note that for even |$l>1$| the contribution from the B configuration in the subsystem, in other words the sum from |$I+1$| to |$I+J$|, vanishes. The last system, where |$l=\nu$|, consists of only one equation. For |$\rm{l} \not{\,=\,}2$| the particular subsystem consists of |$n-4$| equations in the unknown coefficients |$A_i(1-\cos(8t_i))$| and the powered |$z_i^2$| for |$i=1, \ldots, I$|. From Prony’s method (Hildebrand, 1956, pp. 378–382), we then know that |$I$| cannot be less than or equal to |$(n-5)/2$|. Otherwise, following Henrici (1974, p. 603) and Kaltofen et al. (2000), the |$(I+1)\times (I+1)$| matrix $$\left( \begin{matrix} \gamma_8 & \gamma_{10} & \ldots & \gamma_{8+2I} \\ \gamma_{10} & & & \vdots \\ \vdots & & & \\ \gamma_{8+2I} & \ldots & & \gamma_{8+4I} \end{matrix} \right)$$ should be singular, which is a contradiction (for all the weight functions we considered, including the Gegenbauer weight). Let us now take a look at the top subsystem, which consists of |$2\nu+1$| equations. Its first equation can always serve to compute |$A_{I+J+K+1}$|, the weight at the origin. Again using results on Prony-structured systems, we see that for |$J=\nu$| (or for |$K=2\nu$|), the remaining |$2\nu$| equations can be used to compute the |$A_i$| and |$z_i^2$| for |$i=I+1, \ldots, I+J$| (or |$i=I+J+1, \ldots, I+J+K/2$|). For completeness, we give the system of |$n(n+2)/4$| equations in case |$n$| is even, namely |$n=2\nu$|: \[\begin{align} & \quad \,\,\,\left\{ \begin{array}{*{35}{l}} 8\sum\limits_{i=1}^{I}{{{A}_{i}}}+4\sum\limits_{i=I+1}^{I+J}{{{A}_{i}}}+4\sum\limits_{i=I+J+1}^{I+J+K/2}{{{A}_{i}}}={{\gamma }_{0}}, \\ 8\sum\limits_{i=1}^{I}{{{A}_{i}}}z_{i}^{2}+4\sum\limits_{i=I+1}^{I+J}{{{A}_{i}}}z_{i}^{2}+4\sum\limits_{i=I+J+1}^{I+J+K/2}{{{A}_{i}}}z_{i}^{2}={{\gamma }_{2}}, \\ \quad \,\,\,\vdots \\ 8\sum\limits_{i=1}^{I}{{{A}_{i}}}z_{i}^{2n-2}+4\sum\limits_{i=I+1}^{I+J}{{{A}_{i}}}z_{i}^{2n-2}+4\sum\limits_{i=I+J+1}^{I+J+K/2}{{{A}_{i}}}z_{i}^{2n-2}={{\gamma }_{2n-2}} \\ \end{array} \right. \\ & \\ & \left\{ \begin{array}{*{35}{l}} 8\sum\limits_{i=1}^{I}{{{A}_{i}}}z_{i}^{4l}(1-\cos (4l{{t}_{i}}))+4\sum\limits_{i=I+1}^{I+J}{{{A}_{i}}}z_{i}^{4l}(1-\cos l\pi )={{\gamma }_{4l}}, \\ 8\sum\limits_{i=1}^{I}{{{A}_{i}}}z_{i}^{4l+2}(1-\cos (4l{{t}_{i}}))+4\sum\limits_{i=I+1}^{I+J}{{{A}_{i}}}z_{i}^{4l+2}(1-\cos l\pi )={{\gamma }_{4l+2}}, \\ \quad \,\,\,\vdots \\ 8\sum\limits_{i=1}^{I}{{{A}_{i}}}z_{i}^{2n-2}(1-\cos (4l{{t}_{i}}))+4\sum\limits_{i=I+1}^{I+J}{{{A}_{i}}}z_{i}^{2n-2}(1-\cos l\pi )={{\gamma }_{2n-2}} \\ \end{array} \right. \\ & \quad \quad \quad \quad \quad \quad \quad \quad \quad \quad \quad \quad \quad \quad \quad \quad \quad l=1,\ldots ,\nu -1. \\ \end{align}\] (4.2) The same conclusions for |$I, J$| and |$K$| hold as in the case where |$n$| is odd, especially |$2I > n-5$|. Note that now the last system, where |$l=\nu-1$|, consists of two equations. In Section 5, we take a look at the systems (3.3) and (4.1)–(4.2) for particular |$n$|. All the cubature rules listed in Table 1 are rediscovered (we worked out all the details, but only present three typical examples: the cases |$n=3, 5, 7$|)! This new insight should enable us to develop a sharper lower bound for |$N$| in the future. 5. Numerical illustration Let us for simplicity denote |$\rho_i:= z_i^2$| and |$c_i:=\cos(2t_i)$| with |$0 \le t_i \le \pi/2$| in (3.3) and |$C_i:=\cos(4t_i)$| with |$0 \le t_i \le \pi/4$| in (4.1)–(4.2). We consider |$w(\|(x,y)\|_2)=1$|. Rather than going through the details of the recomputation of every cubature rule in Table 1, we illustrate the principle with one case using (3.3) and one case using (4.1)–(4.2). In addition, a third illustration supports the result that no 17-point rule exists for |$n=5$| (Verlinden & Cools, 1992). We start with a simple illustration of (3.3) for |$n=3$|. The system (3.3) has |$n(n+1)/2=6$| equations. The smallest configuration (meaning smallest |$N$|) having at least 6 unknowns is the one with |$H=1, K_x=1$| and |$K_y=0$| (or vice versa), |$O=1$| and |$N=7$| nodes: \[\left\{ \begin{align} &4A_1 + 2A_2 + A_3 = \pi \\ &4A_1 \rho_1 + 2A_2 \rho_2 = {\pi \over 2} \\ &4A_1 \rho_1^2 + 2A_2 \rho_2^2 = {\pi \over 3}\\ &4A_1 \rho_1 c_1 - 2A_2 \rho_2 = 0\\ &4A_1 \rho_1^2 c_1 - 2A_2 \rho_2^2 = 0\\ &4A_1 \rho_1^2 (2 c_1^2-1) + 2 A_2 \rho_2^2 = 0. \end{align} \right.\] From the 4th and 5th equation, we find that |$2A_2\rho_2=4A_1c_1\rho_1$| and |$2A_2\rho_2^2=4A_1c_1\rho_1^2$|. These can be substituted in the 2nd and 3rd equation to result in: \[\left\{ \begin{align} & 4{{A}_{1}}(1+{{c}_{1}}){{\rho }_{1}}=\frac{\pi }{2} \\ & 4{{A}_{1}}(1+{{c}_{1}})\rho _{1}^{2}=\frac{\pi }{3}. \\ \end{align} \right.\] This Prony-structured system can be solved immediately. We find |$4A_1(1+c_1) = 3\pi/4$| and |$\rho_1=2/3$|. Substituting these values in the 4th and 5th equation give another Prony-structured system in the unknowns |$A_2$| and |$\rho_2$|: \[\left\{ \begin{align} &2A_2 \rho_2 = {\pi \over 2} - {8 \over 3} A_1 \\ &2A_2 \rho_2^2 = {\pi \over 3} - {16 \over 9} A_1. \end{align} \right.\] In terms of |$A_1,$| the solution is \[\begin{align} & {{\rho }_{2}}=\frac{\frac{\pi }{3}-\frac{16}{9}{{A}_{1}}}{\frac{\pi }{2}-\frac{8}{3}{{A}_{1}}} \\ & {{A}_{2}}=\frac{2{{\left( \frac{\pi }{4}-\frac{4}{3}{{A}_{1}} \right)}^{2}}}{\frac{\pi }{3}-\frac{16}{9}{{A}_{1}}}. \\ \end{align}\] For |$c_1$| we have the expression |$c_1=-1+3\pi/(16 A_1)$|. Substituting everything in the 6th equation gives us |$A_1=\pi/8$| from which we can deduce all the other values. The obtained cubature formula is the one published in Radon (1948): \begin{equation} \begin{split} \iint_{B_2} p_5(x,y) \; \,{\rm d}x \,& \,{\rm d}y = {\pi \over 8} \left( p_5 \left( \scriptstyle{1 \over \sqrt{2}}, {1 \over \sqrt{6}} \right) + p_5 \left( -\scriptstyle{1 \over \sqrt{2}}, {1 \over \sqrt{6}} \right) + p_5 \left( \scriptstyle{1 \over \sqrt{2}}, -{1 \over \sqrt{6}} \right) + p_5 \left( -\scriptstyle{1 \over \sqrt{2}}, -{1 \over \sqrt{6}} \right) \right) \\ &\qquad + {\pi \over 8} \left( p_5 \left( \sqrt{2/3}, 0 \right) + p_5 \left( -\sqrt{2/3}, 0 \right) \right) + {\pi \over 4} p_5(0,0). \end{split} \nonumber \end{equation} The solution is unique up to some rotations of the configurations R and Ax, since these do not affect the values |$\rho_{1,2}$| or the weights. For |$n=5,$| the choice |$H=2, K_x=2, K_y=2, O=1$| delivers a 17-point configuration and a system of 15 equations in 15 unknowns, of which one can easily prove that it has no solution. This is completely in line with the literature. It is proved in Verlinden & Cools (1992) that no 17-point rule exists on the disk. For |$n=7,$| two configurations are possible. From the setup of Section 3, we find with |$H=2, K_x=3, K_y=2$| a Prony-like system of 28 equations in 29 unknowns. The 35-point rule cited in Table 1 is a solution of the system. The next setup is that of Section 4 with the systems (4.1)–(4.2). We know that |$I \ge 2$| and we choose |$K=6$|, as suggested in the previous section. The number of equations in (4.1) is |$(n+1)^2/4=16$|. Because of the Prony character of the system of equations, we use it as a lower bound for the number of unknowns |$3I+2J+K+1$|. We remind the reader that Prony systems are square nonlinear (polynomial) systems containing an equal number of linear and nonlinear unknowns, which can be solved for the nonlinear and the linear unknowns separately. With |$I=2$| this gives |$2J \ge 3$| or |$J\ge 2$|. We have at least 17 unknowns in (4.1), and hence are expecting nonuniqueness of the solution. The number of nodes for |$I=2, J=2, K=6,$| including the origin, is 37. The system of 16 equations in 17 unknowns is given by: \[\left\{ \begin{align} &8A_1 + 8A_2 + 4A_3 + 4A_4 + 4A_5 +4A_6 + 4A_7 + A_8 = \pi \\ &8A_1 \rho_1 + 8A_2 \rho_2 + 4 A_3\rho_3 + 4 A_4\rho_4 + 4 A_5\rho_5 + 4 A_6\rho_6 + 4 A_7\rho_7 = {\pi \over 2} \\&8A_1 \rho_1^2 + 8A_2 \rho_2^2 + 4 A_3\rho_3^2 + 4 A_4\rho_4^2 + 4 A_5\rho_5^2 + 4 A_6\rho_6^2 + 4 A_7\rho_7^2 = {\pi \over 3} \\&8A_1 \rho_1^3 + 8A_2 \rho_2^3 + 4 A_3\rho_3^3 + 4 A_4\rho_4^3 + 4 A_5\rho_5^3 + 4 A_6\rho_6^3 + 4 A_7\rho_7^3 = {\pi \over 4} \\\ &8A_1 \rho_1^4 + 8A_2 \rho_2^4 + 4 A_3\rho_3^4 + 4 A_4\rho_4^4 + 4 A_5\rho_5^4 + 4 A_6\rho_6^4 + 4 A_7\rho_7^4 = {\pi \over 5} \\&8A_1 \rho_1^5 + 8A_2 \rho_2^5 + 4 A_3\rho_3^5 + 4 A_4\rho_4^5 + 4 A_5\rho_5^5 + 4 A_6\rho_6^5 + 4 A_7\rho_7^5 = {\pi \over 6} \\&8A_1 \rho_1^6 + 8A_2 \rho_2^6 + 4 A_3\rho_3^6 + 4 A_4\rho_4^6 + 4 A_5\rho_5^6 + 4 A_6\rho_6^6 + 4 A_7\rho_7^6 = {\pi \over 7} \\ &8A_1(1-C_1) \rho_1^2 + 8 A_2(1-C_2)\rho_2^2 + 8A_3 \rho_3^2+ 8A_4 \rho_4^2 = {\pi \over 3} \\ &8A_1(1-C_1) \rho_1^3 + 8 A_2(1-C_2)\rho_2^3 + 8A_3 \rho_3^3+ 8A_4 \rho_4^3 = {\pi \over 4} \\ &8A_1(1-C_1) \rho_1^4 + 8 A_2(1-C_2)\rho_2^4 + 8A_3 \rho_3^4+ 8A_4 \rho_4^4 = {\pi \over 5} \\ &8A_1(1-C_1) \rho_1^5 + 8 A_2(1-C_2)\rho_2^5 + 8A_3 \rho_3^5+ 8A_4 \rho_4^5 = {\pi \over 6} \\ &8A_1(1-C_1) \rho_1^6 + 8 A_2(1-C_2)\rho_2^6 + 8A_3 \rho_3^6 + 8A_4 \rho_4^6= {\pi \over 7} \\ &16A_1(1-C_1^2) \rho_1^4 + 16 A_2(1-C_2^2)\rho_2^4 ={\pi \over 5} \\ &16A_1(1-C_1^2) \rho_1^5 + 16 A_2(1-C_2^2)\rho_2^5 ={\pi \over 6} \\ &16A_1(1-C_1^2) \rho_1^6 + 16 A_2(1-C_2^2)\rho_2^6 ={\pi \over 7} \\ &32 A_1C_1 (1-C_1^2) \rho_1^6 + 32 A_2C_2 (1-C_2^2) \rho_2^6 = 0. \end{align} \right.\] The way to tackle this larger system is to break it up in smaller (square) Prony subsystems involving only a subset of the unknowns. We start with the equations 13–15 and add another equation with a symbolic right-hand side, $$16A_1(1-C_1^2)\rho_1^3 + 16A_2(1-C_2^2) \rho_2^3 = \alpha,$$ to have a full Prony system in the unknowns |$B_i=A_i(1-C_i^2), i=1,2$| and |$\rho_1, \rho_2$|. So we can compute expressions for these unknowns in terms of |$\alpha$|. Introducing the parameter |$\alpha$| will actually allow us to characterize the whole infinite set of cubature rules of the form (C,B, A, O) for |$n=7$|. Equation number 16 provides us with a relation between |$C_1$| and |$C_2$|, $$B_1 C_1 \rho_1^6 + B_2 C_2 \rho_2^6 = 0,$$ and hence we can express |$C_2$| in terms of |$B_1, B_2, C_1$| and |$\alpha$|. From the equations 9–12, we find a Prony system in the unknowns |$A_3, A_4, \rho_3$| and |$\rho_4$|. The right-hand sides are in terms of |$C_1$| and |$\alpha$| and so the computed solutions for the unknowns are too. From equation number 8 we obtain |$C_1$| in terms of |$\alpha$|: $${B_1 \over 1+C_1} \rho_1^2 + {B_2 \over 1+C_2} \rho_2^2 + A_3 \rho_3^2 + A_4 \rho_4^2 = {\pi \over 24}.$$ Substituting this into $$B_i=A_i(1-C_i^2), \qquad i=1,2$$ gives us |$A_1$| and |$A_2$| in terms of |$\alpha$|. From the equations 2–7, we obtain expressions for |$A_5, A_6, A_7$| and |$\rho_5, \rho_6, \rho_7$|, also in terms of |$\alpha$|. Finally, from the 1st equation we obtain |$A_8(\alpha)$|. In a nutshell, solving the Prony-like system of 16 equations delivers the cubature rules published in Rabinowitz & Richter (1969) and Cools & Haegemans (1987). The infinity of solutions is described here in terms of a parameter |$\alpha$|, where a valid range for |$\alpha$| is delimited by the constraints \[\begin{align} & -1\le {{C}_{i}}(\alpha )\le 1,\qquad i=1,2, \\ & 0\le {{\rho }_{i}}(\alpha ),\qquad i=1,2,\ldots ,7. \\ \end{align}\] Within this valid region, some of the |$\rho_i$| may become larger than 1 and some weights may become negative. Without presenting the full symbolic analysis of the set of 37-point rules, the following conclusions for cubature rules of the type where |$I=2, J=2, K=6$| are easy to see: From the formulas $$\rho_{1,2} = {-450\alpha + 105\pi \pm \sqrt{15}\sqrt{13\,500\alpha^2-6650\alpha\pi+819\pi^2} \over 42(-25\alpha+6\pi)},$$ we immediately deduce that the valid region is a subset of |$\alpha > 6\pi/25$|. In this region, we can prove that the |$B_{1,2}>0$|, and hence that the weights |$A_{1,2}>0$|, as can be seen in Fig. 3. Fig. 3. View largeDownload slide Expressions |$B_1(\alpha)$| (left) and |$B_2(\alpha)$| (right). Fig. 3. View largeDownload slide Expressions |$B_1(\alpha)$| (left) and |$B_2(\alpha)$| (right). Furthermore, the 37-point rule can be reduced to a 36-point rule by choosing |$\alpha=\alpha^*$|, where |$A_8(\alpha^*)=0$| for |$\alpha^*=0.7751895643351289$|. 6. Discussion When increasing the number of octuples, rectangles, squares or axes pairs in the configurations, only the number of unknowns grows. The number of equations is independent of |$H, I, J, K_x, K_y$| and |$K$|. So the result facilitates the search for a minimal configuration: match the number of unknowns, namely (in case |$n$| is odd) |$3H+2K_x+2K_y$| for (3.3) and |$3I+2J+K$| for (4.1)–(4.2), because of the Prony character of the equations, as closely as possible to the number of these equations, and match the number of nodes, respectively (again in case |$n$| is odd) |$4H+2K_x+2K_y$| and |$8I+4J+2K$|, as closely as possible to the lower bound for |$N$|. We expect to find more general results on cubature in the near future, using the same unifying theory as explained in this article. Already, for some of the cubature rules listed in Table 1, we have been able to reobtain the nodes as common zeros of some linear combinations of the multivariate spherical orthogonal polynomials satisfying (1.4). References Albrecht J. ( 1960 ) Formeln zur numerischen Integration über Kreisbereiche. ZAMM Z. Angew. Math. Mech. , 40 , 514 – 517 . Google Scholar Crossref Search ADS Benouahmane B. & Cuyt A. ( 2000 ) Multivariate orthogonal polynomials, homogeneous Padé approximants and Gaussian cubature. Numer. Algorithms , 24 , 1 – 15 . Google Scholar Crossref Search ADS Brezinski C. ( 1980 ) Padé Type Approximation and General Orthogonal Polynomials . Basel : ISNM 50, Birkhäuser . Cools R. & Haegemans A. ( 1987 ) Construction of fully symmetric cubature formulae of degree 4k-3 for fully symmetric planar regions. J. Comp. Appl. Math. , 17 , 173 – 180 . Google Scholar Crossref Search ADS Cools R. & Haegemans A. ( 1988 ) Construction of symmetric cubature formulae with the number of knots (almost) equal to Möller’s lower bound. Numerical Integration III, (Oberwolfach, 1987) . Internat. Schriftenreihe Numer. Math ., vol. 85 . Basel : Birkhäuser , pp. 25 – 36 . Google Scholar Crossref Search ADS Cools R. & Kim K. ( 2000 ) A survey of known and new cubature formulas for the unit disk. Korean J. Comput. Appl. Math. , 7 , 477 – 485 . Cuyt A. , Benouahmane B. & Verdonk B. ( 2004 ) Spherical orthogonal polynomials and symbolic-numeric Gaussian cubature formulas. LNCS 3037 ( Bubak M. et al. eds). Berlin : Springer , pp. 557 – 560 . ( see extended version printed as Technical Report ). Cuyt A. , Benouahmane B. , Hamsapriye & Yaman I. ( 2011 ) Symbolic-numeric Gaussian cubature rules. Appl. Numer. Math. , 61 , 929 – 945 . Google Scholar Crossref Search ADS Cuyt A. , Yaman I. , Ibrahimoglu A. & Benouahmane B. ( 2012 ) Radial orthogonality and Lebesgue constants on the disk. Numer. Algorithms , 61 , 291 – 313 . Google Scholar Crossref Search ADS Haegemans A. ( 1975 ) Tables of circularly symmetrical integration formulas of degree 2d-1 for two-dimensional circularly symmetrical regions. Technical Report . Belgium : K.U. Leuven . Hammer P. & Stroud A. ( 1958 ) Numerical evaluation of multiple integrals II. Math. Tables Other Aids Computation , 12 , 272 – 280 . Google Scholar Crossref Search ADS Henrici P. ( 1974 ) Applied and Computational Complex Analysis I . New York : John Wiley . Hildebrand F. ( 1956 ) Introduction to Numerical Analysis . New York : McGraw Hill . Kaltofen E. , Lee W.-s. & Lobo A. A. ( 2000 ) Early termination in Ben-Or/Tiwari sparse interpolation and a hybrid of Zippel’s algorithm. Proceedings of the 2000 International Symposium on Symbolic and Algebraic Computation . New York : ACM , pp. 192 – 201 . Kim K.-J. & Song M.-S. ( 1997 ) Symmetric quadrature formulas over a unit disk. Korean J. Comp. Appl. Math. , 4 , 179 – 192 . Li H. , Sun J. & Xu Y. ( 2008 ) Discrete Fourier analysis, cubature, and interpolation on a hexagon and a triangle. SIAM J. Numer. Anal. , 46 , 1653 – 1681 . Google Scholar Crossref Search ADS Möller H. ( 1976 ) Kubaturformeln mit minimaler Knotenzahl. Numer. Math. , 25 , 185 – 200 . Google Scholar Crossref Search ADS Piessens R. & Haegemans A. ( 1975a ) Cubature formulas of degree eleven for symmetric planar regions. J. Comp. Appl. Math. , 1 , 79 – 83 . Google Scholar Crossref Search ADS Piessens R. & Haegemans A. ( 1975b ) Cubature formulas of degree nine for symmetric planar regions. Math. Comp. , 29 , 810 – 815 . Google Scholar Crossref Search ADS Rabinowitz P. & Richter N. ( 1969 ) Perfectly symmetric two-dimensional integration formulas with minimal numbers of points. Math. Comp. , 23 , 765 – 779 . Google Scholar Crossref Search ADS Radon J. ( 1948 ) Zur mechanischen kubatur. Monatsh. Math. , 52 , 286 – 300 . Google Scholar Crossref Search ADS Schmid H. J. & Xu Y. ( 1994 ) On bivariate Gaussian cubature formulae. Proc. Amer. Math. Soc. , 122 , 833 – 842 . Google Scholar Crossref Search ADS Stroud A. ( 1971 ) Approximate Calculation of Multiple Integrals . NJ : Prentice-Hall, Inc . Verlinden P. & Cools R. ( 1992 ) On cubature formulae of degree |$4k+1$| attaining Möller’s lower bound for integrals with circular symmetry. Numer. Math. , 61 , 395 – 407 . Google Scholar Crossref Search ADS Xu Y. ( 2012 ) Minimal cubature rules and polynomial interpolation in two variables. J. Approx. Theory , 164 , 6 – 30 . Google Scholar Crossref Search ADS © The authors 2017. Published by Oxford University Press on behalf of the Institute of Mathematics and its Applications. All rights reserved. This article is published and distributed under the terms of the Oxford University Press, Standard Journals Publication Model (https://academic.oup.com/journals/pages/open_access/funder_policies/chorus/standard_publication_model)
IMA Journal of Numerical Analysis – Oxford University Press
Published: Jan 25, 2019
Read and print from thousands of top scholarly journals.
Already have an account? Log in
Bookmark this article. You can see your Bookmarks on your DeepDyve Library.
To save an article, log in first, or sign up for a DeepDyve account if you don’t already have one.
Copy and paste the desired citation format or use the link below to download a file formatted for EndNote
Access the full text.
Sign up today, get DeepDyve free for 14 days.
All DeepDyve websites use cookies to improve your online experience. They were placed on your computer when you launched this website. You can change your cookie settings through your browser.