[SciPy-dev] MCMC, Kalman Filtering, AI for SciPy?
Fernando.Perez at colorado.edu
Thu Sep 30 15:00:29 CDT 2004
eric jones schrieb:
> Of course, we should get some testing feedback on the main platforms
> with the current CVS before doing this. Besides that, does anyone have
> any specific (small) things to address before we do so? If so, what is
> the timeframe?
I have 2 things:
1. I'd like to have the implementation of Frobenius norms fixed. The current
code, in scipy/linalg/basic.py does this:
elif ord in ['fro','f']:
X = scipy_base.mat(x)
return sqrt(sum(diag(X.H * X)))
This is a straightforward implementation of the formula:
frob_norm(A) = sqrt(Trace(A*hermitian_conjugate(A))), 
where A is a _matrix_, so that * is a true matrix product (not element-wise).
Note, however, that this is equivalent to:
frob_norm(A) = sqrt(sum_i(sum_j(|A_ij|**2))) 
But  is an O(N**3) calculation, while  is O[N**2] for an NxN matrix,
since  requires a matrix multiplication.
The following will produce the right answer for complex matrices:
Here's a quick test of the difference, for a 1024x1024 complex matrix:
In : timing (1,scipy.linalg.norm,a)
In : def fnorm(a):
....: return sqrt(sum((a*conjugate(a)).flat)).real
In : timing (1,fnorm,a)
In : scipy.linalg.norm(a)-fnorm(a)
Untested code valid for real/complex should probably be:
2. The matrixmultiply != bug I recently reported is _very_ serious for people
working with large matrices. I also reported it in the numpy list, but nobody
replied. Since scipy can trivially work around it with the one-liner I
showed, I think this really should be done.
Code for #1 is here (but test it further), and I posted the one-liner for #2
yesterday, so these two (IMHO important) fixes should take no time for a
More information about the Scipy-dev