bug in bases of matroids defined from real matrices

[DaveWitteMorris]

As pointed out in this sage-devel thread, there seems to be a serious bug in the handling of matroids defined from matrices with real entries. Calling the bases method more than once results in incorrect answers (it doesn't even give the correct number of bases), and other methods (such as rank, circuits, ...) also give wrong results. The specific incorrect answers depend on the precision of the real field (and some precisions seem to give correct answers), but the problem persists even at high precision (1000 bits).

sage: for prec in range(990, 1000):
....:     A = matrix(RealField(prec),
....:         [(1.0, 1.0, 0.0, 0.0, 0.0, 5.0),
....:          (1.0, 0.0, 1.0, 0.0, 0.0, 0.0),
....:          (0.0, 1.0, 0.0, 1.0, 0.0, 0.0),
....:          (0.0, 0.0, 0.0, 0.0, 1.0, 1.0)
....:         ])
....:     the_matroid = Matroid(A)
....:     number_of_bases = [len(list(the_matroid.bases())) for _ in range(10)]
....:     # the following should all be `True`
....:     print(prec, [number_of_bases[i] == number_of_bases[0] for i in range(10)])
990 [True, True, True, True, True, True, True, True, True, True]
991 [True, False, False, False, False, False, False, False, False, False]
992 [True, False, False, False, False, False, False, False, False, False]
993 [True, False, False, False, False, False, False, False, False, False]
994 [True, True, True, True, True, True, True, True, True, True]
995 [True, False, False, False, False, False, False, False, False, False]
996 [True, False, False, False, False, False, False, False, False, False]
997 [True, False, False, False, False, False, False, False, False, False]
998 [True, True, True, True, True, True, True, True, True, True]
999 [True, False, False, False, False, False, False, False, False, False]

CC: @sagetrac-Rudi @sagetrac-Stefan

Component: combinatorial designs

Issue created by migration from https://trac.sagemath.org/ticket/33657

[DaveWitteMorris]

comment:1

It seems that the first call gives the correct answer (even when the precision is fairly small), but something goes wrong after that:

sage: actual_bases = [
....:     [0, 1, 2, 4], [0, 1, 2, 5], [0, 1, 3, 4], [0, 1, 3, 5],
....:     [0, 1, 4, 5], [0, 2, 3, 4], [0, 2, 3, 5], [0, 3, 4, 5],
....:     [1, 2, 3, 4], [1, 2, 3, 5], [1, 2, 4, 5], [2, 3, 4, 5]]
....: for prec in range(50, 60):
....:     A = matrix(RealField(prec),
....:         [(1.0, 1.0, 0.0, 0.0, 0.0, 5.0),
....:          (1.0, 0.0, 1.0, 0.0, 0.0, 0.0),
....:          (0.0, 1.0, 0.0, 1.0, 0.0, 0.0),
....:          (0.0, 0.0, 0.0, 0.0, 1.0, 1.0)
....:         ])
....:     the_matroid = Matroid(A)
....:     the_bases = [sorted([sorted(list(B)) for B in the_matroid.bases()]) for _ in range(10)]
....:     # the following should all be `True`
....:     print(prec, [the_bases[i] == actual_bases for i in range(10)])
50 [True, True, True, True, True, True, True, True, True, True]
51 [True, False, False, False, False, False, False, False, False, False]
52 [True, False, False, False, False, False, False, False, False, False]
53 [True, False, False, False, False, False, False, False, False, False]
54 [True, True, True, True, True, True, True, True, True, True]
55 [True, False, False, False, False, False, False, False, False, False]
56 [True, False, False, False, False, False, False, False, False, False]
57 [True, False, False, False, False, False, False, False, False, False]
58 [True, True, True, True, True, True, True, True, True, True]
59 [True, False, False, False, False, False, False, False, False, False]

[DaveWitteMorris]

comment:2

The original report used sagemath 9.5 on Arch linux, but I reproduced the problem with 9.1 and some other sagemath versions on ubuntu 20.04 (cocalc), and also with 9.6b7 on MacOS 11.5.2.

[156bf7a1-cdbf-4d2b-a578-797e3af0730c]

comment:4

The subclass LinearMatroid of BasisExchangeMatroid cannot work with matrices with real entries. The implementation assumes exact arithmetic throughout. Assuming this, we (Stephan and I) could compute the set of bases more efficiently than by computing det(B) for each rxr submatrix of the input matrix, among other things.

I do not know of any method that will consistently extract a matroid from a real matrix and which does not involve making some extra assumptions or choices. When I was writing these classes, I did not want to implement a method which silently makes such choices. E.g. approximating the input matrix A by a rational matrix A' and then creating the matroid of A' will do the job, but when is A' good enough as an approximation? The person who created A is in a much better place to make that call.

Perhaps the documentation does not make it clear enough that if the input to Matroid is a matrix A over F, then F needs to have an implementation with exact arithmetic. A solution could be that Matroid test the input matrix for this and if necessary gives a clear error message.

[DaveWitteMorris]

comment:5

Thanks for the explanation. That makes sense.

Maybe I missed it, but I read the linear matroid docstring, and didn't find anything that says the field needs to be exact. So I think the documentation needs to be fixed. (As suggested, adding an error for inexact fields might also be a good idea.)

Extending the implementation to allow inexact fields (at least, in simple cases like the example in the ticket description, where it is clear what was intended) would be a good enhancement, but that would be a separate ticket.

[mkoeppe]

comment:6

+1 on rejecting inexact fields. We do the same in the Polyhedron code for exactly the same reasons. https://github.com/sagemath/sage-prod/blob/develop/src/sage/geometry/polyhedron/constructor.py#L660