[SciPy-user] how to crash scipy in 3s! (memory leak?)

[David Cournapeau]

Nicolas Chopin wrote:
>     Hi,
> the program below eats all the memory of my computer
> in approx. 3 seconds, making it pretty unstable.
>
> from scipy import *
> def f(a):
>     return rand(1)+a
>
> n = 1e8
> for i in range(n):
>     x += f(0.1)
>
> It looks like a problem with rand()? e.g. rand() forgets
> to free memory? Or am I doing something "forbidden"?
> like returning an array, which creates a reference
> which is never deleted?
> In this case, what is the correct way to do something like
> def f(x):
>     return some_function_of(x, rand(1) )
>
> I use Ubuntu Feisty, Scipy 0.5.2, Numpy 1.01, Python 2.5.1
> (all 3 from Ubuntu repositories).
If it is crashing, it is a bug. If it is eating all your memory, then it 
is working fine :) First, doing loop as you do in an interpreted 
language (like python) is inherently slow for several reasons: don't do 
it. You could say the whole point of numpy is to give you abstractions 
to avoid looping like your doing.

If I understand correctly, you want to create n random scalars, and sum 
them up, right ? First, instead of calling n times rand(1), you can call 
rand(n), which will returns n random values. Then, instead of 
accumulating in a loop, you should first create the array of all 
intermediate values, and then summing the whole array (eg you create a 
temporary array tmp where tmp[i] contains the i-th call of f(0.1)):

"""
import numpy as N
from scipy import rand

def f(a):
    return rand(len(a)) + a

tmp = f(0.1 * N.ones(n))
x = N.sum(a)
"""

(I didn't check whether this code is running). Basically, in 
scipy/numpy, you should avoid using loop as much as possible, and try to 
"vectorize" as much as possible. Most numpy/scipy functions can be used like

numfunc(sequence)

instead of

for i in sequence:
    numfunc(i)

Concrete example:

from scipy import rand

def numfunc(n):
    a = rand(n)
    return numpy.sum(a)

def pyfunc(n):
    x = 0
    for i in xrange(n):
        x += rand(1)
    return x

numfunc is already 200 times faster for n = 1e4 on my computer, and the 
difference will grow for bigger size. It requires some thinking if you 
are not used to it. But if you really need looping with around 1e8 
elements, you won't be able to do it efficiently otherwise (that is if 
you stay in python).

Technically, the difference is that in the first case, the loop is done 
in numpy, which is optimized for this kind of things, and the second one 
is done entirely in python, and loops with functions calls are extremely 
slow compared to compiled language (several order of magnitudes most of 
the cases; this is actually not always true, there are some techniques 
to optimize this kind of thing, but this would take use way beyond the 
point of this thread).

Using xrange as Matthieu suggested would solve the memory part: range(n) 
allocates a list of n integers -> 1e8 * 4 bytes minimum if python is 
using 32 bits integers, I don't remember, and whereas xrange(n) creates 
the new iteration value at each iteration.

David