ExampleB.mw

> with(linalg);

[BlockDiagonal, GramSchmidt, JordanBlock, LUdecomp, QRdecomp, Wronskian, addcol, addrow, adj, adjoint, angle, augment, backsub, band, basis, bezout, blockmatrix, charmat, charpoly, cholesky, col, cold...[BlockDiagonal, GramSchmidt, JordanBlock, LUdecomp, QRdecomp, Wronskian, addcol, addrow, adj, adjoint, angle, augment, backsub, band, basis, bezout, blockmatrix, charmat, charpoly, cholesky, col, cold...[BlockDiagonal, GramSchmidt, JordanBlock, LUdecomp, QRdecomp, Wronskian, addcol, addrow, adj, adjoint, angle, augment, backsub, band, basis, bezout, blockmatrix, charmat, charpoly, cholesky, col, cold...[BlockDiagonal, GramSchmidt, JordanBlock, LUdecomp, QRdecomp, Wronskian, addcol, addrow, adj, adjoint, angle, augment, backsub, band, basis, bezout, blockmatrix, charmat, charpoly, cholesky, col, cold...[BlockDiagonal, GramSchmidt, JordanBlock, LUdecomp, QRdecomp, Wronskian, addcol, addrow, adj, adjoint, angle, augment, backsub, band, basis, bezout, blockmatrix, charmat, charpoly, cholesky, col, cold...[BlockDiagonal, GramSchmidt, JordanBlock, LUdecomp, QRdecomp, Wronskian, addcol, addrow, adj, adjoint, angle, augment, backsub, band, basis, bezout, blockmatrix, charmat, charpoly, cholesky, col, cold...[BlockDiagonal, GramSchmidt, JordanBlock, LUdecomp, QRdecomp, Wronskian, addcol, addrow, adj, adjoint, angle, augment, backsub, band, basis, bezout, blockmatrix, charmat, charpoly, cholesky, col, cold...

> A:=matrix([[1, -1, 0, 0, 0, 0], [4, 5, 0, 0, 0, 0], [0, 0, 0, -1, 0, 0], [0, 0, 1, 0, 0, 0], [0, 0, 0, 0, 0, 1], [0, 0, 0, 0, -1, 0]]);

A := matrix([[1, -1, 0, 0, 0, 0], [4, 5, 0, 0, 0, 0], [0, 0, 0, -1, 0, 0], [0, 0, 1, 0, 0, 0], [0, 0, 0, 0, 0, 1], [0, 0, 0, 0, -1, 0]])

> factor(charpoly(A,t));factor(minpoly(A, t));

(t-3)^2*(t^2+1)^2

(t^2+1)*(t-3)^2

> #The matrix is not diagonalizable. Let us do the primary decomposition theorem in this case. We have C^6 = W1+W2+W3 (direct sum) where:

> W1:=kernel(A - I);

W1 := {vector([0, 0, I, 1, 0, 0]), vector([0, 0, 0, 0, -I, 1])}

> #Obviously, on W1 the map T is diagonal, equal to multiplication by I

> W2:=kernel(A+I);

W2 := {vector([0, 0, 0, 0, I, 1]), vector([0, 0, -I, 1, 0, 0])}

> #Obviously, on W2 the map T is diagonal, equal to multiplication by -I

> W3:=kernel((A-3)^2);

W3 := {vector([1, 0, 0, 0, 0, 0]), vector([0, 1, 0, 0, 0, 0])}

> w:=multiply(A-3, W3[1]);

w := vector([-2, 4, 0, 0, 0, 0])

> #Take as basis for W3 the vectors {w, W3[1]}. (The choice of basis is motivated by the theory of the Jordan canonical form, to be studies soon.) In this basis the transformation is given by

> matrix([[3, 1], [0, 3]]);

matrix([[3, 1], [0, 3]])

> #Let us check that.

> M:=transpose(matrix([W1[1], W1[2], W2[1], W2[2],w, W3[1]]));

M := matrix([[0, 0, 0, 0, -2, 1], [0, 0, 0, 0, 4, 0], [I, 0, 0, -I, 0, 0], [1, 0, 0, 1, 0, 0], [0, -I, I, 0, 0, 0], [0, 1, 1, 0, 0, 0]])

> multiply(inverse(M), A, M);

matrix([[I, 0, 0, 0, 0, 0], [0, I, 0, 0, 0, 0], [0, 0, -I, 0, 0, 0], [0, 0, 0, -I, 0, 0], [0, 0, 0, 0, 3, 1], [0, 0, 0, 0, 0, 3]])

>