In search of beautiful matrices

It is often desirable to have a non-singular  2×2 matrix with small integer entries and with distinct integer eigenvalues. Note that no rational number can occur as an eigenvalue of a matrix with integer entries.

In this problem we want to find all nontriangular (that is a*b ≠ 0) matrices

(({{a, b}, {c, d}}))

such that

a, b, c, 
      d ∈ {-9, -8,  ... , -1, 0, 1, 2,  ...  , 9}

has either real eigenvalues in the set  {-9, -8,  ... , -1, 0, 1, 2,  ...  , 
        9} or the real and imaginary parts of the nonreal eigenvalues are in the set  {-9, -8,  ... , -1, 0, 1, 2,  ...  , 
        9}.

Consider three different cases: (1) the eigenvalues are real and distinct; (2) there is only one real eigenvalue; (3) the eigenvalues are nonreal.  For each case build a table showing the matrix, its  eigenvalues and eigenvectors.  For repeated eigenvalue include a root vector. For nonreal eigenvalues, include only one eigenvalue and a corresponding eigenvector.  Explore  Options[TableForm]  to get your table to display nicely.

There are many similar names in this notebook, so I will set the spellcheck  Off.

Off[General::"spell1"];
Off[General::"spell"];

All matrices with real and distinct eigenvalues

The idea is to look at a matrix as a pair of two pairs of numbers: {{a, b}, {c, d}}.  This is Mathematica's representation for the matrix  (({{a, b}, {c, d}})) .  We will consider only matrices with a nonnegative entry a and relatively prime entries a, b, c, d.   Since all triangular matrices belong to this class we will not include them in the Book.

The set of all possible pairs  {a, b}  is the Cartesian product of the sets  {0, 1, 2,  ... , 9}  and   {-9,  -8,   ... , -1, 0, 1,  ... , 8, 
        9}.   The set of all possible pairs {c, d}  is the Cartesian square of the last set.  The set of matrices that we will start with is the Cartesian product of the previous two Cartesian products.

To get a Cartesian product of two sets I use  Outer[].  This is a complicated command.  For this purpose it is enough to observe that the integer at position  4  in  Outer  indicates to  Outer  to treat as separate elements only sublists at level 1.  Otherwise undesirable output will result. As an example consider

Outer[List, {a, b, c}, {x, y, z}, 1]

allabs = 
        Flatten[Outer[List, Range[0, 9], Range[-9, 9]], 1];

Length[allabs]

allcds = 
        Flatten[Outer[List, Range[-9, 9], Range[-9, 9]], 
          1];

Length[allcds]

Now I produce the starting set of matrices of interest. This set will be further reduced in several steps.

allmats = Flatten[Outer[List, allabs, allcds, 1], 1];

Length[allmats]

allmats[[21249]]

At this step we select non-triangular matrices whose entries are relatively prime.

mats1 = 
        Select[allmats, ((Apply[GCD, Flatten[#]] == 
                    1) ∧ (#[[1, 
                        2]] #[[2, 
                        1]] ≠ 
                    0)) &];

Length[mats1]

Looking at the formula for eigenvalues of a general matrix  (({{a, b}, {c, d}})) , doing it by hand or in Mathematica:

FullSimplify[Eigenvalues[{{a, b}, {c, d}}]]

We see that the eigenvalues will be distinct integers if  (Sqrt[4 b c + (((a - d))^2)]) is a positive integer.  Now I select matrices from  mats1 that satisfy this condition:

Timing[ mats2 = 
          Select[mats1, (IntegerQ[(Sqrt[4 #[[1, 2]]*#[[2, 1]] + (((#[[1, 1]] - #[[2, 2]]))^2)])] ∧ 
                  Positive[(Sqrt[4 #[[1, 2]]*#[[2, 1]] + (((#[[1, 1]] - #[[2, 2]]))^2)])]) &];]

Length[mats2]

The same selection can be done using  Cases[]

Timing[mats2a = 
          Cases[mats1, {{a_, b_}, {c_, d_}} /; 
              IntegerQ[(Sqrt[4 b c + (((a - d))^2)])] ∧ 
                Positive[(Sqrt[4 b c + (((a - d))^2)])]];]

Length[mats2a]

At the final step we drop matrices with too large eigenvalues and too large coordinates in eigenvectors.  These are all desired matrices.

To calculate the eigenvalues and the eigenvectors of a given matrix we use  Eigenystem[]

Eigensystem[mats2[[2356]]]

We are interested only in the absolute value of the numbers that appear in the last object:

Abs[Flatten[Eigensystem[mats2[[2356]]]]]

In fact only the  Max[] of the last set is relevant

Max[Abs[Flatten[Eigensystem[mats2[[2356]]]]]]

Timing[ mats3 = 
Select[ 
mats2, 
(Max[ 
                    Abs[Flatten[Eigensystem[#]]]] < 
                  10) &
              ];
              ]

Length[mats3]

Or using  Cases[]

Timing[mats3a = 
          Cases[mats2, 
            x_ /; Max[Abs[Flatten[Eigensystem[x]]]] < 10];]

Length[mats3a]

mats3[[1]]

Next few cells I will look into the structure of  mats3.

Table[{k, 
        Length[Select[
            mats3, #[[1, 1]] == 
                k &]]}, {k, 0, 9}]

Length[Select[mats3, Det[#] == 0 &]]

Length[Select[mats3, Det[#] > 0 &]]

Length[Select[mats3, Det[#] < 0 &]]

The lack of symmetry in the last two numbers comes from the fact that we are considering only nonnegative top-left eintry.

Recall that our final table has 4292 matrices. (It would have 2*4292-486, that is 8098 matrices if we also considered negative  top-left eintries.)

Length[mats3]

Since we are dealing with a large number of matrices, before making the final table I experiment with a short form of the table.  I do that by using  Take[]

?Take

Options[TableForm]

\(TableForm[
      Partition[Take[mats3, 24], 3, 
          3, {1, 1}, " "] /. {{a_Integer, b_Integer}, {c_Integer, 
              d_Integer}} :> \ \ {MatrixForm[{{a, b}, {c, d}}], 
            TableForm[Map[MatrixForm, Eigensystem[{{a, b}, {c, d}}], {2}], 
              TableSpacing -> \ {1}]}, 
      TableDirections -> {Column, Row, Row, Column}, 
      TableSpacing -> \ {4, 5, 2, 4}, 
      TableHeadings -> \ {None, {"Matrix      Eigensystem", "Matrix      Eigensystem", "Matrix      Eigensystem"}}]\)

?Partition

Partition[Take[mats3, 24], 3, 3, {1, 1}, " "]

\(TableForm[
      Partition[Take[mats3, Length[mats3]], 4, 
          4, {1, 2}, " "] /. {{a_Integer, b_Integer}, {c_Integer, 
              d_Integer}} :> \ \ {MatrixForm[{{a, b}, {c, d}}], 
            TableForm[Map[MatrixForm, Eigensystem[{{a, b}, {c, d}}], {2}], 
              TableSpacing -> \ { .5}]}, 
      TableDirections -> {Column, Row, Row, Column}, 
      TableSpacing -> \ {4, 3, 2, 1}, 
      TableHeadings -> \ {None, {"Matrix      Eigensystem", "Matrix      Eigensystem", "Matrix      Eigensystem", "Matrix      Eigensystem"}}]\)

\(TableForm[
      Partition[\ mats3, 3, 
          3, {1, 1}, " "] /. {{a_Integer, b_Integer}, {c_Integer, 
              d_Integer}} :> \ \ {MatrixForm[{{a, b}, {c, d}}], 
            TableForm[Map[MatrixForm, Eigensystem[{{a, b}, {c, d}}], {2}], 
              TableSpacing -> \ {1}]}, 
      TableDirections -> {Column, Row, Row, Column}, 
      TableSpacing -> \ {4, 5, 2, 4}, 
      TableHeadings -> \ {None, {"Matrix      Eigensystem", "Matrix      Eigensystem", "Matrix      Eigensystem"}}]\)

If we are interested only in matrices with a 0 eigenvalue we select them with

Select[ mats3, Det[#] == 0 &];

\(TableForm[
      Partition[Select[\ mats3, Det[#] == 0 &], 3, 
          3, {1, 1}, " "] /. {{a_Integer, b_Integer}, {c_Integer, 
              d_Integer}} :> \ \ {MatrixForm[{{a, b}, {c, d}}], 
            TableForm[Map[MatrixForm, Eigensystem[{{a, b}, {c, d}}], {2}], 
              TableSpacing -> \ {1}]}, 
      TableDirections -> {Column, Row, Row, Column}, 
      TableSpacing -> \ {4, 5, 2, 4}, 
      TableHeadings -> \ {None, {"Matrix      Eigensystem", "Matrix      Eigensystem", "Matrix      Eigensystem"}}]\)

Matrices ordered by eigenvalues

mats4e = 
        Cases[mats3, x_ :> {x, Eigensystem[x]}];

Length[mats4e]

mats5e = 
        Sort[ mats4e, ((#1[[2, 1, 1]] < #2[[2, 1, 
                        1]]) || (#1[[2, 1, 1]] == #2[[2, 1, 
                          1]] && #1[[2, 1, 2]] < #2[[2, 1, 
                          2]])) &];

Take[ mats5e, 6]

Length[mats5e]

mats6e = 
        Cases[mats5e, 
          x_ :>  
            x[[1]]];

mats6e[[1]]

Length[ mats6e]

\(TableForm[
      Partition[Take[mats6e, Length[mats6e]], 4, 
          4, {1, 2}, " "] /. {{a_Integer, b_Integer}, {c_Integer, 
              d_Integer}} :> \ \ {MatrixForm[{{a, b}, {c, d}}], 
            TableForm[Map[MatrixForm, Eigensystem[{{a, b}, {c, d}}], {2}], 
              TableSpacing -> \ { .5}]}, 
      TableDirections -> {Column, Row, Row, Column}, 
      TableSpacing -> \ {4, 3, 2, 1}, 
      TableHeadings -> \ {None, {"Matrix      Eigensystem", "Matrix      Eigensystem", "Matrix      Eigensystem", "Matrix      Eigensystem"}}]\)

\(TableForm[Partition[\ mats6e, 3, 3, {1, 1}, Null], 
      TableDirections -> {Column, Row, Row, Column}, 
      TableSpacing -> {4, 5, 2, 4}, 
      TableHeadings -> {None, {"Matrix           Eigensystem", "Matrix          Eigensystem", "Matrix      Eigensystem"}}]\)

A way to reduce the number of matrices using "standard rules"

mats4 = 
        Select[mats3, ((Apply[
                      Plus, #[[1]]] ≠ 
                    Apply[Plus, #[[2]]]) ∧ (Apply[
                      Plus, Transpose[#][[1]]] ≠ 
                    Apply[Plus, Transpose[#][[2]]])) &];

Length[mats4]

Table[{k, 
        Length[Select[
            mats4, #[[1, 1]] == 
                k &]]}, {k, 0, 9}]

Length[Select[mats4, Det[#] == 0 &]]

Length[ Cases[mats4, x_ /; Det[x] > 0]]

Length[ Cases[mats4, x_ /; Det[x] < 0]]

Stochastic matrices

Length[mats3]

mats3[[234]]

mats3[[234]][[2, 1]]

matsSt1 = 
        Select[mats3, (#[[1, 
                    1]] + #[[2, 
                    1]] == #[[1, 
                    2]] + #[[2, 
                    2]]) &];

Length[matsSt1]

matsSt2 = 
        Select[matsSt1, 
          And[#[[1, 1]] > 
                0, #[[2, 1]] > 
                0, #[[1, 2]] > 
                0, #[[2, 2]] > 
                0] &];

Length[matsSt2]

matsSte = 
        Cases[matsSt2, x_ :> {x, Eigensystem[x]}];

Length[matsSte]

matsSte1 = 
        Sort[ matsSte, ((#1[[2, 1, 1]] < #2[[2, 1, 
                        1]]) || (#1[[2, 1, 1]] == #2[[2, 1, 
                          1]] && #1[[2, 1, 2]] < #2[[2, 1, 
                          2]])) &];

Take[ matsSte1, 6]

Length[matsSte1]

matsSte2 = 
        Cases[matsSte1, 
          x_ :>  
            x[[1]]];

matsSte2[[1]]

Length[matsSte2]

\(TableForm[
      Partition[Take[matsSte2, Length[matsSte2]], 4, 
          4, {1, 2}, " "] /. {{a_Integer, b_Integer}, {c_Integer, 
              d_Integer}} :> \ \ {MatrixForm[{{a, b}, {c, d}}], 
            TableForm[Map[MatrixForm, Eigensystem[{{a, b}, {c, d}}], {2}], 
              TableSpacing -> \ { .3}]}, 
      TableDirections -> {Column, Row, Row, Column}, 
      TableSpacing -> \ {3, 2, 2, 1}, 
      TableHeadings -> \ {None, {"Matrix    Eigensystem", "Matrix    Eigensystem", "Matrix    Eigensystem", "Matrix    Eigensystem"}}]\)

Complex Integer Eigenvalues

Length[mats1]

Now we select all the matrices with non-real integer eigenvalues:

Timing[mats2c = 
          Cases[mats1, {{a_, b_}, {c_, d_}} /; 
              IntegerQ[(Sqrt[-((((a - d))^2) + 4 b c)])] ∧ 
                Positive[(Sqrt[-((((a - d))^2) + 4 b c)])]];]

Length[mats2c]

mats2c[[60]]

Max[
      Abs[Union[Re[Flatten[Eigensystem[allmats2c[[60]]]]], 
          Re[Flatten[Eigensystem[allmats2c[[60]]]]]]]]

Disregard too large numbers in the solutions: These are all relevant matrices.

Timing[mats3c = 
          Cases[mats2c, 
            x_ /; Max[
                  Abs[Union[Re[Flatten[Eigensystem[x]]], 
                      Re[Flatten[Eigensystem[x]]]]]] < 10];]

Length[mats3c]

TableForm[
      Take[mats3c, {250, 
            260}] /. {{a_, b_}, {c_, 
              d_}} :>   {MatrixForm[{{a, b}, {c, d}}], 
            Map[MatrixForm, Eigensystem[{{a, b}, {c, d}}], {2}]}]

TableForm[
      Take[mats3c, {115, 
            125}] /. {{a_, b_}, {c_, 
              d_}} :>   {MatrixForm[{{a, b}, {c, d}}], 
            Map[MatrixForm, {{Eigensystem[{{a, b}, {c, 
                          d}}][[1, 
                    1]]}, {LCM[
                      Denominator[
                        Re[Eigensystem[{{a, b}, {c, 
                                    d}}][[2, 
                              1]][[1]]]], 
                      Denominator[Eigensystem[{{a, b}, {c, 
                                  d}}][[2, 
                            1]][[1]]]] Eigensystem[{{a, b}, {c, d}}][[2, 
                      1]]}}, {2}]}]

\(TableForm[
      Partition[Take[mats3c, Length[mats2c]], 4, 
          4, {1, 2}, " "] /. {{a_Integer, b_Integer}, {c_Integer, 
              d_Integer}} :> \ \ {MatrixForm[{{a, b}, {c, d}}], 
            TableForm[
              Map[MatrixForm, {{\(Eigensystem[{{a, b}, {c, 
                            d}}]\)[[1, 
                      1]]}, {LCM[
                        Denominator[
                          Re[\(\(Eigensystem[{{a, b}, {c, 
                                      d}}]\)[[2, 
                                1]]\)[\([1]\)]]], 
                        Denominator[\(\(Eigensystem[{{a, b}, {c, 
                                    d}}]\)[[2, 
                              1]]\)[\([1]\)]]]
\(Eigensystem[{{a, b}, {c, d}}]\)[[2, 
                        1]]}}, {2}], 
              TableSpacing -> \ {1}]}, 
      TableDirections -> {Column, Row, Row, Column}, 
      TableSpacing -> \ {4, 3, 2, 1}, 
      TableHeadings -> \ {None, {"Matrix       Eigensystem", "Matrix      Eigensystem", "Matrix      Eigensystem", "Matrix      Eigensystem"}}]\)

Length[mats3c]

Re[Eigensystem[mats3c[[47]]][[1, 1]]]

mats4c = 
        Cases[
mats3c, 
          x_ :> {x, {Re[Eigensystem[x][[1, 1]]], 
                Im[Eigensystem[x][[1, 
                      1]]]}}
];

Length[ mats4c]

mats4c[[354]][[2, 2]]

mats5c = 
        Sort[mats4c, #1[[2, 1]] < #2[[2, 
                    1]] ∨ (#1[[2, 1]] == #2[[2, 
                        1]] ∧ #1[[2, 2]] < #2[[2, 
                        2]]) &];

Length[mats5c]

mats6c = 
        Cases[
mats5c, 
          x_ :> 
            x[[1]]
];

mats6c[[12]]

\(TableForm[
      Partition[Take[mats6c, Length[mats6c]], 4, 
          4, {1, 2}, " "] /. {{a_Integer, b_Integer}, {c_Integer, 
              d_Integer}} :> \ \ {MatrixForm[{{a, b}, {c, d}}], 
            TableForm[
              Map[MatrixForm, {{\(Eigensystem[{{a, b}, {c, 
                            d}}]\)[[1, 
                      1]]}, {LCM[
                        Denominator[
                          Re[\(\(Eigensystem[{{a, b}, {c, 
                                      d}}]\)[[2, 
                                1]]\)[\([1]\)]]], 
                        Denominator[\(\(Eigensystem[{{a, b}, {c, 
                                    d}}]\)[[2, 
                              1]]\)[\([1]\)]]]
\(Eigensystem[{{a, b}, {c, d}}]\)[[2, 
                        1]]}}, {2}], 
              TableSpacing -> \ {1}]}, 
      TableDirections -> {Column, Row, Row, Column}, 
      TableSpacing -> \ {4, 3, 2, 1}, 
      TableHeadings -> \ {None, {"Matrix       Eigensystem", "Matrix      Eigensystem", "Matrix      Eigensystem", "Matrix      Eigensystem"}}]\)

\(TableForm[
      Partition[\ mats6c, 3, 
          3, {1, 1}, " "] /. {{a_Integer, b_Integer}, {c_Integer, 
              d_Integer}} :> \ \ {MatrixForm[{{a, b}, {c, d}}], 
            TableForm[
              Map[MatrixForm, {{\(Eigensystem[{{a, b}, {c, 
                            d}}]\)[[1, 
                      1]]}, {LCM[
                        Denominator[
                          Re[\(\(Eigensystem[{{a, b}, {c, 
                                      d}}]\)[[2, 
                                1]]\)[\([1]\)]]], 
                        Denominator[\(\(Eigensystem[{{a, b}, {c, 
                                    d}}]\)[[2, 
                              1]]\)[\([1]\)]]]
\(Eigensystem[{{a, b}, {c, d}}]\)[[2, 
                        1]]}}, {2}], 
              TableSpacing -> \ {1}]}, 
      TableDirections -> {Column, Row, Row, Column}, 
      TableSpacing -> \ {4, 5, 2, 4}, 
      TableHeadings -> \ {None, {"Matrix      Eigensystem", "Matrix      Eigensystem", "Matrix      Eigensystem"}}]\)

Repeated Eigenvalues

Length[mats1]

Disregard the constant multiple and triangular matrices (b c  = 0):

Collect all repeated eigenvalues solutions:

Timing[mats2r = 
          Cases[mats1, {{a_, b_}, {c_, 
                  d_}} /; ¬ (a d - b c == 0) ∧ ((a^2) + 
                      4 b c - 2 a d + (d^2) == 0)];]

Length[mats2r]

Disregard too large numbers in the solutions: These are all relevant matrices.

Timing[mats3r = 
          Cases[mats2r, 
            x_ /; Max[Abs[Flatten[Eigensystem[x]]]] < 10];]

Length[mats3r]

Off[Solve::"svars"]

Clear[RootVector]; 
    RootVector[M_] := 
      Module[{evl, evc}, 
        evl = Eigensystem[M][[1, 
            1]]; 
        evc = Eigensystem[M][[2, 
            1]]; ({x, 
                y} /. Solve[(M - evl IdentityMatrix[2]) . {x, 
                        y} == evc, {x, 
                    y}][[1]]) /. {y -> 0, x -> 0}]

Length[mats3r]

RootVector[mats3r[[39]]]

\(TableForm[
Take[
          mats3r, {50, 
            60}] /. {{a_, b_}, {c_, 
              d_}} :> \ 
\ {
MatrixForm[{{a, b}, {c, d}}], 
            Map[MatrixForm, 
{
{\(Eigensystem[{{a, b}, {c, 
                          d}}]\)[\([1, 
                      1]\)], "r.v."}, 
Apply[LCM, 
                    Denominator[
                      RootVector[{{a, b}, {c, 
                            d}}]]]
{\(Eigensystem[{{a, 
                            b}, {c, d}}]\)[\([2, 1]\)], \ 
                    RootVector[{{a, b}, {c, 
                          d}}]}
}, 
{2}]
\ \ \ \ \ \ \ \ \ }
\ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ ]\)

\(TableForm[
Take[
          mats3r, {50, 
            60}] /. {{a_, b_}, {c_, 
              d_}} :> \ 
\ {
MatrixForm[{{a, b}, {c, d}}], 
            Map[MatrixForm, 
              ReplacePart[
                  Eigensystem[{{a, b}, {c, d}}], "r.v.", {1, 
                    2}] /. {{0, 0} -> 
                    RootVector[{{a, b}, {c, 
                          d}}]}, {2}]
\ \ \ \ \ \ \ \ \ }
\ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ ]\)

TableForm[
Partition[Take[mats3r, Length[mats3r]], 4, 
          4, {1, 1}, " "] /. {{a_Integer, b_Integer}, {c_Integer, 
              d_Integer}} :>   {MatrixForm[{{a, b}, {c, d}}], 
            TableForm[
              Map[MatrixForm, 
                ReplacePart[
                    Eigensystem[{{a, b}, {c, d}}], "r.v.", {1, 
                      2}] /. {{0, 0} -> 
                      RootVector[{{a, b}, {c, d}}]}, 
{2}], 
              TableSpacing ->  {1}]}, 
      TableDirections -> {Column, Row, Row, Column}, 
      TableSpacing ->  {4, 3, 1, 2}, 
      TableHeadings ->  {None, {"Matrix     Eigensystem", "Matrix \
    Eigensystem", "Matrix     Eigensystem", "Matrix     \
Eigensystem"}}]

Matrices ordered by repeated eigenvalues

Eigensystem[mats3r[[23]]]

Eigensystem[mats3r[[23]]][[1, 1]]

mats4r = 
        Cases[
mats3r, 
          x_ :> {x,  Eigensystem[x][[1, 
                  1]]}
];

Length[mats4r]

mats4r[[54]][[2]]

mats5r = Sort[mats4r, #1[[2]] < #2[[2]] &];

Length[mats5r]

mats6r = 
        Cases[
mats5r, 
          x_ :> 
            x[[1]]
];

TableForm[
Partition[Take[mats6r, Length[mats6r]], 4, 
          4, {1, 1}, " "] /. {{a_Integer, b_Integer}, {c_Integer, 
              d_Integer}} :>   {MatrixForm[{{a, b}, {c, d}}], 
            TableForm[
              Map[MatrixForm, 
                ReplacePart[
                    Eigensystem[{{a, b}, {c, d}}], "r.v.", {1, 
                      2}] /. {{0, 0} -> 
                      RootVector[{{a, b}, {c, d}}]}, 
{2}], 
              TableSpacing ->  {1}]}, 
      TableDirections -> {Column, Row, Row, Column}, 
      TableSpacing ->  {4, 3, 1, 2}, 
      TableHeadings ->  {None, {"Matrix     Eigensystem", "Matrix \
    Eigensystem", "Matrix     Eigensystem", "Matrix     \
Eigensystem"}}]

\(TableForm[
Partition[allmats6r, 3, 
          3, {1, 1}, " "] /. {{a_Integer, b_Integer}, {c_Integer, 
              d_Integer}} :> \ \ {MatrixForm[{{a, b}, {c, d}}], 
            TableForm[
              Map[MatrixForm, 
                ReplacePart[
                    Eigensystem[{{a, b}, {c, d}}], "r.v.", {1, 
                      2}] /. {{0, 0} -> RootVector[{{a, b}, {c, d}}]}, {2}], 
              TableSpacing -> \ {1}]}, 
      TableDirections -> {Column, Row, Row, Column}, 
      TableSpacing -> \ {4, 5, 2, 4}, 
      TableHeadings -> \ {None, {"Matrix      Eigensystem", "Matrix      Eigensystem", "Matrix      Eigensystem"}}]\)

Probabilitiy of a beautiful matrix