In the previous post I said that Jonathan and Peter Borwein figured out how to use the rapid convergence of the AGM to compute various functions, including logarithms. This post will show how to compute logarithms using the AGM.

First I need to define the Jacobi theta functions

Note that if *q* is very small, the series above converge very quickly. We will pick *q* so small that θ_{2}(*q*) and θ_{3}(*q*) are trivial to compute with sufficient accuracy.

Suppose we want to compute the natural log of 10 to a thousand digits. We can use the following equation from [1].

We want *q* to be small, so we want 1/*q* to be large. We set *q* = 10^{–n} to compute log 10^{n} and divide the result by *n*.

We can show from the definition of the theta functions that

Since we want 1000 digit accuracy, we can set *q* = 10^{-250}. Then θ_{2}^{2}(*q*^{4}) = 4×10^{-500} and θ_{3}^{2}(*q*^{4}) = 1 to about 2000 digits of accuracy.

We can calculate the result we want with just 20 iterations of the AGM as the following bc code shows.

define agm(a, b, n) { auto i, c, d for (i = 1; i <= n; i++) { c = (a + b)/2 d = sqrt(a*b) a = c b = d } return a } scale=2000 # a(1) = arctan(1) = pi/4 x = a(1)/(agm(4*10^-500, 1, 20)*250) - l(10) # error l(x)/l(10) # log_10(error)

This returns -999.03…, so we were basically able to get 1000 digits, though the last digit or two might be dodgy. Let’s try again with a smaller value of *q*, setting 10^{-300}.

x = a(1)/(agm(4*10^-600,1, 20)*300) - l(10); l(x)/l(10)

This returns -1199.03… and so we have about 1200 digits. We were using 2000 digit precision in our calculations, but our approximation θ_{3}^{2}(*q*^{4}) = 1 wasn’t good to quite 2000 digits. Using a smaller value of *q* fixed that.

## Related posts

[1] The Borwein brothers, Pi and the AGM. arXiv

Did you just calculated pi/4 by calculating 4*arctan(1) and then dividing by 4?

Yes, force of habit. Whenever I’m writing bc and I see I need π I immediately type 4*a(1).