Mailman 3 February 2006 - SciPy-User

savemat not working in SciPy 0.4.6?
by Albert Strasheim Feb. 22, 2006

Feb. 22, 2006

Hello all I'm trying to use the savemat function in SciPy 0.4.6. The following code doesn't work: from numpy import array from scipy.io import savemat savemat('foo.mat', {'x' : array([0])}) The traceback: Traceback (most recent call last): File "frags.py", line 3, in ? savemat('foo.mat', {'x' : array([0])}) File "C:\Python24\Lib\site-packages\scipy\io\mio.py", line 875, in savemat fid.fwrite(variable+'\x00','char') File "C:\Python24\Lib\site-packages\scipy\io\mio.py", line 223, in write numpyio.fwrite(self,count,data,mtype,bs) numpyio.error: Does not support extended types. Is anybody else experiencing this problem? Regards Albert

2 1

signaltools missing imports
by Albert Strasheim Feb. 22, 2006

Feb. 22, 2006

Hello all It seems signaltools.py is missing some imports in SciPy 0.4.6. The following code: from scipy.signal import triang print triang(3) causes the following traceback: Traceback (most recent call last): File "frags.py", line 10, in ? print triang(3) File "C:\Python24\Lib\site-packages\scipy\signal\signaltools.py", line 579, in triang w = numpy.r_[w, w[-2::-1]] NameError: global name 'numpy' is not defined At the top of signaltools.py we have: import numpy as Numeric This seems incorrect. I've filed this issue in Trac: http://projects.scipy.org/scipy/scipy/ticket/10 Regards Albert

1 0

SciPy Install failure: dfftpack not found
by David Treadwell Feb. 22, 2006

Feb. 22, 2006

Help, Please! SciPy installation is giving me fits. After a recent HD crash, I'm trying to reinstall my Python system. Installing using he bash command in OS X 10.4.5 Foo:~/Desktop/scipy-0.4.6 bar$ sudo python setup.py bdist_egg proceeds just fine until it tries to compile the fftpack. The dfftpack can't be found, even tho the directory ./Lib/fftpack/ dfftpack exists and contains the Fortran source files. See ERROR below. The same error occurs for Foo:~/Desktop/scipy-0.4.6 bar$ sudo python setup.py install I've read and followed, I believe, the instructions from the Installation web page. The info for my Python system as it stands is given below. If there's any more diagnostic info needed, I'll be happy to pass it along. I'd appreciate any help in solving this problem. Thanks, --David Treadwell My Python system information: * ActivePython 2.4.2 Build 10 (ActiveState Corp.) based on Python 2.4.2 (#1, Jan 17 2006, 17:03:20) * GCC 3.3 20030304 (Apple Computer, Inc. build 1666) * fftw2.1.5 * g77 installation: Reading specs from /usr/local/lib/gcc/powerpc-apple- darwin7.9.0/3.4.4/specs Configured with: ../gcc/configure --enable-threads=posix -- enable-languages=f77 gcc version 3.4.4 *NumPy installation: <module 'numpy' from '~/Library/Python2.4/site-packages/numpy-0.9.5- py2.4-macosx-10.4-ppc.egg/numpy/__init__.pyc'> ERROR: /usr/bin/ld: can't locate file for: -ldfftpack collect2: ld returned 1 exit status /usr/bin/ld: can't locate file for: -ldfftpack collect2: ld returned 1 exit status error: Command "/usr/bin/g77 -undefined dynamic_lookup -bundle build/temp.darwin-8.5.0-Power_Macintosh-2.4/build/src/Lib/fftpack/ _fftpackmodule.o build/temp.darwin-8.5.0-Power_Macintosh-2.4/Lib/ fftpack/src/zfft.o build/temp.darwin-8.5.0-Power_Macintosh-2.4/Lib/ fftpack/src/drfft.o build/temp.darwin-8.5.0-Power_Macintosh-2.4/Lib/ fftpack/src/zrfft.o build/temp.darwin-8.5.0-Power_Macintosh-2.4/Lib/ fftpack/src/zfftnd.o build/temp.darwin-8.5.0-Power_Macintosh-2.4/ build/src/fortranobject.o -L/usr/local/lib -L/usr/local/lib/gcc/ powerpc-apple-darwin7.9.0/3.4.4 -Lbuild/temp.darwin-8.5.0- Power_Macintosh-2.4 -ldfftpack -lfftw3 -lg2c -lcc_dynamic -o build/ lib.darwin-8.5.0-Power_Macintosh-2.4/scipy/fftpack/_fftpack.so" failed with exit status 1

2 2

Help on special.hankel1
by Nils Wagner Feb. 22, 2006

Feb. 22, 2006

Hi all, The output of help (special.hankel1) is not very informative. Help on ufunc: hankel1 = <ufunc 'hankel1'> I need the Hankel function of the first kind and zero order H_0^{(1)} (arg), where arg is a complex number. How can I compute this Hankel function with scipy ? Nils

2 1

Re: [SciPy-user] [Numpy-discussion] algorithm, optimization, or other problem?
by Nadav Horesh Feb. 22, 2006

Feb. 22, 2006

You may get a significant boost by replacing the line: w=w+ eta * (y*x - y**2*w) with w *= 1.0 - eta*y*y w += eta*y*x I ran a test on a similar expression and got 5 fold speed increase. The dot() function runs faster if you compile with dotblas. Nadav. -----Original Message----- From: numpy-discussion-admin(a)lists.sourceforge.net on behalf of Bruce Southey Sent: Tue 21-Feb-06 17:15 To: Brian Blais Cc: python-list(a)python.org; numpy-discussion(a)lists.sourceforge.net; scipy-user(a)scipy.net Subject: Re: [Numpy-discussion] algorithm, optimization, or other problem? Hi, In the current version, note that Y is scalar so replace the squaring (Y**2) with Y*Y as you do in the dohebb function. On my system without blas etc removing the squaring removes a few seconds (16.28 to 12.4). It did not seem to help factorizing Y. Also, eta and tau are constants so define them only once as scalars outside the loops and do the division outside the loop. It only saves about 0.2 seconds but these add up. The inner loop probably can be vectorized because it is just vector operations on a matrix. You are just computing over the ith dimension of X. I think that you could be able to find the matrix version on the net. Regards Bruce On 2/21/06, Brian Blais <bblais(a)bryant.edu> wrote: > Hello, > > I am trying to translate some Matlab/mex code to Python, for doing neural > simulations. This application is definitely computing-time limited, and I need to > optimize at least one inner loop of the code, or perhaps even rethink the algorithm. > The procedure is very simple, after initializing any variables: > > 1) select a random input vector, which I will call "x". right now I have it as an > array, and I choose columns from that array randomly. in other cases, I may need to > take an image, select a patch, and then make that a column vector. > > 2) calculate an output value, which is the dot product of the "x" and a weight > vector, "w", so > > y=dot(x,w) > > 3) modify the weight vector based on a matrix equation, like: > > w=w+ eta * (y*x - y**2*w) > ^ > | > +---- learning rate constant > > 4) repeat steps 1-3 many times > > I've organized it like: > > for e in 100: # outer loop > for i in 1000: # inner loop > (steps 1-3) > > display things. > > so that the bulk of the computation is in the inner loop, and is amenable to > converting to a faster language. This is my issue: > > straight python, in the example posted below for 250000 inner-loop steps, takes 20 > seconds for each outer-loop step. I tried Pyrex, which should work very fast on such > a problem, takes about 8.5 seconds per outer-loop step. The same code as a C-mex > file in matlab takes 1.5 seconds per outer-loop step. > > Given the huge difference between the Pyrex and the Mex, I feel that there is > something I am doing wrong, because the C-code for both should run comparably. > Perhaps the approach is wrong? I'm willing to take any suggestions! I don't mind > coding some in C, but the Python API seemed a bit challenging to me. > > One note: I am using the Numeric package, not numpy, only because I want to be able > to use the Enthought version for Windows. I develop on Linux, and haven't had a > chance to see if I can compile numpy using the Enthought Python for Windows. > > If there is anything else anyone needs to know, I'll post it. I put the main script, > and a dohebb.pyx code below. > > > thanks! > > Brian Blais > > -- > ----------------- > > bblais(a)bryant.edu > http://web.bryant.edu/~bblais > > > > > # Main script: > > from dohebb import * > import pylab as p > from Numeric import * > from RandomArray import * > import time > > x=random((100,1000)) # 1000 input vectors > > numpats=x.shape[0] > w=random((numpats,1)); > > th=random((1,1)) > > params={} > params['eta']=0.001; > params['tau']=100.0; > old_mx=0; > for e in range(100): > > rnd=randint(0,numpats,250000) > t1=time.time() > if 0: # straight python > for i in range(len(rnd)): > pat=rnd[i] > xx=reshape(x[:,pat],(1,-1)) > y=matrixmultiply(xx,w) > w=w+params['eta']*(y*transpose(xx)-y**2*w); > th=th+(1.0/params['tau'])*(y**2-th); > else: # pyrex > dohebb(params,w,th,x,rnd) > print time.time()-t1 > > > p.plot(w,'o-') > p.xlabel('weights') > p.show() > > > #============================================= > > # dohebb.pyx > > cdef extern from "Numeric/arrayobject.h": > > struct PyArray_Descr: > int type_num, elsize > char type > > ctypedef class Numeric.ArrayType [object PyArrayObject]: > cdef char *data > cdef int nd > cdef int *dimensions, *strides > cdef object base > cdef PyArray_Descr *descr > cdef int flags > > > def dohebb(params,ArrayType w,ArrayType th,ArrayType X,ArrayType rnd): > > > cdef int num_iterations > cdef int num_inputs > cdef int offset > cdef double *wp,*xp,*thp > cdef int *rndp > cdef double eta,tau > > eta=params['eta'] # learning rate > tau=params['tau'] # used for variance estimate > > cdef double y > num_iterations=rnd.dimensions[0] > num_inputs=w.dimensions[0] > > # get the pointers > wp=<double *>w.data > xp=<double *>X.data > rndp=<int *>rnd.data > thp=<double *>th.data > > for it from 0 <= it < num_iterations: > > offset=rndp[it]*num_inputs > > # calculate the output > y=0.0 > for i from 0 <= i < num_inputs: > y=y+wp[i]*xp[i+offset] > > # change in the weights > for i from 0 <= i < num_inputs: > wp[i]=wp[i]+eta*(y*xp[i+offset] - y*y*wp[i]) > > # estimate the variance > thp[0]=thp[0]+(1.0/tau)*(y**2-thp[0]) > > > > > > > > ------------------------------------------------------- > This SF.net email is sponsored by: Splunk Inc. Do you grep through log files > for problems? Stop! Download the new AJAX search engine that makes > searching your log files as easy as surfing the web. DOWNLOAD SPLUNK! > http://sel.as-us.falkag.net/sel?cmd=lnk&kid=103432&bid=230486&dat=121642 > _______________________________________________ > Numpy-discussion mailing list > Numpy-discussion(a)lists.sourceforge.net > https://lists.sourceforge.net/lists/listinfo/numpy-discussion > ------------------------------------------------------- This SF.net email is sponsored by: Splunk Inc. Do you grep through log files for problems? Stop! Download the new AJAX search engine that makes searching your log files as easy as surfing the web. DOWNLOAD SPLUNK! http://sel.as-us.falkag.net/sel?cmd=k&kid3432&bid#0486&dat1642 _______________________________________________ Numpy-discussion mailing list Numpy-discussion(a)lists.sourceforge.net https://lists.sourceforge.net/lists/listinfo/numpy-discussion

1 0

Re: [SciPy-user] low-discrepancy sequences
by DE MENTEN Sebastien Feb. 22, 2006

Feb. 22, 2006

> > DE MENTEN Sebastien wrote: > > Hi all, > > > > Has anyone worked on low-discrepancy sequences in Scipy? > > Which sequences would you like to see? Halton? Hammersley? Sobol? > Niederreiter-Xing? > > I have a small Halton sequence generator, now, since it was so easy to > implement. The others may take some more time. Apparently, there are a couple of algorithms in fortran already written for this at http://www.csit.fsu.edu/~burkardt/f_src/f_src.html and also in numerical recipes I think. What could be the step for integration in Scipy? Best, Seb

2 2

low-discrepancy sequences
by DE MENTEN Sebastien Feb. 22, 2006

Feb. 22, 2006

Hi all, Has anyone worked on low-discrepancy sequences in Scipy? Sebastien

2 2

inverse of column_stack
by Ryan Krauss Feb. 21, 2006

Feb. 21, 2006

Is there an inverse function to column_stack? i.e. I want to iterate over the columns of a matrix. I could first transpose it and then call tolist, but I was wondering if there was some clean, fast, built in way that I wasn't aware of. Thanks, Ryan

4 7

algorithm, optimization, or other problem?
by Brian Blais Feb. 21, 2006

Feb. 21, 2006

Hello, I am trying to translate some Matlab/mex code to Python, for doing neural simulations. This application is definitely computing-time limited, and I need to optimize at least one inner loop of the code, or perhaps even rethink the algorithm. The procedure is very simple, after initializing any variables: 1) select a random input vector, which I will call "x". right now I have it as an array, and I choose columns from that array randomly. in other cases, I may need to take an image, select a patch, and then make that a column vector. 2) calculate an output value, which is the dot product of the "x" and a weight vector, "w", so y=dot(x,w) 3) modify the weight vector based on a matrix equation, like: w=w+ eta * (y*x - y**2*w) ^ | +---- learning rate constant 4) repeat steps 1-3 many times I've organized it like: for e in 100: # outer loop for i in 1000: # inner loop (steps 1-3) display things. so that the bulk of the computation is in the inner loop, and is amenable to converting to a faster language. This is my issue: straight python, in the example posted below for 250000 inner-loop steps, takes 20 seconds for each outer-loop step. I tried Pyrex, which should work very fast on such a problem, takes about 8.5 seconds per outer-loop step. The same code as a C-mex file in matlab takes 1.5 seconds per outer-loop step. Given the huge difference between the Pyrex and the Mex, I feel that there is something I am doing wrong, because the C-code for both should run comparably. Perhaps the approach is wrong? I'm willing to take any suggestions! I don't mind coding some in C, but the Python API seemed a bit challenging to me. One note: I am using the Numeric package, not numpy, only because I want to be able to use the Enthought version for Windows. I develop on Linux, and haven't had a chance to see if I can compile numpy using the Enthought Python for Windows. If there is anything else anyone needs to know, I'll post it. I put the main script, and a dohebb.pyx code below. thanks! Brian Blais -- ----------------- bblais(a)bryant.edu http://web.bryant.edu/~bblais # Main script: from dohebb import * import pylab as p from Numeric import * from RandomArray import * import time x=random((100,1000)) # 1000 input vectors numpats=x.shape[0] w=random((numpats,1)); th=random((1,1)) params={} params['eta']=0.001; params['tau']=100.0; old_mx=0; for e in range(100): rnd=randint(0,numpats,250000) t1=time.time() if 0: # straight python for i in range(len(rnd)): pat=rnd[i] xx=reshape(x[:,pat],(1,-1)) y=matrixmultiply(xx,w) w=w+params['eta']*(y*transpose(xx)-y**2*w); th=th+(1.0/params['tau'])*(y**2-th); else: # pyrex dohebb(params,w,th,x,rnd) print time.time()-t1 p.plot(w,'o-') p.xlabel('weights') p.show() #============================================= # dohebb.pyx cdef extern from "Numeric/arrayobject.h": struct PyArray_Descr: int type_num, elsize char type ctypedef class Numeric.ArrayType [object PyArrayObject]: cdef char *data cdef int nd cdef int *dimensions, *strides cdef object base cdef PyArray_Descr *descr cdef int flags def dohebb(params,ArrayType w,ArrayType th,ArrayType X,ArrayType rnd): cdef int num_iterations cdef int num_inputs cdef int offset cdef double *wp,*xp,*thp cdef int *rndp cdef double eta,tau eta=params['eta'] # learning rate tau=params['tau'] # used for variance estimate cdef double y num_iterations=rnd.dimensions[0] num_inputs=w.dimensions[0] # get the pointers wp=<double *>w.data xp=<double *>X.data rndp=<int *>rnd.data thp=<double *>th.data for it from 0 <= it < num_iterations: offset=rndp[it]*num_inputs # calculate the output y=0.0 for i from 0 <= i < num_inputs: y=y+wp[i]*xp[i+offset] # change in the weights for i from 0 <= i < num_inputs: wp[i]=wp[i]+eta*(y*xp[i+offset] - y*y*wp[i]) # estimate the variance thp[0]=thp[0]+(1.0/tau)*(y**2-thp[0])

5 4

Interaction of sparse and dense matrices in scipy
by Nils Wagner Feb. 21, 2006

Feb. 21, 2006

Hi all, I need some advice on operations between dense and sparse matrix "objects". Which operations are possible in scipy ? Assume that A and B are sparse while C and D are dense matrices of appropriate dimension. dense sparse dot(C,D) A*B linalg.norm(C) linalg.norm(C.todense()) linalg.kron(C,D) ? identity(n), eye(n) How about interactions of sparse and dense matrices - multiplication of a sparse matrix A with a dense matrix C - Kronecker product of sparse matrices ? - Kronecker product of a dense and a sparse matrix ? - sparse matrix functions expm(A),... ? - sparse eigensolver - to be continued Can I use different sparse formats like csc and csr w.r.t. to multiplication, addition,... And what is the output of A = io.mmread("A.mtx") ? I mean sparse or dense. Any comments ? Nils

4 9