On Tue, Nov 10, 2009 at 2:51 PM, Maggie Xiong < [email protected]> wrote:
> While using PDL::Fit::LM, which uses PDL::Slatec for matrix inversion, > I ran into these numbers, > > perldl> p $a =pdl( [qw( 3.5077326e-32 1.1692442e-32 )], [ qw( > 1.1692442e-32 3.8974806e-33 ) ]) > > [ > [3.5077326e-32 1.1692442e-32] > [1.1692442e-32 3.8974806e-33] > ] > > perldl> p PDL::Slatec::matinv $a > > [ > [-1.6666666e+39 5e+39] > [ 5e+39 -1.5e+40] > ] > > But these numbers are different from PDL::MatrixOps' matrix inversion > function > > perldl> p PDL::MatrixOps::inv $a > > [ > [ 2.9544797e+31 -8.8634392e+31] > [-3.1090575e+30 2.6590318e+32] > ] > > Any idea why they are different? I'm using 32bit Ubuntu 8.04 with > PDL-2.4.4_08. I saw the same thing with 64bit Ubuntu 8.04 and > PDL-2.4.4, and in fact, the fit function using Slatec matinv croaks > with this error on the 64bit system > "uninvertible matrix > [ > [3.5077326e-32 -0.33333333] > [1.1692442e-32 0] > ] > " > I don't even know why the numbers in the error printout is different > from $a, because when I check the pdl that was inverted, it was as in > $a. > > > Thanks! > Maggie > > _______________________________________________ > Perldl mailing list > [email protected] > http://mailman.jach.hawaii.edu/mailman/listinfo/perldl > It looks to me like you're trying to diagonalize a matrix that's not invertible. The determinant of the matrix is -2.34e-72. Even if you rescale the matrix values, you get -2.34e-8. This is a pathological case, and unless you've got a really smart matrix inversion routine, you'll run into trouble where the floating-point roundoff error blows up and gives you meaningless results. Notice that neither the Slatec results nor the MatrixOps results give a matrix that is actually the inverse of the one your originally provided. David
_______________________________________________ Perldl mailing list [email protected] http://mailman.jach.hawaii.edu/mailman/listinfo/perldl
