Trouble getting loop to break

[John Machin]

On Nov 25, 7:00 pm, Dick Moores <r... at rcblue.com> wrote:
> At 01:32 PM 11/20/2007, Fredrik Johansson wrote:
>
>
>
>
>
>
>
> >On Nov 20, 2007 10:00 PM, Dick Moores <r... at rcblue.com> wrote:
> > > And also with the amazing Chudnovsky algorithm for pi. See
> > > <http://python.pastebin.com/f4410f3dc>
>
> >Nice! I'd like to suggest two improvements for speed.
>
> >First, the Chudnovsky algorithm uses lots of factorials, and it's
> >rather inefficient to call mpmath's factorial function from scratch
> >each time. You could instead write a custom factorial function that
> >only uses multiplications and caches results, something like this:
>
> >cached_factorials = [mpf(1)]
>
> >def f(n):
> >     n = int(n)
> >     if n < len(cached_factorials):
> >         return cached_factorials[n]
> >     p = cached_factorials[-1]
> >     for i in range(len(cached_factorials), n+1):
> >         p *= i
> >         cached_factorials.append(p)
> >     return p
>
> >(In some future version of mpmath, the factorial function might be
> >optimized so that you won't have to do this.)
>
> >Second, to avoid unnecessary work, factor out the fractional power of
> >640320 that occurs in each term. That is, change the "denom =" line to
>
> >     denom = (f(3*k) * ((f(k))**3) * (640320**(3*k)))
>
> >and then multiply it back in at the end:
>
> >     print 1/(12*sum/640320**(mpf(3)/2))
>
> >With these changes, the time to compute 1,000 digits drops to only
> >0.05 seconds!
>
> >Further improvements are possible.
>
> >Fredrik
>
> Fredrik,
>
> I'm afraid I'm just too dumb to see how to use your first suggestion
> of cached_factorials. Where do I put it and def()? Could you show me,
> even on-line, what to do? <http://py77.python.pastebin.com/f48e4151c>

(1) Replace line 8 (from mpmath import factorial as f) with Fredrik's
code. (2) Test it to see that it gave the same answers as before ...
[time passes]
Wow! [w/o psyco] Pi to 1000 decimal places takes 13 seconds with
original code and 0.2 seconds with Fredrik's suggestion. 2000: 99
seconds -> 0.5 seconds.