0

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((anp.cos(x))2+(bnp.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)

  1. 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.

minmax
  • 877

1 Answers1

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.

enter image description here

PierreCarre
  • 23,092
  • 1
  • 20
  • 36
  • I think I showed the iteration converges in my OP. The problem is to determine the limit $L(x_0, y_0)$ for $x_0>y_0>0$. I can only do that numerically. I believe there is no analytical expression, since it is directly related to the complete elliptic integral of the first kind. – minmax Nov 18 '21 at 14:02
  • @minmax I attached a figure to my post. It seems that the attracting set for each fixed point is a hyperbolic type curve. – PierreCarre Nov 18 '21 at 14:08
  • thanks! but I think don't understand your figure. Are these curves level sets? The function $L(a,b)$ is homogenous, that is $L(\lambda a, \lambda b)=\lambda L(a, b)$. I think this could help. I just added the simple code I used to verify the theory. – minmax Nov 18 '21 at 14:32
  • @minmax If you pick $(x_0,y_0)$ and use it as initial approximation for the fixed point method, you will observe convergence to some point $(x_(x_0,y_0), x_(x_0,y_0))$. The picture represents the contours of $f(x_0,y_0)= x_*(x_0,y_0)$. For instance, points on the contour level 2 produce sequences converging to $(2,2)$. – PierreCarre Nov 18 '21 at 14:37