Nathan -

On Mon, Jun 21, 2010 at 11:16 PM, Nathan Jahnke <[email protected]> wrote:

> Hi,
>
> PDL is rocking. It's great to do in a real language what I could
> previously do only in MATLAB.
>

Yes, that is nice. :-)


> That said, I've run into what appears to be a bug with svd() in PDL
> 2.4.6 (built using the included cpan client) under Strawberry 5.10.0
> (Windows XP SP3). Searching in the mailing list archive appears to be
> broken and Google is not turning up anything recent, so I am seeking
> advice on this problem here.
>
> I constructed a minimal pair. svd() succeeds here:
>
>
> use PDL;
>
> my $A = pdl(
> [1,2,3],
> [4,5,6]
> );
> print $A;
>
> my ($u, $s_diag, $v) = svd($A);
> my $s = zeroes($v);
> $s->diagonal(0,1) .= $s_diag;
> my $vt = transpose($v);
>
> my $matcopy = $u x $s x $vt;
> print $matcopy;
>
>
> Output:
> [
>  [1 2 3]
>  [4 5 6]
> ]
>
> [
>  [1 2 3]
>  [4 5 6]
> ]
>
> The original matrix has been restored (cf.
> http://www.mail-archive.com/[email protected]/msg00302.html
> )
>
> But svd() fails if the matrix's second row is also [1,2,3]:
>
>
> (snip)
> my $A = pdl(
> [1,2,3],
> [1,2,3]
> );
> (snip)
>
> Output:
> [
>  [1 2 3]
>  [1 2 3]
> ]
>
> [
>  [-1.#IND -1.#IND -1.#IND]
>  [-1.#IND -1.#IND -1.#IND]
> ]
>
> What?? :)
>
> It also fails with matrices with only a single row (apparently, any
> single row EXCEPT a 1x1 matrix):
>
> (snip)
>
> my $A = pdl(
> [4,5,6],
> );
> (snip)
>
> Output:
> [4 5 6]
> [
>  [-1.#IND -1.#IND -1.#IND]
> ]
>
>
> It seems like the single row case (but not the multi-row case) is also
> broken under the same version of PDL under Perl 5.10.0 on my i7 Mac
> laptop (running OS X 10.6.4). Output:
>
> [4 5 6]
> [
>  [nan nan nan]
> ]
>
>
> Singular value decomposition is a great tool and it would be great to
> have it working at full strength in PDL even on the "less commonly
> used" platforms. ;)
>
>
> Thanks for your time,
>
>
> Nathan
>

I was able to reproduce this bug on my Linux machine. I have not figured out
exactly what is causing it, but the code for the algorithm used for SVD
dates to the turn of the century, so it's been around for a while.

That raises a very interesting question: If it's been around so long, why
has nobody noticed it, or fixed it? The answer, I believe, is that when you
are working with real-world data (or good simulations), you will never have
two rows that are strictly identical unless you have a hardware malfunction
with your experiment or a problem with the simulation that produced the
matrix. Under such circumstances, would the broken output be a bug or a
feature? :-)

Anyway, I just wanted to let you know that I'm looking into it. I'll get
back if I have more details to report.

David

-- 
Sent via my carrier pigeon.
_______________________________________________
Perldl mailing list
[email protected]
http://mailman.jach.hawaii.edu/mailman/listinfo/perldl

Reply via email to