Changes between Initial Version and Version 2 of Ticket #21579


Ignore:
Timestamp:
09/24/16 23:46:57 (3 years ago)
Author:
bober
Comment:

Legend:

Unmodified
Added
Removed
Modified
  • Ticket #21579

    • Property Component changed from modular forms to linear algebra
    • Property Summary changed from Modular symbols computations are inconsistent and occasionally wrong to Errors calculating characteristic polynomials of rational matrices
  • Ticket #21579 – Description

    initial v2  
     1I'm leaving the description below because that's where I first ran into problems, but this seems like a much more basic issue. I suspect some sort of memory bug, because I'm getting different results based on the order in which I execute functions.
     2
     3This could just be a problem with polynomials. I haven't quite diagnosed it yet.
     4
     5{{{#!python
     6A = Matrix([[ 10/7,   3/4,     1,   8/7,     3,  -1/7,    -1,    -1,   7/4,    -2,     3,   1/2,  -8/3,   3/8, -10/7,    -3,  -1/2,  -4/7,   9/8,    -4],
     7            [  3/2,   3/5, -3/10,  -4/3,  -2/7,  -3/2,    -1,     3,     2,     0,   4/3,    -2,   4/9,  -1/2,  -1/5,     1,   7/8,     1,  -5/2,  10/9],
     8            [ 9/10,     7,     3,     0,  -3/2,   1/3,    -5,   8/3,   6/7,     0,   7/4,  -4/3,  -1/2,  -1/6,  -1/4,  -3/7,  -6/5,  -1/9,  -1/3,  -8/5],
     9            [    8,   1/2,  -8/7,     1,  -3/8,    -7, -10/7,     1,   4/3,     0,     0,     4,    -3,  -1/4,     5,    -1,    -1,     0,  -2/7, -10/3],
     10            [  8/3,     5,  -2/3,   9/4,  -9/8,  -1/3,  -2/3,    -1,    -8,   4/3,  -3/4,  -5/4,  -2/3,  -1/2,     5,  -1/3,    10,   1/2,  -5/3,   9/2],
     11            [-10/9,    -3,    -7,   1/2,   5/9,     0,     3,    -1,   9/4,  -3/8,   -10,     0,     6,  -1/6,    -2,  -9/8,   2/3,     6,     0,     0],
     12            [ -2/5,   6/5,     0,     0,  -6/7, -9/10,  -2/7, -1/10,  -1/3,   7/2,   1/3,  -3/2,   7/2,   3/2,   1/3,     0,   7/2,  -2/9,   2/3,  -1/8],
     13            [  1/3,     3,  10/3,     1,  -1/2,  -3/8,    -5,   3/2,   1/3,     0,  -1/2,   4/9,  -7/5,     0,  -3/2,  9/10,    -1,   1/2,    -1, -9/10],
     14            [   -1,   3/8,     1,  -2/5,   1/3,     1,  -3/2,    -4,    -1,    -4,     1,     1, -9/10,     0,     1,    -1,  -7/6,     0, -10/9,     0],
     15            [   -6,     0,  -8/5,  -9/7,     7,  -1/3,  -8/9,    -8,   7/9,   1/2,     1,    -1,     1,  -2/7,   8/3,  -1/3,  -4/5,   2/3,    -5,  -1/7],
     16            [  5/3,   1/3,  -9/7,    -2,   4/5,  -5/7,   1/4,  -5/8,    -1,   5/7,   3/2,  -3/2,   2/5,     0,  -1/2,   1/4,     0,     3,     1,     0],
     17            [ -7/8,    -1,     6,  -9/7,   3/5,   1/2,     2,   7/2,     1,    -6,  -7/2,   1/9,    -1,  -7/5,    -5,    -5,  -2/7,  -5/4,    -3,    -3],
     18            [ -7/3,   1/5,    -2,    -8,   1/2,     2,     0,   2/7,  -1/4,     0,     0,  -4/7,   1/3,    -1,    -5,   3/2,    10,  3/10,   1/3,   3/4],
     19            [  8/5,     0,  -4/5,     0,   1/4,  -8/5,  -4/9,     5,  -7/9,   1/9,  -7/8,  -5/8,  -1/4,    10,    -1,  -8/9,   1/2,   1/6,  -8/7,   1/7],
     20            [ -2/5,  -5/4,     2,   1/2,  -3/5,  -5/4, -10/9,   1/2,   3/5,   9/7,  -5/6,  -3/7,   7/6,  -4/5,  -9/5,  -7/3,  -1/3,     1,     5,  -2/5],
     21            [ -1/2,     2,     0,     3,   5/6,     0,    -3,  -1/3,   1/2,     1,   2/9,  -5/7,    -1,    -8,   1/2,    10,    -5,   1/2,     1,     5],
     22            [  8/3,    -2,    -2,  -4/3,   9/2,  -4/5,  10/7,  -1/2,     9,  -1/3, -10/7,    -4,     7, -1/10,  10/9,  -7/6,   2/7,   2/7,     0,     8],
     23            [ -2/3,   1/2,     4,  -3/5,  -2/3,  -2/7,  9/10,  -1/4,  -5/3,     6,  -1/3,    -1,  -5/4,  -4/3,    -6,  -8/7,    -2,     1,   5/6,   1/3],
     24            [-10/3,     5,     0,   2/5,  -4/3,  -1/2,  -1/2,  -9/4,     6,   1/5,  -1/9,    -2,  -6/5,    -2,     4,  -3/4,    -1,   1/4,  -3/8,   1/2],
     25            [ -7/6, -7/10,     0,   1/5,   3/4,  -3/5,    -2,  -3/2,    -1,   6/7,  -5/7,     0,   9/5,  -1/2,    -1,    -2,    -7,     0,  -4/5,  -9/2]])
     26
     27k = int(sys.argv[1])
     28
     29if k == 0:
     30    sage.matrix.matrix_integer_dense.USE_LINBOX_POLY = True # True is the default
     31    A._clear_cache()
     32    f = A.charpoly() # should be using linbox
     33
     34    sage.matrix.matrix_integer_dense.USE_LINBOX_POLY = False # Setting this False should force
     35                                                             # integer charpoly to use flint when
     36                                                             # algorithm=='linbox'. (Rational
     37                                                             # matrices don't have a flint option.)
     38    A._clear_cache()
     39    g = A.charpoly() # should use flint, I think
     40
     41    A._clear_cache()
     42    h = A.charpoly(algorithm='generic')
     43elif k == 1:
     44    sage.matrix.matrix_integer_dense.USE_LINBOX_POLY = False
     45    A._clear_cache()
     46    g = A.charpoly() # should use flint, I think
     47
     48    sage.matrix.matrix_integer_dense.USE_LINBOX_POLY = True
     49    A._clear_cache()
     50    f = A.charpoly() # should be using linbox
     51
     52    A._clear_cache()
     53    h = A.charpoly(algorithm='generic')
     54
     55
     56if k == 0 or k == 1:
     57    print 'flint == linbox:', f == g
     58    print 'flint == generic:', g == h
     59    print 'linbox == generic:', f == h
     60
     61    print 'linbox (maybe) correct? ', f(A) == 0
     62    print 'flint (maybe) correct?', g(A) == 0
     63    print 'generic (maybe) correct?', h(A) == 0
     64
     65sage.matrix.matrix_integer_dense.USE_LINBOX_POLY = True
     66
     67if k in [0,1,2]:
     68    # let's try computing the charpoly the way that
     69    # B.charpoly() does it to see what goes wrong.
     70    B, denom = A._clear_denom()
     71    B._clear_cache()
     72    f1 = B.charpoly(algorithm='linbox')
     73    B._clear_cache()
     74    g1 = B.charpoly(algorithm='flint')
     75    print 'flint == linbox', f1 == g1
     76    print 'linbox (maybe) correct?', f1(B) == 0
     77    print 'flint (maybe) correct?', g1(B) == 0
     78
     79    x = f1.parent().gen()
     80    F = f1(x * denom) * (1 / (denom**f1.degree()))
     81
     82    x = g1.parent().gen()
     83    G = g1(x * denom) * (1 / (denom**g1.degree()))
     84
     85    print 'flint == linbox', F == G
     86    if x == 0 or x == 1:
     87        print 'linbox == linbox', f == F
     88        print 'flint == flint?', g == G
     89    print '2nd linbox (maybe) correct?', F(A) == 0
     90    print '2nd flint (maybe) correct?', G(A) == 0
     91
     92if k == 3:
     93    # let's try computing the charpoly the way that
     94    # B.charpoly() does it to see what goes wrong.
     95    B, denom = A._clear_denom()
     96    f1 = B.charpoly(algorithm='linbox')
     97    print 'linbox (maybe) correct?', f1(B) == 0
     98
     99    x = f1.parent().gen()
     100    F = f1(x * denom) * (1 / (denom**f1.degree()))
     101
     102    print '2nd linbox (maybe) correct?', F(A) == 0
     103
     104if k == 4:
     105    # let's try computing the charpoly the way that
     106    # B.charpoly() does it to see what goes wrong.
     107    B, denom = A._clear_denom()
     108    g1 = B.charpoly(algorithm='flint')
     109    print 'flint (maybe) correct?', g1(B) == 0
     110
     111    x = g1.parent().gen()
     112    G = g1(x * denom) * (1 / (denom**g1.degree()))
     113
     114    print '2nd flint (maybe) correct?', G(A) == 0
     115
     116}}}
     117
     118On the most recent version of the develop branch, here's output that I get consistently:
     119
     120{{{
     121jb12407@lmfdb5:~/sage-bug$ sage bug.sage 0
     122flint == linbox: False
     123flint == generic: False
     124linbox == generic: True
     125linbox (maybe) correct?  True
     126flint (maybe) correct? False
     127generic (maybe) correct? True
     128flint == linbox True
     129linbox (maybe) correct? True
     130flint (maybe) correct? True
     131flint == linbox True
     1322nd linbox (maybe) correct? False
     1332nd flint (maybe) correct? False
     134
     135jb12407@lmfdb5:~/sage-bug$ sage bug.sage 1
     136flint == linbox: True
     137flint == generic: False
     138linbox == generic: False
     139linbox (maybe) correct?  False
     140flint (maybe) correct? False
     141generic (maybe) correct? True
     142flint == linbox True
     143linbox (maybe) correct? True
     144flint (maybe) correct? True
     145flint == linbox True
     1462nd linbox (maybe) correct? False
     1472nd flint (maybe) correct? False
     148
     149jb12407@lmfdb5:~/sage-bug$ sage bug.sage 2
     150flint == linbox True
     151linbox (maybe) correct? True
     152flint (maybe) correct? True
     153flint == linbox True
     1542nd linbox (maybe) correct? False
     1552nd flint (maybe) correct? False
     156
     157jb12407@lmfdb5:~/sage-bug$ sage bug.sage 3
     158linbox (maybe) correct? True
     1592nd linbox (maybe) correct? True
     160
     161jb12407@lmfdb5:~/sage-bug$ sage bug.sage 4
     162flint (maybe) correct? True
     1632nd flint (maybe) correct? False
     164}}}
     165
     166All those `False`s ought to be `True`. Option 3 does get the right answer, but I was getting (less directly repeatable) bugs just by calling charpoly() repeatedly on the same (rational) matrix, which should have been using only linbox.
     167
     168In another (nonrepeated, interactive) test, I wound up in a situation that was something like
     169{{{
     170sage: f = A.charpoly(algorithm='flint')
     171sage: g = A.charpoly(algorithm='linbox')
     172sage: f == g
     173True
     174sage: f(A) == 0
     175True
     176sage: g(A) == 0
     177False
     178sage: print f
     179...
     180sage: print g
     181...
     182# yes, the certainly looked the same!
     183}}}
     184(I might have the True and False mixed up, though.)
     185
     186---
     187
     188Old Description:
     189
    1190I ran the following code for a bit
    2191