Oops, luckily, Julian was eavesdropping! ---------- Forwarded message ---------- From: Julian Ziegler Hunts <julianj.zh@gmail.com> Date: Sat, Apr 23, 2011 at 8:40 PM Subject: Re: [math-fun] How can I accurately calculate (p/q)^n mod 1, for large n To: billgosper@gmail.com
Timing[Do[10^16*Floor[Mod[1.5^k, 1]], {k, 99999}]]
What? That's always 0. Mod[x,1] is ≥0 and <1, so Floor[Mod[x,1]]=0 for all x. The two methods should never differ, since even the numerics know that 0≤Mod[x,1]<1. Maybe you want Floor[10^16*Mod[1.5^k,1]], which diverges from the exact computation at k=21, or, with 10000 (decimal) digits of precision, at k=56758. Julian On Sat, Apr 23, 2011 at 7:32 PM, Bill Gosper <billgosper@gmail.com> wrote:
*From: *
At 05:20 PM 4/23/2011, Dan Asimov wrote:
Maybe I should be more specific.
Mathematica can do a loop of length 10000, say calculating
100*Floor[Mod[(3/2)^k,1]]
for 10001 <= k < 20000,
but when I try to do something like that for
a loop of size 100,000 or 1,000,000, it seems to take forever (or hang, I'm not sure which).
--Dan
It's computing the numerators and denominators, e.g., 2^100000, exactly.
In[622]:= Timing[Do[If[10^16*Floor[Mod[(3/2)^k, 1]] != 10^16*Floor[Mod[1.5^k, 1]],Print[k]], {k, 99999}]]
Out[622]= {59.0911, Null}
About 50 times faster:
In[623]:= Timing[Do[10^16*Floor[Mod[1.5^k, 1]], {k, 99999}]]
Out[623]= {1.21972, Null}
Or 60 this way: In[625]:= Block[{t = 10^16},Timing[Do[t*Floor[Mod[1.5^k, 1]], {k, 99999}]]]
Out[625]= {0.997455, Null}
At some point, the two methods will probably differ by 1 (in 10^16). --rwg (Boy, Neil must be busy, not to have beaten me to this one.)