Take the 2-minute tour ×
Mathematica Stack Exchange is a question and answer site for users of Mathematica. It's 100% free, no registration required.

Is there a way to calculate the Rational Canonical Form of an $n\times n$ integer matrix using Mathematica?

I have been perusing the documentation and web, but nothing so far.

share|improve this question

1 Answer 1

up vote 6 down vote accepted

I found the problem. What is needed is to factor the CharacteristicPolynomial by the MatrixMinimalPolynomial continously in order to obtain all minimal polynomials.

I did not know how to do that using one of the polynomial functions in Mathematica, so I wrote a function to do it.

Result

This shows test matrices. The left side is the matrix, and the right side is the matrix rational form. This function called rationalMatrixForm requires 2 helper functions: companionMatrix and MatrixMinimalPolynomial

This is not optimal function by any means (uses AppendTo for example), and I am sure it can be improved, but first I wanted to make sure it is correct. I tested it on number of matrices and the results agrees with wikipdia example and mupad results. If you find a bug in it, please let me know.

rationalMatrixForm[{{2, -2, 14}, {0, 3, -7}, {0, 0, 2}}]

Mathematica graphics

Grid[{MatrixForm[#], MatrixForm@rationalMatrixForm[#]} & /@ tests, Frame -> All]

Mathematica graphics

Code

CompanionMatrix[p_, x_] := Module[
  {n, w = CoefficientList[p, x]},
  w = -w/Last[w];
  n = Length[w] - 1;
  SparseArray[{{i_, n} :> w[[i]], {i_, j_} /; i == j + 1 -> 1}, {n, 
    n}]]

MatrixMinimalPolynomial[a_List?MatrixQ, x_] := 
 Module[{i, n = 1, qu = {}, 
   mnm = {Flatten[IdentityMatrix[Length[a]]]}}, 
  While[Length[qu] == 0, AppendTo[mnm, Flatten[MatrixPower[a, n]]];
   qu = NullSpace[Transpose[mnm]];
   n++];
  First[qu].Table[x^i, {i, 0, n - 1}]]

and

rationalMatrixForm[a_?(MatrixQ[#, NumericQ] &)] := Module[(*version 8/24/13 2PM*)
  {p, q, min, c = {}, moreFactors = True, z, x},

  p = CharacteristicPolynomial[a, x];
  min = MatrixMinimalPolynomial[a, x];

  While[moreFactors,
   q = PolynomialQuotient[p, min, x];

   If[q === 0,

    moreFactors = False;
    If[Not[FreeQ[p, x]],
     z = CompanionMatrix[p, x];
     AppendTo[c, z]
     ],

    z = CompanionMatrix[min, x];
    AppendTo[c, z];
    p = q

    ] (* if *)

   ]; (*end WHILE more factorization needed*)

  SparseArray[Band[{1, 1}] -> c]
  ]

test matrices

 tests = {
   {{2, -2, 14}, {0, 3, -7}, {0, 0, 2}},
   {{3, 4, 0}, {-1, -3, -2}, {1, 2, 1}},
   {{-2, -1, -2, -1, 1, 0}, {-2, -1, -2, -1, 1, 1}, {2, 1, 2, 1, 0, 
     0}, {2, 1, 0, 1, -3, -1}, {-2, 0, -2, 0, 0, 0}, {2, -2, 0, 0, 0, 
     0}},
   {{0, -4, 85}, {1, 4, -30}, {0, 0, 3}},
   {{2, -2, 14, 5, 6, 7}, {0, 3, -7, 9, 20, 33}, {0, 0, 2, 9, 0, 
     3}, {2, -2, 14, 5, -8, 7}, {2, 2, 14, 23, 6, 7}, {2, 2, 14, 23, 
     6, 70}}
   };
share|improve this answer
    
Very nice, however produces an incorrect result for mat = {{2, -2, 14}, {0, 3, -7}, {0, 0, 2}}, but correct for mat = {{0, -4, 85}, {1, 4, -30}, {0, 0, 3}}. These have same characteristic polynomial, but different rational canonical forms. Regards –  Amzoti Aug 23 '13 at 5:08
    
@Nasser: Are you saying no way to easily fix? –  Amzoti Aug 23 '13 at 5:19
    
@Nasser: It might have to do with the fact that the first matrix in my comment has a minimal polynomial $(x-2)(x-3)$ while the second has $(x-2)^2(x-3)$, even though they both have the same characteristic polynomial. –  Amzoti Aug 23 '13 at 5:23
    
@Nasser: Is there a way to print $P$, as in $P^{-1}AP=r$. In other words, you found $r$, but is it possible to also provide $P$? –  Amzoti Aug 23 '13 at 5:25
    
@Amzoti I get this !Mathematica graphics –  Nasser Aug 23 '13 at 5:45

Your Answer

 
discard

By posting your answer, you agree to the privacy policy and terms of service.

Not the answer you're looking for? Browse other questions tagged or ask your own question.