[Numpy-discussion] when does numpy create temporaries ?

[Arnd Baecker]

On Tue, 30 May 2006, Simon Burton wrote:

> Consider these two operations:
>
> >>> a=numpy.empty( 1024 )
> >>> b=numpy.empty( 1024 )
> >>> a[1:] = b[:-1]
> >>> a[1:] = a[:-1]
> >>>
>
> It seems that in the second operation
> we need to copy the view a[:-1]
> but in the first operation we don't need
> to copy b[:-1].
>
> How does numpy detect this, or does it
> always copy the source when assigning to a slice ?
>
> I've poked around the (numpy) code a bit and
> tried some benchmarks, but it's still not so
> clear to me.

Hi,

not being able to give an answer to this question,
I would like to emphasize that this can be a very important issue:
Firstly, I don't know how to monitor the memory usage
*during* the execution of a line of code.
(Putting numpy.testing.memusage() before and after
that line does not help, if I remember things correctly).

With Numeric I ran into a memory problem with the code appended below.
It turned out, that internally a copy had been made
which for huge arrays brought my system into swapping.
(numpy behaves the same as Numeric. Moreover,  it
seems to consume around 8.5 MB more memory than Numeric?!)

So I strongly agree that it would be nice to know in advance
when temporaries are created.
In addition it would be good to be able to debug memory allocation.
(For example, with f2py there is the option
-DF2PY_REPORT_ON_ARRAY_COPY=1
Note that this only works when generating the wrapper library,
i.e., there is no switch to turn this on or off afterwards,
at least as far as I know).

Best, Arnd



##########################################################
from Numeric import *
#from numpy import *
import os

pid=os.getpid()
print "Process id: ",pid

N=200
NB=30             # number of wavefunctions
NT=20             # time steps
print "Expected size of `wfk`      (in KB):", N*N*NB*8/1024.0
print "Expected size of `time_arr` (in KB):", N*N*NT*16/1024.0

wfk=zeros( (N,N,NB),Float)
phase=ones(NB,Complex)
time_arr=zeros( (N,N,NT),Complex)

print "press enter and watch the memory"
raw_input("(eg. with pmap %d | grep total)" % (pid)  )

# this one does a full copy of wfk, because it is complex !!!
#while 1:
#    for tn in range(NT):
#        time_arr[:,:,tn]+=dot(wfk, phase)
#
#
# memory usage: varies around:
# - 38524K/57276K with Numeric
# - 46980K/66360K with numpy

while 1:
    for tn in range(NT):
        time_arr[:,:,tn]+=dot(wfk, phase.real)+1j*dot(wfk, phase.imag)
#
# memory usage: varies around:
# -  38524K/40412K with Numeric
# -  46984K/47616K with numpy

################################################################