[math-fun] Shanks weirdness, DanA's puzzle sum
Trying for several digits of DanA's Out[197]= Sum[Product[1 - 1/Prime[k], {k, 1, -1 + n}]/Prime[n], {n, 1, Infinity}] (No need to special-case n=1 -- Mma does the right thing on empty products. And the *wrong* thing when they're of negative length!) In[198]:= shanks[L_List] := (#1*#3 - #2^2)/(#1 - 2*#2 + #3) &[Drop[L, -2], L[[2 ;; -2]], Drop[L, 2]] Take exponentially longer segments of the series: In[187]:= numMProd[ Table[MProd[ MatrixForm[{{1 - 1./Prime[i], 1/Prime[i]}, {0, 1}}], {i, 2^n + 1, 2^(n + 1)}], {n, 0, 11}]] Out[187]= {{{0.666667, 1/3}, {0, 1}}, {{0.685714, 0.314286}, {0., 1.}}, {{0.74823, 0.25177}, {0., 1.}}, {{0.795719, 0.204281}, {0., 1.}}, {{0.83033, 0.16967}, {0., 1.}}, {{0.856913, 0.143087}, {0., 1.}}, {{0.876387, 0.123613}, {0., 1.}}, {{0.891033, 0.108967}, {0., 1.}}, {{0.902836, 0.0971637}, {0., 1.}}, {{0.912323, 0.0876769}, {0., 1.}}, {{0.920145, 0.079855}, {0., 1.}}, {{0.926663, 0.0733374}, {0., 1.}}} Restore the first term, and accumulate: In[188]:= FoldList[Dot, {{1, 1}/2, {0, 1}}, %] Out[188]= {{{1/2, 1/2}, {0, 1}}, {{0.333333, 0.666667}, {0., 1.}}, {{0.228571, 0.771429}, {0., 1.}}, {{0.171024, 0.828976}, {0., 1.}}, {{0.136087, 0.863913}, {0., 1.}}, {{0.112997, 0.887003}, {0., 1.}}, {{0.0968287, 0.903171}, {0., 1.}}, {{0.0848595, 0.915141}, {0., 1.}}, {{0.0756126, 0.924387}, {0., 1.}}, {{0.0682658, 0.931734}, {0., 1.}}, {{0.0622804, 0.93772}, {0., 1.}}, {{0.057307, 0.942693}, {0., 1.}}, {{0.0531043, 0.946896}, {0., 1.}}} The partial sums are in the upper right [[1,2]]: In[189]:= shanks[#[[1, 2]] & /@ %] Out[189]= {0.948718, 0.899118, 0.917897, 0.932005, 0.94094, 0.949258, \ 0.955796, 0.96014, 0.964033, 0.967136, 0.969815} Second tier: In[190]:= shanks[%] Out[190]= {0.912739, 0.974621, 0.956372, 1.06137, 0.97982, 0.968741, \ 0.997558, 0.979348, 0.986716} Third (losing both ends each time.) In[191]:= shanks[%] Out[191]= {0.960528, 0.971919, 1.01547, 0.966999, 0.976743, 0.986399, \ 0.984594} Fourth. In[192]:= shanks[%] Out[192]= {0.956494, 0.99253, 0.975112, 2.0492, 0.984878} A tall pole in the tent! The immediate explanation is that 0.966999, 0.976743, 0.986399 is nearly in arithmetic progression, provoking division by nearly 0. In fact, the irregularity of the primes makes the process so noisy that the second differences change sign wantonly, and this whole computation is meaningless. But suppose we smooth it out, replacing Prime[n] with n*Log[n]. That's an interesting sum, too. Suppose we sum 10^6 terms directly, and Shanks the tail: numMProd[Table[ MProd[MatrixForm[{{1 - 1`33/(n*Log[n]), 1/(n*Log[n])}, {0, 1}}], {n, 2^k + 1000000, 999999 + 2^(k + 1)}], {k, 0, 11}]]; (Here, numMProd just performs the matrix products, to 33 places due to the 1`33.) Form the sequence of exponentially longer partial sums: FoldList[Dot, IdentityMatrix[2], %]; Run one Shanks: shanks[#[[1, 2]] & /@ %] {-7.23825846114180023820963784936630936956*10^-8, \ -7.2383081775534954986936726671956993755*10^-8, \ -7.238507042592270807744603740319705184*10^-8, \ -7.239302497883386891719990656592430349*10^-8, \ -7.24248428013694135359940241455532980*10^-8, \ -7.25521109787942075472466273513777568*10^-8, \ -7.3061158789240489622769357685310492*10^-8, \ -7.509715087677159821674775767630926*10^-8, \ -8.323952662885734390860079101633870*10^-8, \ -1.157962990187933870110393101592588*10^-7, \ -2.459217059666940024454702903981525*10^-7} How can these conceivably be <0?? Try the 2nd order Shanks: shanks2[#[[1, 2]] & /@ %198] {-7.2382418887852421433705793226704*10^-8, \ -7.238241887452633043864365057607*10^-8, \ -7.238241876791840693390633496226*10^-8, \ -7.238241791506788978747515379514*10^-8, \ -7.23824110924696740587213638954*10^-8, \ -7.23823565149782811295872404833*10^-8, \ -7.23819199477433435746658414992*10^-8, \ -7.2378428252583618557144316956*10^-8, \ -7.2350508161419247635566250253*10^-8} Still negative! The only glimmer of sanity comes from summing bursts of length 3^n, shanks[#[[1, 2]] & /@ %] {-7.23827917631153038724969102081745002169*10^-8, \ -7.2385774734354090466183845664296486922*10^-8, \ -7.241262106902162871015196793342498037*10^-8, \ -7.26542271070426051438560418211001551*10^-8, \ -7.4828385233749802945842107930570293*10^-8, \ -9.438781720404053555835539434832801*10^-8, \ -2.702074800500521053964558787965094*10^-7, \ -1.846816420328610065496655139531113*10^-6, \ -0.00001588394699143661453072982592933523, \ -0.0001383570856200631342075355905173333, \ -0.001152755525647733026707277401399129} and the using Shanks2: shanks2[#[[1, 2]] & /@ %%] {-7.23824188680716125090267195281569*10^-8, \ -7.2382418304287692662035666190882*10^-8, \ -7.238240308325929783457067074808*10^-8, \ -7.238199220757133558198054009671*10^-8, \ -7.237090600928513738314351741322*10^-8, \ -7.20721785747766461202604502914*10^-8, \ -6.40543875690906529608871726514*10^-8, 1.48716801605978891479361949759*10^-7, 5.62589299046691128953938459235*10^-6} I have never before seen this wrong sign behavior on a smooth sequence. Unnerving. --rwg Does anyone have a nice little Richardson accelerator to share? Is it time to give up and beg Dan for his answer?
Okay, this has gone on entirely too long! Shall I explain the sum (1/2) + (1/3)(1-1/2) + (1/5)(1-1/2)(1-1/3) + (1/7)(1-1/2)(1-1/3)(1-1/5) + ... ? Or give a hint? --Dan On 2014-01-10, at 2:21 PM, Bill Gosper wrote:
Trying for several digits of DanA's
Out[197]= Sum[Product[1 - 1/Prime[k], {k, 1, -1 + n}]/Prime[n], {n, 1, Infinity}]
(No need to special-case n=1 -- Mma does the right thing on empty products. ...
A tall pole in the tent!
Bill! This is a family list!
... Is it time to give up and beg Dan for his answer?
[SPOILER] . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 1! I actually have two solutions for this, one somewhat more general (though not very rigorous) and one which is probably Dan's intended solution. The first one rearranges the sum into: (1/2)+(1-1/2)(1/3+(1-1/3)(1/5+(1-1/5)( ... ))) which can be seen as the composition of an infinite number of functions A_n[z]=1/Prime[n]+(1-1/Prime[n])*z (Bill points out this is essentially a matrix product, and that I'm just obscuring things by writing it out like this.) And now the risky bit: looking at the tail end of the infinite sum, we'll have something like A_1[A_2[ ... A_inf[A_inf[A_inf[ ... ]]]...]] Of course the only fixed point of A, and thus the fixed point of A_inf (whatever that is) is 1, and then the nested expression simplifies to 1. This can also be generalized- just replace 1/Prime[n] with any real, non-one sequence in (0,2) which converges to a non-negative number, and the sum will still come out to 1. (Caution: some of these constraints are only a few minutes old!) Second solution: Expand the entire sum (parentheses added for clarity): (1/2) + (1/3 - (1/3)*(1/2)) + (1/5 - (1/5)*(1/2)-(1/5)*(1/3) + (1/5)*(1/2)*(1/3)) +... and it's easy to see that this looks like an application of the inclusion-exclusion principle, except with multiplication instead of intersection, and asymptotic density instead of cardinality. One can carry the analogy further and claim this is the asymptotic density of the union of <composites divisible by 2>, <composites divisible by 3>, and so forth- in other words, the asymptotic density of the set of composites. Since PrimePi[x]<x/(Log[x]-4) (for x>55) and is thus less than linear, the sum must be equal to 1, although it will converge very, /very/ slowly. This could also be combined with the first solution for some set-theoretical* results which I'm almost certain have been discovered before. --Neil *number-theoretical? measure-theoretical? On Wed, Jan 15, 2014 at 10:46 AM, Dan Asimov <dasimov@earthlink.net> wrote:
Okay, this has gone on entirely too long!
Shall I explain the sum
(1/2) + (1/3)(1-1/2) + (1/5)(1-1/2)(1-1/3) + (1/7)(1-1/2)(1-1/3)(1-1/5) + ...
?
Or give a hint?
--Dan
On 2014-01-10, at 2:21 PM, Bill Gosper wrote:
Trying for several digits of DanA's
Out[197]= Sum[Product[1 - 1/Prime[k], {k, 1, -1 + n}]/Prime[n], {n, 1, Infinity}]
(No need to special-case n=1 -- Mma does the right thing on empty products. ...
A tall pole in the tent!
Bill! This is a family list!
... Is it time to give up and beg Dan for his answer?
_______________________________________________ math-fun mailing list math-fun@mailman.xmission.com http://mailman.xmission.com/cgi-bin/mailman/listinfo/math-fun
participants (3)
-
Bill Gosper -
Dan Asimov -
Neil Bickford