Re: Re: A problem in Pi digits as Lattice space filling

*To*: mathgroup at smc.vnet.net*Subject*: [mg93976] Re: [mg93843] Re: A problem in Pi digits as Lattice space filling*From*: danl at wolfram.com*Date*: Sun, 30 Nov 2008 07:00:02 -0500 (EST)*References*: <ggj7co$j2i$1@smc.vnet.net>

> A friend and long time Mathematica programmer came up with this: > Clear[f, n] > f[0] = 1; > f[n_Integer /; n > 0] := Module[{cnt, m, p}, > cnt = 0; m = Table[0, Evaluate[Apply[Sequence, Table[{10}, {n}]]]]; > p = Partition[RealDigits[Pi, 10, 20*10^n][[1]], n, 1]; > For[i = 1, i <= 20*10^n - n + 1, i++, > If[m[[Apply[Sequence, p[[i]] + 1]]]++ == 0, > If[++cnt == 10^n, Print[i]; Break[]]]; > ] > ]; > Table[f[n], {n, 0, 4}] > From In[9]:= > 33 > > From In[9]:= > 606 > > From In[9]:= > 8554 > > From In[9]:= > 99847 > There my machine gives up in real time... > I pretty much trust this and my older clunky programs. > When run on a faster machine the answers are: > starting at f[0]=1; > 1,33,606,8554,99847,1369561, 14118308 > the hyperbolic probability limit: > > Limit[10^n/f[n],n->Infinity] > > appears to be zero from these numbers. > That result would suggest that the Pi digits are really not normal at all. I do not see that implication. > But with just these few numbers it is hard to say. > Since the BBP approach is actually based on the normal character of the > Pi digits The existance of the BBP formula probably has no bearing on normality of Pi. That at least is what I heard Helaman Fergusen say about it, 11 years ago. > ( one of the assumptions necessary for a PSLQ digit solve to actually > work! I'm not certain > of the theory, but that is what the papers say.). If there are integer relations to be found amongst the digits, PSLQ can find them. I am not aware of a connection here to normality, or lack thereof. > I'm not willing to argue the theory of the point, as I really don't much > understand it. > The idea of the space filling character > of the Pi digits came from them being normal. > If the lattice diverges as it seems to toward Infinity, > that would imply > that the digits of Pi aren't really normal, > I think. No. That is to say, if the lattice behaves as this one is behaving, there is nothing to indicate non-normality. I'll explain further below. > As a result the BBP Pi digits base 16 forms are not derived > correctly? > Or at least there is a question, > now, involved. The derivation was proven, many years ago. It's not conjectural. > Me, I'm an experimenter and not an axiomatic theoretical type. > But in any case the question of the space filling limit seems important > to my mind. > > > Respectfully, Roger L. Bagula > 11759Waterhill Road, Lakeside,Ca 92040-2905,tel: 619-5610814 > :http://www.geocities.com/rlbagulatftn/Index.html > alternative email: rlbagula at sbcglobal.net There are ways to investigate this. A relatively simple approach is to use random large numbers, and check the behavior of the same lattice filling as you are doing with digits of Pi. I will use my version of the filling time code, from an earlier post (your friend's is similar but slightly slower due to more work in the lookup phase). latticeFillMark[n_?NumericQ, d_Integer] /; d > 0 := Catch[Module[{digits, vals, len = 10^d, hitlist, count = 0, elem}, digits = First[RealDigits[n, 10, 10*Ceiling[d/2]*len]]; vals = Map[FromDigits, Partition[digits, d, 1]]; hitlist = ConstantArray[0, len]; Do[elem = vals[[j]] + 1; If[hitlist[[elem]] == 1, Continue[]]; hitlist[[elem]] = 1; count++; If[count == len, Throw[j]];, {j, Length[vals]}]; $Failed]] Here we generate random reals, of sufficient length to enable the experimental analysis. randomFillMarks[digits_, len_] := Table[ latticeFillMark[ RandomReal[WorkingPrecision -> 10*Ceiling[digits/2]*10^digits], digits] , {len}] We create a sample of 100 cases wherein we check how long to fill the 10x10 lattice of 100 values. r2sample = N[randomFillMarks[2, 100]]; Now I show that the mean and standard deviation make a value of 606 quite plausible. Mean[r2sample] Out[77]= 517.09 In[86]:= StandardDeviation[r2sample] Out[86]= 130.531 Here I perform a similar experiment on 20 random values, this time going after the 10x10x10 lattice of 1000 values. r3sample = N[randomFillMarks[3, 20]] Out[80]= {7145., 6530., 7058., 5903., 6476., 7655., 7303., 6636., \ 11536., 7789., 6639., 8176., 6445., 6553., 6606., 5620., 6604., \ 10084., 6408., 10454.} In[82]:= Mean[r3sample] Out[82]= 7381. In[84]:= StandardDeviation[r3sample] Out[84]= 1567.77 Again, the value for Pi (8554) is within one (estimated) standard deviation of the experimental mean. And we see similarly for the 10x10x10x10 case as well. r4sample = N[randomFillMarks[4, 20]]; In[182]:= Mean[r4sample] Out[182]= 92447.8 In[183]:= StandardDeviation[r4sample] Out[183]= 10790.9 There is an important caveat here. It is imperative that the random number generator be extremely good, that is, not suffer from defects of correlation between generated digits. The default RNG is based on a cellular automaton and has, I believe, performed well on the l'Ecuyer test suite. So I trust it. But if you do not, you could first do SeedRandom[Method -> "Rule30CA"] This uses the Rule 30 cellular automaton RNG. It's slower than the default one, but has always seemed quite robust in terms of behavior. Finally I will indicate an exact approach to these computations of expected value and standard deviation. We assume the digits of our input (Pi, in this case) are uniformly independently distributed (uid), and try to analyze the expected value of our lattice filling thresholds. The code below does this. In brief, it works as follows. We let p[j,n,d] denote the quantity of cases of n uid values, each in range 1-d, filling exactly j positions. This is seen to be the number filling at most that many, minus the number that exactly fit 1, 2, ..., j-1 positions. Then we let q[j,n,d] denote the probability of this happening (the number of cases, divided by total number of such ensembles, which is d^n). Why p for quantity and q for probability? I think my typing fingers got away from me (both of them). Clear[p, q] p[0, n_, d_] = 0; q[0, n_, d_] = 0; p[j_, n_, d_] := p[j, n, d] = Together[Binomial[d, j]*(j^n - Sum[Binomial[j, k] p[k, n, k], {k, j - 1}])] q[j_, n_, d_] := q[j, n, d] = p[j, n, d]/d^n I have not been able to come up with any closed form. But this suffices for computing at least the 10 and 10x10 cases exactly. The values of interest are q[d,n,d] with d=10 or 100, and n>=d. In[13]:= q10[n_] = q[10, n, 10] Out[13]= 10^-n (-10 + 45 2^n + 45 2^(3 n) + 105 2^(1 + 2 n) - 10 3^(2 n) - 40 3^(1 + n) - 252 5^n + 35 6^(1 + n) - 120 7^n + 10^ n) In[29]:= d10[n_] = Together[q10[n] - q10[n - 1]] Out[29]= -(1/63) 2^(-2 - n) 5^-n (-22680 + 2835 2^(3 n) + 2835 2^(4 + n) + 19845 2^(2 + 2 n) - 280 3^(2 n) - 7840 3^(2 + n) + 245 2^(4 + n) 3^(2 + n) - 63504 5^n - 12960 7^n) q100[n_] = q[100, n, 100]; d100[n_] = Together[q100[n] - q100[n - 1]]; We can now compute expected values. In[192]:= mean10 = Sum[n*d10[n], {n, 10, Infinity}] Out[192]= 7381/252 In[193]:= N[mean10] Out[193]= 29.2897 Recall that this is where Pi gave 33. Now we'll compute the standard deviation. In[194]:= var10 = Sum[n^2*d10[n], {n, 10, Infinity}] - mean10^2; In[195]:= std10 = Sqrt[var10]; In[196]:= N[std10] Out[196]= 11.211 So 33 is well within a standard deviation of the mean. We'll repeat this exact computation for the 10x10 case. mean100 = Sum[n*d100[n], {n, 100, Infinity}]; In[204]:= N[mean100] Out[204]= 518.738 In[205]:= var100 = Sum[n^2*d100[n], {n, 100, Infinity}] - mean100^2; In[206]:= std100 = Sqrt[var100]; N[std100] Out[207]= 125.822 Again we see that the computed value for Pi, 8554, is within a standard deviation of the expected value. This method runs out of steam around here due to computational complexity. I have not as yet found a more tractable formulation. In particular I have not found an exact form for the general case (that is, one for all d,n pairs). I would guess this has been studied in the literature, and that if one exists it is well known, just not to me. Daniel Lichtblau Wolfram Research