Mailman 3 August 2002 - SciPy-Dev

RE: [SciPy-dev] Reliability of linalg.eig for complex matrices A,B
by Chuck Harris Aug. 29, 2002

Aug. 29, 2002

from scipy import * # # Transcendental eigenvalue problem # # A(x) y = 0 # You are looking for nontrivial vectors in the null space? This isn't an eigenvalue problem, really, since you already know the eigenvalue --- 0. necessary & sufficient condition: determinant(A(x)) = 0 which is a zero finding problem. Solutions then follow by reduction to row echelon form and back substitution. # # B(x) = -d A / d x # You doing linearization of an ODE? Or what? > > > Hi, > > I tried to solve a transcendental eigenvalue problem with scipy. > However It works more or less for r e a l arguments x, but > for complex arguments I cannot observe convergency to any > eigenvalue. I guess its a problem of the eigenvalue solver. > Please find attached my small example. > > Any ideas or suggestions ? > > Nils > Chuck

1 0

Reliability of linalg.eig for complex matrices A,B
by My VDI Freemail Aug. 29, 2002

Aug. 29, 2002

Hi, I tried to solve a transcendental eigenvalue problem with scipy. However It works more or less for r e a l arguments x, but for complex arguments I cannot observe convergency to any eigenvalue. I guess its a problem of the eigenvalue solver. Please find attached my small example. Any ideas or suggestions ? Nils

1 0

Re: [SciPy-dev] Re: scipy / SGI MIPSpro compiler (part 3)
by Pearu Aug. 27, 2002

Aug. 27, 2002

Steve, On Sat, 24 Aug 2002, Steve M. Robbins wrote: > On Sat, Aug 24, 2002 at 05:03:19PM -0500, Pearu wrote: > > Actually there is no exception raised. Only the exception message is shown > > as a warning of missing cblas. cblas (and clapack) is build only if ATLAS > > is available. > > Oh. So you're saying that both "No module named clapack" and "No > module named cblas" came from a single warning message, and that I > simply mis-parsed them as two? No. When running tests then there were also the same warnings for clapack. I see no problem here, though. > > Could you get the latest > > f2py from CVS and run the following tests: > > > > cd f2py2e && python setup.py install > > This bombed with a message about scipy_distutils missing. > I made a symlink scipy_distutils -> ../scipy-cvs/scipy_distutils > where the scipy-cvs directory holds scipy checked out from > CVS within the last day or so. Is this OK? That is OK. Normally scipy_distutils is shipped also with f2py. To install scipy, one needs f2py first. f2py installs f2py2e and scipy_distutils packages. And then scipy re-installs scipy_distutils. > Incidentally, it looks like f2py installs scipy_distutils itself. > Am I right? If so, is there a way to inhibit this? Yes, you are right. f2py and scipy_distutils are developed very closely, so usually there is no problem if either scipy or f2py installs scipy_distutils, just make sure that you have the latest scipy_distutils from scipy CVS. I also use symlink from scipy to f2py2e directory while working with the codes. However, if you really want to stop f2py to re-installing scipy_distutils, then you must edit f2py2e/setup.py. > > > cd f2py2e/tests/f77 > > python return_integer.py --quiet > > Sadly, this failed. Here is the non-quiet output from > "python return_integer.py 2>&1": > > > cc-3303 cc: WARNING File = /tmp/(a)425198.0/ext_return_integermodule.c, Line = 179 > A type qualifier on a return type is meaningless. > > const void (*f2py_func)(int*,int*)) { > ^ These warnings are now fixed. > Traceback (most recent call last): > File "return_integer.py", line 99, in ? > runtest(t) > File "return_integer.py", line 56, in runtest > assert t(-123)==-123 > AssertionError > > > > I ran the test "t(-123)" by hand, and it is only the "t1" > version that fails: > > >>> t1(-123) > 133 > > So the integer*1 is unsigned, perhaps? The SGI cc manual says "The > default is to treat values of type char as if they had type unsigned > char". I guess the generated code should use "signed char". Yes, fixed in f2py CVS. > > python return_real.py --quiet > > This time I got a bus error from python. :-( > Running the tests manually, I find that > > abs(t([234])-234)<=err > > with t=t0 is the first test that coredumps. This is known bug with SGI compilers (when using optimization flags) and Python 2.1. See http://sourceforge.net/tracker/index.php?func=detail&aid=435026&group_id=54… I have worked around this bug for f2py and SGI. > > python return_complex.py --quiet > > try: raise RuntimeError,`t(10l**400)` > RuntimeError: (inf+0j) Fixed. > > python return_character.py --quiet > Traceback (most recent call last): > File "return_character.py", line 62, in ? > runtest(t) > File "return_character.py", line 50, in runtest > assert t(array(77))=='M' > AssertionError Fixed. > Both t0(array(77)) and t1(array(77)) return ''. > > There may be some other weirdness going on. While running some of > the "return_character.py" tests by hand, I had the python process > killed by the system, with a note in the syslog: > > Aug 24 20:48:14 1A:wart unix: |$(0x6dd)ALERT: Process [python] 428471 generated trap, but has signal 11 held or ignored > Aug 24 20:48:14 6A:wart unix: epc 0xfa3a38c ra 0xfa3abd4 badvaddr 0x4973f130 > Aug 24 20:48:14 6A:wart unix: Process has been killed to prevent infinite loop > Aug 24 20:49:10 1A:wart unix: |$(0x6dd)ALERT: Process [python] 428406 generated trap, but has signal 11 held or ignored > Aug 24 20:49:10 6A:wart unix: epc 0xfa3a38c ra 0xfa3abd4 badvaddr 0x4973f180 > Aug 24 20:49:10 6A:wart unix: Process has been killed to prevent infinite loop > > This diagnostic is a new one to me. Hmm, can you reproduce this one and send the details to me? Thanks for the bug reports. I believe they are all fixed now in f2py CVS. Could you verify that? Thanks, Pearu

6 21

Re: scipy / SGI MIPSpro compiler (part 3)
by Steve M. Robbins Aug. 24, 2002

Aug. 24, 2002

Hi again, I've managed to build and install scipy (using the two patches sent previously) on an IRIX 6.5 system with SGI's MIPSpro C compiler and GNU Fortran (GCC 3.2). I'm using CVS scipy. I built *without* ATLAS, using the BLAS and LAPACK sources from netlib. (Many thanks to whoever made it possible for scipy to build blas & lapack from source.) I've moved on to running the scipy tests. I'm not sure which failures are worth worrying about. For instance, a simple "import scipy" results in a warning: >>> import scipy /home/bic/stever/irix-6/lib/python2.1/site-packages/scipy/linalg/lapack.py:24: UserWarning: exceptions.ImportError: No module named clapack warnings.warn(clapack.__doc__) exceptions.ImportError: No module named cblas I don't think this is a problem -- the comments in lapack.py suggest that either Fortran or C lapack is okay. I'm not sure why the exception is raised. (??) I'm not sure where the exception about cblas comes from -- is cblas required? (I'm using the Fortran BLAS sources, compiled during the scipy build) Now, when I run the test suite itself, python -c 'import scipy; scipy.test(level=1)' 2>err | tee out I get some more worrying messages. Here are the standard output and standard error, respectively. I'd appreciate any comments if I've done something wrong in building scipy. Are these messages about "No test suite found for ..." normal? --------------------- out ---------------------------------- exceptions.ImportError: No module named cblas No test suite found for scipy.helper No test suite found for scipy.__cvs_version__ No test suite found for scipy.helpmod creating test suite for: scipy.common No test suite found for scipy.scipy_tempfile No test suite found for scipy.sync No test suite found for scipy.pilutil No test suite found for scipy.proc No test suite found for scipy.ga No test suite found for scipy.io No test suite found for scipy.cow No test suite found for scipy.integrate creating test suite for: scipy.stats.distributions No test suite found for scipy.stats.pstat creating test suite for: scipy.stats.morestats No test suite found for scipy.stats.rv creating test suite for: scipy.stats.stats No test suite found for scipy.stats.rv2 creating test suite for: scipy_base.limits creating test suite for: scipy_base.index_tricks creating test suite for: scipy_base.shape_base creating test suite for: scipy_base.function_base creating test suite for: scipy_base.matrix_base No test suite found for scipy_base.scimath creating test suite for: scipy_base.type_check No test suite found for scipy_base.polynomial No test suite found for scipy.special No test suite found for scipy.cluster FAILURE to import scipy.weave.accelerate_tools <string>:0: AttributeError: 'scipy.weave' module has no attribute 'accelerate_tools' (in ?) No test suite found for scipy.weave.scalar_info No test suite found for scipy.weave.blitz_info No test suite found for scipy.weave.cxx_info creating test suite for: scipy.weave.catalog No test suite found for scipy.weave.blitz_spec No test suite found for scipy.weave.dumbdbm_patched creating test suite for: scipy.weave.scalar_spec No test suite found for scipy.weave.lib2def creating test suite for: scipy.weave.sequence_spec No test suite found for scipy.weave.misc No test suite found for scipy.weave.code_blocks No test suite found for scipy.weave.inline_info No test suite found for scipy.weave.converters No test suite found for scipy.weave.unicode_spec creating test suite for: scipy.weave.inline_tools No test suite found for scipy.weave.bytecodecompiler creating test suite for: scipy.weave.slice_handler No test suite found for scipy.weave.unicode_info creating test suite for: scipy.weave.build_tools creating test suite for: scipy.weave.blitz_tools No test suite found for scipy.weave.base_info No test suite found for scipy.weave.swig_info creating test suite for: scipy.weave.ast_tools creating test suite for: scipy.weave.ext_tools building extensions here: /home/bic/stever/.python21_compiled/4190480 building extensions here: /home/bic/stever/.python21_compiled/4190481 No test suite found for scipy.weave.wx_info No test suite found for scipy.weave.common_info No test suite found for scipy.weave.wx_spec creating test suite for: scipy.weave.common_spec No test suite found for scipy.weave.dumb_shelve No test suite found for scipy.weave.base_spec No test suite found for scipy.weave.conversion_code No test suite found for scipy.weave.standard_array_info creating test suite for: scipy.weave.size_check No test suite found for scipy.weave.conversion_code_old creating test suite for: scipy.weave.standard_array_spec No test suite found for scipy.interpolate No test suite found for scipy.signal No test suite found for anneal No test suite found for scipy.optimize.minpack No test suite found for scipy.optimize.optimize No test suite found for scipy.optimize.common_routines creating test suite for: scipy.optimize.zeros No test suite found for scipy.linalg.interface_gen creating test suite for: scipy.linalg.blas No test suite found for scipy.linalg.lapack creating test suite for: scipy.linalg.basic exceptions.ImportError: No module named cblas No test suite found for scipy.linalg.flinalg creating test suite for: scipy.linalg.decomp creating test suite for: scipy.linalg.matfuncs !! FAILURE building test for scipy.linalg.matfuncs <string>:1: ImportError: No module named test_matfuncs (in ?) Testing uniform Testing stnorm Testing norm Testing lognorm Testing expon Testing beta Testing power Testing bradford Testing burr Testing fisk Testing cauchy Testing halfcauchy Testing foldcauchy Testing gamma Testing gengamma Testing loggamma Testing alpha Testing anglit Testing arcsine Testing betaprime Testing erlang Testing dgamma Testing extreme3 Testing exponweib Testing exponpow Testing frechet Testing gilbrat Testing f Testing ncf Testing chi2 Testing chi Testing nakagami Testing genpareto Testing genextreme Testing genhalflogistic Testing pareto Testing lomax Testing halfnorm Testing halflogistic Testing fatiguelife Testing foldnorm Testing ncx2 Testing t Testing nct Testing weibull Testing dweibull Testing maxwell Testing rayleigh Testing genlogistic Testing logistic Testing gumbel Testing gompertz Testing hypsecant Testing laplace Testing reciprocal Testing triang Testing tukeylambda Ties preclude use of exact statistic. Ties preclude use of exact statistic. warning: specified build_dir '_bad_path_' does not exist or is or is not writable. Trying default locations warning: specified build_dir '..' does not exist or is or is not writable. Trying default locations warning: specified build_dir '_bad_path_' does not exist or is or is not writable. Trying default locations warning: specified build_dir '..' does not exist or is or is not writable. Trying default locations TESTING CONVERGENCE zero should be 1 function f2 cc.bisect : 1.0000000000001952 cc.ridder : 1.0000000000004661 cc.brenth : 0.9999999999999997 cc.brentq : 0.9999999999999577 function f3 cc.bisect : 1.0000000000001952 cc.ridder : 1.0000000000000000 cc.brenth : 1.0000000000000009 cc.brentq : 1.0000000000000011 function f4 cc.bisect : 1.0000000000001952 cc.ridder : 1.0000000000001454 cc.brenth : 0.9999999999993339 cc.brentq : 0.9999999999993339 function f5 cc.bisect : 1.0000000000001952 cc.ridder : 1.0000000000004574 cc.brenth : 0.9999999999991444 cc.brentq : 0.9999999999991444 function f6 cc.bisect : 1.0000000000001952 cc.ridder : 1.0000000000004949 cc.brenth : 0.9999999999990239 cc.brentq : 1.0000000000001117 ----------------------------- err ----------------------------------- /home/bic/stever/irix-6/lib/python2.1/site-packages/scipy/linalg/lapack.py:24: UserWarning: exceptions.ImportError: No module named clapack warnings.warn(clapack.__doc__) /home/bic/stever/irix-6/lib/python2.1/site-packages/scipy/linalg/lapack.py:24: UserWarning: exceptions.ImportError: No module named clapack warnings.warn(clapack.__doc__) ...................................................................F..................................................................................................................................................................................................................EE..........................................................................................................................F..FF......scopy:n=4 ..scopy:n=3 ....dcopy:n=4 ..dcopy:n=3 ....ccopy:n=4 ..ccopy:n=3 ....zcopy:n=4 ..zcopy:n=3 .....saxpy:n=4 ..saxpy:n=3 .....daxpy:n=4 ..daxpy:n=3 .....caxpy:n=4 ..caxpy:n=3 .....zaxpy:n=4 ..zaxpy:n=3 ...sscal:n=4 ...dscal:n=4 ...cscal:n=4 ...zscal:n=4 ....sswap:n=4 ..sswap:n=3 ....dswap:n=4 ..dswap:n=3 ....cswap:n=4 ..cswap:n=3 ....zswap:n=4 ..zswap:n=3 .............................................................. ====================================================================== ERROR: check_add_function_ordered (test_catalog.test_catalog) ---------------------------------------------------------------------- Traceback (most recent call last): File "/usr/local/unstable/lib/python2.1/site-packages/scipy/weave/tests/test_catalog.py", line 286, in check_add_function_ordered q.add_function('f',string.upper) File "/usr/local/unstable/lib/python2.1/site-packages/scipy/weave/catalog.py", line 575, in add_function self.add_function_persistent(code,function) File "/usr/local/unstable/lib/python2.1/site-packages/scipy/weave/catalog.py", line 615, in add_function_persistent module = getmodule(function) File "/usr/local/unstable/lib/python2.1/site-packages/scipy/weave/catalog.py", line 66, in getmodule if mod and object in mod.__dict__.values(): File "/usr/local/unstable/lib/python2.1/site-packages/scipy_distutils/misc_util.py", line 43, in __getattr__ raise self._info[0],self._info[1] ImportError: No module named cblas ====================================================================== ERROR: Test persisting a function in the default catalog ---------------------------------------------------------------------- Traceback (most recent call last): File "/usr/local/unstable/lib/python2.1/site-packages/scipy/weave/tests/test_catalog.py", line 274, in check_add_function_persistent1 q.add_function_persistent('code',i) File "/usr/local/unstable/lib/python2.1/site-packages/scipy/weave/catalog.py", line 615, in add_function_persistent module = getmodule(function) File "/usr/local/unstable/lib/python2.1/site-packages/scipy/weave/catalog.py", line 66, in getmodule if mod and object in mod.__dict__.values(): File "/usr/local/unstable/lib/python2.1/site-packages/scipy_distutils/misc_util.py", line 43, in __getattr__ raise self._info[0],self._info[1] ImportError: No module named cblas ====================================================================== FAIL: check_basic (test_morestats.test_shapiro) ---------------------------------------------------------------------- Traceback (most recent call last): File "/usr/local/unstable/lib/python2.1/site-packages/scipy/stats/tests/test_morestats.py", line 30, in check_basic assert_almost_equal(w,0.90047299861907959,8) File "/usr/local/unstable/lib/python2.1/site-packages/scipy_base/testing.py", line 283, in assert_almost_equal assert round(abs(desired - actual),decimal) == 0, msg AssertionError: Items are not equal: DESIRED: 0.900472998619 ACTUAL: 0.900472760201 ====================================================================== FAIL: check_asum (test_blas.test_fblas1_simple) ---------------------------------------------------------------------- Traceback (most recent call last): File "/usr/local/unstable/lib/python2.1/site-packages/scipy/linalg/tests/test_blas.py", line 59, in check_asum assert_almost_equal(f([3,-4,5]),12) File "/usr/local/unstable/lib/python2.1/site-packages/scipy_base/testing.py", line 283, in assert_almost_equal assert round(abs(desired - actual),decimal) == 0, msg AssertionError: Items are not equal: DESIRED: 12 ACTUAL: 0.0 ====================================================================== FAIL: check_dot (test_blas.test_fblas1_simple) ---------------------------------------------------------------------- Traceback (most recent call last): File "/usr/local/unstable/lib/python2.1/site-packages/scipy/linalg/tests/test_blas.py", line 66, in check_dot assert_almost_equal(f([3,-4,5],[2,5,1]),-9) File "/usr/local/unstable/lib/python2.1/site-packages/scipy_base/testing.py", line 283, in assert_almost_equal assert round(abs(desired - actual),decimal) == 0, msg AssertionError: Items are not equal: DESIRED: -9 ACTUAL: 0.0 ====================================================================== FAIL: check_nrm2 (test_blas.test_fblas1_simple) ---------------------------------------------------------------------- Traceback (most recent call last): File "/usr/local/unstable/lib/python2.1/site-packages/scipy/linalg/tests/test_blas.py", line 75, in check_nrm2 assert_almost_equal(f([3,-4,5]),math.sqrt(50)) File "/usr/local/unstable/lib/python2.1/site-packages/scipy_base/testing.py", line 283, in assert_almost_equal assert round(abs(desired - actual),decimal) == 0, msg AssertionError: Items are not equal: DESIRED: 7.07106781187 ACTUAL: 3.73560645334e+270 ---------------------------------------------------------------------- Ran 559 tests in 11.403s FAILED (failures=4, errors=2) ---------------------------------------------------------------------- Thanks, -Steve

2 2

scipy / SGI MIPSpro compiler (part 2)
by Steve M. Robbins Aug. 23, 2002

Aug. 23, 2002

Hi again, Got stuck in scipy_base/_compiled_base.c with three errors and a warning. cc -O -OPT:Olimit=0 -I/usr/local/unstable/include/python2.1 -c /tmp/smr/scipy/scipy_base/_compiled_base.c -o build/temp.irix64-6.5-2.1/_compiled_base.o error: command 'cc' failed with exit status 2 swmgr@wart{scipy}cc -O -OPT:Olimit=0 -I/usr/local/unstable/include/python2.1 -c /tmp/smr/scipy/scipy_base/_compiled_base.c cc-3316 cc: ERROR File = /tmp/smr/scipy/scipy_base/_compiled_base.c, Line = 45 The expression must be a pointer to a complete object type. for (k=0; k < asize; k++,ip+=instride) { ^ cc-3316 cc: ERROR File = /tmp/smr/scipy/scipy_base/_compiled_base.c, Line = 48 The expression must be a pointer to a complete object type. for (j=0; j < copied; j++,op2+=elsize) { ^ cc-3316 cc: ERROR File = /tmp/smr/scipy/scipy_base/_compiled_base.c, Line = 58 The expression must be a pointer to a complete object type. op += elsize; /* Get ready to put next match */ ^ cc-1174 cc: WARNING File = /tmp/smr/scipy/scipy_base/_compiled_base.c, Line = 87 The variable "f1" was declared but never referenced. PyObject *m, *d, *s, *f1; ^ My understanding is that you can't do pointer arithmetic on void* types. I changed them to char*, assuming that the strides are given in bytes. Index: scipy_base/_compiled_base.c =================================================================== RCS file: /home/cvsroot/world/scipy/scipy_base/_compiled_base.c,v retrieving revision 1.1 diff -u -b -B -r1.1 _compiled_base.c --- scipy_base/_compiled_base.c 16 Aug 2002 13:09:42 -0000 1.1 +++ scipy_base/_compiled_base.c 23 Aug 2002 20:29:44 -0000 @@ -15,8 +15,8 @@ int asize, abytes, new; int copied=0, nd; int instride=0, elsize, k, j, dims[1]; - void *ip, *op; /* Current memory buffer */ - void *op2; + char *ip, *op; /* Current memory buffer */ + char *op2; static char *kwlist[] = {"input", NULL}; @@ -84,7 +84,7 @@ }; DL_EXPORT(void) init_compiled_base(void) { - PyObject *m, *d, *s, *f1; + PyObject *m, *d, *s; /* Create the module and add the functions */ m = Py_InitModule("_compiled_base", methods);

1 0

license for djbfft?
by eric jones Aug. 23, 2002

Aug. 23, 2002

Dr. Bernstein, I've copied the scipy-dev list on this mail. We develop an Open Source (BSD style) package of scientific computing modules for Python at www.scipy.org. We've been using fftpack and, optionally fftw, for our fft module. fftw's license is GPL, otherwise we would be using it as the standard package. One of the community developers (Pearu Peterson) has been working hard on improving the fft module lately and ran across your fft package. Congrats on making such a nice tool. Combining your algorithms for power of 2 and fftpack's generic length ffts provides pretty much equivalent or better performance than fftw on most fft lengths. Pearu posted benchmarks below: http://www.scipy.org/site_content/mailman?fn=scipy-dev/2002-August/00110 1.html I searched for a license in the djbfft source but didn't see one. I also didn't see any mention of djbfft on http://cr.yp.to/distributors.html. So, the question is, can we include djbfft source code in SciPy and distribute it under its BSD license? The SciPy community is certainly willing to contribute every change we make to the source code. Our changes are typically made to help with compatibility across platforms and only occasionally to the core algorithms. Thanks for making your tools available on the web, and I look forward to your reply regarding the licensing issue. Regards, eric

1 0

scipy / SGI MIPSpro compiler
by Steve M. Robbins Aug. 23, 2002

Aug. 23, 2002

Hello, I'm trying to build scipy on an IRIX 6.5 system using the vendor compiler. (That's what I built python with). The compiler complained bitterly about signal/S_bspline_util.c because there are function prototypes mixed with old-style function declarations. Here's a patch that changes them all to "new style" (I forget the technical term) declarations. Cheers, -Steve Index: signal/S_bspline_util.c =================================================================== RCS file: /home/cvsroot/world/scipy/signal/S_bspline_util.c,v retrieving revision 1.4 diff -u -b -B -r1.4 S_bspline_util.c --- signal/S_bspline_util.c 30 Jun 2002 23:55:05 -0000 1.4 +++ signal/S_bspline_util.c 23 Aug 2002 20:11:36 -0000 @@ -28,12 +28,8 @@ /* with a given starting value loaded into the array */ void -S_IIR_order1 (a1, a2, x, y, N, stridex, stridey) - float a1; - float a2; - float *x; - float *y; - int N, stridex, stridey; +S_IIR_order1 ( float a1, float a2, float *x, float *y, + int N, int stridex, int stridey ) { float *yvec = y+stridey; float *xvec = x+stridex; @@ -51,13 +47,8 @@ /* y[n] = a1 * x[n] + a2 * y[n-1] + a3 * y[n-2] */ /* with two starting values loaded into the array */ void -S_IIR_order2 (a1, a2, a3, x, y, N, stridex, stridey) - float a1; - float a2; - float a3; - float *x; - float *y; - int N, stridex, stridey; +S_IIR_order2 ( float a1, float a2, float a3, float *x, float *y, + int N, int stridex, int stridey ) { float *yvec = y+2*stridey; float *xvec = x+2*stridex; @@ -86,14 +77,8 @@ */ void -S_IIR_order2_cascade (cs, z1, z2, y1_0, x, yp, N, stridex, stridey) - float cs; - float z1; - float z2; - float y1_0; - float *x; - float *yp; - int N, stridex, stridey; +S_IIR_order2_cascade ( float cs, float z1, float z2, float y1_0, float *x, + float *yp, int N, int stridex, int stridey) { float *yvec = yp+stridey; float *xvec = x+stridex; @@ -137,13 +122,8 @@ */ int -S_IIR_forback1 (c0, z1, x, y, N, stridex, stridey, precision) - float c0; - float z1; - float *x; - float *y; - int N, stridex, stridey; - float precision; +S_IIR_forback1 ( float c0, float z1, float *x, float *y, + int N, int stridex, int stridey, float precision ) { float *yp = NULL; float *xptr = x; @@ -189,12 +169,8 @@ /* h must be odd length */ /* strides in units of sizeof(float) bytes */ void -S_FIR_mirror_symmetric (in, out, N, h, Nh, instride, outstride) - float *in; - float *out; - int N, Nh; - float *h; - int instride, outstride; +S_FIR_mirror_symmetric ( float *in, float *out, int N, float *h, int Nh, + int instride, int outstride ) { int n, k; int Nhdiv2 = Nh >> 1; @@ -254,14 +230,9 @@ } int -S_separable_2Dconvolve_mirror(in, out, M, N, hr, hc, Nhr, - Nhc, instrides, outstrides) - float *in; - float *out; - int M, N; - float *hr, *hc; - int Nhr, Nhc; - int *instrides, *outstrides; +S_separable_2Dconvolve_mirror( float *in, float *out, int M, int N, + float *hr, float *hc, int Nhr, int Nhc, + int *instrides, int *outstrides ) { int m, n; float *tmpmem; @@ -305,10 +276,7 @@ static float S_hs(int,float,double,double); float -S_hc(k, cs, r, omega) - int k; - float cs; - double r, omega; +S_hc( int k, float cs, double r, double omega ) { if (k < 0) return 0.0; if (omega == 0.0) @@ -319,10 +287,7 @@ } float -S_hs(k, cs, rsq, omega) - int k; - float cs; - double rsq, omega; +S_hs( int k, float cs, double rsq, double omega ) { float cssq; float c0; @@ -383,12 +348,8 @@ */ int -S_IIR_forback2 (r, omega, x, y, N, stridex, stridey, precision) - double r,omega; - float *x; - float *y; - int N, stridex, stridey; - float precision; +S_IIR_forback2 ( double r, double omega, float *x, float *y, + int N, int stridex, int stridey, float precision ) { float cs; float *yp = NULL; @@ -497,13 +458,8 @@ */ int -S_cubic_spline2D(image, coeffs, M, N, lambda, strides, cstrides, precision) - float *image; - float *coeffs; - int M, N; - double lambda; - int *strides, *cstrides; - float precision; +S_cubic_spline2D( float *image, float *coeffs, int M, int N, double lambda, + int *strides, int *cstrides, float precision ) { double r, omega; float *inptr; @@ -590,13 +546,9 @@ */ int -S_quadratic_spline2D(image, coeffs, M, N, lambda, strides, cstrides, precision) - float *image; - float *coeffs; - int M, N; - double lambda; - int *strides, *cstrides; - float precision; +S_quadratic_spline2D( float *image, float *coeffs, int M, int N, + double lambda, int *strides, int *cstrides, + float precision ) { double r; float *inptr;

1 0

djbfft - The fastest FFT in the universe ;)
by Pearu Aug. 23, 2002

Aug. 23, 2002

Hi, I have found a fft library djbfft by D. J. Bernstein that is even faster than FFTW. I have implemented its (currently optional) support to fftpack2. djbfft provides power-of-2 FFT for single/double real/complex data. Therefore, djbfft must live with some other fft library that would provide non-power-of-2 FFT. However, there is upto 1.7 times speed improvement for power-of-2 FFTs when using djbfft. djbfft home page is at http://cr.yp.to/djbfft.html Eric, can you check if there are any license issues with djbfft that would prevent us including it to scipy? At least to me it is not clear if djbfft is in public domain or not? The author notes that "djbfft can be used in my own code" but he does not specify if it can be distributed in a package like scipy. See also relevant notes in http://cr.yp.to/distributors.html but the link above does not mention djbfft. Here are the latest testing results (only with complex input): fftpack2 = f2py interface to fftpack FFT = f2c version of fftpack (from NumPy) scipy = current scipy interface to Fortran fftpack Fast Fourier Transform ================================== size |fftpack2 | FFT | scipy | fftw2 100 | 0.36 | 0.39 | 0.37 | skipped (secs for 7000 calls) 1000 | 0.56 | 0.63 | 0.53 | skipped (secs for 2000 calls) 256 | 0.60 | 0.74 | 0.61 | skipped (secs for 10000 calls) 512 | 0.96 | 1.29 | 0.98 | skipped (secs for 10000 calls) 1024 | 0.17 | 0.25 | 0.17 | skipped (secs for 1000 calls) 2048 | 0.35 | 0.49 | 0.35 | skipped (secs for 1000 calls) 4096 | 0.52 | 0.70 | 0.52 | skipped (secs for 500 calls) 8192 | 2.50 | 2.83 | 2.59 | skipped (secs for 500 calls) fftpack2 = djbfft & fftpack Fast Fourier Transform ================================== size |fftpack2 | FFT | scipy | fftw2 100 | 0.37 | 0.39 | 0.37 | skipped (secs for 7000 calls) 1000 | 0.57 | 0.63 | 0.53 | skipped (secs for 2000 calls) 256 | 0.57 | 0.74 | 0.61 | skipped (secs for 10000 calls) 512 | 0.81 | 1.28 | 0.97 | skipped (secs for 10000 calls) 1024 | 0.14 | 0.25 | 0.17 | skipped (secs for 1000 calls) 2048 | 0.27 | 0.58 | 0.35 | skipped (secs for 1000 calls) 4096 | 0.42 | 0.70 | 0.51 | skipped (secs for 500 calls) 8192 | 1.40 | 2.99 | 2.66 | skipped (secs for 500 calls) . fftpack2 = djbfft & fftw Fast Fourier Transform ================================== size |fftpack2 | FFT | scipy | fftw2 100 | 0.35 | 0.39 | 0.37 | skipped (secs for 7000 calls) 1000 | 0.47 | 0.63 | 0.53 | skipped (secs for 2000 calls) 256 | 0.57 | 0.74 | 0.61 | skipped (secs for 10000 calls) 512 | 0.82 | 1.29 | 0.98 | skipped (secs for 10000 calls) 1024 | 0.14 | 0.25 | 0.17 | skipped (secs for 1000 calls) 2048 | 0.27 | 0.48 | 0.35 | skipped (secs for 1000 calls) 4096 | 0.42 | 0.71 | 0.52 | skipped (secs for 500 calls) 8192 | 1.49 | 2.95 | 2.57 | skipped (secs for 500 calls) Regards, Pearu

1 0

Resultant of two polynomials - Sylvester matrix.
by Nils Wagner Aug. 21, 2002

Aug. 21, 2002

Hi, I am looking for a function which computes the resultant of two polynomials. It is defined to be the determinant of their Sylvester matrix. A small example would be appreciated. Thanks in advance. Nils

1 0

fftpack preliminary testing results
by Pearu Aug. 20, 2002

Aug. 20, 2002

Hi, I have started to create a testing site for fftpack and here are some preliminary benchmark results: Fast Fourier Transform ================================== size |fftpack2 | FFT | scipy | fftw2 100 | 1.55 | 1.07 | 1.65 | 0.42 (secs for 7000 calls) 1000 | 1.77 | 1.82 | 3.01 | 0.65 (secs for 2000 calls) 512 | 3.97 | 4.56 | 3.20 | 1.50 (secs for 10000 calls) 2048 | 1.30 | 1.82 | 1.26 | 0.68 (secs for 1000 calls) 8192 | 5.65 | 6.61 | 5.68 | 2.39 (secs for 500 calls) fftpack2 refers to f2py'd version of fftpack FFT refers to fftpack (f2c version of fftpack) shipped with Numeric scipy refers to scipy.fftpack fftw2 refers to quickly-put-together interface to FFTW All tests are carried out on complex 1D contiguous data. All fft routines pass correctness check (i.e. they all return the same result). Notes: * Timings for fftpack2 and scipy fluctuate strongly. In average, their performance is the same. * fftpack2/scipy perform better than FFT for large data sets. * FFTW beats them all (which makes more sense than my earlier tests). Pearu

1 0