On Sun, 16 Apr 2023 at 15:52, Trevor Karn <trevor.k.k...@gmail.com> wrote: > > I don't have much understanding of how floating point arithmetic works, but > the argument > > >If you're writing python code, you should expect 2*0.5 to return a > >float. But if you're learning linear algebra for the first time and > >typing a matrix into the Sage notebook, typing 0.5 instead of 1/2 > >should not ruin the entire assignment without so much as a warning. > > seems important to me. I also sometimes worry that Sage seems to be written > for users who are both experts in [insert mathematical field here] and > Python. This prevents (a) mathematical experts who don't know python from > using Sage and (b) more junior students from using sage. Maybe one can make > the argument that mathematical experts should just learn python (although I > would disagree), but I think (b) is a genuine problem. Mathematica doesn't > seem to require the user to be expert in either math or Wolfram language (or > maybe I haven't used it enough). I don't have an actionable complaint here, > but it is something I encounter when helping folks around my math department > with their code. I'm trying to keep it in mind as I develop code.
This hardly requires expert knowledge of Python but rather just the understanding that 0.5 means a floating point number and then some understanding of the implications of floating point numbers on exactness. This is not uniquely a feature of Python by any stretch. The use of a decimal point to denote approximate or floating point numbers is the same in Mathematica, Maple, Maxima, Matlab, and not just Python but most general purpose programming languages. You can see the Mathematica and Maple docs explaining this here: https://www.maplesoft.com/support/help/maple/view.aspx?path=float https://reference.wolfram.com/language/howto/ChangeTheFormatOfNumbers.html Select quotes from there: > > > Enter an approximate real number by explicitly including a decimal point: > > > The presence of a floating-point number in an expression generally > > > implies that the computation will use floating-point evaluation. > > > any calculation mixing exact and approximate numbers gives an approximate > > > result The fact that writing 0.5 (rather than 1/2) implies the use of floating point which in turn implies rounding errors and inexactness is certainly something that can surprise beginners. In the context of Sage it could have been a valid choice not to do that and to have 0.5 mean an exact rational number when it appears as a literal in the user's code. There is no difficulty in making that work any more than making 1/2 return a rational number (since in Python itself 1/2 is a float). This is a design choice though and using decimal literals for floats is useful for anyone who actually wants to do calculations with floats since it means that there is a succinct syntax for floats as well as a separate syntax for exact rationals. > I see that > > sage: A = matrix([[-3, 2, 1 ],[ 2,-4, 4 ], [ 1, 2,-5 ]]) > sage: B = (2*0.5*A) > sage: rank(B) > 3 > > What about floating point arithmetic makes this break? I print (B) and see > > [-3.00000000000000 2.00000000000000 1.00000000000000] > [ 2.00000000000000 -4.00000000000000 4.00000000000000] > [ 1.00000000000000 2.00000000000000 -5.00000000000000] > > which doesn't have any obvious reason to be not rank-2. Could someone ELI5? Arithmetic is needed to be able to compute the rank. With floats that arithmetic is going to be approximate and will have some rounding error. Knowing that the rank is 2 here requires knowing that the third row is *exactly* equal to some linear combination of the first and second but rounding errors will make the linear combinations not be exactly equal. It is basically impossible to compute the rank of a matrix using inexact arithmetic since the tiniest rounding errors will make most matrices appear to be full rank. That does not stop there being floating point routines for computing the rank but they work differently, are inherently heuristic and need thresholds as discussed here: https://numpy.org/doc/stable/reference/generated/numpy.linalg.matrix_rank.html NumPy's matrix_rank function does give the expected answer here and is probably a better implementation of this for machine precision: sage: np.linalg.matrix_rank(B) 2 The Mathematica docs have a caveat about computing the rank of a matrix of floats: https://reference.wolfram.com/language/ref/MatrixRank.html Under the "possible issues" section it says: > > > Use machine arithmetic. Machine numbers cannot distinguish between 9 and > > > 9 + 1e-20 It then shows an example in which the rank of a matrix of floats is computed incorrectly. With that example NumPy fails for the same reason: sage: A = np.array([[1, 2, 3], [4, 5, 6], [7, 8, 9 + 1e-20]]) sage: A array([[1., 2., 3.], [4., 5., 6.], [7., 8., 9.]]) sage: np.linalg.matrix_rank(A) 2 It is not that the matrix_rank function gets anything wrong here but the use of floats even means that 9 + 1e-20 is equal to 9 so the matrix A already has the wrong values before matrix_rank is called. The matrix_rank function returned the exact correct answer for the inputs that it was given. This also shows why converting all floats to rational and then computing the rank does not work in general. If the matrix arrived at the rank function with floats in it then most likely rounding errors have already happened. In that case computing the exact rank based on the exact values of the floats in the matrix is probably not what the user intended. -- Oscar -- You received this message because you are subscribed to the Google Groups "sage-devel" group. To unsubscribe from this group and stop receiving emails from it, send an email to sage-devel+unsubscr...@googlegroups.com. To view this discussion on the web visit https://groups.google.com/d/msgid/sage-devel/CAHVvXxRADFToFBuOQJqbme1aLwrrAn2GQNFGBvvRHbVTL0V5jQ%40mail.gmail.com.
