Arithmetic mean geometric mean 2d map

Question

Consider the 2d iteration \begin{aligned} x_{n+1} &=\frac{x_n+y_n}2,\\ y_{n+1} &=\sqrt{x_ny_n}. \end{aligned} We assume that $x_0>y_0>0$. Hence, $x_n>x_{n+1}>y_{n+1}>y_n$ due to the AM-GM inequality. Since the sequence $x_n$ is bounded from below and decreases monotonically, there is a limit $x_\infty$. Also, $y_n$ grows monotonically and is bounded from above. Hence, $\lim_{n\rightarrow\infty}y_n=y_\infty$. We also know that $$ x_{n+1}-y_{n+1}= \frac{(\sqrt{x_n}-\sqrt{y_n})^2}2=\frac12\frac{\sqrt{x_n}-\sqrt{y_n}}{\sqrt{x_n}+\sqrt{y_n}}(x_n-y_n)<\frac12(x_n-y_n)<2^{-n-1}(x_0-y_0). $$ Therefore, we find that both sequences have the same limit $L(x_0, y_0)=x_\infty=y_\infty$.

1. I cannot find the fixed-point $(x^*, y^*)$ of the iteration above, since both equations yield only $x^*=y^*$. I wonder if there is an analytic expression for $L(x_0, y_0)$ or at least an aproximate expression.

2. In the paper by Borwein and Borwein "The arithmetic-geometric mean and fast computation of elementary functions", SIAM V26, 351 (1984), they claim that elliptic integral of the second kind $$ I(a,b)= \int_0^{\pi/2}\frac{d\theta}{\sqrt{a^2\cos^2\theta+b^2\sin^2\theta}}, $$ where $I(a,b)=I(L(a,b))=\frac{\pi}{2L(a,b)}$. By numerical integration, I got something else though. Edit: I fixed this numerical problem, it was a bug in my code. I added the code below. Unfortunately, it seems the code formatting tool of math stackexchange is not working properly.

import numpy as np

import scipy.integrate

from scipy.special import ellipk

a = 2.0

b = 1.0

x = a

y = b

###########################################################

def fun(x, a, b):

  f = 1.0/np.sqrt((a*np.cos(x))**2+(b*np.sin(x))**2)

  return f

###########################################################

def fun2(x, m):

  f = 1.0/np.sqrt(1-m*(np.sin(x))**2)

  return f

##########################################################

integral, err = scipy.integrate.quad(fun,0., np.pi/2,args=(a,b,))

print('Original Gaussian quadrature integral', integral)

for i in range(100):

  plt.plot(x, y, 'ro', ms=2.0)

  x_n = (x+y)/2

  y_n = np.sqrt(x*y)

  x = x_n

L = x

print('iteration result', np.pi/(2*L))

m = 1-(b/a)**2

Im = ellipk(m)/a

print('ellipk', Im)

integral, err = scipy.integrate.quad(fun2, 0., np.pi/2,args=(m,))

print('Gaussian quadrature integral', integral/a)

3. How to show that $I(a, b)=I\left(\frac{a+b}2, \sqrt{ab}\right)$?

According to the authors, this problem has a long history, dating back to Gauss and Legendre.

Answer 1

Not sure what is your question in 1... it is clear that all points $(x^*,x^*)$ are fixed points of the underlying map. Another issue is to decide, for any particular $(x_0,y_0)$, if the sequence is converging, and to which of these fixed points.

In the picture bellow you can find the contour plots of $(x_0,y_0)\mapsto L(x_0,y_0)$, i.e. initial approximations in the same contour lead to convergence to the same fixed point.