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: |
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
comment:2 Changed 11 years ago by
- Cc spancratz added
comment:3 Changed 11 years ago by
- 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
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
- Cc nthiery added
- Status changed from new to needs_review
comment:6 Changed 11 years ago by
- 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
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.
comment:8 Changed 11 years ago by
- 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
- 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: ↓ 14 Changed 11 years ago by
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
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
- 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
- Cc jason added
comment:14 in reply to: ↑ 10 Changed 11 years ago by
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
- 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
- Merged in set to sage-4.6.1.alpha0
- Resolution set to fixed
- Status changed from positive_review to closed
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,