Pairwise combining the terms of A&S 6.5.31 (upper Gamma CF), then adding last week's lower gamma CF: In[185]:= Equal[Gamma[ z], ((y^z*((1/(cfk[n*(z - n), -z + y + 2*n + 1, List[n, Infinity]] - z + y + 1)) + (1/(z + cfk[n*y, z - y + n, List[n, Infinity]] - y))))/(E^y))] z 1 Gamma[z] == (y (------------------------------------------- + -y + z + cfk[n y, n - y + z, {n, Infinity}] 1 y -----------------------------------------------------------)) / E 1 + y - z + cfk[n (-n + z), 1 + 2 n + y - z, {n, Infinity}] where cfk:=ContinuedFractionK and y>0 can be chosen to optimize convergence. Specializing, In[186]:= % /. z -> 1/2 /. y -> 1 Out[186]= 1 Sqrt[Pi] == (-------------------------------------- + 1 1 -(-) + cfk[n, -(-) + n, {n, Infinity}] 2 2 1 ------------------------------------------) / E 3 1 3 - + cfk[(- - n) n, - + 2 n, {n, Infinity}] 2 2 2 In[192]:= N[List @@ %186 /. cfk[num_, den_, _] :> ContinuedFractionK[num, den, {n, 9}]] Out[192]= {1.77245, 1.77244} This is the first non-asymptotic sum or CF I've seen for sqrt(pi). But the sum behind the lowergamma cf has been there all along, so I must not have been paying attention. --rwg
Julian has found that a convergence estimator I reported last year is less accurate than I thought, but hardly useless. With his correction, I just found a surprise in the convergence rate of the slower CF in the Γ formula below. For this quadratic/linear with z:=1/2, the formula predicts an error proportional to E^(-4 √(n y)). While this is asymptotically 0 bits/term, it is enormously better than any (q+1)Fq[1], whose error after n terms is 1/polynomial(n). Roughly speaking, instead of 10^n terms, we need only n^2. And note the surprising sensitivity to the free parameter y, a presumably small, fixed number added to the growing variable 2n. In fact, the formula predicts nonconvergence for y=0. It's almost right. Another oddity: simply changing the sign of the CF numerators, i.e. cfk[n*(n - z), -z + y + 2*n + 1,{n,∞}], almost completely desensitizes to y the convergence. The formula says the second, faster CF converges like y^(-1 + n)/(-1 + n)! Note that, unlike the first cf, this error increases with y. Suppose we are willing to compute 69 terms of each CF. Equating the two error terms and solving, Mathematica gave an exact large Lambert-W expression, approximately y=7.1, yielding about 36 digits of √π. At "0 bits/term"! But this remains a non-competitive Γ method unless we find a (presumably matrix) speedup of the siower CF. --rwg On Tue, Aug 3, 2010 at 4:09 PM, Bill Gosper <billgosper@gmail.com> wrote:
Pairwise combining the terms of A&S 6.5.31 (upper Gamma CF), then adding last week's lower gamma CF:
In[185]:= Equal[Gamma[ z], ((y^z*((1/(cfk[n*(z - n), -z + y + 2*n + 1, List[n, Infinity]] - z + y + 1)) + (1/(z + cfk[n*y, z - y + n, List[n, Infinity]] - y))))/(E^y))]
z 1 Gamma[z] == (y (------------------------------------------- + -y + z + cfk[n y, n - y + z, {n, Infinity}]
1 y -----------------------------------------------------------)) / E 1 + y - z + cfk[n (-n + z), 1 + 2 n + y - z, {n, Infinity}]
where cfk:=ContinuedFractionK and y>0 can be chosen to optimize convergence. Specializing,
In[186]:= % /. z -> 1/2 /. y -> 1
Out[186]= 1 Sqrt[Pi] == (-------------------------------------- + 1 1 -(-) + cfk[n, -(-) + n, {n, Infinity}] 2 2
1 ------------------------------------------) / E 3 1 3 - + cfk[(- - n) n, - + 2 n, {n, Infinity}] 2 2 2
In[192]:= N[List @@ %186 /. cfk[num_, den_, _] :> ContinuedFractionK[num, den, {n, 9}]]
Out[192]= {1.77245, 1.77244}
This is the first non-asymptotic sum or CF I've seen for sqrt(pi). But the sum behind the lowergamma cf has been there all along, so I must not have been paying attention. --rwg
Traffic signs can lie. I saw one that said No turn on red, and then the red turned on.
The slower converging "half" of this formula is the [w,oo) integral, w^z/(E^w*Gamma[z, w]) == MProd[MatrixForm[{{-1 + 2*n + w - z, n*(-n + z)}, {1, 0}}], {n, 1, Infinity}] where MProd is matrix product, and the matrices a1 b1 . a2 b2 . . . 1 0 1 0 compute the continued fraction a1+b1/(a2+b2/(... . Unfortunately, this K quadratic/linear is not a special case of the path invariant system at the bottom of p13 of http://www.tweedledum.com/rwg/stanfordn3.pdf, so I couldn't accelerate it. Finally, I found the {n,j} pair n: {{a + j + 2 n, -n (j + n)}, {1, 0}}, j: {{-a - n, (-1 + n) (j + n)}, {-1, j + n}} but the first few contours proved disappointing, so I sought a third dimension k. Usually this would be no easier (nor likelier to happen) than dredging up the original j matrix, but I found a cheap trick which seems too good to be true (and probably is, in most cases): The system is free of k, so the k matrix is obviously {1,0;0,1} : In[330]:= updatemat[k, IdentityMatrix[2]] Out[330]= {{n, {{a + j + 2 n, -n (j + n)}, {1, 0}}}, {k, {{1, 0}, {0, 1}}}, {j, {{-a - n, (-1 + n) (j + n)}, {-1, j + n}}}} Then In[331]:= recoord[n -> n + k] Out[331]= {{n, {{a + j + 2 (k + n), -(k + n) (j + k + n)}, {1, 0}}}, {k, {{a + j + 2 (k + n), -(k + n) (j + k + n)}, {1, 0}}}, {j, {{-a - k - n, (-1 + k + n) (j + k + n)}, {-1, j + k + n}}}} making the k matrix the same as the n matrix. Then In[332]:= recoord[j -> j - k] Out[332]= {{n, {{a + j + k + 2 n, -(j + n) (k + n)}, {1, 0}}}, {k, {{-((a + j + n)/(1 + a)), ((-1 + j + n) (k + n))/(1 + a)}, {-(1/(1 + a)), (k + n)/(1 + a)}}}, {j, {{-a - k - n, (j + n) (-1 + k + n)}, {-1, j + n}}}} Now the n matrix is symmetric in j and k, so k is highly legitimate! But why aren't the j and k matrices analogous? They are, but all the elements of k are divided by 1+a. Scaling the whole matrix (even by a function of k) is irrelevant to path invariance, so In[333]:= updatemat[k, (1 + a)*mats[k]] Out[333]= {{n, {{a + j + k + 2 n, -(j + n) (k + n)}, {1, 0}}}, {k, {{-a - j - n, (-1 + j + n) (k + n)}, {-1, k + n}}}, {j, {{-a - k - n, (j + n) (-1 + k + n)}, {-1, j + n}}}} In[334]:= picheck[%] Out[334]= True (Path invariance check.) This matrix functions package is from Corey Ziegler Hunts. Early results from these {j,k,n} matrices are inconclusive. I'm hitting a lot of numerical instability, which seems inherent: When a CF matrix a b c d converges, a/c approaches b/d. And thus the matrix approaches singularity. --rwg On Tue, Aug 3, 2010 at 4:09 PM, Bill Gosper <billgosper@gmail.com> wrote:
Pairwise combining the terms of A&S 6.5.31 (upper Gamma CF), then adding last week's lower gamma CF:
In[185]:= Equal[Gamma[ z], ((y^z*((1/(cfk[n*(z - n), -z + y + 2*n + 1, List[n, Infinity]] - z + y + 1)) + (1/(z + cfk[n*y, z - y + n, List[n, Infinity]] - y))))/(E^y))]
z 1 Gamma[z] == (y (------------------------------------------- + -y + z + cfk[n y, n - y + z, {n, Infinity}]
1 y -----------------------------------------------------------)) / E 1 + y - z + cfk[n (-n + z), 1 + 2 n + y - z, {n, Infinity}]
where cfk:=ContinuedFractionK and y>0 can be chosen to optimize convergence. Specializing,
In[186]:= % /. z -> 1/2 /. y -> 1
Out[186]= 1 Sqrt[Pi] == (-------------------------------------- + 1 1 -(-) + cfk[n, -(-) + n, {n, Infinity}] 2 2
1 ------------------------------------------) / E 3 1 3 - + cfk[(- - n) n, - + 2 n, {n, Infinity}] 2 2 2
In[192]:= N[List @@ %186 /. cfk[num_, den_, _] :> ContinuedFractionK[num, den, {n, 9}]]
Out[192]= {1.77245, 1.77244}
This is the first non-asymptotic sum or CF I've seen for sqrt(pi). But the sum behind the lowergamma cf has been there all along, so I must not have been paying attention. --rwg
It's easy to use calculus to prove the most efficient can (cylinder + top + bottom) has a square cross-section, where efficient means it has the least area for a given volume. The well-known answer is that it has a square cross-section: height = diameter. This is easy to prove with standard calculus optimization methods. But what is the fastest way to see this? --Dan
participants (2)
-
Bill Gosper -
Dan Asimov