jordan_form with transformation=true returns non-invertible transformation

[syazdani77]

The following code returns an incorrect result:

mm=Matrix(GF(2),[[1,0,1,0,0,0,1],[1,0,0,1,1,1,0],[1,1,0,1,1,1,1],[1,1,1,0,1,1,1],[1,1,1,0,0,1,0],[1,1,1,0,1,0,0],[1,1,1,1,1,1,0]])
_,S = mm.jordan_form(transformation=True)
S.rank()

S should be invertible, so the rank should be 7, but the rank of the above is 5.

CC: @jasongrout

Component: linear algebra

Keywords: jordan_form, transformation

Author: Sebastian Pancratz

Reviewer: Rob Beezer, Minh Van Nguyen

Merged: sage-4.3.3.alpha0

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

[jasongrout]

comment:3

According to Magma:

u,v:=JordanForm(Matrix(GF(2),7,7,StringToIntegerSequence("1 0 1 0 0 0 1 1 0 0 1 1 1 0 1 1 0 1 1 1 1 1 1 1 0 1 1 1 1 1 1 0 0 1 0 1 1 1 0 1 0 0 1 1 1 1 1 1 0")));

print "jordan form:";
u;
print "transformation";
v;

gives

Magma V2.15-14    Tue Nov 24 2009 10:36:58    [Seed = 4156006409]
   -------------------------------------

jordan form:
[1 0 0 0 0 0 0]
[0 1 1 0 0 0 0]
[0 0 1 0 0 0 0]
[0 0 0 1 1 0 0]
[0 0 0 0 1 0 0]
[0 0 0 0 0 1 1]
[0 0 0 0 0 0 1]
transformation
[0 0 0 0 1 1 0]
[0 0 0 0 0 0 1]
[1 1 1 1 1 1 1]
[0 0 0 0 0 1 0]
[1 1 1 0 1 1 0]
[1 0 1 0 0 0 0]
[1 1 0 1 1 1 0]

Total time: 0.180 seconds, Total memory usage: 7.27MB

[jasongrout]

comment:4

Or more verbosely:

m:=Matrix(GF(2),7,7,StringToIntegerSequence("1 0 1 0 0 0 1 1 0 0 1 1 1 0 1 1 0 1 1 1 1 1 1 1 0 1 1 1 1 1 1 0 0 1 0 1 1 1 0 1 0 0 1 1 1 1 1 1 0"));
jor,trans:=JordanForm(m);

print "jordan form:";
jor;
print "transformation";
trans;
print "m*transformation";
trans*m;
print "transformation*jordan form";
jor*trans; 

gives

Magma V2.15-14    Tue Nov 24 2009 10:41:01    [Seed = 3482294566]
   -------------------------------------

jordan form:
[1 0 0 0 0 0 0]
[0 1 1 0 0 0 0]
[0 0 1 0 0 0 0]
[0 0 0 1 1 0 0]
[0 0 0 0 1 0 0]
[0 0 0 0 0 1 1]
[0 0 0 0 0 0 1]
transformation
[0 0 0 0 1 1 0]
[0 0 0 0 0 0 1]
[1 1 1 1 1 1 1]
[0 0 0 0 0 1 0]
[1 1 1 0 1 1 0]
[1 0 1 0 0 0 0]
[1 1 0 1 1 1 0]
m*transformation
[0 0 0 0 1 1 0]
[1 1 1 1 1 1 0]
[1 1 1 1 1 1 1]
[1 1 1 0 1 0 0]
[1 1 1 0 1 1 0]
[0 1 1 1 1 1 0]
[1 1 0 1 1 1 0]
transformation*jordan form
[0 0 0 0 1 1 0]
[1 1 1 1 1 1 0]
[1 1 1 1 1 1 1]
[1 1 1 0 1 0 0]
[1 1 1 0 1 1 0]
[0 1 1 1 1 1 0]
[1 1 0 1 1 1 0]

Total time: 0.170 seconds, Total memory usage: 7.27MB

[sagetrac-spancratz]

comment:5

Attachment: trac6942_jordan.patch.gz

The above patch fixes this problem, and also resolves ticket #6932.

[sagetrac-spancratz]

Additional doctests

[sagetrac-spancratz]

comment:6

Attachment: trac6942_jordan_tests.patch.gz

The second patch adds additional doctests. There are three of them, all for 10 by 10 matrices over the rationals and with only one eigenvalue. The Jordan blocks are of sizes (a) 3,3,3,1, (b) 3,3,2,2, and (c) 3,2,2,2,1.

Sebastian

[rbeezer]

comment:7

This is looking pretty good. But I'll have to spend some more time with it. Until then, here's another 10x10 matrix with a nice Jordan form and nearly no fractions in the transformation matrix.

matrix(QQ, [
[-6,  9,  -7,  -5,  5,  12,  -22,  14,  8,  21 ],
[ -3,  5,  -3,  -1,  2,  7,  -12,  9,  1,  12 ],
[8,  -9,  8,  6,  0,  -14,  25,  -13,  -4,  -26 ],
[-7,  9,  -7,  -5,  0,  13,  -23,  13,  2,  24 ],
[0,  -1,  0,  -1,  -3,  -2,  3,  -4,  -2,  -3 ],
[3,  2,  1,  2,  9,  -1,  1,  5,  5,  -5 ],
[-1,  3,  -3,  -2,  4,  3,  -6,  4,  4,  3 ],
[3,  -4,  3,  2,  1,  -5,  9,  -5,  1,  -9 ],
[0,  2,  0,  0,  2,  2,  -4,  4,  2,  4 ],
[-4,  4,  -5,  -4,  -1,  6,  -11,  4,  1, 10]
])

[rbeezer]

Attachment: trac_6942-reviewer.patch.gz

One-character reviewer doctest fix

[rbeezer]

comment:8

Hi Sebastian,

Very nice!

1. One line in a list indented one-too-many characters. Fix in reviewer patch for your convenience. Accept or not.

2. Most of the linear algebra matrix decompositions return all the pieces, all the time. I've always thought it inconsistent Jordan form has this transformation switch, sometimes returning one matrix, other times two. Now would be a good time to change this behavior. But I see that the form is almost a combinatorial routine (once you have the spectrum), while the basis vectors require a whole lot more work, so maybe best not to compute it unless it is asked for?

3. The list Vsmall+Y ends up in a span() in the "vector_in_difference" routine. The span has an already_echelonized switch, so I spent a lot of time trying to decide how (or if it was possible) to somehow "keep" Vsmall+Y echelonized, since Vsmall will start that way. Couldn't see a reasonable way to do anything too clever, though.

4. Do you want the 10x10 matrix above for the doctests? I'd be happy to package it up with the one-character patch. ;-)

Other than the doctest fix, this is ready to go. If you want to accept the doctest fix, then go ahead and mark this as "positive review" - everything else is just for your consideration.

Great to see all your good work since the summer, including this.

Rob