Hi,
one hint for doing *fast* computations:
Use Cython, rather than pure Python!
>From my point of view, Cython is one of the main strengths of Sage:
You can do your computations with C data types, and you should do so,
because it is faster than using "complicated" data types such as Sage-
Integer. But you still have the nice feel of Python.
A few months ago, I made up a Mandelbrot creator myself, just for fun
and for practicing. It consists of a class "Mandelbrot". In the
initiation, you can pass various parameters (the x- and y-range, the
number of pixel on both axes, the maximal number of iterations).
If i understand the meaning of your parameters correctly, the
following produces the image that you want:
sage: attach Mandelbrot.pyx
# the code is below
sage: MB=Mandelbrot(xmin=-2,xmax=1,xsteps=200,
ymin=-1,ymax=1,ysteps=200, maxIter=200)
# i think this is your example
sage: time MB.make()
CPU times: user 0.34 s, sys: 0.00 s, total: 0.34 s
Wall time: 0.34 s
sage: show(MB.img)
# Objects of class <Mandelbrot> have an attribute img, which is
computed with the method make()
# and shown with show()
Main reasons why the code is fast:
1. The computations are done with data of type <double>. For example,
cdef double cr
cdef double ci
describes real and imaginary part of the complex number for which we
iterate
2. All loop variables are of type <int>. You do
for i in range(n):
...
which is slow. Faster is
for i in xrange(n):
...
Using C types is still faster:
cdef int i
cdef int N=n # hence, the <Integer> n is converted into a C
<int>
for i from 0 <= i < N:
...
One very useful feature of Sage for developing fast code is
"annotation". By typing
> sage -cython -a Mandelbrot.pyx
from the command line, i get a file Mandelbrot.html, which i can view
with a browser. That file shows me my own code -- with hints on how
efficient the code is (the more yellow a line of code the more slow it
is executed).
The code below also relies on an experimental package for image
manipulation, that you can get (from the command line) with
> sage -i PIL-1.1.5
Cheers
Simon
###################################################
## Here is the code of Mandelbrot.pyx
import sage
import sage.all
import Image
import colorsys
# requires PIL
class Mandelbrot:
def __init__(self,xmin=-2.4,xmax=1.3,xsteps=600,
ymin=-1.7,ymax=1.7,ysteps=600, maxIter=120):
"""
Class for creating a picture of the Mandelbrot set.
MB = Mandelbrot()
The picture can be created with MB.make() and shown with
show(MB.img)
"""
self.colors = {} # used for cacheing
self.img = Image.new("RGB",(xsteps,ysteps),(0,0,0))
self.d2 = float(4) # sequence is unbounded if it comes outside
a circle of radius 2
self.N = int(maxIter) # maximal number of iterations
if xmin>=xmax:
raise ValueError, "xmax<=xmin"
self.xmin=float(xmin)
self.xmax=float(xmax)
if xsteps==0:
raise ValueError, "xsteps==0"
self.dx=(self.xmax-self.xmin)/float(xsteps)
self.xsteps = int(xsteps)
if ymin>=ymax:
raise ValueError, "ymax<=ymin"
self.ymin=float(ymin)
self.ymax=float(ymax)
if ysteps==0:
raise ValueError, "ysteps==0"
self.dy=(self.ymax-self.ymin)/float(ysteps)
self.ysteps = int(ysteps)
def color(self,i):
"computes the color of a pixel if the iteration leaves the
circle of radius 2 after i steps"
# the result is cached in self.colors
if i==self.N:
return (0,0,0)
CL = self.colors.get(i)
if CL:
return CL
(r,g,b) = colorsys.hsv_to_rgb(1.0-float(i)/self.N,1.0,1.0)
CL = (int(256*r)-1,int(256*g)-1,int(256*b)-1)
self.colors[i] = CL
return CL
def point(self,x,y):
"""
return the number of iterations after which point (x,y)
leaves the circle of radius 2.
Maximum number of iterations is self.N
"""
# Morally, this is C code!
cdef double cr = self.xmin+x*self.dx
cdef double ci = self.ymin+y*self.dy
cdef double zr = 0
cdef double zi = 0
cdef double zrneu,zineu
cdef double zr2=0
cdef double zi2=0
cdef double d2 = self.d2
cdef int N = self.N
cdef int i
for i from 0 <= i <= N:
zrneu = zr2-zi2 + cr
zineu = 2*zr*zi + ci
zr=zrneu
zi=zineu
zr2=zr**2
zi2=zi**2
if (zr2+zi2)>d2:
self.img.putpixel((x,y),self.color(i))
return i
return self.N
def make (self):
"Create a picture of the Mandelbrot set"
cdef int i
cdef int j
cdef int xm = self.xsteps
cdef int ym = self.ysteps
for i from 0 <= i < xm:
for j from 0 <= j < ym:
self.point(i,j)
--~--~---------~--~----~------------~-------~--~----~
To post to this group, send email to sage-support@googlegroups.com
To unsubscribe from this group, send email to [EMAIL PROTECTED]
For more options, visit this group at
http://groups.google.com/group/sage-support
URLs: http://www.sagemath.org
-~----------~----~----~----~------~----~------~--~---
