Model each three consecutive terms of the numerical sequence a_n as three points on the curve c+f*d^-n and form a new sequence b_{1,n} from these values of c. Likewise form b_{2,n} from b_{1,n}, etc. Then b_{k,n} is an array of estimates of lim_{n->oo} a_n:
(d31) b := if k = 0 then a k, n n 2 b b - b k - 1, n - 1 k - 1, n + 1 k - 1, n else ------------------------------------------- b - 2 b + b k - 1, n + 1 k - 1, n k - 1, n - 1
Note that this is 0 for all n if the preceding sequence progressed geometrically, and infinity if it progressed arithmetically.
These estimates tend to improve gradually with n and dramatically with k up to some point where noise and ripple in the a_n sequence gets amplified to the point of swamping everything. But let a_n be a sequence of (noise-free) continued fractions, e.g.,
(c34) a[n] := (cfexpand(d24[n]),%%[1]/%%[2])[1]
%% 1 (d34) a := (cfexpand(d24 ), ---) n n %% 1 2
(c24) apply('matrix,makelist([3,7,15,1,k],k,2,8))
[ 3 7 15 1 2 ] [ ] [ 3 7 15 1 3 ] [ ] [ 3 7 15 1 4 ] [ ] (d24) [ 3 7 15 1 5 ] [ ] [ 3 7 15 1 6 ] [ ] [ 3 7 15 1 7 ] [ ] [ 3 7 15 1 8 ]
with some term progressing arithmetically. Then the first acceleration (c35) apply(matrix,makelist(cf(b[1,n]),n,2,6))
[ 3 7 15 1 6 1 15 7 ] [ ] [ 3 7 15 1 8 1 15 7 ] [ ] (d35) [ 3 7 15 1 10 1 15 7 ] [ ] [ 3 7 15 1 12 1 15 7 ] [ ] [ 3 7 15 1 14 1 15 7 ]
produces a sequence with that same term progressing twice as fast.
(c36) apply(matrix,makelist(cf(b[2,n]),n,3,5))
[ 3 7 15 1 18 1 4 2 1 1 1 2 ] [ ] (d36) [ 3 7 15 1 22 1 4 2 1 1 1 2 ] [ ] [ 3 7 15 1 26 1 4 2 1 1 1 2 ]
(c37) apply(matrix,makelist(cf(b[3,n]),n,4,4))
(d37) [ 3 7 15 1 46 1 1 3 3 1 3 ]
What's going on? Each of these progressions is a ratio of arithmetic progressions. Given
(c39) a[n] := (a*n+b)/(c*n+d)
a n + b (d39) a := ------- n c n + d
then the first "accelerant" is (c40) ratsimp(b[1,n])
2 a c n + a d + b c (d40) ------------------- 2 2 c n + 2 c d
or, dividing through by 2 c,
(c41) opsubst("+" = lambda([[l]],2*c,%%*multthru(1/%%,apply("+",l))),%)
a d b a n + --- + - 2 c 2 (d41) ------------- . c n + d
Remarkably, only one coefficient has changed. After n accelerations,
(c43) subst(b(0),b0,differenceq(b(n) = a*d/2/c+b(n-1)/2,b(n),b(0) = b0))
a d a d - b(0) c (d43) b(n) = --- - ------------ c n c 2
which clearly approaches a*d/c, whereat a_n approaches (c44) factor(subst(a*d/c,b,a[n]))
a (d44) -, c
i.e., we progressed geometrically through the sequence of arithmetic progressions. --rwg
See http://en.wikipedia.org/wiki/Shanks_transformation , but note there are usually better ways to accelerate sums. Also, the Gregory's series for pi=4 atan 1 is discussed in HAKMEM--it is unclear whether this process converges in the long run. Things get *very* shaky after a few dozen iterations. many iterations.
I decided it was time to repeat this experiment on modern hardware (= an old laptop), and found the method to be even more wildly impractical than I remembered. If you use exact rational arithmetic, the integers grow like a double exponential (making the simple behavior on continued fractions the more surprising). If you use approximate arithmetic, ferocious subtractive cancellation requires kilodigits of significance to get 100 digits of pi. Here is the end of a run starting with 13333 digit approximations to 4/1,4/3,4/5,... . Note that it divides by 0 on the 126th step. Note also the crazily irregular convergence(?), where it makes three or four consecutive corrections(?) many orders of magnitude smaller than needed, followed by a (sometimes too) big one and then more tiny ones. k b_{k,2k+1} 96 3.14159265358979323846264338327950288419716939937510582097494459230781639134930408761923749865559466 97 3.14159265358979323846264338327950288419716939937510582097494459230781639134930408761386651551062205 98 3.14159265358979323846264338327950288419716939937510582097494459230781640628565979016487075572426161 99 3.14159265358979323846264338327950288419716939937510582097494459230781640628565979016487075572426161 100 3.14159265358979323846264338327950288419716939937510582097494459230781640628565979016487075572426161 101 3.14159265358979323846264338327950288419716939937510582097494459230781640628565979016487075572426161 102 3.14159265358979323846264338327950288419716939937510582097494459230781640628565987264008565988274034 103 3.14159265358979323846264338327950288419716939937510582097494459230781640628565987264008565988274034 104 3.14159265358979323846264338327950288419716939937510582097494459230781640628565987264008565988274034 105 3.14159265358979323846264338327950288419716939937510582097494459230781640628557017500807143906315824 125 3.141592653589793238462643383279502884197169399375105820974944592307816406286208998687220540166 Power::infy: Infinite expression 1/0.*10^-4 encountered. >> Pi= 3.14159265358979323846264338327950288419716939937510582097494459230781640628620899862803482534211707