Some determinants can not be computed

[edd8e884-f507-429a-b577-5d554626c0fe]

To use some fast algorithms that only work for fields, the determinant method uses the is_field method for rings which uses the is_maximal method for ideals. Unfortunately, this method is not always implemented.

The determinant is a very fundamental function, and there are many division-free algorithms that work in every ring, hence this method should never raise an error.

Here is an example of the bug i encountered while testing a conjecture in combinatorics:

sage: A = GF(2)['x,y,z']                                                                                
sage: A.inject_variables()                                                                              
Defining x, y, z
sage: R = A.quotient(x^2 + 1).quotient(y^2 + 1).quotient(z^2 + 1)                                       
sage: R.inject_variables()                                                                              
Defining xbarbarbar, ybarbarbar, zbarbarbar
sage: M = matrix([[1,1,1,1],[xbarbarbar,ybarbarbar,1,1],[0,1,zbarbarbar,1],[xbarbarbar,zbarbarbar,1,1]])
sage: M.determinant()
---------------------------------------------------------------------------
NotImplementedError                       Traceback (most recent call last)

/tmp/<ipython console> in <module>()

/opt/sagemath/sage-4.5.3/local/lib/python2.6/site-packages/sage/matrix/matrix2.so in sage.matrix.matrix2.Matrix.determinant (sage/matrix/matrix2.c:6942)()
   1092         # TODO: Find a reasonable cutoff point.  (This is field specific, but
   1093         # seems to be quite large for Q[x].)
-> 1094         if (R.is_field() and R.is_exact() and algorithm is None) or (algorithm == "hessenberg"):
   1095             try:
   1096                 c = self.charpoly('x', algorithm="hessenberg")[0]

/opt/sagemath/sage-4.5.3/local/lib/python2.6/site-packages/sage/rings/quotient_ring.pyc in is_field(self, proof)
    420         """
    421         if proof:
--> 422             return self.defining_ideal().is_maximal()
    423         else:
    424             try:

/opt/sagemath/sage-4.5.3/local/lib/python2.6/site-packages/sage/rings/ideal.pyc in is_maximal(self)
    559             NotImplementedError
    560         """
--> 561         raise NotImplementedError
    562 
    563     def is_primary(self, P=None):

NotImplementedError: 

A simple (maybe dirty) solution might be to add a try at the beginning of the method and an except that use a basic division-free formula in case of error.

I am not an algebraist, so i prefer not to fix this bug by myself.

CC: @sagetrac-sage-combinat @sagetrac-spancratz @aghitza @nthiery @jasongrout

Component: commutative algebra

Keywords: determinant, ring, ideal

Author: Thierry Monteil

Reviewer: Mike Hansen, Sébastien Labbé

Merged: sage-4.6.1.alpha0

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

[nthiery]

comment:1

Hi Thierry!

I already stumbled into similar annoyances. Adding an exception handling in order to fall back to the default division-free algorithm when Sage can't determine if the base ring is a field sounds very reasonable to me!

Cheers, Nicolas

[edd8e884-f507-429a-b577-5d554626c0fe]

comment:3

ok, i will make a very basic (purely logical (a simple try: followed by except:pass works at least for my example)) patch for it, and add my example as a test for the determinant method.

[edd8e884-f507-429a-b577-5d554626c0fe]

comment:4

Here is a patch that fixes the problem. I also took the opportunity of this bugfix to accelerate the method a bit when the division-free algorithm df is explicitly chosen. For this, i changed the ordering in the and boolean expression: the method first test if the chosen algorithm in None before checking if the ring is a field.

[mwhansen]

comment:6

We don't want to put a bare "except:" clause since that will catch things like KeyboardInterrupt exceptions.  We should be explicit about the type of exceptions that we want to catch.

[mwhansen]

Reviewer: Mike Hansen

[edd8e884-f507-429a-b577-5d554626c0fe]

comment:7

Ok, since a general robust (division-free) algorithm follows, i thought any error could be corrected, i didn't thought to the KeyboardInterrupt exception. Thanks for the hint! I will be explicit about the type of exception.

[edd8e884-f507-429a-b577-5d554626c0fe]

Attachment: trac_10063-determinant_not_computed_in_some_rings_bugfix-tm.patch.gz

Tested on 4.5.3

[edd8e884-f507-429a-b577-5d554626c0fe]

comment:8

Only the NotImplementedError is catched now, so that the bug i found is still fixed.

Other types of errors (possibly) arising from what is behind is_field or is_exact methods will not be catched (which is bad for the user, but good for debugging).

Since both the patch and the modification are small, i replaced the existing patch by another one with the same name.

[edd8e884-f507-429a-b577-5d554626c0fe]

Author: Thierry Monteil

[seblabbe]

Changed reviewer from Mike Hansen to Mike Hansen, Sébastien Labbé

[seblabbe]

comment:9

Dear Thierry,

You included too much code inside of the try statement. You are "trying" a lot more than what is necessary, that is, R.is_field() in this case. From the PEP0008, cited in the Python coding conventions of the Sage Developpers Guide, one can read :

    - Additionally, for all try/except clauses, limit the 'try' clause
      to the absolute minimum amount of code necessary.  Again, this
      avoids masking bugs.

      Yes:

          try:
              value = collection[key]
          except KeyError:
              return key_not_found(key)
          else:
              return handle_value(value)

      No:

          try:
              # Too broad!
              return handle_value(collection[key])
          except KeyError:
              # Will also catch KeyError raised by handle_value()
              return key_not_found(key)

Hence, you could do something like :

    try:
        R_is_field = R.is_field()
    except NotImplementedError:
        ...

    if (algorithm is None and R_is_field and blablablal) or ...:
        ...

Finally, you need an empty line after the :: for the documentation to build properly.

Hence, I change the status of this ticket to needs work.

Sébastien

[edd8e884-f507-429a-b577-5d554626c0fe]

comment:10

Hi Sebastien,

thanks for your comments. As you know, i am not a computer scientist nor a programmer, so i do not know much about design principles, but if i had to code a determinant function from scratch, knowing that:

• the determinant is an fundamental function
• there exist a division-free algorithm that works in any case
• there are some fast/fancy/pretty/lounge/fieldic... optimizations that are very good in some particular cases but not completely finished, i will have written a code like (even before meeting any bug):

try:
    all the kind of optimized algorithms
except StandardError:
    pass #or send an automatic bug report if the error seems new.
the division-free algorithm that always works

I do not feel it is a dirty code, since its structure shows that one algorithm works everywhere and the other are fragile optimizations.

Making patches that only repair the bugs found will make a code whose architecture depend on the history of the bugs discovery, not necessarily a readable code (why is there a test for the is_field method and not for is_exact,...).

The aim of the Python convention is to avoid masking bugs, but if the exception is NotImplementedError, then the "bug" is already known, since it means that someone wrote an empty method which raises this NotImplementedError.

Note also that the Python convention is a bit slower since i have to test the is_field method even if the algorithm is not None. Should i replace:

try:
    R_is_field_attempt = R.is_field()
except NotImplementedError:
    R_is_field_attempt = False

by:

if algorithm is None:
    try:
        R_is_field_attempt = R.is_field()
    except NotImplementedError:
        R_is_field_attempt = False

in order to skip the test when another algorithm is is called?

Anyway, here is a patch that implements your recommendations (i rename it to keep the previous version on the trac server and follow the discussion).

I fixed the problem of the documentation as well (i took the opportunity to fix the existing test, which i copied).

Should we think to clean the whole determinant code?

[edd8e884-f507-429a-b577-5d554626c0fe]

Tested on 4.5.3

[edd8e884-f507-429a-b577-5d554626c0fe]

comment:11

Attachment: trac_10063-determinant_not_computed_in_some_rings_bugfix_attempt_3-tm.patch.gz

In this version, the is_field method is tested only if algorithm is not None. If another algorithm is chosen, the test is skipped to save some time.

[edd8e884-f507-429a-b577-5d554626c0fe]

Tested on 4.5.3

[edd8e884-f507-429a-b577-5d554626c0fe]

comment:12

Attachment: trac_10063-determinant_not_computed_in_some_rings_bugfix_attempt_4-tm.patch.gz

A typo in the previous comment: the is_field method is tested only if algorithm is None.

[seblabbe]

comment:14

Should we think to clean the whole determinant code?

Maybe. I don't know. But, I don't feel I am the determinant expert that could rethink/refactor/improve that code.

Although, I reviewed your most recent patch (trac_10063-determinant_not_computed_in_some_rings_bugfix_attempt_4-tm.patch).

I was able to reproduce the problem on my computer. The problem is indeed fixed by the patch. All test passed in the repository sage/matrix. Documentation builds fine.

The problem originated from a NotImplementedError when computing R.is_field() in the case where algorithm=None (see below). So your new code instead tries to compute R.is_field() and if R.is_field() raises a NotImplementedError you consider this as R.is_field() == False. This computation is done only if algorithm is None so that the method doesn't get slower. You also followed my earlier advise, that is, the try statement only tries what is necessary.

Hence, from the knowledge I have, I am OK with giving a positive review to this ticket. Maybe Mike Hanson or Jason have comments, so I wait 24 hours, and then I will change the status to positive review.

Sébastien

PS : The is_field method is not implemented for the following ring (should we open a ticket for it?):

sage: A = GF(2)['x,y,z']
sage: A.inject_variables()
Defining x, y, z
sage: R = A.quotient(x^2 + 1).quotient(y^2 + 1).quotient(z^2 + 1) 
sage: R.is_field()
Traceback (most recent call last):

NotImplementedError: 

[jdemeyer]

comment:15

I did a full ptestlong against sage-4.6.alpha3, works fine. Positive review considering the comments of Sébastien Labbé.

[jdemeyer]

Merged: sage-4.6.1.alpha0