Michael Orlitzky wrote: > > You're going to have a rather hard time convincing me that, upon > encountering a decimal point, any other math software is going to flip > out and silently return incorrect results for basic matrix operations. > > Recall: > > sage: m = matrix([ [-0.3, 0.2, 0.1], > [0.2, -0.4, 0.4], > [0.1, 0.2, -0.5] ]) > > sage: m.echelon_form() > > [ 1.00000000000000 0.000000000000000 0.000000000000000] > [0.000000000000000 1.00000000000000 0.000000000000000] > [0.000000000000000 0.000000000000000 1.00000000000000] >
The problem here is that no one has really implemented a numerically stable version of echelon_form for RR. I thought we called scipy for rank calculations over RDF, but apparently we don't even do that! Scipy answers correctly: sage: import scipy sage: scipy.rank(m.numpy()) 2 Very little attention has been paid to numerically stable linear algebra for non-exact rings. We do get some numerically stable things for matrices over RDF and CDF because those are based on numpy and scipy (which call blas and lapack routines). However, apparently we don't call numpy/scipy to calculate the rank, but instead rely on our unstable echelon form computation! Scipy and numpy don't have facilities for calculating the echelon form because it does not make sense to calculate the echelon form in numerical settings (paraphrasing their words to the best of my memory). That said, there is plenty of work for someone who would like to implement numerically stable, arbitrary precision linear algebra in Sage. I've talked with at least one person who has worked on such a C library (arbitrary precision, numerically stable linear algebra) about contributing to Sage in the last year, but politics have held them up from contributing any code. We don't even need huge experts in the area; even someone implementing basic textbook algorithms would be a huge improvement, I think. If someone wanted to take make a good library of numerically stable linear algebra routines based on mpfr, I think it would be absolutely awesome. Thanks, Jason -- Jason Grout --~--~---------~--~----~------------~-------~--~----~ To post to this group, send email to sage-support@googlegroups.com To unsubscribe from this group, send email to sage-support-unsubscr...@googlegroups.com For more options, visit this group at http://groups.google.com/group/sage-support URL: http://www.sagemath.org -~----------~----~----~----~------~----~------~--~---