Convergence acceleration of a recursively defined sequence $a_{n+1} = (1-a_n)^{\frac1p},\ p>1$

Question

On this question of recursive sequence, I have proved its convergence. As the sequence oscillates around its limit, the convergence rate can be accelerated. Here is the definition of sequence (series) convergence acceleration. How would one accelerate the convergence rate of this sequence? Please provide proof and references.

Answer 1

The limit of the sequence $$ a_{n+1}=(1-a_n)^{\frac1p}\tag1 $$ if it exists, would be the solution to $$ x^p+x-1=0\tag2 $$ Since $$ \overbrace{1-p^{-\frac{p}{p-1}}}^{x}\gt\overbrace{\ \ p^{-\frac1{p-1}}\ \ }^{(1-x)^{\frac1p}}\tag3 $$ The root of $(2)$ must be less than $1-p^{-\frac{p}{p-1}}$.

For some $\xi$ between $a_{n-1}$ and $a_n$, we have $$ \begin{align} a_{n+1}-a_n &=(1-a_n)^{\frac1p}-(1-a_{n-1})^{\frac1p}\\ &=\frac{(1-\xi)^{\frac1p-1}}p\,(a_{n-1}-a_n)\tag4 \end{align} $$ If $\xi\lt1-p^{-\frac{p}{p-1}}$, $(4)$ represents a contraction map and will converge. This is precisely the condition satisfied by the root of $(2)$, as shown in $(3)$. Thus, in a neighborhood of the limit, the sequence $(1)$ converges geometrically with a ratio close to $\frac{(1-\xi_0)^{\frac1p-1}}p$ (the smaller the neighborhood, the closer the ratio), where $\xi_0$ is the root of $(2)$ in $(0,1)$.

That is, in a neighborhood of the limit, the error in $\boldsymbol{(1)}$ decreases by a factor of $\boldsymbol{\frac{(1-\xi)^{\frac1p-1}}p}$ with each iteration.

---

Newton's Method gives the solution to $(2)$ as the limit of the sequence $$ \begin{align} b_{n+1} &=b_n-\frac{b_n^p+b_n-1}{pb_n^{p-1}+1}\\ &=\frac{(p-1)b_n^p+1}{pb_n^{p-1}+1}\tag{5} \end{align} $$ It can be shown that Newton's Method converges quadratically.

That is, in a neighborhood of the limit, the error in $\boldsymbol{(5)}$ is squared with each iteration.

---

Example: $\boldsymbol{p=2}$

Iterating $(1)$ with $a_0=0.5$ yields $$ \begin{array}{l} 0.5\\ 0.70710678118654752440\\ 0.54119610014619698440\\ 0.67735064763666020695\\ 0.56802231678283538462\\ 0.65725009183503703439\\ 0.58544846755710529028\\ 0.64385676391794992909\\ 0.59677737564526528464\\ 0.63499812941042174015\\ 0.60415384678869525497\\ 0.62916305772931769061\\ 0.60896382673413558793\\ \end{array} $$ This is converging, but quite slowly (about $1$ digit every $4.5$ iterations).

Iterating $(5)$ with $b_0=0.5$ yields $$ \begin{array}{l} 0.5\\ 0.625\\ 0.61805555555555555556\\ 0.61803398895790200138\\ 0.61803398874989484822\\ 0.61803398874989484820\\ 0.61803398874989484820 \end{array} $$ After $5$ iterations, we have $20$ digits of convergence (doubling the number of digits every iteration).

---

Proof of $\boldsymbol{(3)}$

Since $p\gt1$, Bernoulli's Inequality says that $$ (1+(p-1))^{\frac{p}{p-1}}\gt1+p\tag6 $$ Multiplying both sides by $p^{-\frac{p}{p-1}}$ gives $$ \begin{align} 1 &\gt(p+1)\,p^{-\frac{p}{p-1}}\\ &=p^{-\frac1{p-1}}+p^{-\frac{p}{p-1}}\tag7 \end{align} $$ which is equivalent to $(3)$.

Answer 2

From the previous question, we know that

$$a_{2n}<a<a_{2n+1}$$

And so, for any $m>0$, the weighted average of consecutive terms is closer to $a$.

$$a_{2n}<\frac{a_{2n}+ma_{2n+1}}{1+m}<a_{2n+1}$$

So $\frac{a_{2n}+ma_{2n+1}}{1+m}$ is closer to $a$ than $a_{2n}$ and $a_{2n+1}$, yielding faster convergence. Specifically,

\begin{align}\small\lim_{n\to\infty}\frac{a-\frac{a_n+ma_{n+1}}{1+m}}{a-a_n}&=\frac1{1+m}+\frac m{1+m}\lim_{n\to\infty}\frac{a-(1-a_n)^{1/p}}{a-a_n}\\&=\frac1{1+m}-\frac{ma}{p(1+m)(1-a)}\end{align}

For convergence to be fastest, we're interested in

$$0=\frac1{1+m}-\frac{ma}{p(1+m)(1-a)}$$

Which has the solution

$$m=\frac{p(1-a)}a$$

However, in order for this to be useful, you'll need to know what $a$ is. Replacing it with $a_n$ is good enough,

$$a=\lim_{n\to\infty}\frac{a_n^2+p(1-a_n)^{(p+1)/p}}{a+p(1-a)}$$

which converges quadratically. (as fast as Newton's method)

The technique of weighting terms such that the resulting limit is zero can easily be extended to more generally things such as

$$f(x)=\frac{x+m_1(1-x)^{1/p}+m_2(1-(1-x)^{1/p})^{1/p}+\dots}{1+m_1+m_2+\dots}$$

where you'll want

$$f^{(k)}(a)=0$$

for as many $k$ as you are interested in. (It results in quadratic, cubic, quartic, etc. convergence depending on how many terms you use.)