Mailman 3 July 2007 - SciPy-Dev

gcc bzonfig.h
by John Hunter July 16, 2007

July 16, 2007

On our x86 solaris platform with gcc johnh@flag:svn> gcc --version gcc (GCC) 3.4.1 johnh@flag:svn> uname -a SunOS flag 5.10 Generic_118855-15 i86pc i386 i86pc we are getting an error with blitz inline (eg in Prabhu's laplace.py example) because we are getting the gnu/bzconfig.h macro defined /* define if the compiler has isnan function in namespace std */ #ifndef BZ_ISNAN_IN_NAMESPACE_STD #define BZ_ISNAN_IN_NAMESPACE_STD #endif and then a compiler error that isnan is not in the std::namespace. We tried to "undec_macros" this variable, but that was no help, so we patched the file with the scipy svn patch below. This shouldn't adversely affect anyone, because the default behavior will be the same, but will allow us to compile weave blitz code with the define_macros definte_macros=[('NO_BZ_ISNAN_IN_NAMESPACE_STD',1)] Although the double negative syntax is a bit irksome, I didn't know a better way to turn off this macro. If this, or something like it, could be included in scipy, that would be great. Thanks. JDH Index: Lib/weave/blitz/blitz/gnu/bzconfig.h =================================================================== --- Lib/weave/blitz/blitz/gnu/bzconfig.h (revision 3166) +++ Lib/weave/blitz/blitz/gnu/bzconfig.h (working copy) @@ -298,9 +298,11 @@ #endif /* define if the compiler has isnan function in namespace std */ +#ifndef NO_BZ_ISNAN_IN_NAMESPACE_STD #ifndef BZ_ISNAN_IN_NAMESPACE_STD #define BZ_ISNAN_IN_NAMESPACE_STD #endif +#endif /* define if the compiler has C math functions in namespace std */ #ifndef BZ_MATH_FN_IN_NAMESPACE_STD

1 0

ANN: SciPy 2007 Conference Updates
by Travis Vaught July 15, 2007

July 15, 2007

Greetings, We're excited to have *Ivan Krstić*, the director of security architecture for the One Laptop Per Child project as our Keynote Speaker this year. The planning for the SciPy 2007 Conference is moving along. Please see below for some important updates. Schedule Available ------------------ The full schedule of talks has been posted here: http://www.scipy.org/SciPy2007/ConferenceSchedule Early Registration Extended --------------------------- If you haven't yet registered for the conference, the early registration deadline has been extended to Wednesday, July 18th, 2007. For more information on the conference see: http://www.scipy.org/SciPy2007 Student Sponsorship ------------------- Enthought, Inc. (http://www.enthought.com) is sponsoring the registration fees for up to 5 college or graduate students to attend the conference. To apply, please send a short description of what you are studying and why you'd like to attend to info(a)enthought.com. Please include telephone contact information. BOFs & Sprints -------------- If you're planning to attend and are interested in selecting BOF or Sprint Session topics, please weigh in at: BOFs: http://www.scipy.org/SciPy2007/BoFs Sprints: http://www.scipy.org/SciPy2007/Sprints We're looking forward to a great conference this year! Best, Travis

1 0

sparse comments and questions
by Nathan Bell July 14, 2007

July 14, 2007

Currently the sparse matrix storage formats (e.g. sparse.csr_matrix ) only support integer indices ( numpy.intc ) and floating point data (float32, float64, complex64, and complex128). While this addresses the common usage, it would be nice to support other index and data formats also. For example, when using a csr_matrix to represent a graph the dtype is mostly irrelevant, so smaller types like uint8 would be desirable. Likewise, on a 64-bit machine, one might like to use 32-bit indices (instead of the default C-int size) to economize on space. The sparsetools routines are fully templated, so it's a simple matter to support the other types there. However, I don't know whether other subpackages that use sparse matrices would require changes. Any thoughts? I have a few other questions and comments regarding the sparse package. 1) Does anyone depend on the .ftype attribute of the sparse formats? Can it be eliminated? If not, can it be trapped by __getattr__() so that it doesn't become unsynchronized from data.dtype? 2) When operating on multiple types, what's the best way to get the correct upcast type? For example int32 + float32 -> float64. I know this can be done experimentally ( http://www.scipy.org/Tentative_NumPy_Tutorial/Discussion ). Is there a better way? 3) Should we allow modifications to CSR and CSC matrices that force O(N) updates/reallocations. For instance, adding a new value into a CSR matrix is essentially as expensive as building a new one from scratch. Could we instead prohibit such modifications and raise an exception informing the user that lil_matrix is the proper format for such operations? Note that changing an existing matrix value is not so costly (typically O(1)), so not all calls to __setitem__ would be banned. 4) lil_matrix now supports more sophisticated slicing, but the implementation is pretty hairy and not especially fast. Is anyone looking at this? Does the slicing work correctly? I fixed a bug related to negative indices, but I'm not completely confident in the robustness of my fix. 5) Scipy's spdiags() works differently than the MATLAB version. Should it be changed to coincide? If we later support a sparse diagonal format, should the output of spdiags change to dia_matrix instead of the present csc_matrix? 6) When working with sparse matrices, do people expect the result to be of the same type? For example, it's far more efficient to make transpose(csr_matrix) output a csc_matrix (as is currently done in SciPy) since that's a O(1) operation (as opposed to an O(n) cost otherwise). So in general, should a method to return the result in most efficient output format, relying on the user to convert it to another format if so desired? I.e. if you need a csr_matrix at a given point, then you should call .tocsr() to be sure. If it's already a csr_matix, then no conversion will be done. 7) Where do functions that operate on sparse matrices belong? For example, suppose I wrote a function to extract the lower diagonal entries of a matrix. Should it be added to sparse.py (like spdiags) or somewhere else? I'd like to add some functionality to allow the extraction of submatrices from CSR/CSC matrices. This could be used when removing the DOFs corresponding to boundary conditions from a bilinear form. Does anyone have a suggestion for an interface for this functionality? It's possible to implement this via sparse matrix multiplication, albeit with some additional overhead. -- Nathan Bell wnbell(a)gmail.com

4 8

proposed patches: reviewers welcome
by saintmlx July 12, 2007

July 12, 2007

Hello SciPy Developers, I recently submitted patches that, I hope, will help closing two opened tickets: the conversion of scalars using int() (#549, opened by myself), and... the notorious "agitation" around the "standard" standard deviation (#502 et al.). Any reviews would be welcome. If these could make it to the svn repository, I'd be more than happy :) I've attached diff files to both tickets. Thanks for any help, Xavier Saint-Mleux

1 0

Sorted Eigenvalues/vectors
by Jeremy Conlin July 10, 2007

July 10, 2007

I have a need to have the eigenvalues/vectors calculated in my code returned sorted by eigenvalue. Of course this necessitates sorting the eigenvectors accordingly. If I understand correctly, a method for this doesn't exist in the current distribution so I wrote a simple one myself. Is this something that the community is interested in? Should a method like this be included in the scipy distribution? It is small, but I think very useful. Thanks, Jeremy Conlin

2 3

scikits infrastructure: ticket Mailing list, list of projects
by David Cournapeau July 10, 2007

July 10, 2007

Hi, Since several projects are now under the scikits umbrella (mlabwrap, learn, optimization toolbox and audio IO tools), I would like to know if some infrastructure could be set. In particular, having tickets/commit Mailing lists would be useful. cheers, David

4 5

Current SVN seems to be broken
by Matthieu Brucher July 9, 2007

July 9, 2007

Hi, I just updated scipy to the HEAD and when building, I got this : Lib/sparse/sparsetools/sparsetools_wrap.cxx: In function 'int SWIG_Python_ConvertPtr(PyObject*, void**, swig_type_info*, int)': Lib/sparse/sparsetools/sparsetools_wrap.cxx:1231: error: invalid conversion from 'const char*' to 'char*' Lib/sparse/sparsetools/sparsetools_wrap.cxx: In function 'void SWIG_Python_FixMethods(PyMethodDef*, swig_const_info*, swig_type_info**, swig_type_info**)': Lib/sparse/sparsetools/sparsetools_wrap.cxx:20206: error: invalid conversion from 'const char*' to 'char*' Lib/sparse/sparsetools/sparsetools_wrap.cxx: In function 'int SWIG_Python_ConvertPtr(PyObject*, void**, swig_type_info*, int)': Lib/sparse/sparsetools/sparsetools_wrap.cxx:1231: error: invalid conversion from 'const char*' to 'char*' Lib/sparse/sparsetools/sparsetools_wrap.cxx: In function 'void SWIG_Python_FixMethods(PyMethodDef*, swig_const_info*, swig_type_info**, swig_type_info**)': Lib/sparse/sparsetools/sparsetools_wrap.cxx:20206: error: invalid conversion from 'const char*' to 'char*' error: Command "g++ -pthread -fno-strict-aliasing -DNDEBUG -g -O3 -Wall -fPIC -ILib/sparse/sparsetools -I/home/brucher/local/lib/python2.5/site-packages/numpy/core/include -I/home/brucher/local/include/python2.5 -c Lib/sparse/sparsetools/sparsetools_wrap.cxx -o build/temp.linux-i686-2.5/Lib/sparse/sparsetools/sparsetools_wrap.o" failed with exit status 1 With gcc 4.1.1 on FC6 This file was updated since my last checkout, can it be repared ? Matthieu

2 2

restructuredtext en docformat
by dmitrey July 9, 2007

July 9, 2007

I noticed that http://twistedmatrix.com/pipermail/twisted-bugs/2006-October/001609.html says: epydoc has a epytext variant that's valid restructured text. You can either set `__docformat__ = "restructuredtext en"` in your module, or pass `--docformat restructuredtext` to `epydoc` to use it. So, maybe we can omit putting the line (just do latter proposition), + taking into account that in future more appropriate format can appear and then we should make all those changes in all source files? D.

2 1

GSoC weekly report
by dmitrey July 7, 2007

July 7, 2007

Hi all, today's topics: 1. Automatic gradient checks had been implemented 2. f,c,h patterns: need more time & efforts (I describe the problems here according to Alan Isaac's request) 3. Letters from ALGENCAN developers. 4. My optimization department staff have explained me their own opinion about OpenOpt renaming, that is insisted by Jarrod Millman and Alan Isaac (and maybe some other scipy developers?). 1. This week automatic gradient checks for f->min, c[i](x)<=0, h[j](x)=0 had been implemented. I guess it's done in much more convenient way that MATLAB fmincon do, isn't it? See http://openopt.blogspot.com/2007/07/automatic-gradient-check-ready.html for more details, and/or see the example below: ... OpenOpt checks user-supplied gradient df (size: (15,)) according to: prob.diffInt = 1e-07# step for numerical gradient obtaining prob.check.maxViolation = 1e-07# lines where difference is less than the number will not be shown, default 1e-5 df num user-supplied numerical difference 0 +6.000e+00 +6.000e+00 -2.122e-07 2 -1.666e+01 -1.666e+01 +1.048e-06 3 -2.584e+01 -2.584e+01 +1.400e-06 4 -2.046e+01 -2.046e+01 +9.434e-07 5 -5.461e+00 -5.461e+00 +1.512e-06 6 +5.363e+00 +5.363e+00 +4.905e-07 9 -2.458e+01 -2.458e+01 +5.861e-07 10 -2.343e+01 -2.343e+01 +1.162e-06 11 -9.929e+00 -9.929e+00 +1.588e-07 12 +3.502e+00 +3.502e+00 +7.354e-07 14 -7.812e+00 -7.812e+00 -8.602e-07 max(abs(df_user - df_numerical)) = 1.51219384126e-06 (is registered in df number 5) sum(abs(df_user - df_numerical)) = 9.37081698371e-06 ======================== OpenOpt checks user-supplied gradient dc (size: (15, 2)) according to: prob.diffInt = 1e-07 prob.check.maxViolation = 1e-07 dc num i,j:dc[i]/dx[j] user-supplied numerical difference 0 0 / 0 +4.096e+03 +4.096e+03 -4.787e-05 2 1 / 0 +0.000e+00 -3.638e-05 +3.638e-05 3 1 / 1 +8.645e+00 +8.645e+00 -1.301e-07 4 2 / 0 +0.000e+00 -3.638e-05 +3.638e-05 5 2 / 1 -6.658e+00 -6.658e+00 -1.350e-07 6 3 / 0 +0.000e+00 -3.638e-05 +3.638e-05 8 4 / 0 +0.000e+00 -3.638e-05 +3.638e-05 10 5 / 0 +0.000e+00 -3.638e-05 +3.638e-05 12 6 / 0 +0.000e+00 -3.638e-05 +3.638e-05 14 7 / 0 +0.000e+00 -3.638e-05 +3.638e-05 16 8 / 0 +0.000e+00 -3.638e-05 +3.638e-05 18 9 / 0 +0.000e+00 -3.638e-05 +3.638e-05 20 10 / 0 +0.000e+00 -3.638e-05 +3.638e-05 22 11 / 0 +0.000e+00 -3.638e-05 +3.638e-05 24 12 / 0 +0.000e+00 -3.638e-05 +3.638e-05 26 13 / 0 +0.000e+00 -3.638e-05 +3.638e-05 28 14 / 0 +0.000e+00 -3.638e-05 +3.638e-05 max(abs(dc_user - dc_numerical)) = 4.78662564092e-05 (is registered in dc number 0) sum(abs(dc_user - dc_numerical)) = 0.000557448428236 ======================== OpenOpt checks user-supplied gradient dh (size: (15, 2)) according to: prob.diffInt = 1e-07 prob.check.maxViolation = 1e-07 dh num i,j:dh[i]/dx[j] user-supplied numerical difference 27 13 / 1 +7.642e+02 +7.642e+02 -2.108e-05 28 14 / 0 +3.312e+03 +3.312e+03 -5.292e-03 max(abs(dh_user - dh_numerical)) = 0.00529207172031 (is registered in dh number 28) sum(abs(dh_user - dh_numerical)) = 0.00531315227033 ======================== (to see the messages you need to turn prob.check.df=1, prob.check.dh=1, prob.check.dc=1) 2. Patterns: I should decide which way of implementing those ones is the best. For example, consider the problem n,m,s = 10000,1000,1000 ((x+15)^2).sum() -> min# size(x)=n acoording to constraints x1^2 + x2^2 + ... + xm^2 <=c1 x2^2 + x3^2 + ... + x[m+1]^2 <=c2 ... xs^2+x[s+1]^2+... + x[m+s]^2 <= c[s] so cPattern will look like 1 1 1...1 0 0 0 0 0 1 1 ..1 1 0 0 0 0 0 1 ..1 1 1 0 0 ... 0 0 0 ..0 0 0 0 0 0 0 0 ..0 0 0 0 0 if I will store dc as 1) dense matrix 2) sparse matrix it will require 1) O(n x s) ~= 10000*1000*sizeof(float) memory bytes 2) O(m x s) = 1000*1000*(sizeof(float) + 2*sizeof(int)) memory bytes However, lots of solvers require just dL/dx = df/dx + dc/dx + dh/dx = grad(f) + sum(dot(mu,c(x))) + sum(dot(lambda,h(x))) so, if dc/dx will be obtained sequentially as dc[i]/dx or dc/dx[i], we can just summarize these ones step by step and amount of memory needed will be O(n) , approximately sizeof(float)*n However, unfortunately not all solvers need just dL/dx, or L is not just f + <mu,c> + <lambda, h>, like for example in Augmented Lagrangian - related solvers. On the other hand, seems like it should always work with fPattern, because f = Fwhole is always just mere sum of f[j] One more problem is how to implement data storing. The most obvious approach is of course double cycle for i in xrange(n) for j in xrange(len(c))#number of c inequalities dc[indexes related to the vector (and its' length) returned from c[j](x) evaluation] = (c[j](x)-c0[j])/p.diffInt ('cause c[j] can return several numbers, if we consider the most general case) However, double cycle is too slow for Python. I intend to write single cycle, based on for i in xrange(n) j = where(cPattern[i]) C = <somehow evaluate all c[j](x)> dc[some ind] = (C[some ind] - C0[some ind])/p.diffInt here some ind could be something like [1,2,3] + [4,5] + [8,9,10] + ... + [100, 101] - values obtained from those j: cPattern[j]=1 of couse in code I will divide dcmatrix only once, using ufunc /, when whole matrix dc is obtained However, some problems still exist. 3.This week I've received (and answered to) some letters from ALGENCAN developers, they are interested in connecting their solvers to Python and OO. Their solvers often works better than IPOPT, some results are attached in their articles http://www.ime.usp.br/~egbirgin/tango/publications.php (Their software is Augmented Lagrangian - based) + They informed me that they had changed lots of their solvers licenses from "free for noncommercial usage" (that are not OSI-approved) to GPL (some their software remain commercial). 4. So, I have informed my department about your intention to change OpenOpt name. They have answered me: Dmitrey, we had agreed to make OpenOpt our department's vendor, like TomOpt, PENOPT, IPOPT, SolvOpt, CVXOPT are, as well as construct (in future) a website like openopt.net and create a banner like "member of OpenOpt project" for link exchange and/or our partners. So you propose to turn it into something like "member of scikits.optimize project"? So it's your own decision, but please, don't annoy with your questions anymore - we have our own work, and we do not see any benefits of spending our time with your one, moreover, for free. Please, just finish your education as quickly as it can be and go away - we need people working in our projects. Particularly, it means that I will not be able to implement the QP/QPQC solver that I told about (this one is needed as default lincher QP subproblem solver with BSD license). Also, unfortunately it means problems with consulting about global GRASP-based solver (and implementing this one into scikits.optimization), because I still misunderstood some dr. Shilo's algorithm steps. In brief that's all for the week. BTW it's last one before 1st GSoC milestone (July 9th) and I had done all mentioned in http://projects.scipy.org/scipy/scipy/wiki/OptimizationProblems except patterns. Regards, D.

5 9

PEP 368: Standard image protocol and class
by Lino Mastrodomenico July 6, 2007

July 6, 2007

[Sorry for the cross-posting, but I think this may be relevant for both NumPy and ndimage.] Hello everyone, I have submitted to the Python core developers a new PEP (Python Enhancement Proposal): http://www.python.org/dev/peps/pep-0368/ It proposes two things: * the creation of a standard image protocol/interface that can be hopefully implemented interoperably by most Python libraries that manipulate images; * the addition to the Python standard library of a basic implementation of the new protocol. The new image protocol is heavily inspired by a subset of the NumPy array interface, with a few image-specific additions and changes (e.g. the "size" attribute of an image is a tuple (width, height)). Of course it would be wonderful if these new image objects could interoperate out-of-the-box with numpy arrays and ndimage functions. There is another proposal that would be very useful for that, PEP 3118 by Travis Oliphant and Carl Banks: http://www.python.org/dev/peps/pep-3118/ The image PEP (368) currently lists only modes based on uint8/16/32 numbers, but the final version will probably also include modes based on float32 and float16 (converted in software to/from float32/64 when necessary). A discussion about it is currently going on in the python-3000 mailing list: <http://mail.python.org/pipermail/python-3000/2007-July/thread.html#8648> Any suggestion, comment or criticism from the NumPy/SciPy people would be very useful, but IMHO keeping the discussion only on the python-3000 ML may be a good idea, to avoid duplicating answers on different mailing lists. Thanks in advance. -- Lino Mastrodomenico E-mail: l.mastrodomenico(a)gmail.com

3 3