[wxPython] Mandelbrot - further improvements

classic Classic list List threaded Threaded
10 messages Options
Reply | Threaded
Open this post in threaded view
|

[wxPython] Mandelbrot - further improvements

Markus Gritsch
Hi!

Some time ago (Tue, 02 Nov 1999), there was an interesting thread on this
mailing-list about a mandelbrot optimization.  Well, I took a look at the
result,
which was posted by John (guy who started the thread), and I found, that there
is the possibility for some improvements, which he left out:

*) I modified the while loop of the mandelbrot function in such a way, that the

calculations of z_ = z_ * z_ and z_ = z_ + c_ are done in-place.  This reduced
the computation time of this function on my box from 30 seconds to 26 seconds.

*) I found it annoying, that the preparation of the mandelbrot data after the
first OnPaint event takes about half as long as the mandelbrot function itself.

Ionel suggested a better solution, but is wasn't included in the final code.
Now it is used, and it takes 0.076 seconds on my box compared to 13.8 seconds
before it was used.  Thats about 180 times faster :-))

*) I also fixed the size of the wxFrame after the Fit().  Now it has exactly
the size of the wxImage.

So the overall computation time is now reduced from 44 seconds to 26 seconds.

If MAX_ITERATIONS is greater than 255, the output behaves a little bit strange,

but in John's version it was even worse, because the mandelbrot function was
called on and on without displaying anything.  Maybe this could be improved
further.

Because the application is blocking during calculation Robin suggested to use
a wxIdle-handler.  Has someone a code snippet which shows how to uses such
a handler?  Or is someone able to include it in the mandelbrot-program?

Have fun, and here's the code

---
Markus


"""
Mandelbrot program, for no particular reason.
John Farrell 1999/11/02

CREDITS:
Thanks to Robin Dunn for showing me how to use images and bitmaps.
Thanks to Ionel Simionescu for showing me how to use the Numeric package,
and consequently writing all of the computation code.
Thanks to Alexander Smishlajev for suggestions on optimising loops
and drawing.
"""

from wxPython.wx import *
import Numeric

WIDTH = 200
HEIGHT = 200

#coords = ((-1.0, -1.5), (0.0, -0.5))
coords = ((-2.0, -1.5), (1.0, 1.5))
WEIGHTS = (16, 1, 32)
MAX_ITERATIONS = 255
zero = complex(0.0, 0.0)

def setup():
    # we use a single array for the real values corresponding to the x
coordinates
    diff = coords[1][0] - coords[0][0]
    xs = 0j + (coords[0][0] + Numeric.arange(WIDTH).astype(Numeric.Float32) *
diff / WIDTH )
    # we use a single array for the imaginary values corresponding to the
ycoordinates
    diff = coords[1][1] - coords[0][1]
    ys = 1j * (coords[0][1] + Numeric.arange(HEIGHT).astype(Numeric.Float32) *
diff / HEIGHT )
    # we build <c> in direct correpondence with the pixels in the image
    c = Numeric.add.outer(ys, xs)
    z = Numeric.zeros((HEIGHT, WIDTH)).astype(Numeric.Complex32)
    c.shape = (HEIGHT * WIDTH,)
    z.shape = (HEIGHT * WIDTH,)
    return (c, z)

def mandelbrot(progressMonitor):
    import time
    t = time.clock()
    (c, z) = setup()
    iterations = 0
    i_no = Numeric.arange(WIDTH * HEIGHT)       # non-overflow indices
    data = MAX_ITERATIONS + Numeric.zeros(WIDTH * HEIGHT)
    # initialize the "todo" arrays;
    # they will contain just the spots where we still need to iterate
    c_ = Numeric.array(c)
    z_ = Numeric.array(z)
    while iterations <= MAX_ITERATIONS and len(i_no):
        progressMonitor.progress(iterations, MAX_ITERATIONS)
        z_ = z_ * z_ + c_
        overflow = Numeric.greater_equal(abs(z_), 2.0)
        not_overflow = -(overflow - 1)
        # get the indices where overflow occured
        overflowIndices = Numeric.compress(overflow, i_no)
        # set the pixel indices there
        for idx in overflowIndices:
            data[idx] = iterations
        # compute the new array of non-overflow indices
        i_no = Numeric.compress(not_overflow, i_no)
        # update the todo arrays
        c_ = Numeric.compress(not_overflow, c_)
        z_ = Numeric.compress(not_overflow, z_)
        iterations = iterations + 1
    progressMonitor.progress(MAX_ITERATIONS, MAX_ITERATIONS)
    print time.clock() - t
    return data

class MandelCanvas(wxWindow):
    def __init__(self, parent, id = -1):
        wx.wxWindow.__init__(self, parent, id)
        self.parent = parent
        self.border = (1,1)
        self.SetSize(wxSize(WIDTH, HEIGHT))
        self.SetBackgroundColour(wx.wxNamedColour("black"))
        self.bitmap = None
        self.colours = Numeric.zeros( (MAX_ITERATIONS + 1, 3) )
        arangeMax = Numeric.arange(0, MAX_ITERATIONS + 1)
        self.colours[:,0] = Numeric.clip(arangeMax * WEIGHTS[0], 0,
MAX_ITERATIONS)
        self.colours[:,1] = Numeric.clip(arangeMax * WEIGHTS[1], 0,
MAX_ITERATIONS)
        self.colours[:,2] = Numeric.clip(arangeMax * WEIGHTS[2], 0,
MAX_ITERATIONS)

    def OnPaint(self, event):
        if not self.bitmap:
            data = mandelbrot(self.parent)
            import time
            t = time.clock()
            pixels = [chr(0)] * len(data) * 3
            for i in range(len(data)):
                di = data[i]
                i3 = i * 3
                for j in range(3):
                   pixels[i3 + j] = chr(self.colours[di, j])
            img = wxEmptyImage(WIDTH, HEIGHT)
            d = string.join(pixels, '')
            img.SetData(d)
            self.bitmap = img.ConvertToBitmap()
            print time.clock() - t
        dc = wxPaintDC(self)
        dc.BeginDrawing()
        dc.DrawBitmap(self.bitmap, 0, 0, false)
        dc.EndDrawing()

class MyFrame(wxFrame):
    def __init__(self, parent, ID, title):
        wxFrame.__init__(self, parent, ID, title)
        self.CreateStatusBar()
        self.SetStatusText("")
        self.Centre(wxBOTH)
        mdb = MandelCanvas(self)

    def progress(self, done, outof):
        self.SetStatusText("%d/%d" % (done, outof))

class MyApp(wxApp):
    def OnInit(self):
        frame = MyFrame(NULL, -1, "Mandelbrot")
        frame.Fit()
        frame.Show(true)
        self.SetTopWindow(frame)
        return true

app = MyApp(0)
app.MainLoop()



_______________________________________________
wxPython-users maillist  -  [hidden email]
http://starship.python.net/mailman/listinfo/wxpython-users



Reply | Threaded
Open this post in threaded view
|

Re: [wxPython] Mandelbrot - further improvements

Markus Gritsch
> Have fun, and here's the code

Well, somehow I managed to attach the original code in my last posting ?-)

Here is the improved one (now really :-)

#! /usr/bin/env python
#######################################################################################94x78##

"""
Mandelbrot program, for no particular reason.
John Farrell 1999/11/02

CREDITS:
Thanks to Robin Dunn for showing me how to use images and bitmaps.
Thanks to Ionel Simionescu for showing me how to use the Numeric package,
and consequently writing all of the computation code.
Thanks to Alexander Smishlajev for suggestions on optimising loops
and drawing.
"""

from wxPython.wx import *
import Numeric

WIDTH = 200
HEIGHT = 200

#coords = ((-1.0, -1.5), (0.0, -0.5))
coords = ((-2.0, -1.5), (1.0, 1.5))
WEIGHTS = (16, 1, 32)
MAX_ITERATIONS = 255
zero = complex(0.0, 0.0)

def setup():
    # we use a single array for the real values corresponding to the x coordinates

    diff = coords[1][0] - coords[0][0]
    xs = 0j + (coords[0][0] + Numeric.arange(WIDTH).astype(Numeric.Float32) * diff
/ WIDTH )
    # we use a single array for the imaginary values corresponding to the
ycoordinates
    diff = coords[1][1] - coords[0][1]
    ys = 1j * (coords[0][1] + Numeric.arange(HEIGHT).astype(Numeric.Float32) *
diff / HEIGHT )
    # we build <c> in direct correpondence with the pixels in the image
    c = Numeric.add.outer(ys, xs)
    z = Numeric.zeros((HEIGHT, WIDTH)).astype(Numeric.Complex32)
    c.shape = (HEIGHT * WIDTH,)
    z.shape = (HEIGHT * WIDTH,)
    return (c, z)

def mandelbrot(progressMonitor):
    import time
    t = time.clock()
    (c, z) = setup()
    iterations = 0
    i_no = Numeric.arange(WIDTH * HEIGHT)       # non-overflow indices
    data = MAX_ITERATIONS + Numeric.zeros(WIDTH * HEIGHT)
    # initialize the "todo" arrays;
    # they will contain just the spots where we still need to iterate
    c_ = Numeric.array(c).astype(Numeric.Complex32)
    z_ = Numeric.array(z).astype(Numeric.Complex32)
    while iterations <= MAX_ITERATIONS and len(i_no):
        progressMonitor.progress(iterations, MAX_ITERATIONS)
        # do the calculations in-place
        # On my box the mandelbrot calculation now takes 26 seconds
        # compared to 30 seconds which were needed without in-place notation.
        Numeric.multiply(z_, z_, z_)
        Numeric.add(z_, c_, z_)
        overflow = Numeric.greater_equal(abs(z_), 2.0)
        not_overflow = -(overflow - 1)
        # get the indices where overflow occured
        overflowIndices = Numeric.compress(overflow, i_no)
        # set the pixel indices there
        for idx in overflowIndices:
            data[idx] = iterations
        # compute the new array of non-overflow indices
        i_no = Numeric.compress(not_overflow, i_no)
        # update the todo arrays
        c_ = Numeric.compress(not_overflow, c_)
        z_ = Numeric.compress(not_overflow, z_)
        iterations = iterations + 1
    progressMonitor.progress(MAX_ITERATIONS, MAX_ITERATIONS)
    print time.clock() - t
    return data

class MandelCanvas(wxWindow):
    def __init__(self, parent, id = -1):
        wx.wxWindow.__init__(self, parent, id)
        self.parent = parent
        self.border = (1,1)
        self.SetSize(wxSize(WIDTH, HEIGHT))
        self.SetBackgroundColour(wx.wxNamedColour("black"))
        self.bitmap = None
        self.colours = Numeric.zeros( (MAX_ITERATIONS + 1, 3) )
        arangeMax = Numeric.arange(0, MAX_ITERATIONS + 1)
        self.colours[:,0] = Numeric.clip(arangeMax * WEIGHTS[0], 0,
MAX_ITERATIONS)
        self.colours[:,1] = Numeric.clip(arangeMax * WEIGHTS[1], 0,
MAX_ITERATIONS)
        self.colours[:,2] = Numeric.clip(arangeMax * WEIGHTS[2], 0,
MAX_ITERATIONS)

    def OnPaint(self, event):
        if not self.bitmap:
            data = mandelbrot(self.parent)
            import time
            t = time.clock()

            #-- John Farrell's version: Takes 13.8 seconds on my box,
            #-- that's about half as long as the mandelbrot calculation,
            #-- which takes 26 seconds on my box.
            ##pixels = [chr(0)] * len(data) * 3
            ##for i in range(len(data)):
            ##    di = data[i]
            ##    i3 = i * 3
            ##    for j in range(3):
            ##       pixels[i3 + j] = chr(self.colours[di, j])
            ##img = wxEmptyImage(WIDTH, HEIGHT)
            ##d = string.join(pixels, '')
            ##img.SetData(d)
            ##self.bitmap = img.ConvertToBitmap()

            #-- Ionel Simionescu's version: takes 0.076 seconds on my box,
            #-- 180 times faster than John's version :-)
            data.shape = (HEIGHT, WIDTH)
            # build the pixel values
            pixels = Numeric.take( self.colours, data )
            # create the image data
            self.bitmap = pixels.astype(Numeric.UnsignedInt8).tostring()
            # create the image itself
            image = wxEmptyImage(WIDTH, HEIGHT)
            image.SetData(self.bitmap)
            self.bitmap = image.ConvertToBitmap()

            print time.clock() - t
        dc = wxPaintDC(self)
        dc.BeginDrawing()
        dc.DrawBitmap(self.bitmap, 0, 0, false)
        dc.EndDrawing()

class MyFrame(wxFrame):
    def __init__(self, parent, ID, title):
        wxFrame.__init__(self, parent, ID, title)
        self.CreateStatusBar()
        self.SetStatusText("")
        self.Centre(wxBOTH)
        mdb = MandelCanvas(self)

        self.SetAutoLayout(true)
        self.Layout()
        self.Fit()
        sz = self.GetClientSize()
        self.SetClientSize(wxSize(sz.width-7, sz.height-14))

    def progress(self, done, outof):
        self.SetStatusText("%d/%d" % (done, outof))

class MyApp(wxApp):
    def OnInit(self):
        frame = MyFrame(NULL, -1, "Mandelbrot")
        frame.Show(true)
        self.SetTopWindow(frame)
        return true

app = MyApp(0)
app.MainLoop()



_______________________________________________
wxPython-users maillist  -  [hidden email]
http://starship.python.net/mailman/listinfo/wxpython-users



Reply | Threaded
Open this post in threaded view
|

Re: [wxPython] Mandelbrot - further improvements

Markus Gritsch
Ionel Simionescu wrote:

> Hi Markus,
>
> I do not understand how this snippet intends to work:
>
>         Numeric.multiply(z_, z_, z_)
>         Numeric.add(z_, c_, z_)
>
> In my tests, Numeric.add(a, b, c) does not write-in any of the arrays passed
> as arguments, but simply returns the sum.

Well, for me it works - maybe this was a bug which is now fixed?  I use NumPy
1.0 which is available as a Windows installer.

The only thing I want to claim about is, that the NumPy reference pages don't
give any information about in-place operations.  The only source of information
I could find was the NumPy tutorial, but even there it isn't clear if the array
in which the result is stored must be the first or the last argument.  But this
can be found out easily by doing something like you suggested (Numeric.add(a,
b, c)), I did it myself, and as already sayed, for me it worked fine.

Have fun, Markus



_______________________________________________
wxPython-users maillist  -  [hidden email]
http://starship.python.net/mailman/listinfo/wxpython-users



Reply | Threaded
Open this post in threaded view
|

Re: [wxPython] Mandelbrot - further improvements

Markus Gritsch
Ionel Simionescu wrote:

> Hi,
>
> It seems to work indeed: if there are three arguments and the first and the
> last are the same, this argument is both source and target.
>
> Did not know that.

I didn't know that too.  I thought, I saw this notation working for three
arbitrary specified arrays.  o.k. if it only works if the first and last
argument is the same, there is no possibility of specifying an input and an
output array ;-)

Be well,
Markus


_______________________________________________
wxPython-users maillist  -  [hidden email]
http://starship.python.net/mailman/listinfo/wxpython-users



Reply | Threaded
Open this post in threaded view
|

Re: [wxPython] Mandelbrot - further improvements

Markus Gritsch
In reply to this post by Markus Gritsch
Ionel Simionescu wrote:

> Here's a little comparison between the standard way to write in place and
> the one just being revealed.
>
> The difference is worth noting.
>
> ---
>
> from Numeric import *
> import time
>
> N = 16000*16
> T = 40
>
> a = ones((N,), Float)
> b = ones((N,), Float)
>
> t0 = time.clock()
> for t in range(T): a[:] = a+b
> t1 = time.clock()
> print ' ready! Total running time:  %f ' %  (t1 - t0)
>
> t0 = time.clock()
> for t in range(T): add(a,b,a)
> t1 = time.clock()
> print ' ready! Total running time:  %f ' %  (t1 - t0)


_______________________________________________
wxPython-users maillist  -  [hidden email]
http://starship.python.net/mailman/listinfo/wxpython-users



Reply | Threaded
Open this post in threaded view
|

Re: [wxPython] Mandelbrot - further improvements

Markus Gritsch
In reply to this post by Markus Gritsch
Ionel Simionescu wrote:

> As a side note: the standard write-in-place (a[:] = a+b) is significantly
> slower and does not seem to save any memory in comparison to the normal
> "name-rebinding" (a = a+b).

That's clear because writing a[:] = a+b doesn't do anything "in-place" at all.
I even wouldn't call this "in-place".  Because when the right hand side is
evaluated, a new object is generated, and therefore a new memory block is
allocated, and thats quite costly.  So it's clear, that no memory is saved
compared to the normal name rebinding.  And the fact that it is slower may be
caused by the additional slicing you used in a[:] = a+b.

---
Markus

_______________________________________________
wxPython-users maillist  -  [hidden email]
http://starship.python.net/mailman/listinfo/wxpython-users



Reply | Threaded
Open this post in threaded view
|

RE: [wxPython] Mandelbrot - further improvements

Mike Fletcher-4
In reply to this post by Markus Gritsch
Okay, following version is about 2x the speed on my PIII-500 (~7.3 seconds
down to ~3.25 seconds).  Mostly is a rewriting of compress to just return
indices, and elimination of the expensive arange creations on every
iteration (which looks weird, but keeps things fast). Did a few
trivial-impact tweaks as I played. Finally, inlined the rewritten compress
function (it's only two lines), doesn't do a lot, but does reduce some
overhead.

Optimisation has got down to the point where some knowledge of the algo
and/or performance of the C functions is required (I have neither). (Why
doesn't profile report C function time anyway?  Doesn't need to trace
through sub-calls in the DLL, but it would be nice to know how fast the
various functions are...)  Did try localising variables, wasn't a
spectacular speedup.

I did add a timer-based callback, but the system is still non-responsive
during calculation.

Enjoy all,
Mike



mandelbrot.py (7K) Download Attachment
Reply | Threaded
Open this post in threaded view
|

RE: [wxPython] Mandelbrot - further improvements

Mike Fletcher-4
Summary:
        Mandelbrot Module, calculates (and displays) Mandelbrot set using Numeric
extensions.  Has been optimised slightly on the wxPython list, module
attached...

Another few 1/10ths of a second improvement (~3.3s down to ~2.9s), as my
part has ceased to have any relevance to wxPython, I should probably switch
to the main list :) , this will be my last post to the wxPython list on this
topic :) .  I added code to the module so that it will run without wxPython
installed (it will use wxPython if it's available, otherwise just prints
status to command line).

Enjoy,
Mike



_______________________________________________
wxPython-users maillist  -  [hidden email]
http://starship.python.net/mailman/listinfo/wxpython-users



Reply | Threaded
Open this post in threaded view
|

[wxPython] last argument of ufuncs

Markus Gritsch
In reply to this post by Markus Gritsch
Ionel Simionescu wrote:

> Hi,
>
> It seems to work indeed: if there are three arguments and the first and the
> last are the same, this argument is both source and target.
>
> Did not know that.

Yesterday I gave Numeric another try and I was able to do the following without
any problem:

a=arange(10)
b=arange(10)
c=arange(10)

multiply(a,b,c)

and the square of each number was stored in c ...

works for me too if the first and last argument are different arrays

Be well,
Markus

_______________________________________________
wxPython-users maillist  -  [hidden email]
http://starship.python.net/mailman/listinfo/wxpython-users



Reply | Threaded
Open this post in threaded view
|

Re: [wxPython] last argument of ufuncs

Ionel Simionescu

Hi Markus,

I attach a post I sent to comp.lang.python few days ago, and some follow-ups
from David Ascher.

----------------------------------------------------------------------------
-------------


Hi,

In a post on the wxPython list, Markus Gritsch suggested that array ops such
as
a = a+b
can be done faster by using write-in place:
Numric.add(a,b,a)

Playing around I see now that other ufuncs allow a write-in target
specification as their last argument. For example,

Numeric.sqrt(a,b)

will put the result in b (if its dimension and type are Ok).

The docs do not tell much about this.
Can anyone confirm/clarify/comment ?

Thanks,
ionel

---
>
> The docs do not tell much about this.
> Can anyone confirm/clarify/comment ?

1) Confirm: yes
2) Clarify: well, I think you pretty much summed it up
3) Comment: Watch out for 'slices'.  If you do

a = arange(10)
b = a[::-1]
add(a,b,b)
print b

you might be suprised =)

--david ascher

----------------------------------------------------------------------------
-------------


> I am [surprised] indeed. :-)
> How do the slices work under the hood to trigger that ?

  a = arange(10)

that's a contiguous, no-fancy-anything array with a contiguous block of
increasing numbers in its data space.

  b = a[::-1]

this makes no copies (unlike lists!).  Instead it creates an array which
shares a's data space, but notes that the data space for b starts at the
end of b and goes backwards through a's data space.

  add(a,b,a)

this goes through a's data space forwards (to read a's data) and
backwards (to read b's data) and writes them in a's data space (forwards).
A mess ensues.

=)

--david

PS: IMHO, it's a design bug.  Too late to fix.


----------------------------------------------------------------------------
-------------


For some reason, David did not post his messages.


ionel


_______________________________________________
wxPython-users maillist  -  [hidden email]
http://starship.python.net/mailman/listinfo/wxpython-users