By using this site, you agree to our updated Privacy Policy and our Terms of Use. Manage your Cookies Settings.
439,941 Members | 1,753 Online
Bytes IT Community
+ Ask a Question
Need help? Post your question and get tips & solutions from a community of 439,941 IT Pros & Developers. It's quick & easy.

code optimization (calc PI)

P: n/a
mm
(Yes, I konw whats an object is...)
BTW. I did a translation of a pi callculation programm in C to Python.
(Do it by your own... ;-)
--------
Calc PI for 800 digs(?). (german: Stellen)
------
int a=10000,b,c=2800,d,e,f[2801],g;main(){for(;b-c;)f[b++]=a/5;
for(;d=0,g=c*2;c-=14,printf("%.4d",e+d/a),e=d%a)for(b=c;d+=f[b]*a,
f[b]=d%--g,d/=g--,--b;d*=b);}

$ ./a.exe
31415926535897932384626433832795028841971693993751 058209749445923078164062862089
98628034825342117067982148086513282306647093844609 550582231725359408128481117450
28410270193852110555964462294895493038196442881097 566593344612847564823378678316
52712019091456485669234603486104543266482133936072 602491412737245870066063155881
74881520920962829254091715364367892590360011330530 548820466521384146951941511609
43305727036575959195309218611738193261179310511854 807446237996274956735188575272
48912279381830119491298336733624406566430860213949 463952247371907021798609437027
70539217176293176752384674818467669405132000568127 145263560827785771342757789609
17363717872146844090122495343014654958537105079227 968925892354201995611212902196
08640344181598136297747713099605187072113499999983 729780499510597317328160963185
--------

But Python is much slower here then C here. I used a while-loop within a
while-loop. Not that much calculation here.
Jan 3 '07 #1
Share this Question
Share on Google+
21 Replies


P: n/a
mm wrote:
(Yes, I konw whats an object is...)
BTW. I did a translation of a pi callculation programm in C to Python.
(Do it by your own... ;-)
Is that a question on how to optimize code you won't show us? If yes, I'm
sorry to tell you that crystal balls are short these days. Too much
new-year-outlooks.

Diez
Jan 3 '07 #2

P: n/a
mm

Hmm... it's a question. It was not that easy to translate this #@*?%!
C-Program into readable code and then to Python. But it works.

There are only two while-loops (a while within an other while).

I konw, that for example while-loops in Perl are very slow. Maybe this
is also known in Pyhton. Then, I can translate the while-loops in to
for-loops, for example.
More general, maybe there is a speed optimazation docu out there.

(Anyway, this code for C-calc is complex. And I realy don't understand
it right now... 8-)
But anyway, there is not that much calculation, it has more something to
do with the too while-loops. Here is some code:

(In general, it has basically nothing to do with PI-calc.)
c=2800 ## a counter
while c*2:
## do some calc. did not change c here.
...
b=c ## number of elements

while (b-1):
## so some calc.
...
b=b-1 ## just for while-loop condition.
c = c-14; ## this is code vor the 1st while-loop, BUT bust run after
the 2nd-while-loop.
pi = pi + str("%04d" % int(e + d/a)) ## this should be fast?! I dont
know.
There are a output string, a list, and integers. No complex data
structures; and no complex calculations.

Diez B. Roggisch wrote:
mm wrote:

>>(Yes, I konw whats an object is...)
BTW. I did a translation of a pi callculation programm in C to Python.
(Do it by your own... ;-)


Is that a question on how to optimize code you won't show us? If yes, I'm
sorry to tell you that crystal balls are short these days. Too much
new-year-outlooks.

Diez
Jan 3 '07 #3

P: n/a

Ah, no. It was a HASH (assoziative Array or somethings like that).
mm wrote:
>
I konw, that for example while-loops in Perl are very slow. Maybe this
is also known in Pyhton. Then, I can translate the while-loops in to
for-loops, for example.
More general, maybe there is a speed optimazation docu out there.
Jan 3 '07 #4

P: n/a
Using the '+' operator for string concatonation can be slow, especially
when done many times in a loop.
pi = pi + str("%04d" % int(e + d/a)) ## this should be fast?! I dont
The accepted solution would be to make pi an array and append to the
end...

pi = [] #create the array (empty)
....
....
pi.append(str("%04d"%int(e+d/a))) # append to it

And when it is time to print the results do the following:

print "".join(pi)

It may look strange, but it is a common Python idiom. You might also
look into an optimizer such as psycho: http://psyco.sourceforge.net/.

Jan 3 '07 #5

P: n/a
On 2007-01-03 16:50, mm wrote:
More general, maybe there is a speed optimazation docu out there.
At least Alex Martellis "Python in a Nutshell" has a section on
optimization.

I presented this at the last EuroPython conference:
http://sschwarzer.com/download/optim...python2006.pdf

Stefan

--
Dr.-Ing. Stefan Schwarzer
SSchwarzer.com - Softwareentwicklung für Technik und Wissenschaft
http://sschwarzer.com
Jan 3 '07 #6

P: n/a
Hmm.. thanks. I did this changes, but without any performance profits.
Matimus wrote:
Using the '+' operator for string concatonation can be slow, especially
when done many times in a loop.

> pi = pi + str("%04d" % int(e + d/a)) ## this should be fast?! I dont


The accepted solution would be to make pi an array and append to the
end...

pi = [] #create the array (empty)
...
...
pi.append(str("%04d"%int(e+d/a))) # append to it

And when it is time to print the results do the following:

print "".join(pi)

It may look strange, but it is a common Python idiom. You might also
look into an optimizer such as psycho: http://psyco.sourceforge.net/.
Jan 3 '07 #7

P: n/a

Yes. But it still runns very slowly.

If someone is really interested in speed optimization, I can publish my
PI-calc code.

Maybe for some Python compiler/interpreter hackers... ;-)
(There are only 2 while-loops, one within another, and some simple basic
calculations. Nothing special.)
Stefan Schwarzer wrote:
On 2007-01-03 16:50, mm wrote:
>>More general, maybe there is a speed optimazation docu out there.


At least Alex Martellis "Python in a Nutshell" has a section on
optimization.

I presented this at the last EuroPython conference:
http://sschwarzer.com/download/optim...python2006.pdf

Stefan
Jan 4 '07 #8

P: n/a
At Wednesday 3/1/2007 12:50, mm wrote:
>Hmm... it's a question. It was not that easy to translate this #@*?%!
C-Program into readable code and then to Python. But it works.
Because that code was written with C in mind, and uses some C
subtlecies (uhm, got it right?) obscuring the original algorithm.
BTW, do you know where the algorithm was taken from?
>while c*2:
`while c:` does *exactly* the same without computing (and discarding)
an intermediate object.
while (b-1):
If b>=1 initially (as it should, else this will be a loooooong loop),
this is the same as `while b>1:` but without computing (and
discarding) an intermediate object.
b=b-1 ## just for while-loop condition.
I'm not sure if the compiler is smart enough to translate this into `b -= 1`
c = c-14; ## this is code vor the 1st while-loop, BUT bust run after
Same as above
pi = pi + str("%04d" % int(e + d/a)) ## this should be fast?! I dont
Someone else already pointed out how to improve this.

You don't show the rest of the code...
--
Gabriel Genellina
Softlab SRL


__________________________________________________
Preguntá. Respondé. Descubrí.
Todo lo que querías saber, y lo que ni imaginabas,
está en Yahoo! Respuestas (Beta).
¡Probalo ya!
http://www.yahoo.com.ar/respuestas

Jan 4 '07 #9

P: n/a
If someone is really interested in speed optimization, I can publish my
PI-calc code.
I wouldn't mind seeing it. Chances are you will get much better help if
you post your code anyway.

-Matt

Jan 4 '07 #10

P: n/a

Ok, here is the code. It is a translation of the following code, found
on the internet.

* The C is very fast, Python not.
* Target: Do optimization, that Python runs nearly like C.
Auf 800 Stellen in 160 Zeichen...
------
int a=10000,b,c=2800,d,e,f[2801],g;main(){for(;b-c;)f[b++]=a/5;
for(;d=0,g=c*2;c-=14,printf("%.4d",e+d/a),e=d%a)for(b=c;d+=f[b]*a,
f[b]=d%--g,d/=g--,--b;d*=b);}

$ ./a.exe
31415926535897932384626433832795028841971693993751 058209749445923078164062862089
98628034825342117067982148086513282306647093844609 550582231725359408128481117450
28410270193852110555964462294895493038196442881097 566593344612847564823378678316
52712019091456485669234603486104543266482133936072 602491412737245870066063155881
74881520920962829254091715364367892590360011330530 548820466521384146951941511609
43305727036575959195309218611738193261179310511854 807446237996274956735188575272
48912279381830119491298336733624406566430860213949 463952247371907021798609437027
70539217176293176752384674818467669405132000568127 145263560827785771342757789609
17363717872146844090122495343014654958537105079227 968925892354201995611212902196
08640344181598136297747713099605187072113499999983 729780499510597317328160963185

Here the Python:
----------------------------------------------------
#!/usr/bin/python
# -*- coding: utf-8 -*-

## http://de.wikipedia.org/wiki/Pi_(Kreiszahl)
from __future__ import division
from array import array
import decimal
import time
#int a=10000,b,c=2800,d,e,f[2801],g;main(){for(;b-c;)f[b++]=a/5;
# for(;d=0,g=c*2;c-=14,printf("%.4d",e+d/a),e=d%a)for(b=c;d+=f[b]*a,
# f[b]=d%--g,d/=g--,--b;d*=b);}
#
print "\nTiming a 1 million loop 'for loop' ..."
start = time.clock()
for x in range(1000000):
y = x # do something
end = time.clock()
print "Time elapsed = ", end - start, "seconds"
start_prg = time.clock()
#pi=[] ## create an empty array
pi=''

a=10000
b=0
c=5600 ## c=2800
d=0
e=0
f=[] # f[2801]
g=0
counter=c
while 0<=counter+4000:
f.append(2000) # f.append( int(a/5) )
counter=counter-1
# b=b+1

#print "DEBUG: b counter: ", b, counter
while (c*2):
d=0

g=c*2 ## see while condition
## ---------------------------
b=c ## anzahl elemente
#print "DEBUG: before 3 while loop..."
#print "DEBUG: b=", b
while (b-1):
d = d + f[b]*a
g=g-1
f[b] = d%g
d = int(d/g) ## needs cast to int
g=g-1

d=d*b

b=b-1 ## see while condition

#print "DEBUG: d=", d
c = c-14;
#pi.append(str("%04d" % int(e + d/a))) # append to it
pi = pi + str("%04d" % int(e + d/a))
#print "%04d" % int(e + d/a),
## need cast to int
#print int(e + int(d/a))
e = d%a;

#print "".join(pi)
print pi

end_prg = time.clock()
print "Total Time elapsed = ", end_prg - start_prg, "seconds"
## EOF.
-------------------------------

Matimus wrote:
>>If someone is really interested in speed optimization, I can publish my
PI-calc code.


I wouldn't mind seeing it. Chances are you will get much better help if
you post your code anyway.

-Matt
Jan 4 '07 #11

P: n/a
At Wednesday 3/1/2007 22:10, Michael M. wrote:
>Ok, here is the code. It is a translation of the following code, found
on the internet.

* The C is very fast, Python not.
* Target: Do optimization, that Python runs nearly like C.
Why? Python is strong in other aspects, *not* on computation speed.

Anyway, for a nice and rather fast Python implementation (yielding
unlimited digits *without* unlimited storage), see
http://mail.python.org/pipermail/pyt...er/414347.html
The code, for reference:

def pi():
k, a, b, a1, b1 = 2, 4, 1, 12, 4
while 1:
# Next approximation
p, q, k = k*k, 2*k+1, k+1
a, b, a1, b1 = a1, b1, p*a+q*a1, p*b+q*b1
# Yield common digits
d, d1 = a/b, a1/b1
while d == d1:
yield str(d)
a, a1 = 10*(a%b), 10*(a1%b1)
d, d1 = a/b, a1/b1

import sys
sys.stdout.writelines(pi())

(The converted C code does a lot of nonsense in Python terms - maybe
you should try to interpret *what* it does and then reimplement that
using Python)
--
Gabriel Genellina
Softlab SRL


__________________________________________________
Preguntá. Respondé. Descubrí.
Todo lo que querías saber, y lo que ni imaginabas,
está en Yahoo! Respuestas (Beta).
¡Probalo ya!
http://www.yahoo.com.ar/respuestas

Jan 4 '07 #12

P: n/a
Michael M.:
* The C is very fast, Python not.
* Target: Do optimization, that Python runs nearly like C.
Python can't be fast as C for that kind of programs.
Note that your original C program gives less digits than the Python
program.
Your original takes about ~15.2 s on my PC. The following version (I
have just cleaned it up a bit, probably it can be improved) needs about
~0.7 seconds with Psyco and ~5.7 s without (Psyco compilation time and
Python startup times aren't included, if you include them the timings
are ~0.94 s and ~5.96 s).
Compiled with ShedSkin this program needs ~0.29 s (redirecting the
output to a file, because ShedSkin printing is slow still).
from time import clock

def compute_pi():
pi = []
a = 10000
b = d = e = g = 0
c = 5600
f = [2000] * (c + 4000 + 1)
while c:
d = 0
g = c * 2
b = c
while b 1:
d += f[b] * a
g -= 1
f[b] = d % g
d = d // g
g -= 1
d *= b
b -= 1
c -= 14
pi.append("%04d" % int(e + d // a))
e = d % a
return "".join(pi)

import psyco; psyco.bind(compute_pi)
start_time = clock()
print compute_pi()
print "Total time elapsed:", round(clock() - start_time, 2), "s"

Bye,
bearophile

Jan 4 '07 #13

P: n/a
On Thu, 04 Jan 2007 02:10:44 +0100, Michael M. wrote:
print "\nTiming a 1 million loop 'for loop' ..."
start = time.clock()
for x in range(1000000):
y = x # do something
Why not "pass # do nothing"?

end = time.clock()
print "Time elapsed = ", end - start, "seconds"
Are you aware that you're timing the CREATION of a rather large list, as
well as the loop?

[snip code]
#print "".join(pi)
print pi

end_prg = time.clock()
and also the time taken for I/O printing the digits? That can be quite
slow.

For timing functions and small snippets of code, you should investigate
the timeit module. It is very flexible, and will often give more accurate
results than roll-you-own timing.

For example:
>>import timeit
t = timeit.Timer("for i in range(1000000): pass", "pass")
t.timeit(100) # time the loop one hundred times
29.367985010147095

Now compare it to the speed where you only create the list once.
>>t = timeit.Timer("for i in L: pass", "L = range(1000000)")
t.timeit(100)
16.155359983444214
--
Steven D'Aprano

Jan 4 '07 #14

P: n/a

Michael M. wrote:
Ok, here is the code. It is a translation of the following code, found
on the internet.

* The C is very fast, Python not.
* Target: Do optimization, that Python runs nearly like C.
There is an error in the translated code. It returns 1600 digits
instead of 800 digits.
>
counter=c
while 0<=counter+4000:
f.append(2000) # f.append( int(a/5) )
counter=counter-1
# b=b+1
This creates a list f with length 9601. It should have a length of
2801.

I found an explanation of the original C program at
http://rooster.stanford.edu/~ben/maths/pi/code.html

Using the variable names from the above explanation and editing the
code posted by bearophile to match the explanation, I have the
following:

from time import clock

def compute_pi():
pi = []
a = 10000
i = k = b = d = c = 0
k = 2800
r = [2000] * 2801
while k:
d = 0
i = k
while True:
d += r[i] * a
b = 2 * i - 1
r[i] = d % b
d //= b
i -= 1
if i == 0: break
d *= i
k -= 14
pi.append("%04d" % int(c + d // a))
c = d % a
return "".join(pi)

start_time = clock()
pi = compute_pi()
print pi
print "Total time elapsed:", round(clock() - start_time, 2), "s"
print len(pi)

You're original version takes 2.8 seconds on my computer. The above
version takes .36 seconds. I tried a couple of optimizations but
couldn't make any more improvements.

casevh

Jan 4 '07 #15

P: n/a
In <cu***************@nntpserver.swip.net>, Michael M. wrote:
* The C is very fast, Python not.
* Target: Do optimization, that Python runs nearly like C.
As someone else already asked: Why? You can't beat a compiled to machine
code language with an interpreted one when doing integer arithmetic.
counter=c
while 0<=counter+4000:
f.append(2000) # f.append( int(a/5) )
counter=counter-1
# b=b+1
Why do you count so complicated here? It's hard to see that this creates
a list with 2801 elements.
d = int(d/g) ## needs cast to int
Instead of converting the numbers to float by "real" division and then
converting back the result you should do integer division:

d = d // g

or

d //= g
pi = pi + str("%04d" % int(e + d/a))
The `str()` call is unnecessary. And again: try to stick to integer
arithmetic with ``//``.

Here is my attempt to convert the C code, not written with speed in mind
and I was too lazy too time it. :-)

from itertools import izip

def pi():
result = list()
d = 0
e = 0
f = [2000] * 2801
for c in xrange(2800, 0, -14):
for b, g in izip(xrange(c, 1, -1), xrange((c * 2) - 1, 0, -2)):
d += f[b] * 10000
h, f[b] = divmod(d, g)
d = h * b
h, i = divmod(d, 10000)
result.append('%.4d' % (e + h))
e = i
return ''.join(result)

Ciao,
Marc 'BlackJack' Rintsch
Jan 4 '07 #16

P: n/a
Here is my attempt to convert the C code, not written with speed in mind
and I was too lazy too time it. :-)

from itertools import izip

def pi():
result = list()
d = 0
e = 0
f = [2000] * 2801
for c in xrange(2800, 0, -14):
for b, g in izip(xrange(c, 1, -1), xrange((c * 2) - 1, 0, -2)):
d += f[b] * 10000
h, f[b] = divmod(d, g)
d = h * b
h, i = divmod(d, 10000)
result.append('%.4d' % (e + h))
e = i
return ''.join(result)
Your version: .36 seconds

It's ugly, but unrolling the divmod into two separate statements is
faster for small operands. The following takes .28 seconds:

from time import clock
from itertools import izip

def pi():
result = list()
d = 0
e = 0
f = [2000] * 2801
for c in xrange(2800, 0, -14):
for b, g in izip(xrange(c, 1, -1), xrange((c * 2) - 1, 0, -2)):
d += f[b] * 10000
f[b] = d % g
h = d // g
d = h * b
h, i = divmod(d, 10000)
result.append('%.4d' % (e + h))
e = i
return ''.join(result)

start_time = clock()
pi = pi()
print pi
print "Total time elapsed:", round(clock() - start_time, 2), "s"
print len(pi)

casevh

Jan 4 '07 #17

P: n/a
mm <mi*****@mustun.chwrote:
(Yes, I konw whats an object is...)
BTW. I did a translation of a pi callculation programm in C to Python.
(Do it by your own... ;-)
Calc PI for 800 digs(?). (german: Stellen)
int a=10000,b,c=2800,d,e,f[2801],g;main(){for(;b-c;)f[b++]=a/5;
for(;d=0,g=c*2;c-=14,printf("%.4d",e+d/a),e=d%a)for(b=c;d+=f[b]*a,
f[b]=d%--g,d/=g--,--b;d*=b);}
Below you'll find my translation. Note that I improved the algorithm
greatly. Note the nice easy to read algorithm also! This calculated
800 digits of pi in 0.1 seconds and 10,000 digits in 20 seconds.

------------------------------------------------------------
"""
Standalone Program to calculate PI using python only

Nick Craig-Wood <ni**@craig-wood.com>
"""

import sys
from time import time

class FixedPoint(object):
"""A minimal immutable fixed point number class"""
__slots__ = ['value'] # __weakref__
digits = 25 # default number of digits
base = 10**digits

def __init__(self, value=0):
"""Initialise the FixedPoint object from a numeric type"""
try:
self.value = long(value * self.base)
except (TypeError, ValueError), e:
raise TypeError("Can't convert %r into long", value)

def set_digits(cls, digits):
"""Set the default number of digits for new FixedPoint numbers"""
cls.digits = digits
cls.base = 10**digits
set_digits = classmethod(set_digits)

def copy(self):
"Copy self into a new object"
new = FixedPoint.__new__(FixedPoint)
new.value = self.value
return new
__copy__ = __deepcopy__ = copy

def __add__(self, other):
"""Add self to other returning the result in a new object"""
new = self.copy()
new.value = self.value + other.value
return new
__radd__ = __add__

def __sub__(self, other):
"""Subtract self from other returning the result in a new object"""
new = self.copy()
new.value = self.value - other.value
return new
def __rsub__(self, other):
new = self.copy()
new.value = other.value - self.value
return new

def __mul__(self, other):
"""
Multiply self by other returning the result in a new object. The
number of digits is that of self.

Note that FixedPoint(x) * int(y) is much more efficient than
FixedPoint(x) * FixedPoint(y)
"""
new = self.copy()
if isinstance(other, (int, long)):
# Multiply by scalar
new.value = self.value * other
else:
new.value = (self.value * other.value) // other.base
return new
__rmul__ = __mul__

def __div__(self, other):
"""
Divide self by other returning the result in a new object. The
number of digits is that of self.

Note that FixedPoint(x) / int(y) is much more efficient than
FixedPoint(x) / FixedPoint(y)
"""
new = self.copy()
if isinstance(other, (int, long)):
# Divide by scalar
new.value = self.value // other
else:
new.value = (self.value * other.base) // other.value
return new
def __rdiv__(self, other):
if not isinstance(other, FixedPoint):
other = FixedPoint(other, self.digits)
return other.__div__(self)

def __str__(self):
"""
Convert self into a string. This will be the integral part of
the number possibly preceded by a - sign followed by a . then
exactly n digits where n is the number of digits specified on
creation.
"""
if self.value < 0:
s = str(-self.value)
else:
s = str(self.value)
s = s.zfill(self.digits + 1)
s = s[:-self.digits] + "." + s[-self.digits:]
if self.value < 0:
s = "-" + s
return s

def __abs__(self):
"""Return the absolute value of self"""
new = self.copy()
if new.value < 0:
new.value = - new.value
return new

def __cmp__(self, other):
"""Compare self with other"""
return cmp(self.value, other.value)

def sqrt(n):
"""
Return the square root of n. It uses a second order
Newton-Raphson convgence. This doubles the number of significant
figures on each iteration. This will work for any finite
precision numeric type.
"""
x = n / 2 # FIXME find a better guess!
old_delta = None
while 1:
x_old = x
x = (x + n / x) / 2
delta = abs(x - x_old)
if old_delta is not None and delta >= old_delta:
break
old_delta = delta
return x

def pi_agm(one = 1.0, sqrt=sqrt):
"""
Brent-Salamin Arithmetic-Geometric Mean formula for pi. This
converges quadratically.

FIXME can increase the number of digits in each iteration to speed up?
"""
_1 = one
_2 = _1 + _1
a = _1
b = _1/sqrt(_2)
t = _1/2
k = 1
two_to_the_k = 2L
old_delta = None

while 1:
s = (a + b ) / 2
d = a - s
d = d * d
a = s
s = s * s
t = t - two_to_the_k * d
b = sqrt(s - d)
delta = abs(a - b)
# print k, delta, old_delta
if old_delta is not None and delta >= old_delta:
break
old_delta = delta
k = k + 1
two_to_the_k *= 2

a = a + b
return ( a * a ) / ( _2 * t)

if __name__ == "__main__":
try:
digits = int(sys.argv[1])
except:
digits = 100
print "Calculating",digits,"digits of pi"
start = time()
FixedPoint.set_digits(digits)
_1 = FixedPoint(1)
print pi_agm(_1, sqrt=sqrt)
dt = time() - start
print "That took",dt,"seconds"
------------------------------------------------------------
But Python is much slower here then C here. I used a while-loop within a
while-loop. Not that much calculation here.
There is more than one way to skin this cat. A better algorithm helps
a lot.

If I really wanted lots of digits of pi from python I'd do it like
this. This is an example of the right tool for the job. 0.5ms looks
pretty good to me!
>>import math, gmpy
bits = int(800/math.log10(2))
start = time(); pi=gmpy.pi(bits); dt = time()-start
print "That took",dt,"seconds"
That took 0.00055193901062 seconds
>>pi
mpf('3.1415926535897932384626433832795028841971693 99375105820974944592307816406286208998628034825342 11706798214808651328230664709384460955058223172535 94081284811174502841027019385211055596446229489549 30381964428810975665933446128475648233786783165271 20190914564856692346034861045432664821339360726024 91412737245870066063155881748815209209628292540917 15364367892590360011330530548820466521384146951941 51160943305727036575959195309218611738193261179310 51185480744623799627495673518857527248912279381830 11949129833673362440656643086021394946395224737190 70217986094370277053921717629317675238467481846766 94051320005681271452635608277857713427577896091736 37178721468440901224953430146549585371050792279689 25892354201995611212902196086403441815981362977477 13099605187072113499999983729780499510597317328160 96318595024459455e0',2657)
>>>
--
Nick Craig-Wood <ni**@craig-wood.com-- http://www.craig-wood.com/nick
Jan 4 '07 #18

P: n/a

Mainly, it was fload-div. Changed to int-div (python //) runs faster.

Yes, this "gmpy" sounds good for calc things like that.
But not available on my machine.
ImportError: No module named gmpy

Anyway, thanks for posting. This gmpy module can be very intersting.
But right now, the focus was, if it's possible to translate a strange C
code into Python. And it is. Sure ;-)

Maybe, someone can also translate a very simple algorithm for calc a
range of PI in Python. (Available Code for C.)
http://crd.lbl.gov/~dhbailey/

http://crd.lbl.gov/~dhbailey/bbp-codes/piqpr8.c

But seems to be removed? It was working last days. (Maybe removed from
server.)
An other well known is: (but someone implemented this already in Python,
slow, and, depends on the previouse calc values, like most other
algorithmes.)

pi 01 01 01 01 01 01 01 01 01 01
__ = __ - __ + __ - __ + __-__+__-__+__-__
04 01 03 05 07 09 11 13 15 17 19 ...

pi/4 = 1/1 - 1/3 + 1/5 - 1/7 + 1/9 ...
Or things like that:
http://numbers.computation.free.fr/C...am/pifast.html

There is more than one way to skin this cat. A better algorithm helps
a lot.

If I really wanted lots of digits of pi from python I'd do it like
this. This is an example of the right tool for the job. 0.5ms looks
pretty good to me!
>>import math, gmpy
>>bits = int(800/math.log10(2))
>>start = time(); pi=gmpy.pi(bits); dt = time()-start
>>print "That took",dt,"seconds"
That took 0.00055193901062 seconds
>>pi
mpf('3.1415926535897932384626433832795028841971693 99375105820974944592307816406286208998628034825342 11706798214808651328230664709384460955058223172535 94081284811174502841027019385211055596446229489549 30381964428810975665933446128475648233786783165271 20190914564856692346034861045432664821339360726024 91412737245870066063155881748815209209628292540917 15364367892590360011330530548820466521384146951941 51160943305727036575959195309218611738193261179310 51185480744623799627495673518857527248912279381830 11949129833673362440656643086021394946395224737190 70217986094370277053921717629317675238467481846766 94051320005681271452635608277857713427577896091736 37178721468440901224953430146549585371050792279689 25892354201995611212902196086403441815981362977477 13099605187072113499999983729780499510597317328160 96318595024459455e0',2657)
>>>
Jan 4 '07 #19

P: n/a
Yes, this "gmpy" sounds good for calc things like that.
But not available on my machine.
ImportError: No module named gmpy
What type of machine?

The home page for gmpy is http://sourceforge.net/projects/gmpy/

I have Windows versions available at http://home.comcast.net/~casevh/

casevh

Jan 4 '07 #20

P: n/a
ca****@comcast.net wrote:
>>Yes, this "gmpy" sounds good for calc things like that.
But not available on my machine.
ImportError: No module named gmpy


What type of machine?
Windoof with Cygwin.
>
The home page for gmpy is http://sourceforge.net/projects/gmpy/

I have Windows versions available at http://home.comcast.net/~casevh/

casevh
If the binary works, I don't know. Maybe there is also a gmpy module
already in Cygwin.
Jan 4 '07 #21

P: n/a
Michael M. <mi*****@mustun.chwrote:
Yes, this "gmpy" sounds good for calc things like that.
But not available on my machine.
ImportError: No module named gmpy
Sorry, I should have said - you'll need to download that from

http://gmpy.sourceforge.net/
Anyway, thanks for posting. This gmpy module can be very intersting.
But right now, the focus was, if it's possible to translate a strange C
code into Python. And it is. Sure ;-)

Maybe, someone can also translate a very simple algorithm for calc a
range of PI in Python. (Available Code for C.)
http://crd.lbl.gov/~dhbailey/
Here are some more pi calculation methods. You'll need the FixedPoint
class in my last posting, or you cah use gmpy mpf~s also.

_1 = FixedPoint(1)
# or
_1 = gmpy.mpf(1, bits)
# or (for testing)
_1 = 1.0

_3 = _1 * 3

def arctan(x):
"""
Calculate arctan(x)

arctan(x) = x - x**3/3 + x**5/5 - ... (-1 <= x <= 1)
"""
total = x
power = x
divisor = 1
old_delta = None
while 1:
power *= x
power *= x
power = -power
divisor += 2
old_total = total
total += power / divisor
delta = abs(total - old_total)
if old_delta is not None and delta >= old_delta:
break
old_delta = delta
return total

def pi_machin():
return 4*(4*arctan(_1/5) - arctan(_1/239))

def pi_ferguson():
return 4*(3*arctan(_1/4) + arctan(_1/20) + arctan(_1/1985))

def pi_hutton():
return 4*(2*arctan(_1/3) + arctan(_1/7))

def pi_gauss():
return 4*(12*arctan(_1/18) + 8*arctan(_1/57) - 5*arctan(_1/239))

def pi_euler():
return 4*(5*arctan(_1/7) + 2*arctan(_3/79))
--
Nick Craig-Wood <ni**@craig-wood.com-- http://www.craig-wood.com/nick
Jan 5 '07 #22

This discussion thread is closed

Replies have been disabled for this discussion.