MSC. Primary 65C20, 60F05; Secondary 05C80, 60D05.

Keywords. Random polygons; Nonlinear extrapolations; Dirichlet distribution; Central limit theorems; Cramér’s theorem.

1 Introduction

Given a convex set \( K \subseteq \mathbb {R}^d \), the stochastic properties of the convex hull \( K_n \) generated by \(n\) independent random points on \(K\), such as the area, volume and number of vertices of \( K_n \), their probability distributions and asymptotic behavior have attracted extensive attention (see, e.g., [ 7 , 10 , 11 , 16 , 17 , 20 ] ). In the case of \(n\) points randomly selected on a unit circle in \( \mathbb {R}^2 \), the resulting convex hull is a random \(n\)-gon inscribed in the circle which is obtained by connecting all adjacent vertices on the circle. Using the same set of random points, one may also construct (w.p.\(1\)) a random circumscribing \(n\)-gon which is tangent to the circle at each of the \(n\) random points. In the simplest case when the vertices are independent and uniformly distributed on the circle, it is known that the semiperimeter \( S_n \) and area \( A_n \) of such random inscribed polygons, and the semiperimeter (or area) \( S'_n \) of the random circumscribing polygons all converge to \(\pi \) w.p.\(1\) as \( n \to \infty \) and their distributions are also asymptotically Gaussian [ 4 , 24 ] . Furthermore, by using extrapolation techniques [ 12 , 15 ] originating exactly from the famous Archimedean approximations of \( \pi \) based on regular polygons [ 3 , 13 , 19 ] , it has been shown [ 22 , 23 , 25 ] that simple weighted averages such as \( \frac{4}{3} S_n - \frac{1}{3} A_n \), \( \frac{2}{3} S_n + \frac{1}{3} S_n' \) and \( \frac{16}{15} S_n - \frac{1}{5} A_n + \frac{2}{15} S_n' \), etc., provide much more accurate approximations of \( \pi \), and at the same time also satisfy similar central limit theorems as \( n \to \infty \). We note that extrapolation methods are useful in many important applications such as numerical evaluation of integrals, numerical solution of differential equations, and polynomial interpolations, etc. To accelerate the convergence associated with existing low-precision approximations, extrapolation seeks to combine them in a way such that the leading order error terms are cancelled out as much as possible. For example, in the simpler case of the Archimedean approximation of \( \pi \), while both \( S_n = n \sin (\pi /n) \) and \( A_n = \frac{1}{2} n \sin (2 \pi /n) = S_{n/2} \) converge to \( \pi \) with errors of order \( {\mathcal O}(n^{-2}) \), a closer look reveals that \( S_n - \pi \approx - \frac{\pi ^3}{6n^2} \) and \( A_n - \pi \approx - \frac{2 \pi ^3}{3n^2} \). This implies the error of \( A_n \) roughly quadruples that of \( S_n \). Consequently, the weighted average \( \frac{4}{3} S_n - \frac{1}{3} A_n \), or equivalently, \( \frac{4}{3} S_n - \frac{1}{3} S_{n/2} \), exactly cancels out the leading order error terms in \( S_n \) and \( A_n \) to yield an improved estimate of \( \pi \) with a reduced error now of order \( {\mathcal O}(n^{-4}) \).

More recently, in [ 26 ] , the authors have initiated the study of novel nonlinear extrapolation estimates of \( \pi \) in the forms \( \mathcal{X}_n = S_n^{\alpha } A_n^{\beta } \) and \( \mathcal{Y}_n (p) = (\alpha S_n^p + \beta A_n^p)^{1/p} \) where \( \alpha + \beta = 1 \) and \( p \neq 0 \). By deriving probabilistic asymptotic expansions with carefully controlled error estimates, it is shown that, for both \( \mathcal{X}_n \) and \( \mathcal{Y}_n (p) \), the same choice \( \alpha = 4/3 \), \( \beta = -1/3 \) minimizes the approximation error with \( \mathcal{X}_n = \pi + n^{-3+\delta } o(1) \), \( \mathcal{Y}_n (p) = \pi + n^{-3+\delta } o(1) \) where \( \delta {\gt} 0 \) is any positive number and \( o(1) \) represents a random variable which converges to \(0\) w.p.\(1\) as \( n \to \infty \). Furthermore, \( \mathcal{X}_n \) and \( \mathcal{Y}_n (p) \) are also asymptotically normal with \( \mathcal{X}_n \sim {\rm AN} (\pi - 2 \pi ^5 / n^4, 2496 \pi ^{10} / n^9) \), \( \mathcal{Y}_n(p) \sim {\rm AN} (\pi - 2(p+1) \pi ^5 / n^4, (160p^2 + 960p + 2496 ) \pi ^{10} / n^9 ) \) where for a sequence of random variables \( \{ Z_n \} \) and \( \mu _n \in \mathbb {R} \), \( \sigma _n {\gt} 0 \), the notation \( Z_n \sim {\rm AN} (\mu _n, \sigma _n^2 ) \) means \( \left( Z_n - \mu _n \right) / \sigma _n \stackrel{\mathcal{L}}{\longrightarrow } \mathcal{N} (0, 1) \). In particular, for \( p=1 \), \( \mathcal{Y}_n (1) \) reduces to the optimal linear extrapolation estimate \( \frac{4}{3} S_n - \frac{1}{3} A_n \). Moreover, \( \mathcal{X}_n \) may be viewed as the limit of \( \mathcal{Y}_n (p) \) when \( p \to 0 \), a reflection of the relation \( \lim _{ p \to 0 } (\alpha x^p + \beta y^p )^{1/p} = x^{\alpha } y^{\beta } \) for any \( x, y {\gt} 0 \) and \( \alpha + \beta = 1 \).

In this paper, we aim to further develop nonlinear random extrapolation methods for approximating \( \pi \). On the one hand, it would be natural to include also random circumscribing polygons in the approximation process. Motivated by the work in [ 26 ] , we study nonlinear functions of \( S_n \), \( A_n \) and \( S_n' \) in the forms \( \mathcal{W}_n = S_n^\alpha A_n^\beta S_n'^{\gamma } \), \( \mathcal{W}_n(p) = (\alpha S_n^p + \beta A_n^p + \gamma S_n'^{\, p})^{1/p} \) where \( \alpha + \beta + \gamma = 1 \) and \( p \neq 0 \). On the other hand, we are also interested in extending the theory to more general random cyclic polygons whose vertices are not independently and uniformly distributed on the circle. While this is a very challenging problem in general, as a first step, we focus on the particular case of random cyclic polygons generated from symmetric Dirichlet distributions with an arbitrary concentration parameter \( a{\gt}0 \). We note that in such cases, it has been proved [ 21 ] that the respective semiperimeters and areas, again denoted by \( S_n \), \( A_n \) and \( S_n' \), satisfy similar convergence estimates and central limit theorems as in [ 4 , 24 ] for the uniform case, which in fact corresponds to the special case \( a=1 \) of the Dirichlet distribution.

Clearly, as in [ 26 ] , the case \( p = 1 \) reduces to linear extrapolations based on \( S_n \), \( A_n \) and \( S_n' \), and due to the relation \( \lim _{p \to 0} \, (\alpha x^p + \beta y^p + \gamma z^p)^{1/p} = x^\alpha y^\beta z^\gamma \), we expect to recover \( \mathcal{W}_n \) from \( \mathcal{W}_n (p) \) in the limit \( p \to 0 \). More importantly, based on similar asymptotic expansion results established in [ 21 ] for \( S_n \), \( A_n \) and \( S_n' \) in terms of various power sums of the underlying Dirichlet distribution, we derive rigorous probabilistic asymptotic expansions with carefully controlled error estimates for various nonlinear functions of \( S_n \), \( A_n \) and \( S_n' \), particularly \( \mathcal{W}_n \) and \( \mathcal{W}_n (p) \). Such probabilistic asymptotic expansions resemble the well-known Taylor series expansions in many deterministic approximation problems and provide a cornerstone for establishing the corresponding probability convergence estimates and central limit theorems.

It turns out that, for both \( \mathcal{W}_n \) and \( \mathcal{W}_n(p) \), the optimal approximation occurs when \( \alpha + 4 \beta - 2 \gamma = 0 \) with \( \mathcal{W}_n = \pi + n^{-3 + \delta } o(1)\) and \( \mathcal{W}_n (p) = \pi + n^{-3 + \delta } o(1) \). Such results are comparable with those obtained in [ 26 ] for the case \( \gamma = 0 \) and are actually weaker than the corresponding optimal linear extrapolation estimate \( \mathcal{W}_n(1) = \frac{16}{15} S_n - \frac{1}{5} A_n + \frac{2}{15} S_n' = \pi + n^{-5 + \delta } o(1) \), see Theorems 11–13 below for details. Note that together with \( \alpha + \beta + \gamma = 1 \), the condition \( \alpha + 4 \beta - 2 \gamma = 0 \) implies that \( \alpha = 4/3 - 2 \gamma \), \( \beta = - 1/3 + \gamma \) where \( \gamma \) is an arbitrary constant. Due to complicated nonlinear effects, however, the extra “free" parameter \( \gamma \) can no longer be used to further improve the approximation associated with \( \mathcal{W}_n \) and \( \mathcal{W}_n (p) \). Nevertheless, by further combining such nonlinear extrapolation estimates with different values of \( \gamma \), it is possible achieve additional improvements better than the linear case.

Finally, it is interesting to note that in the case of the classical Archimedean polygons, for both \( \mathcal{W}_n \) and \( \mathcal{W}_n (p) \), the optimal estimates also occur when \( \alpha + 4 \beta - 2 \gamma = 0 \) with \( \mathcal{W}_n = \pi + \frac{1}{45} \frac{\pi ^5}{n^4} + {\mathcal O}(n^{-6}) \), \( \mathcal{W}_n (p) = \pi + \frac{1}{180} \frac{\pi ^5}{n^4} [ 45 \gamma p - 10 p + 4] + {\mathcal O}(n^{-6}) \). In fact, for \( \mathcal{W}_n (p) \), by choosing \( \gamma = \frac{10 p - 4}{ 45 p } \) and \( p = p_{\pm } = \frac{-21 \pm \sqrt{721}}{70} \), we can further obtain \( \mathcal{W}_n (p_{\pm }) = \pi -\frac{119 \pm \sqrt{721}}{661500} \frac{\pi ^9}{n^8} + {\mathcal O}(n^{-10}) \), which is two orders of magnitude higher than the optimal linear extrapolation estimate \( \mathcal{W}_n (1) = \frac{16}{15} S_n - \frac{1}{5} A_n + \frac{2}{15} S_n' = \pi + \frac{\pi ^7}{105 n^6} + {\mathcal O}(n^{-8}) \). However, for \( \mathcal{W}_n \), the result turns out to be completely independent of \( \gamma \). This is due to the relation \( A_n S_n' = S_n^2 \), a variant of Archimedes’s celebrated geometric mean relation, which implies \( \mathcal{W}_n = S_n^\alpha A_n^\beta S_n'^{\, \gamma } = S_n^{\alpha + 2\gamma } A_n^{\beta - \gamma } \).

The remainder of the paper is organized as follows. In section 2, we present some useful preliminary results related to the Dirichlet distribution and its various power sums. section 3 is devoted to the study of nonlinear extrapolation estimates of \( \pi \) based on random inscribed and circumscribing polygons generated from symmetric Dirichlet distributions. Finally, we offer several additional remarks in section 4 and some concrete numerical simulation results in section 5 to conclude our study on nonlinear random extrapolation approximations.

2 PRELIMINARIES

2.1 Basic properties of Dirichlet distributions

Recall that a random vector \( {\bf Y}'=(Y_1, \cdots , Y_{n-1}) \in \mathbb {R}^{n-1} \), \( n \geq 2 \), is said to have Dirichlet distribution [ 2 ] with parameters \( {\bf a}=(a_1, \cdots , a_{n-1}; a_n) \in \mathbb {R}^{n+} \) if it has joint probability density function

where \( y_i{\gt}0 \), \( \sum _{i=1}^{n-1} y_i {\lt} 1 \), \( y_n = 1 - \sum _{i=1}^{n-1} y_i \), and \( \Gamma (a) = \int _0^{\infty } v^{a-1} e^{-v} \, dv \) is the gamma function defined for all \( a{\gt}0 \). Let \( Y_n = 1 - \sum _{i=1}^{n-1} Y_i \). With slight abuse of notation, we also refer to \( {\bf Y}=({\bf Y}', Y_n) \in \mathbb {R}^n \) as Dirichlet distribution and write \( {\bf Y} \sim \mathbb {D}{\rm ir}({\bf a}) \), \( {\bf Y}' \sim \mathbb {D}{\rm ir}' ({\bf a}) \).

In this paper, we focus on symmetric Dirichlet distributions, that is, \( a_i = a {\gt} 0 \) for all \( 1 \leq i \leq n \). In such cases, all \( Y_i \sim \mathbb {B}{\rm eta}(a,(n-1)a) \) have identical Beta distribution.

Let \( \mathcal{D}_{n,k} = \sum ^n_{i=1} Y_i^k \), \( k \in \mathbb {N} \). Then \( n^{-(k-1)} \leq \mathcal{D}_{n,k} \leq 1 \) w.p.1. Furthermore, for large \(n\), by using the above Dirichlet integrals, it is easy to compute \( \mathbb {E}(\mathcal{D}_{n,k}) = n \gamma _k (a) / \gamma _k (na) \approx n^{-(k-1)} m_k \) and \( \mathbb {V}{\rm ar}(\mathcal{D}_{n,k}) \approx n^{-(2k-1)} \sigma _k^2 \) where \( \gamma _k(a) = \Gamma (a)^{-1} \Gamma (a+k) = a (a + 1) \cdots (a + k - 1) \), \( m_k = a^{-k} \gamma _k(a) \), and \( \sigma _k^2 = m_{2k} - (1 + k^2/a) m_k^2 \).

The following two lemmas provide additional asymptotic convergence results for various nonlinear expressions of \( \mathcal{D}_{n,k} \). Their proofs are slightly lengthy and are deferred to the Appendix.

2.2 Random cyclic polygons under symmetric Dirichlet distributions

The Dirichlet distribution \( {\bf Y} \sim \mathbb {D}{\rm ir}({\bf a}) \) is naturally associated with the (non-uniform) random division \( 0 = X_0 {\lt} X_1 {\lt} \cdots {\lt} X_{n-1} {\lt} X_n = 1 \) of the unit interval where \( X_0 = 0 \) and \( X_i = \sum _{j=1}^i Y_j \) for \( 1 \leq i \leq n \). In the special case \( a=1 \), this corresponds to the classical uniform random division [ 6 , 14 ] generated by \( n-1 \) independent and uniformly distributed random points on \( (0, 1) \). With the rescaling \( X_i \mapsto \theta _i = 2 \pi X_i \), this can be further mapped to a random division of the unit circle, separated by points \( P_i (\cos \theta _i, \sin \theta _i) \), \( 0 \leq i\leq n \), in counterclockwise direction where \( P_n \) represents the same point as \( P_0 \). By connecting these points consecutively, we obtain an inscribed random \(n\)-gon with its semiperimeter \( S_n \) and area \( A_n \) given by \( S_n = \sum ^n_{i=1} \sin \pi (X_i-X_{i-1}) = \sum _{i=1}^n \sin \pi Y_i \), \( A_n = \frac{1}{2} \sum ^n_{i=1} \sin 2 \pi (X_i-X_{i-1}) = \frac{1}{2} \sum ^n_{i=1} \sin 2 \pi Y_i \). Similarly, using the same random vertices \( P_i \), we can also construct w.p.\(1\) a circumscribing random \(n\)-gon which is tangent to the circle at each point \( P_i \) with its semiperimeter and area both given by \( S_n' = \sum _{i=1}^n \tan \pi (X_i-X_{i-1}) = \sum _{i=1}^n \tan \pi Y_i \). Note that in the event (which has probability \(0\)) all vertices are equally spaced, that is, \( Y_i =1/n \) for all \( 1 \leq i \leq n \), these random \(n\)-gons happen to be regular \(n\)-gons inscribed in or circumscribed about the circle with \( S_n = n \sin (\pi /n) \), \( A_n = \frac{1}{2} n \sin (2\pi /n) = S_{n/2} \) and \( S_n' = n \tan (\pi /n) \). Additionally, such random cyclic polygons generated from the symmetric Dirichlet distribution \( {\bf Y} \sim \mathbb {D}{\rm ir}({\bf a}) \) also degenerate to regular \(n\)-gons in the limit as \( a \to \infty \).

By using the Taylor series expansion of the sine and tangent functions, it is easy to obtain, at least formally, the following probabilistic asymptotic expansions for \( S_n \), \( A_n \) and \( S_n' \):

where \( B_j \) is the \( j \)th Bernoulli number. Note that by Lemma 3, the random infinitesimal terms \( \mathcal{D}_{n,k} \) in the above expansions decrease progressively in order of magnitude. The validity of these asymptotic expansions is rigorously justified by the following lemma.

3 Nonlinear extrapolation estimates

3.1 Probabilistic asymptotic expansions for nonlinear functions of \( S_n \), \( A_n \) and \( S_n' \)

In this section, we study nonlinear random extrapolation estimates of \(\pi \) based on the semiperimeters and areas of both inscribed and circumscribed random polygons with an aim to construct more accurate nonlinear extrapolation estimates than in [ 26 ] . To facilitate the derivation of asymptotic expansions for nonlinear functions of \( S_n \), \( A_n \) and \( S_n' \) in the forms \( \mathcal{W}_n = S_n^\alpha A_n^\beta S_n'^{\gamma } \) and \( \mathcal{W}_n(p) = (\alpha S_n^p + \beta A_n^p + \gamma S_n'^{\, p})^{1/p} \) where \( \alpha + \beta + \gamma = 1 \) and \( p \neq 0 \), we follow the development in [ 26 ] and write

We mention that while the analysis in [ 26 ] is carried out for uniform random divisions only, with ??, it is straightforward to verify that exactly the same asymptotic expansion results in fact extend to the case of symmetric Dirichlet distributions \( {\bf Y} \sim \mathbb {D} {\rm ir} ({\bf a}) \) for arbitrary \( a{\gt}0 \). Note also that, in view of ??, we have, for any \( \delta {\gt} 0 \), \( \mathcal{W}_{n,1} = n^{-1+\delta } o(1) \), \( \mathcal{U}_{n,p,1} = n^{-1+\delta } o(1) \), \( \mathcal{W}_{n,2} = n^{-3+\delta } o(1) \), \( \mathcal{U}_{n,p,2} = n^{-3+\delta } o(1) \), \( \mathcal{W}_{n,3} = n^{-5+\delta } o(1) \), \( \mathcal{U}_{n,p,3} = n^{-5+\delta } o(1) \). Thus, the above probabilistic asymptotic expansions in lemma 8 imply, \( \log \left( S_n / \pi \right) = - \pi ^2 \mathcal{W}_{n,1} + n^{-3 + \delta } o(1) = - \pi ^2 \mathcal{W}_{n,1} - \pi ^4 \mathcal{W}_{n,2} + n^{-5 + \delta } o(1) \) and \( ( S_n / \pi )^p = 1 - \mathcal{U}_{n,p,1} \pi ^2 + n^{-3 + \delta } o(1) = 1 - \mathcal{U}_{n,p,1} \pi ^2 + \mathcal{U}_{n,p,2} \pi ^4 + n^{-5 + \delta } o(1) \).

Next, to derive asymptotic expansions for nonlinear functions of \( S_n' \) such as \( \log (S_n' / \pi ) \) and \( S_n'^{\, p} \), we apply the same Taylor series expansions of \( \log (1+x) \) and \( (1+x)^p \) on \( S_n' / \pi - 1 = \frac{1}{3} \pi ^2 \mathcal{D}_{n,3} + \frac{2}{15} \pi ^4 \mathcal{D}_{n,5} + \frac{17}{315} \pi ^6 \mathcal{D}_{n,7} + n^{-7 + \delta } o(1) \). Similar to lemma 8, we now obtain

where

With \( \mathcal{T}_{n,1} = \pi \) and \( \mathcal{D}_{n,k} \leq \Delta _n^{k-1} \leq t^{k-1} \), it is clear that we may choose \( t \) suitably small such that \( 0 {\lt} S_n'/\pi - 1 \leq 1/2 \) and \( 0 \leq \mathcal{T}_{n,2m-1}/\pi - 1 \leq 1/2 \). By the mean value theorem, we thus obtain

We now take \( m = 4 \). By inserting \( \mathcal{T}_{n,7}/\pi -1 = \tfrac {1}{3} \pi ^2 \mathcal{D}_{n,3} + \tfrac {2}{15} \pi ^4 \mathcal{D}_{n,5} + \tfrac {17}{315} \pi ^6 \mathcal{D}_{n,7} = n^{-1+\delta } o(1) \) into the Taylor series approximation \( \log (1+x) = x - \frac{1}{2} x^2 + \frac{1}{3} x^3 + {\mathcal O}(1) x^4 \) for \( |x| \leq 1/2 \) and keeping only terms at order \( n^{-7+\delta } \) and below, we obtain \( \log (\mathcal{T}_{n,7}/\pi ) 1_{ \{ \Delta _n \leq t \} } = \frac{1}{3} \pi ^2 \mathcal{D}_{n,3} + \left(-\frac{1}{18} \mathcal{D}_{n,3}^2 + \frac{2}{15} \mathcal{D}_{n,5}\right) \pi ^4 + \left(\frac{1}{81} \mathcal{D}_{n,3}^3 - \frac{2}{45} \mathcal{D}_{n,3} \mathcal{D}_{n,5} + \frac{17}{315} \mathcal{D}_{n,7}\right) \pi ^6 + n^{-7 + \delta } o(1) = \mathcal{M}_{n,1} \pi ^2 + \mathcal{M}_{n,2} \pi ^4 + \mathcal{M}_{n,3} \pi ^6 \) \(+ n^{-7+\delta } o(1). \) Consequently, we have \( \log (S_n'/\pi ) = \mathcal{M}_{n,1} \pi ^2 + \mathcal{M}_{n,2} \pi ^4 + \mathcal{M}_{n,3} \pi ^6 + n^{-7+\delta } o(1) \). Similarly, we can verify \(( S_n' / \pi )^p = 1 + \mathcal{V}_{n,p,1} \pi ^2 + \mathcal{V}_{n,p,2} \pi ^4 + \mathcal{V}_{n,p,3} \pi ^6 + n^{-7+\delta } o(1)\).

3.2 Nonlinear extrapolations of the form \( \mathcal{W}_n = S_n^\alpha A_n^\beta S_n'^{\, \gamma } \)

With the above preparations, we are now ready to derive probabilistic asymptotic expansions for the nonlinear extrapolation estimates \( \mathcal{W}_n = S_n^\alpha A_n^\beta S_n'^{\, \gamma } \) and \( \mathcal{W}_n(p) = ( \alpha S_n^p + \beta A_n^p + \gamma S_n'^{\, p} )^{1/p}\) where \( \alpha + \beta + \gamma = 1 \) and \( p \neq 0 \). By taking the exponential of the linear combination of \( \log (S_n / \pi ) \), \( \log (A_n / \pi ) \) and \( \log (S_n' / \pi ) \) in 0, 1 and 4, or by multiplying \( S_n^\alpha \), \( A_n^\beta \) and \( S_n'^{\, \gamma } \) directly from 2, 3, and 5, we obtain

Let \( \eta = \alpha + 4\beta - 2\gamma \). Then with \( \alpha + \beta + \gamma = 1 \), we may write \( \alpha = \frac{4}{3} - 2\gamma -\frac{1}{3} \eta \), \( \beta = -\frac{1}{3} + \gamma + \frac{1}{3} \eta \). Clearly if \( \eta \neq 0 \), we have \( \mathcal{W}_n = \pi + n^{-1 + \delta } o(1) \, = \pi - \frac{1}{6} \eta \pi ^3 \mathcal{D}_{n,3} + n^{-3 + \delta } o(1). \) Then by Slutsky’s theorem, it follows that \( \mathcal{W}_n(p) \sim {\rm AN} (\pi - \frac{1}{6} n^{-2} \eta \pi ^3 m_3, \frac{1}{36} n^{-5} \eta ^2 \pi ^6 \sigma _3^2) \).

However, if \( \eta = 0 \), that is, \( \alpha = 4/3 - 2 \gamma \), \( \beta = -1/3 + \gamma \) and \( \gamma \) is an arbitrary constant, it is then possible to eliminate the leading error term involving \( \mathcal{D}_{n,3} \) in 5 to obtain

In particular, this implies \( \mathcal{W}_n = \pi - \tfrac {1}{4} \pi ^5 \left[ \left( \gamma - \tfrac {2}{9} \right) \mathcal{D}_{n,3}^2 +\left( \frac{2}{15} - \gamma \right) \mathcal{D}_{n,5} \right] + n^{-5 + \delta } o(1) \). By lemma 5 and remark 6, it is clear that \( T_n=T_n(\gamma - \tfrac {2}{9}, \tfrac {2}{15} - \gamma )= \left( \gamma - \tfrac {2}{9} \right) \mathcal{D}_{n,3}^2 + \left( \tfrac {2}{15} - \gamma \right) \mathcal{D}_{n,5} \) is nondegenerate for all \( \gamma \in \mathbb {R} \) and is asymptotically normal with \( T_n \sim {\rm AN} \left( n^{-4} \mu _T(\gamma - \tfrac {2}{9}, \tfrac {2}{15} - \gamma ), n^{-9} \sigma _T^2 (\gamma - \tfrac {2}{9}, \tfrac {2}{15} - \gamma ) \right). \)

By further applying Slutsky’s theorem, a related central limit theorem can be established for the optimal nonlinear extrapolation estimate \( \mathcal{W}_n = \) \(S_n^{4/3 - 2 \gamma } A_n^{-1/3 + \gamma } S_n'^{\, \gamma } \). Note that in such cases, it is impossible to take advantage of \( \gamma \) to further eliminate the leading order error term in 6 to achieve \( \mathcal{W}_n = \pi + n^{-5 + \delta } o(1) \).

3.3 Nonlinear extrapolations of the form \( \mathcal{W}_n(p) = (\alpha S_n^p + \beta A_n^p + \gamma S_n'^{\, p})^{1/p}\)

Next, we consider \( \mathcal{W}_n(p) = (\alpha S_n^p + \beta A_n^p + \gamma S_n'^{\, p})^{1/p} \). By taking the linear combination of \( (S_n / \pi )^p \), \( (A_n / \pi )^p \) and \( (S_n' / \pi )^p \) in 2, 3 and 5, and applying Newton’s generalized binomial formula for \( (1+x)^{1/p} \), we derive

Again, let \( \eta = \alpha + 4\beta - 2\gamma \) so that \( \alpha = \frac{4}{3} - 2\gamma - \frac{1}{3} \eta \), \( \beta = -\frac{1}{3} + \gamma + \frac{1}{3} \eta \). Thus as in the case of \( \mathcal{W}_n \), if \( \eta \neq 0 \), then \( \mathcal{W}_n(p) = \pi - \frac{1}{6} \eta \pi ^3 \mathcal{D}_{n,3} + n^{-3 + \delta } o(1) \) and \( \mathcal{W}_n(p) \sim {\rm AN} (\pi - \frac{1}{6} n^{-2} \eta \pi ^3 m_3, \frac{1}{36} n^{-5} \eta ^2 \pi ^6 \sigma _3^2) \). However, if \( \eta = 0 \), that is, \( \alpha = 4/3 - 2 \gamma \), \( \beta = -1/3 + \gamma \), a further improvement as in theorem 11 is possible. In such cases, we have

4 Additional Remarks

We offer some additional remarks to conclude our study on the nonlinear random extrapolation estimates \( \mathcal{W}_n = S_n^{\alpha } A_n^{\beta } S_n'^{\, \gamma } \) and \( \mathcal{W}_n (p) = (\alpha S_n^p + \beta A_n^p + \gamma S_n'^{\, p})^{1/p} \). For brevity, we address for both \( \mathcal{W}_n \) and \( \mathcal{W}_n (p) \) the optimal case only with \( \alpha + 4 \beta - 2 \gamma = 0 \), that is, \( \alpha = 4/3 - 2 \gamma \), \( \beta = -1/3 + \gamma \).

4.1 Special cases of \( \alpha = 0 \), \( \beta = 0 \), or \( \gamma = 0 \)

Note that for \( \gamma = 0 \), \( \gamma = 1/3 \) and \( \gamma = 2/3 \), \( \mathcal{W}_n \) reduces to \( \mathcal{X}_n = S_n^{4/3} A_n^{-1/3} \), \( \mathcal{Y}_n = S_n^{2/3} S_n'^{\, 1/3} \) and \( \mathcal{Z}_n = A_n^{1/3} S_n'^{\, 2/3} \) respectively with \( \mathcal{X}_n = \pi + \frac{\pi ^5}{90} \left( 5 \mathcal{D}_{n,3}^2 - 3 \mathcal{D}_{n,5} \right) + n^{-5 + \delta } o(1) \), \( \mathcal{Y}_n = \pi - \frac{\pi ^5}{180} \left( 5 \mathcal{D}_{n,3}^2 - 9 \mathcal{D}_{n,5} \right) + n^{-5 + \delta } o(1) \), \( \mathcal{Z}_n = \pi - \frac{\pi ^5}{45} \left( 5 \mathcal{D}_{n,3}^2 - 6 \mathcal{D}_{n,5} \right) + n^{-5 + \delta } o(1) \).

Similarly, for \( \mathcal{W}_n (p) \), the same choices \( \gamma = 0 \), \( \gamma = 1/3 \), \( \gamma = 2/3 \) yield \( \mathcal{X}_n (p) = \left( \frac{4}{3} S_n^p - \frac{1}{3} A_n^p \right)^{1/p} = \pi - \frac{1}{90} \pi ^5 \left[ 5(p-1) \mathcal{D}_{n,3}^2 + 3 \mathcal{D}_{n,5} \right] + n^{-5 + \delta } o(1)\), \( \mathcal{Y}_n (p) = \left(\frac{2}{3} S_n^p + \frac{1}{3} S_n'^{\, p} \right)^{1/p} = \pi + \frac{1}{180} \pi ^5 \left[ 5(p-1) \mathcal{D}_{n,3}^2 + 9 \mathcal{D}_{n,5} \right] + n^{-5 + \delta } o(1) \), \( \mathcal{Z}_n (p) = \left(\frac{1}{3} A_n^p + \frac{2}{3} S_n'^{\, p}\right)^{1/p} = \pi + \frac{1}{45} \pi ^5 \left[ 5(p-1) \mathcal{D}_{n,3}^2 + 6 \mathcal{D}_{n,5} \right] + n^{-5 + \delta } o(1)\).

4.2 Uniform spacings

In the special case of \( a=1 \), the symmetric Dirichlet distribution corresponds to the uniform spacings generated by \( n-1 \) independent and uniformly distributed random points on the unit interval. Thus, by setting \( a=1 \), we immediately obtain optimal nonlinear extrapolation estimates for random polygons generated by independent and uniformly distributed random points on the unit circle.

4.3 Regular polygons

Finally, we remark that in the case of regular polygons, with \( Y_i = 1/n \) and \( \mathcal{D}_{n,k} = n^{-(k-1)} \), it is straightforward to check that the optimal estimate for both \( \mathcal{W}_n \) and \( \mathcal{W}_n (p) \) occurs also when \( \alpha + 4 \beta - 2 \gamma = 0 \) with

Note that the result for \( \mathcal{W}_n \) actually does not depend on the parameter \( \gamma \) at all. This is due to the fact \( A_n S_n' = S_n^2 \), a variant of Archimedes’s celebrated geometric mean relation, which implies \( \mathcal{W}_n = S_n^\alpha A_n^\beta S_n'^{\, \gamma } = S_n^{\alpha + 2\gamma } A_n^{\beta - \gamma } \). In particular, in such cases, \( \mathcal{X}_n = S_n^{4/3} A_n^{-1/3} \), \( \mathcal{Y}_n = A_n^{2/3} S_n'^{\, 1/3} \) and \( \mathcal{Z}_n = A_n^{1/3} S_n'^{\, 2/3} \) all yield exactly the same result. While this relation no longer holds for random polygons, it helps explain why throwing in an extra term does not always increase the accuracy of the extrapolation estimates.

For \( \mathcal{W}_n(p) \), however, by taking \( \gamma = \frac{10 p - 4}{45 p} \), it is possible to eliminate the leading error term \( {\mathcal O}(n^{-4}) \) to obtain \( \mathcal{W}_n(p) = \pi + \frac{1}{5670} \frac{\pi ^7}{n^6} (35p^2 + 21p-2) - \frac{1}{32400} \frac{\pi ^9}{n^8} (25p^3 - 75p^2-52p+12) + {\mathcal O}(n^{-10}) \). Thus if we choose \( p_{\pm } = \frac{-21 \pm \sqrt{721}}{70} \), we can further obtain \( \mathcal{W}_n (p_{\pm }) = \pi -\frac{119 \pm \sqrt{721}}{661500} \frac{\pi ^9}{n^8} + {\mathcal O}(n^{-10}) \), which is two orders of magnitude higher than the optimal linear estimate \( \mathcal{W}_n (1) = \frac{16}{15} S_n - \frac{1}{5} A_n + \frac{2}{15} S_n' = \pi + \frac{1}{105} \frac{\pi ^7}{n^6} + \frac{1}{360}\frac{\pi ^9}{n^8} + {\mathcal O}(n^{-10}) \).

5 Numerical simulations

In this section, we present numerical simulation results to confirm the main probabilistic convergence estimates obtained in this paper. For this purpose, we use the MATLAB command gamrnd(a,1) to first generate \(n\) independent gamma random variables \( {\bf V} = (V_1, V_2, \cdots V_n) \) with shape parameter \( a{\gt}0 \). The symmetric Dirichlet random vector \( {\bf Y} \sim \mathbb {D}{\rm ir} (a {\bf 1}) \) is then obtained by normalization \( {\bf Y} = {\bf V} / \| {\bf V} \| _{1} \). Next, we compute \( S_n = \sum _{i=1}^n \sin \pi Y_i \), \( A_n = \frac{1}{2} \sum _{i=1}^n \sin 2 \pi Y_i \), \( S_n' = \sum _{i=1}^n \tan \pi Y_i \), and subsequently, \( \mathcal{W}_n = S_n^{\alpha } A_n^{\beta } S_n'^{\gamma } \), \( \mathcal{W}_n (p) = (\alpha S_n^p + \beta A_n^p + \gamma S_n'^p)^{1/p} \) for \( n=96 \times 2^k \), \( k=4,5,6,7 \). For simplicity, we consider two \(p\) values: \( p=2 \), \( p=-1 \) and choose \( \alpha = 16/15 \), \( \beta = -1/5 \), \( \gamma = 2/15 \) which clearly satisfies the optimality condition \( \eta = \alpha + 4 \beta - 2 \gamma = 0 \). We repeat these simulations for \( m=100,000 \) times. For each of these random samples, we compute its empirical mean \( \hat{\mu } (X_n) \) and empirical standard deviation \( \hat{\sigma } (X_n) \) with normalization \( \tilde{X}_n = (X_n - \hat{\mu } (X_n))/\hat{\sigma }(X_n) \). The histograms (with bin size 400) for the empirical PDFs of normalized \( S_n \), \( A_n \), \( S_n' \), \( W_n \), \( W_n (p) \) are displayed in figure 1 and figure 2 below for specified parameter values.

To effectively compare the empirical data with their theoretical asymptotic values, suitable scaling factors are used for \( \hat{\mu } (X_n) \) and \( \hat{\sigma } (X_n) \) in the tables below. For \( \mathcal{W}_n \) and \( \mathcal{W}_n(p) \), by ??, it is clear that \( \hat{\mu } (X_n) = \pi + {\mathcal O}(1) n^{-4} \pi ^5 \), and \( \hat{\sigma } (X_n) = {\mathcal O}(1) n^{-9/2} \pi ^5 \). Thus, for easy numerical comparison, we display \( \hat{\mu }_n = n^4 \pi ^{-5} (\hat{\mu } (X_n) - \pi ) \) and \( \hat{\sigma }_n = n^{9/2} \pi ^{-5} \hat{\sigma } (X_n) \) instead, together with similarly scaled limiting values. For \( S_n \), \( A_n \) and \( S_n' \), the scaling factors for the mean and standard deviation are \( n^2 \pi ^{-3} \) and \( n^{5/2} \pi ^{-3} \) respectively.

Finally, we note that for \( p=1 \), the optimal linear extrapolation \( \mathcal{W}_n (1) = \frac{16}{15} S_n - \frac{1}{5} A_n + \frac{2}{15} S_n' \) converges most rapidly with an asymptotical mean \( \hat{\mu } (\mathcal{W}_n (1))\) \( = \pi + {\mathcal O}(1) n^{-6} \) and standard deviation \( \hat{\sigma } (\mathcal{W}_n (1)) = {\mathcal O}(1) n^{-13/2} \). At such “atomic" scales, however, for \(n\) in the range \( 10^3 \sim 10^4 \), the usual double precision computation may not be enough to prevent severe loss of significant digits. In such cases, the distribution of the rescaled simulated data would appear to be more “discrete" with unusually large variance. As a compromise, partial numerical evidence of the convergence results may be witnessed by using relatively smaller values of \(n\) instead. See figure 3 and table 7. Such phenomena also occur to \( \mathcal{W}_n \) and \( \mathcal{W}_n (p) \), but to a much lesser extent.

Appendix

Proof for lemma 4

We use similar ideas in [ 21 ] to prove ??. First, by using the multinomial expansion formula

where \( p \) is a positive integer, \( {\bf m} = (m_1, m_2, \cdots , m_n) \) is a multi-index with each \( m_j \geq 0 \) and \( | {\bf m} | = \sum _{j=1}^n m_j = p \), and \( \binom {p}{\bf m} \) is given by

we obtain for each \( 1 \leq i \leq l \)

where \( {\bf m}_i = (m_{i,1}, m_{i,2}, \cdots , m_{i,n}) \) such that \( | {\bf m}_i | = \sum _{j=1}^n m_{i,j} = p_i \). Then for all \( 1 \leq i \leq l \), we obtain

By using 1 and \( \sum _{j=1}^n \sum _{i=1}^l k_i m_{i,j} = \sum _{i=1}^l \sum _{j=1}^n k_i m_{i,j} = \sum _{i=1}^l k_i p_i \), we obtain

The key is to estimate \( \prod _{j=1}^n \frac{ {\Gamma } ( a + q_j) }{ {\Gamma } (a) } \) with \( q_j = \sum _{i=1}^l k_i m_{i,j} \). Note that

Since

the number \(r\) of indices \( \# \{ 1 \leq j \leq n: q_j \geq 1 \} \) is at most \( m_* \). Then if \( q_j \geq 1 \), we have

where \( \widehat{a} = \lceil a \rceil \) is the smallest integer that is greater than or equal to \( a \), \( {\mathcal O}(1) \) is some positive constant independent of \(n\). By using Stirling’s formula [ 1 ]

as \( z \to \infty \), we obtain

Substituting 3 and 4 into 0, and using

we have

By applying Markov inequality, for any \( \delta {\gt} 0 \) and \( \varepsilon {\gt} 0 \), we obtain

as \( n \to \infty \), which implies \( n^{ \left( \sum _{i=1}^l (k_i -1) p_i \right) - \delta } \prod _{j=1}^l \mathcal{D}_{n,k_j}^{p_j} \to 0 \) in probability as \( n \to \infty \). In addition, we have

By Borel-Cantelli Lemma, it follows that \( n^{ \left( \sum _{i=1}^l (k_i -1) p_i \right) - 1 -\delta } \prod _{i=1}^l \mathcal{D}_{n,k_i}^{p_i} \to 0 \) with probability \(1\) as \( n \to \infty \). This completes the proof of lemma 4.

Proof for lemma 5

By applying the equivalent representation [ 2 ]

where \( V_1, V_2, \cdots , V_n \) are independent gamma random variables with \( V_i \sim \Gamma (a, 0, 1) \), \( i = 1, 2, \cdots , n \), for \( a {\gt} 0 \), we may rewrite

With the above reformulation of \( T_n \), we may consider the joint asymptotic distribution of the three sums \( \sum _{i=1}^n V_i, \sum _{i=1}^n V_i^3, \sum _{i=1}^n V_i^5 \). By using the multivariate central limit theorem, we obtain, as \(n\rightarrow \infty \),

where \( \mu _1 = \mathbb {E}(V_i) = a \), \( \mu _2 = \mathbb {E}(V_i^3) = \gamma _3(a) \), \( \mu _3 = \mathbb {E}(V_i^5) = \gamma _5(a) \), and \( \Sigma \) is the covariance matrix of the random vector \( (V_i, V_i^3, V_i^5) \) with

Next, we apply Cramér’s theorem [ 9 ] to obtain

where \( f \) is a mapping: \( \mathbb {R}^3 \rightarrow \mathbb {R} \) such that \( \nabla f(\mu ) \) is continuous in a neighborhood of \( \mu = (\mu _1, \mu _2, \mu _3) \in \mathbb {R}^3 \), and \( \sigma _F^2 = \nabla f(\mu ) \Sigma (\nabla f(\mu ))^T \). To do so, we choose \( f(x,y,z) = \alpha x^{-6}y^2 + \beta x^{-5}z \) with \( \frac{\partial f}{\partial x} = -6 \alpha x^{-7}y^2 - 5 \beta x^{-6}z \), \( \frac{\partial f}{\partial y} = 2 \alpha x^{-6}y \), \( \frac{\partial f}{\partial z} = \beta x^{-5} \). Then we have

with \( f(\mu ) = \alpha m_3^2 + \beta m_5 = a^{-4} (a+1) (a+2) \big[ \alpha (a+1)(a+2) + \beta (a+3)(a+4) \big] \), \( \nabla f(\mu ) = (-6 \alpha a^{-1} m_3^2 - 5 \beta a^{-1} m_5, 2 \alpha a^{-6} \gamma _3(a), \beta a^{-5}) \). Then \( \sigma _T^2 = \nabla f(\mu ) \Sigma (\nabla f(\mu ))^T = 4 \alpha ^2 m_3^2 \sigma _3^2 + 4 \alpha \beta (m_3 m_8 - (1 + 15/a) m_3^2 m_5) + \beta ^2 \sigma _5^2 = 8 a^{-9} (a + 1) (a + 2) \big[3 \alpha ^2 (a+1)^2 (a+2)^2 (3a+7) + 5 \beta ^2 (a+3) (a+4) (5a^3 + 60a^2 + 250a + 363) + 30 \alpha \beta (a+1) (a+2) (a+3)^2 (a+4)\big] \). This completes the proof of lemma 5.

1

M. Abramowitz, I.A. Stegun, Handbook of Mathematical Functions, Dover, 1970.

2

N. Balakrishnan, V.B. Nevzorov, A Primer on Statistical Distributions, Wiley, 2003, https:doi.org/10.1002/0471722227.

3

P. Beckmann, A History of \( \pi \), Fifth edition, The Golem Press, 1982.

4

C. Bélisle, On the polygon generated by \(n\) random points on a circle, Statist. Probab. Lett., 81 (2011), pp. 236–242, https://doi.org/10.1016/j.spl.2010.11.012.

5

P. Billingsley, Probability and Measure, Third edition, Wiley, 1995.

6

D.A. Darling, On a class of problems related to the random division of an interval, Ann. Math. Statistics, 24 (1953), pp. 239–253, https:doi.org/10.1214/aoms/1177729030.

7

B. Efron, The convex hull of a random set of points, Biometrika, 52 (1965), pp. 331–343, https://doi.org/10.1093/biomet/52.3-4.331.

8

W. Feller, An Introduction to Probability Theory and its Applications, vol. II, Wiley, 1966.

9

T.S. Ferguson, A Course in Large Sample Theory, Chapman & Hall, New York, 1996.

10

T. Hsing, On the asymptotic distribution of the area outside a random convex hull in a disk, Ann. Appl. Probab, 4 (1994) no. 2, pp. 478–493.

11

I. Hueter, The convex hull of a normal sample, Adv. Appl. Probab, 26 (1994) no. 4, pp. 855–875, https://doi.org/10.2307/1427894.

12

D.C. Joyce, Survey of extrapolation processes in numerical analysis, SIAM Rev., 13 (1971), pp. 435–490, https://doi.org/10.1137/1013092.

13

T. Li, Essays about \(\pi \), Higher Education Press, 2007 (Chinese).

14

R. Pyke, Spacings, J. Roy. Statist. Soc. Ser. B, 27 (1965), pp. 395–449, https://doi.org/10.1111/j.2517-6161.1965.tb00602.x.

15

P. Rabinowitz, Extrapolation methods in numerical integration, Numer. Algor., 3 (1992), pp. 17–28, https://doi.org/10.1007/BF02141912.

16

A.Rényi, R. Sulanke, Über die konvexe Hülle von \(n\) zufällig gewählten Punkten I, Z. Wahrscheinlichkeitstheorie und Verw. Gebiete, 2 (1963), pp. 75–84, https://doi.org/10.1007/BF00535300.

17

A. Rényi, R. Sulanke, Über die konvexe Hülle von \(n\) zufällig gewählten Punkten II, Z. Wahrscheinlichkeitstheorie und Verw. Gebiete, 3 (1964), pp. 138–147.

18

R.J. Serfling, Approximation Theorems of Mathematical Statistics, Wiley, 1980, https://doi.org/10.1002/9780470316481.

19

T. Tang, Computational mathematics originated from the computation of \(\pi \), Higher Education Press, 2018 (Chinese).

20

V. Vu, Central limit theorems for random polytopes in a smooth convex set, Adv. Math., 207 (2006) no. 1, pp. 221–243, https://doi.org/10.1016/j.aim.2005.11.011.

21

S. Wang, W.-Q. Xu, Random cyclic polygons from Dirichlet distributions and approximations of \(\pi \), Statist. Probab. Lett., 140 (2018), pp. 84–90, https://doi.org/10.1016/j.spl.2018.05.007.

22

S. Wang, W.-Q. Xu, J. Liu, Random polygons and optimal extrapolation estimates of \( \pi \), Stat. Optim. Inf. Comput., 9 (2021), pp. 241–249, https://doi.org/10.19139/soic-2310-5070-1193.

23

W.-Q. Xu, Extrapolation methods for random approximations of \( \pi \), J. Numer. Math. Stoch., 5 (2013), pp. 81–92.

24

W.-Q. Xu, Random circumscribing polygons and approximations of \( \pi \), Statist. Probab. Lett., 106 (2015), pp. 52–57, https://doi.org/10.1016/j.spl.2015.06.026.

25

W.-Q. Xu, L. Meng, Y. Li, Random polygons and estimations of \( \pi \), Open Math., 17 (2019), pp. 575–581, https://doi.org/10.1515/math-2019-0049.

26

W.-Q. Xu, S. Wang, D. Xu, Nonlinear extrapolation estimates of \(\pi \), Acta Math. Sin. Appl. Engl. Ser., 40 (2024), no. 1, pp. 91–108.