Opened 11 years ago

Closed 11 years ago

#10063 closed defect (fixed)

Some determinants can not be computed

Reported by: tmonteil Owned by: tmonteil
Priority: critical Milestone: sage-4.6.1
Component: commutative algebra Keywords: determinant, ring, ideal
Cc: sage-combinat, spancratz, AlexGhitza, nthiery, jason Merged in: sage-4.6.1.alpha0
Authors: Thierry Monteil Reviewers: Mike Hansen, Sébastien Labbé
Report Upstream: N/A Work issues:
Branch: Commit:
Dependencies: Stopgaps:

Status badges

Description

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.

Attachments (3)

Change History (19)

comment:1 Changed 11 years ago by nthiery

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

comment:2 Changed 11 years ago by tmonteil

  • Cc spancratz added

comment:3 Changed 11 years ago by tmonteil

  • Cc AlexGhitza added
  • Component changed from algebra to commutative algebra
  • Owner changed from AlexGhitza to tmonteil

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.

comment:4 Changed 11 years ago by tmonteil

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.

comment:5 Changed 11 years ago by tmonteil

  • Cc nthiery added
  • Status changed from new to needs_review

comment:6 Changed 11 years ago by mhansen

  • Reviewers set to Mike Hansen
  • Status changed from needs_review to needs_work

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.

comment:7 Changed 11 years ago by tmonteil

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.

Changed 11 years ago by tmonteil

Tested on 4.5.3

comment:8 Changed 11 years ago by tmonteil

  • Authors set to Thierry Monteil
  • Status changed from needs_work to needs_review

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.

comment:9 Changed 11 years ago by slabbe

  • Reviewers changed from Mike Hansen to Mike Hansen, Sébastien Labbé
  • Status changed from needs_review to needs_work

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

comment:10 follow-up: Changed 11 years ago by tmonteil

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?

comment:11 Changed 11 years ago by tmonteil

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.

comment:12 Changed 11 years ago by tmonteil

  • Status changed from needs_work to needs_review

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

comment:13 Changed 11 years ago by jason

  • Cc jason added

comment:14 in reply to: ↑ 10 Changed 11 years ago by slabbe

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: 

comment:15 Changed 11 years ago by jdemeyer

  • Milestone set to sage-4.6.1
  • Status changed from needs_review to positive_review

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

comment:16 Changed 11 years ago by jdemeyer

  • Merged in set to sage-4.6.1.alpha0
  • Resolution set to fixed
  • Status changed from positive_review to closed
Note: See TracTickets for help on using tickets.