Mailman 3 April 2008 - SciPy-User

Help required regding Eigenvectors
by Sudheer Phani 14 Mar '22

14 Mar '22

Hello All, I have been encountering some strange problems with the calculation of eigen vectors using python. The problem is with the sign of the eigen vectors calculated. On Very few occasions the eigen vectors calculated using python have an opposite sign than the actual eigenvectors, though the magnitude is correct. Can any one of you plz help me in this regard and suggest me some solution if you have encoutered similar problems . Any help in this regd is appreciated Thanks in Advace regds Sudheer

5 5

3D plotting
by Robert Cimrman 10 Sep '09

10 Sep '09

Hi, I need to do some simple 3D plots, so I tried the two Cookbook recipes (xplt, mplot3d) and failed: * xplt: does not work when using a custom installation directory such as python setup.py install --root=/home/share/software apparently the gistdata directory path is wrong: gist.error: failed to create drawing -- Gist work.gs style sheet missing * mplot3d: does not work with my version of matplotlib (0.80). I have made the changes mentioned in the Cookbook to no avail. (Axes.__init__() args apparently changed, as well as some other matplotlib object attributes) Any ideas? mplot3d looks great, I would really like to use it! thanks, r.

6 6

assign to diagonal values?
by Carlos Scheidegger 07 Jul '09

07 Jul '09

Say I have a matrix m whose diagonal values I want to replace with the contents of vector new_diagonal. Right now, I'm using m -= diag(m) m += diag(new_diagonal) but that seems pretty wasteful (especially if diag() creates new arrays, which I suspect it does). Is there a standard (or any efficient) way to assign the diagonal values directly without an explicit for loop? Googling for "site:www.scipy.org assign diagonal" only returned results for sparse matrices, which is not the case here - these are plain numpy arrays. Thanks in advance, -carlos

11 17

Is numpy's argsort lying about its numpy.int32 types?
by Rob Clewley 15 May '09

15 May '09

Hi, I'm having a problem comparing some types when using numpy's argsort. I'm using numpy 1.0.2. I can reproduce it simply: Python 2.4.4 (#71, Oct 18 2006, 08:34:43) [MSC v.1310 32 bit (Intel)] << snip >> In [1]: from numpy import argsort, sort, int32, array In [2]: x=array([1,3,2]) In [3]: aa=argsort(x) In [4]: as=sort(x) In [5]: type(aa[0]) Out[5]: <type 'numpy.int32'> In [6]: type(as[0]) Out[6]: <type 'numpy.int32'> In [7]: int32 Out[7]: <type 'numpy.int32'> In [8]: type(as[0])==int32 Out[8]: True In [9]: type(aa[0])==int32 Out[9]: False Any of the three indices in aa give me the same problem. Can someone explain if I should be doing this a different way, and if this is a bug?! Thanks, Rob

5 6

The IO library and image file formats -- compare with with PIL
by Sebastian Haase 16 Jun '08

16 Jun '08

Hi, Could someone here please summarize the plans for the io library !? I saw that there are parts of it going into numpy already. And more stuff is already in SciPy -- IIRC. What is the status / plans regrading image formats like TIFF. Are you guys planning to duplicate the efforts of the Python Imaging Library ( PIL ) ? Or can you just copy-paste some of the code ? The reason I'm asking now, is that I just submitted a patch to PIL that adds writing capabilities for MultiPage TIFF files. I did not get any response on the PIL list, and I'm actually still not clear on their handling of free vs. commercial. The commercial PIL license sounds like you have to pay $2000 to get access to svn, i.e. the current development version of PIL. Anyone here who could comment on this ? Thanks, Sebastian Haase

7 26

scipy.stats rv objects from data
by Erik Tollerud 02 May '08

02 May '08

I'm finding the scipy.stats documentation somewhat difficult to follow, so maybe the answer to this question is in there... I can't really find it, though. What I have is a sequence of numbers X_i . Two things I'd like to be able to do with this: 1. Create a discrete probability distribution (class rv_discrete) from this data so as to use the utility functions that take rv_discrete objects. The rv_discrete documentation suggests should be easy. I did the following >>>ddist=rv_discrete(values=(x,[1/len(x) for i in x]),name='test') >>>ddist.pmf(50) array(0.0) Any value I try to get of the pmf seems to be 0. Do I have to explicitly subclass rv_discrete with my data and a _pmf method or something? This seems like a very natural thing to want to do, and hence it seems odd to not have some helper like make_dist(x,name='whatever') . I can take a shot at creating such a function, but I don't want to do so if one exists. 2. Create a continuous probability distribution from something like spline fitting or simple linear interpolation of a the data in X_i. Does this require explict subclassing, or is there a straightforward way to do it that's builtin? I'm not sure if this step is strictly necessary - what I really want to do is be able to draw from the discrete distribution in 1 just by sampling the cdf... maybe this is how it's supposed to work with the discrete distribution, but when I tried to sample it using ddist.rvs, I would always get the input values I specified rather random values sampled from the cdf. I'm on scipy 0.6.0 and numpy 1.0.4

4 4

Polynomial interpolation
by Anne Archibald 01 May '08

01 May '08

Hi, It appears that scipy does not have a facility for using the Lagrange polynomial to interpolate data. (Or did I miss it?) I am well aware of the numerical difficulties this poses, but it (and its generalization, the Hermite polynomial) does prove useful on occasion. I have written prototype implementations of two algorithms for evaluating this polynomial, and I'd like comments before submitting them for inclusion in scipy. One implementation is based on Krogh 1970, "Efficient Algorithms for Polynomial Interpolation and Numerical Differentiation"; it allows the construction of Hermite polynomials (where some derivatives at each point may also be specified) and the evaluation of arbitrary derivatives. It is based on a Neville-like algorithm, so that it does O(n^2) work when constructed but only O(n) per point evaluated, or O(n^2) per point for which all derivatives must be evaluated. (n here is the degree of the polynomial.) The other implementation is based on Berrut and Trefethen 2004, "Barycentric Lagrange Interpolation". This implementation uses a rewriting of the Lagrange polynomial as a rational function, and is efficient and numerically stable. It also allows efficient updating of the y values at which interpolation occurs, as well as addition of new x values. It does not support evaluation of derivatives or construction of Hermite polynomials. Finally, I have written a "PiecewisePolynomial" class for constructing splines in which each piece may have arbitrary degree, and for which the function values and some derivatives are specified at each knot. The intent is that this be used to represent solutions obtained from ODE solvers, using one polynomial for each solver step, with the same order as the solver's internal polynomial solution, and with (some) derivatives matching at the ends. Such a Trajectory object would contain additional information about the solution (for example stiffness or error estimates) beyond what is in PiecewisePolynomial. (I tried implementing Trajectory using splines, which are evaluated in compiled code, but their maximum degree is 5 while the solvers will go up to degree 12.) They all need work, in particular, efficiency would be improved by making the y values vectors, error checking needs to be more robust, and documentation is not in reST form. Ultimately, too, the evaluation functions at least should be written in a compiled language (cython?). But I thought I'd solicit comments on the code first - is the object-oriented design cumbersome? Is including the algorithm in the name confusing to users? Is the calling convention for Hermite polynomials too confusing or error-prone? Thanks, Anne

12 52

numpy.fastCopyAndTranspose Segfaults
by Dave 30 Apr '08

30 Apr '08

I tried posting to the numpy mailing list but it didn't seem to make it through so I'm re-posting in case it's of interest to the scipy community. As the subject says the below code will crash the Python interpreter. I'm running 32-bit Win XP on an x64 xeon cpu. Python 2.5.1 (r251:54863, Apr 18 2007, 08:51:08) [MSC v.1310 32 bit (Intel)] on win32 Type "help", "copyright", "credits" or "license" for more information. from datetime import datetime, timedelta from numpy import asarray, fastCopyAndTranspose import numpy numpy.__version__ '1.0.5.dev5008' timestamps = [datetime(2007,1,1)] while timestamps[-1] < datetime(2007,12,31): timestamps.append(timestamps[-1] + timedelta(1)) x = [asarray(timestamps),asarray(range(len(timestamps)))] fastCopyAndTranspose(x) -Dave

2 2

[ANN] EuroSciPy Abstracts Deadline Reminder
by Travis Vaught 30 Apr '08

30 Apr '08

Greetings, Just a reminder: the abstracts for the EuroSciPy Conference in Leipzig are due by midnight tonight (CST, US [UTC -6]) April, 30. If you'd like to present, please submit your abstract as a PDF, MS Word or plain text file to euroabstracts(a)scipy.org. For more information on the EuroSciPy Conference, please see: http://www.scipy.org/EuroSciPy2008

1 0

Linear Interpolation Question
by Andrea Gavana 30 Apr '08

30 Apr '08

Hi All, I have 2 matrices coming from 2 different simulations: the first column of the matrices is a date (time) at which all the other results in the matrix have been reported (simulation step). In these 2 matrices, very often the simulation steps do not coincide, so I just want to interpolate the results in the second matrix using the dates in the first matrix. The problem is, I have close to 13,000 columns in every matrices, and repeating interp1d all over the columns is quite expensive. An example of what I am doing is as follows: # Loop over all the columns for indx in indices: # Set up a linear interpolation with: # x = dates in the second simulation # y = single column in the second matrix simulation function = interp1d(secondaryMatrixDates, secondaryMatrixResults[:, indx], kind='linear') # Interpolate the second matrix results using the first simulation dates interpolationResults = function(mainMatrixDates) # I need the difference between the first simulation and the second newMatrix[:, indx] = mainMatrixResults[:, indx] - interpolationResults This is somehow a costly step, as it's taking up a lot of CPU (increasing at every iteration) and quite a long time (every column has about 350 data). Is there anything I can do to speed up this loop? Or may someone suggest a better approach? Thank you very much for your suggestions. Andrea. "Imagination Is The Only Weapon In The War Against Reality." http://xoomer.alice.it/infinity77/

3 3