Alex Ghitza wrote: > -----BEGIN PGP SIGNED MESSAGE----- > Hash: SHA1 > > Hi, > > There are some inconsistencies in linear algebra over fields like CDF. > Some trouble was already reported in #2256. The following is another issue: > > sage: M = matrix(CDF, 2, 2, [(-1 - 2*I, 5 - 6*I), (-2 - 4*I, 10 - 12*I)]) > sage: M.is_invertible() > True > sage: M.determinant() > 5.3290705182e-15 + 1.7763568394e-15*I > sage: M.inverse() > [ 1.01330991616e+15 - 2.58956978574e+15*I -5.06654958079e+14 + > 1.29478489287e+15*I] > [ 5.62949953421e+14 + 5.62949953421e+14*I -2.81474976711e+14 - > 2.81474976711e+14*I] > > So because of roundoff errors, Sage thinks that we have an invertible > matrix. But the code for echelon_form knows that it's not invertible: > > sage: M.echelon_form() > [ 1.0 1.4 > + 3.2*I] > [-2.22044604925e-16 - 4.4408920985e-16*I > ~ 0] > sage: M.rank() > 1 > > This behavior is inconsistent. Either there is enough precision to > decide that the matrix is invertible, or there isn't. Doing this in two > different ways should not yield two mutually exclusive answers. > > Similar problems occur over CC. > > I hope it's clear by now that I have no idea how CDF and friends > actually work, so I don't know what would be a good solution or even > whether such a solution should exist. My initial guess would be that > finding the determinant requires more computations so there is more > precision loss than in computing the echelon form; in that case, maybe > over such fields we should not base is_invertible() and inverse() on > determinant() but rather on echelon_form(). > > Note also that this problem does not seem to occur over RR or RDF.
Do you know how Matlab and friends treat a similar case? Should we be computing singular values and letting the user decide on a cut-off point for deciding if a singular value really is zero? Jason --~--~---------~--~----~------------~-------~--~----~ To post to this group, send email to sage-devel@googlegroups.com To unsubscribe from this group, send email to [EMAIL PROTECTED] For more options, visit this group at http://groups.google.com/group/sage-devel URLs: http://www.sagemath.org -~----------~----~----~----~------~----~------~--~---
