---
title: "Eigensystem"
language: "en"
type: "Symbol"
summary: "Eigensystem[m] gives a list {values, vectors} of the eigenvalues and eigenvectors of the square matrix m. Eigensystem[{m, a}] gives the generalized eigenvalues and eigenvectors of m with respect to a. Eigensystem[m, k] gives the eigenvalues and eigenvectors for the first k eigenvalues of m. Eigensystem[{m, a}, k] gives the first k generalized eigenvalues and eigenvectors."
keywords: 
- eigen decomposition
- eigen system
- matrix diagonalization
- characteristic polynomial
- eigenvalues
- eigenvectors
- characteristic values
- characteristic vectors
- EISPACK
- generalized eigen decomposition
- generalized eigen system
- generalized similarity transformation
- LAPACK
- linear matrix pencils
- similarity transformation
- ARPACK
- FEAST
- Arnoldi
- Lanczos
- EIGENQL
- TRIQL
- eigsys
- eig
- eigs
canonical_url: "https://reference.wolfram.com/language/ref/Eigensystem.html"
source: "Wolfram Language Documentation"
related_guides: 
  - 
    title: "Matrix Operations"
    link: "https://reference.wolfram.com/language/guide/MatrixOperations.en.md"
  - 
    title: "Matrix Decompositions"
    link: "https://reference.wolfram.com/language/guide/MatrixDecompositions.en.md"
  - 
    title: "Unsupervised Machine Learning"
    link: "https://reference.wolfram.com/language/guide/UnsupervisedMachineLearning.en.md"
related_functions: 
  - 
    title: "EigenvalueDecomposition"
    link: "https://reference.wolfram.com/language/ref/EigenvalueDecomposition.en.md"
  - 
    title: "Eigenvalues"
    link: "https://reference.wolfram.com/language/ref/Eigenvalues.en.md"
  - 
    title: "Eigenvectors"
    link: "https://reference.wolfram.com/language/ref/Eigenvectors.en.md"
  - 
    title: "NDEigensystem"
    link: "https://reference.wolfram.com/language/ref/NDEigensystem.en.md"
  - 
    title: "DEigensystem"
    link: "https://reference.wolfram.com/language/ref/DEigensystem.en.md"
  - 
    title: "NullSpace"
    link: "https://reference.wolfram.com/language/ref/NullSpace.en.md"
  - 
    title: "JordanDecomposition"
    link: "https://reference.wolfram.com/language/ref/JordanDecomposition.en.md"
  - 
    title: "SchurDecomposition"
    link: "https://reference.wolfram.com/language/ref/SchurDecomposition.en.md"
  - 
    title: "SingularValueDecomposition"
    link: "https://reference.wolfram.com/language/ref/SingularValueDecomposition.en.md"
  - 
    title: "QRDecomposition"
    link: "https://reference.wolfram.com/language/ref/QRDecomposition.en.md"
  - 
    title: "MatrixFunction"
    link: "https://reference.wolfram.com/language/ref/MatrixFunction.en.md"
related_tutorials: 
  - 
    title: "Eigenvalues and Eigenvectors"
    link: "https://reference.wolfram.com/language/tutorial/LinearAlgebra.en.md#9501"
---
# Eigensystem

Eigensystem[m] gives a list {values, vectors} of the eigenvalues and eigenvectors of the square matrix m. 

Eigensystem[{m, a}] gives the generalized eigenvalues and eigenvectors of m with respect to a. 

Eigensystem[m, k] gives the eigenvalues and eigenvectors for the first k eigenvalues of m. 

Eigensystem[{m, a}, k] gives the first k generalized eigenvalues and eigenvectors.

## Details and Options

* ``Eigensystem`` finds numerical eigenvalues and eigenvectors if ``m`` contains approximate real or complex numbers.

* For approximate numerical matrices ``m``, the eigenvectors are normalized. »

* For exact or symbolic matrices ``m``, the eigenvectors are not normalized. »

* All the nonzero eigenvectors given are independent. If the number of eigenvectors is equal to the number of nonzero eigenvalues, then corresponding eigenvalues and eigenvectors are given in corresponding positions in their respective lists. The eigenvalues correspond to rows in the eigenvector matrix.

* If there are more eigenvalues than independent eigenvectors, then each extra eigenvalue is paired with a vector of zeros. »

* If they are numeric, eigenvalues are sorted in order of decreasing absolute value.

* The eigenvalues and eigenvectors satisfy the matrix equation ``m.Transpose[vectors] == Transpose[vectors].DiagonalMatrix[values]``. »

* The generalized finite eigenvalues and eigenvectors satisfy ``m.Transpose[vectors] == a.Transpose[vectors].DiagonalMatrix[values]``. »

* Ordinary eigenvalues are always finite; generalized eigenvalues can be infinite. Infinite generalized eigenvalues correspond to those $v$ for which $m.v=\lambda  v \land a.v=0$. »

* When matrices ``m`` and ``a`` have a dimension‐$d$ shared null space, then $d$ of their generalized eigenvalues will be ``Indeterminate`` and will be paired with zero vectors in the generalized eigenvector list. »

* ``{vals, vecs} = Eigensystem[m]`` can be used to set ``vals`` and ``vecs`` to be the eigenvalues and eigenvectors, respectively. »

* For numeric eigenvalues, ``Eigensystem[m, k]`` gives the ``k`` eigenvalues that are largest in absolute value and corresponding eigenvectors.

* ``Eigensystem[m, -k]`` gives the ``k`` eigenvalues that are smallest in absolute value and corresponding eigenvectors.

* ``Eigensystem[m, spec]`` is equivalent to applying ``Take[…, spec]`` to each element of ``Eigensystem[m]``.

* ``Eigensystem[m, UpTo[k]]`` gives ``k`` eigenvalues and corresponding eigenvectors, or as many as are available.

* ``SparseArray`` objects and structured arrays can be used in ``Eigensystem``.

* ``Eigensystem`` has the following options and settings:

|           |           |                                           |
| --------- | --------- | ----------------------------------------- |
| Cubics    | False     | whether to use radicals to solve cubics   |
| Method    | Automatic | select a method to use                    |
| Quartics  | False     | whether to use radicals to solve quartics |
| ZeroTest  | Automatic | test for when expressions are zero        |

* The ``ZeroTest`` option only applies to exact and symbolic matrices.

* Explicit ``Method`` settings for approximate numeric matrices include:

|           |                                                                                                    |
| --------- | -------------------------------------------------------------------------------------------------- |
| "Arnoldi" | Arnoldi iterative method for finding a few eigenvalues                                             |
| "Banded"  | direct banded matrix solver for Hermitian matrices                                                 |
| "Direct"  | direct method for finding all eigenvalues                                                          |
| "FEAST"   | FEAST iterative method for finding eigenvalues in an interval (applies to Hermitian matrices only) |

* The ``"Arnoldi"`` method is also known as a Lanczos method when applied to symmetric or Hermitian matrices.

* The ``"Arnoldi"`` and ``"FEAST"`` methods take suboptions ``Method -> {"name", opt1 -> val1, …}``, which can be found in the ``Method`` subsection.

---

## Examples (74)

### Basic Examples (5)

Eigenvalues and eigenvectors of an exact, two-dimensional matrix:

```wl
In[1]:= Eigensystem[{{-3, 2}, {-15, 8}}]

Out[1]= {{3, 2}, {{1, 3}, {2, 5}}}
```

---

Set ``vals`` and ``vecs`` to be the eigenvalues and eigenvectors, respectively:

```wl
In[1]:=
m = {{1, 2}, {3, 4}};
{vals, vecs} = Eigensystem[m];
```

The eigenvalues are a simple list:

```wl
In[2]:= vals

Out[2]= {(1/2) (5 + Sqrt[33]), (1/2) (5 - Sqrt[33])}
```

The eigenvectors are a list of lists:

```wl
In[3]:= vecs

Out[3]= {{(1/6) (-3 + Sqrt[33]), 1}, {(1/6) (-3 - Sqrt[33]), 1}}
```

The result of multiplying each eigenvector by ``m`` is to scalar multiply by the corresponding eigenvalue:

```wl
In[4]:= {m.Subscript[vecs, [[1]]] == Subscript[vals, [[1]]]Subscript[vecs, [[1]]], m.Subscript[vecs, [[2]]] == Subscript[vals, [[2]]]Subscript[vecs, [[2]]]}

Out[4]= {True, True}
```

---

Eigenvalues and eigenvectors of a three-dimensional matrix:

```wl
In[1]:= Eigensystem[{{1, 2, 3}, {4, 5, 6}, {7, 8, 9}}]

Out[1]= {{(3/2) (5 + Sqrt[33]), -(3/2) (-5 + Sqrt[33]), 0}, {{-(-7 - Sqrt[33]/2 (11 + 2 Sqrt[33])), -(-29 - 5 Sqrt[33]/4 (11 + 2 Sqrt[33])), 1}, {-(7 - Sqrt[33]/2 (-11 + 2 Sqrt[33])), -(29 - 5 Sqrt[33]/4 (-11 + 2 Sqrt[33])), 1}, {1, -2, 1}}}
```

---

Eigenvalues and eigenvectors computed with machine precision:

```wl
In[1]:= {val, vec} = Eigensystem[{{1.1, 2.2, 3.25}, {0.76, 4.6, 5}, {0.1, 0.1, 6.1}}];

In[2]:= val

Out[2]= {6.60674, 4.52536, 0.667901}

In[3]:= vec//MatrixForm

Out[3]//MatrixForm=
(⁠|           |           |            |
| --------- | --------- | ---------- |
| 0.48687   | 0.833694  | 0.260598   |
| -0.479424 | -0.873368 | 0.085911   |
| 0.985096  | -0.171352 | -0.0149803 |⁠)
```

---

Symbolic eigenvalues and eigenvectors:

```wl
In[1]:= Eigensystem[{{a, b}, {c, d}}]

Out[1]= {{(1/2) (a + d - Sqrt[a^2 + 4 b c - 2 a d + d^2]), (1/2) (a + d + Sqrt[a^2 + 4 b c - 2 a d + d^2])}, {{-(-a + d + Sqrt[a^2 + 4 b c - 2 a d + d^2]/2 c), 1}, {-(-a + d - Sqrt[a^2 + 4 b c - 2 a d + d^2]/2 c), 1}}}
```

### Scope (20)

#### Basic Uses (6)

Find the eigensystem of a machine-precision matrix:

```wl
In[1]:= MatrixForm /@ Eigensystem[Table[N[1 / (i + j + 1), 18], {i, 3}, {j, 3}]]

Out[1]=
{(⁠|                         |
| ----------------------- |
| 0.657051428297579245    |
| 0.0189263109740709143   |
| 0.000212736918826030734 |⁠), (⁠|                      |                       |                       |
| -------------------- | -- ... | --------------------- |
| 0.703152857567645584 | 0.549267614369702795  | 0.451531999640191372  |
| 0.668534800678373420 | -0.294440006648162146 | -0.682910171813949300 |
| 0.242151355925135888 | -0.782055094152352984 | 0.574240847148626711  |⁠)}
```

---

Approximate 18-digit precision eigenvalues and eigenvectors:

```wl
In[1]:= Eigensystem[Table[N[1 / (i + j + 1), 18], {i, 3}, {j, 3}]]

Out[1]= {{0.657051428297579245, 0.0189263109740709143, 0.000212736918826030734}, {{0.703152857567645584, 0.549267614369702795, 0.451531999640191372}, {0.668534800678373420, -0.294440006648162146, -0.682910171813949300}, {0.242151355925135888, -0.782055094152352984, 0.574240847148626711}}}
```

---

Eigensystem of a complex matrix:

```wl
In[1]:= Eigensystem[{{1.1 - .2I, 2.2, 3.25}, {0.76, 4.6, 5 - 2I}, {0.1, 0.1 + I, 6.1}}]

Out[1]= {{7.64039  + 1.15405 I, 3.53464  - 1.16549 I, 0.624978  - 0.188568 I}, {{0.429431  + 0.0578862 I, 0.796511  + 0. I, 0.300904  + 0.295406 I}, {0.540025  - 0.095532 I, 0.787863  + 0. I, -0.157159 - 0.231991 I}, {0.981801  + 0. I, -0.179757 - 0.045745 I, -0.0218189 + 0.0344192 I}}}
```

---

A real, exact eigensystem:

```wl
In[1]:= Eigensystem[{{(1/3), (1/2), (3/5)}, {(1/2), (4/5), 1}, {(3/5), 1, (9/7)}}]

Out[1]= {{Root[-1 + 1195*#1 - 25400*#1^2 + 10500*#1^3 & , 3, 0], Root[-1 + 1195*#1 - 25400*#1^2 + 10500*#1^3 & , 2, 0], Root[-1 + 1195*#1 - 25400*#1^2 + 10500*#1^3 & , 1, 0]}, {{Root[3150 - 6603*#1 - 1400*#1^2 + 3528*#1^3 & , 2, 0], Root[-3850 - 6635*#1 +  ... & , 3, 0], 1}, {Root[3150 - 6603*#1 - 1400*#1^2 + 3528*#1^3 & , 1, 0], Root[-3850 - 6635*#1 + 9856*#1^2 + 5880*#1^3 & , 2, 0], 1}, {Root[3150 - 6603*#1 - 1400*#1^2 + 3528*#1^3 & , 3, 0], Root[-3850 - 6635*#1 + 9856*#1^2 + 5880*#1^3 & , 1, 0], 1}}}
```

A complex, exact eigensystem:

```wl
In[2]:= Eigensystem[{{π, (1/3)}, {I, 5}}]

Out[2]= {{(1/6) (15 + 3 π + Sqrt[3 ((75 + 4 I) - 30 π + 3 π^2)]), (1/6) (15 + 3 π - Sqrt[3 ((75 + 4 I) - 30 π + 3 π^2)])}, {{-(1/6) I (-15 + 3 π + Sqrt[3 ((75 + 4 I) - 30 π + 3 π^2)]), 1}, {(1/6) I (15 - 3 π + Sqrt[3 ((75 + 4 I) - 30 π + 3 π^2)]), 1}}}
```

---

Eigensystem of a symbolic matrix:

```wl
In[1]:=
Eigensystem[(|   |   |   |
| - | - | - |
| a | 2 | 0 |
| 2 | 3 | 1 |
| 0 | 1 | 7 |)]

Out[1]= {{Root[28 - 20 a + (16 + 10 a) #1 + (-10 - a) #1^2 + #1^3&, 1], Root[28 - 20 a + (16 + 10 a) #1 + (-10 - a) #1^2 + #1^3&, 2], Root[28 - 20 a + (16 + 10 a) #1 + (-10 - a) #1^2 + #1^3&, 3]}, {{(1/2) (20 - 10 Root[28 - 20 a + (16 + 10 a) #1 + (-10 - a ... (-10 - a) #1^2 + #1^3&, 2], 1}, {(1/2) (20 - 10 Root[28 - 20 a + (16 + 10 a) #1 + (-10 - a) #1^2 + #1^3&, 3] + Root[28 - 20 a + (16 + 10 a) #1 + (-10 - a) #1^2 + #1^3&, 3]^2), -7 + Root[28 - 20 a + (16 + 10 a) #1 + (-10 - a) #1^2 + #1^3&, 3], 1}}}
```

---

The eigenvalues and eigenvectors of large numerical matrices are computed efficiently:

```wl
In[1]:= m = RandomReal[{0, 5}, {100, 100}];

In[2]:= Eigensystem[m] ;//Timing

Out[2]= {0.051753, Null}
```

#### Subsets of Eigenvalues (5)

Compute the three largest eigenvalues and their corresponding eigenvectors:

```wl
In[1]:= {vals, vecs} = Eigensystem[Table[If[Abs[i - j] < 3, 1.0, 0], {i, 100}, {j, 100}] , 3];
```

Visualize the three vectors, using the eigenvalue as a label:

```wl
In[2]:= ListPlot[vecs, PlotLegends -> vals]

Out[2]= [image]
```

---

Eigensystem corresponding to the three smallest eigenvalues:

```wl
In[1]:= Eigensystem[Table[N[1 / (i + j + 1), 12], {i, 6}, {j, 6}], -3]

Out[1]= {{0.0000531681803613, 7.9933981112676432005931990090801`12.*^-7, 5.374160219636231802894455446`12.*^-9}, {{0.156642430099, -0.604475265692, 0.298411689270, 0.498707568405, 0.134069731544, -0.504320666110}, {0.0397884774870, -0.303195059511, 0.647138724596, -0.233751807251, -0.533115450788, 0.385809684708}, {-0.00607971468708, 0.0785483344143, -0.345051984627, 0.677082251733, -0.611143120843, 0.206872024132}}}
```

---

Find the eigensystem corresponding to the four largest eigenvalues, or as many as there are if fewer:

```wl
In[1]:=
Eigensystem[(⁠|       |       |       |
| ----- | ----- | ----- |
| (1/3) | (1/4) | (1/5) |
| (1/4) | (1/5) | (1/6) |
| (1/5) | (1/6) | (1/7) |⁠), UpTo[4]]

Out[1]= {{Root[-1 + 4755*#1 - 255600*#1^2 + 378000*#1^3 & , 3, 0], Root[-1 + 4755*#1 - 255600*#1^2 + 378000*#1^3 & , 2, 0], Root[-1 + 4755*#1 - 255600*#1^2 + 378000*#1^3 & , 1, 0]}, {{Root[126 - 251*#1 - 196*#1^2 + 196*#1^3 & , 3, 0], Root[140 - 337*#1 - 56*#1^2 + 196*#1^3 & , 3, 0], 1}, {Root[126 - 251*#1 - 196*#1^2 + 196*#1^3 & , 1, 0], Root[140 - 337*#1 - 56*#1^2 + 196*#1^3 & , 2, 0], 1}, {Root[126 - 251*#1 - 196*#1^2 + 196*#1^3 & , 2, 0], Root[140 - 337*#1 - 56*#1^2 + 196*#1^3 & , 1, 0], 1}}}
```

---

Repeated eigenvalues are considered when extracting a subset of the eigensystem:

```wl
In[1]:= mat = {{(7/2), 0, (1/2), 0}, {0, 3, 0, 1}, {(1/2), 0, (7/2), 0}, {0, 1, 0, 3}};

In[2]:= Eigensystem[mat]

Out[2]= {{4, 4, 3, 2}, {{0, 1, 0, 1}, {1, 0, 1, 0}, {-1, 0, 1, 0}, {0, -1, 0, 1}}}

In[3]:= Eigensystem[mat, 3]

Out[3]= {{4, 4, 3}, {{0, 1, 0, 1}, {1, 0, 1, 0}, {-1, 0, 1, 0}}}

In[4]:= Eigensystem[mat, -3]

Out[4]= {{4, 3, 2}, {{1, 0, 1, 0}, {-1, 0, 1, 0}, {0, -1, 0, 1}}}
```

---

Zero vectors are used when there are more eigenvalues than independent eigenvectors:

```wl
In[1]:= Eigensystem[{{2, 1, 0}, {0, 2, 0}, {0, 0, 1}}]

Out[1]= {{2, 2, 1}, {{1, 0, 0}, {0, 0, 0}, {0, 0, 1}}}
```

#### Generalized Eigenvalues (4)

Compute machine-precision generalized eigenvalues and eigenvectors:

```wl
In[1]:= a = {{1., 2.}, {3., 4.}};

In[2]:= b = {{1., 4.}, {9., 16.}};

In[3]:= Eigensystem[{a, b}]

Out[3]= {{0.25  + 0.193649 I, 0.25  - 0.193649 I}, {{0.848472  - 0.0858378 I, -0.498043 - 0.157099 I}, {0.848472  + 0.0858378 I, -0.498043 + 0.157099 I}}}
```

---

Generalized exact eigensystem:

```wl
In[1]:= a = {{1, 1, 1}, {1, 0, 1}, {0, 0, 1}};

In[2]:= b = {{0, 1, 1}, {0, 1, 1}, {1, 0, 0}};

In[3]:= Eigensystem[{a, b}]

Out[3]= {{∞, (1/2) (1 + Sqrt[5]), (1/2) (1 - Sqrt[5])}, {{0, -1, 1}, {(1/2) (-1 + Sqrt[5]), 0, 1}, {(1/2) (-1 - Sqrt[5]), 0, 1}}}
```

Compute the result at finite precision:

```wl
In[4]:= Eigensystem[N[{a, b}, 20]]//Chop

Out[4]= {{∞, 1.6180339887498948482, -0.61803398874989484820}, {{0, 1.0000000000000000000, -1.0000000000000000000}, {-0.61803398874989484820, 0, -1.0000000000000000000}, {-1.0000000000000000000, 0, 0.61803398874989484820}}}
```

---

Compute a symbolic generalized eigensystem:

```wl
In[1]:= a = {{x, 1 + x}, {1 - x, x}};

In[2]:= b = {{1, 1}, {1, 1}};

In[3]:= Eigenvectors[{a, b}]

Out[3]= {{-1, 1}, {-(1/-1 + 2 x), 1}}
```

---

Find the two smallest generalized eigenvalues and corresponding generalized eigenvectors:

```wl
In[1]:=
a = {{1, 2, 3}, {4, 5, 6}, {7, 8, 9}};
b = {{11, 12, 13}, {1, 15, 16}, {17, 18, 19}};
Eigensystem[{a, b}, -2]

Out[1]= {{1, 0}, {{0, -1, 1}, {1, -2, 1}}}
```

#### Special Matrices (5)

Eigensystems of sparse matrices:

```wl
In[1]:= SparseArray[{{1, 3} -> 1, {2, 2} -> 2, {3, 1} -> 3, {4, 2} -> 5}, {4, 4}]

Out[1]= SparseArray[Automatic, {4, 4}, 0, {1, {{0, 1, 2, 3, 4}, {{3}, {2}, {1}, {2}}}, {1, 2, 3, 5}}]

In[2]:= Eigensystem[%]

Out[2]= {{2, -Sqrt[3], Sqrt[3], 0}, {{0, 2, 0, 5}, {-(1/Sqrt[3]), 0, 1, 0}, {(1/Sqrt[3]), 0, 1, 0}, {0, 0, 0, 1}}}

In[3]:= SparseArray[Band[{1, 1}, {-1, -1}] -> {{{1, 2}, {2, 1}}}, {4, 4}]

Out[3]=
SparseArray[Automatic, {4, 4}, 0, {1, {{0, 2, 4, 6, 8}, {{2}, {1}, {2}, {1}, {4}, {3}, {4}, {3}}}, 
  {2, 1, 1, 2, 2, 1, 1, 2}}]

In[4]:= Eigensystem[%]

Out[4]= {{3, 3, -1, -1}, {{0, 0, 1, 1}, {1, 1, 0, 0}, {0, 0, -1, 1}, {-1, 1, 0, 0}}}
```

---

Eigensystems of structured matrices:

```wl
In[1]:= SymmetrizedArray[{{1, 1} -> 2.2, {1, 2} -> 1.1, {3, 2} -> 5.7}, {3, 3}, Symmetric[All]]

Out[1]=
SymmetrizedArray[StructuredArray`StructuredData[{3, 3}, 
  {{{1, 1} -> 2.2, {1, 2} -> 1.1, {2, 3} -> 5.7}, Symmetric[{1, 2}]}]]

In[2]:= Eigensystem[%]

Out[2]= {{5.86736, -5.77635, 2.109}, {{0.210326, 0.701218, 0.681217}, {-0.0976925, 0.708391, -0.699027}, {0.972738, -0.0804737, -0.217497}}}

In[3]:= QuantityArray[{{1, 2}, {3, 4}}, "Meters"]

Out[3]= QuantityArray[StructuredArray`StructuredData[{2, 2}, {{{1, 2}, {3, 4}}, "Meters", {{1}, {2}}}]]
```

The units of a ``QuantityArray`` object are in the eigenvalues, leaving the eigenvectors dimensionless:

```wl
In[4]:= Eigensystem[%]

Out[4]= {{Quantity[(1/2)*(5 + Sqrt[33]), "Meters"], Quantity[(1/2)*(5 - Sqrt[33]), "Meters"]}, {{(1/6) (-3 + Sqrt[33]), 1}, {(1/6) (-3 - Sqrt[33]), 1}}}
```

---

All the eigenvalues of ``IdentityMatrix[n]`` are 1, and it eigenvectors form the standard basis for $\mathbb{C}^n$ :

```wl
In[1]:= Eigensystem[IdentityMatrix[4]]

Out[1]= {{1, 1, 1, 1}, {{0, 0, 0, 1}, {0, 0, 1, 0}, {0, 1, 0, 0}, {1, 0, 0, 0}}}
```

---

Eigenvectors of ``HilbertMatrix`` :

```wl
In[1]:= Eigensystem[HilbertMatrix[3]]

Out[1]= {{Root[-1 + 381*#1 - 3312*#1^2 + 2160*#1^3 & , 3, 0], Root[-1 + 381*#1 - 3312*#1^2 + 2160*#1^3 & , 2, 0], Root[-1 + 381*#1 - 3312*#1^2 + 2160*#1^3 & , 1, 0]}, {{Root[20 - 92*#1 - 95*#1^2 + 50*#1^3 & , 3, 0], Root[30 - 29*#1 - 30*#1^2 + 25*#1^3 & , 3, 0], 1}, {Root[20 - 92*#1 - 95*#1^2 + 50*#1^3 & , 1, 0], Root[30 - 29*#1 - 30*#1^2 + 25*#1^3 & , 2, 0], 1}, {Root[20 - 92*#1 - 95*#1^2 + 50*#1^3 & , 2, 0], Root[30 - 29*#1 - 30*#1^2 + 25*#1^3 & , 1, 0], 1}}}
```

If the matrix is first numericized, the eigenvectors (but not eigenvalues) change significantly:

```wl
In[2]:= Eigensystem[N@HilbertMatrix[3]]

Out[2]= {{1.40832, 0.122327, 0.00268734}, {{-0.827045, -0.459864, -0.323298}, {0.547448, -0.52829, -0.649007}, {0.127659, -0.713747, 0.688672}}}
```

This is because the eigenvectors are normalized for numeric input:

```wl
In[3]:= Norm /@ Last[%]

Out[3]= {1., 1., 1.}
```

---

Eigensystem of a ``CenteredInterval`` matrix:

```wl
In[1]:= SeedRandom[777];(m = Map[CenteredInterval, RandomReal[{-10, 10}, {3, 3}, WorkingPrecision -> 10], {2}])//MatrixForm

Out[1]//MatrixForm=
(⁠|                                                            |                                                              |                                                              |
| ------------------------------------------------------- ... edInterval[{{-12098865075, -32, 649552873, -61}, 44}]  |
| CenteredInterval[{{9652099725, -31, 1036386317, -61}, 44}] | CenteredInterval[{{-28987783505, -32, 778134889, -60}, 44}]  | CenteredInterval[{{-36931291055, -32, 991366796, -60}, 44}]  |⁠)

In[2]:= {vals, vecs} = Eigensystem[m]

Out[2]=
{{CenteredInterval[{{{-6160307369087, -40, 627908656, -54}, {-14680229171717, -41, 627628227, -54}}, 
  44}], CenteredInterval[{{{-6972228065791, -39, 627269961, -54}, {623616021747, -82, 622875997, -54}}, 44}], CenteredInterval[{{{-12320614738175, ... 23, -49, 622170095, -55}}, 
  44}], CenteredInterval[{{{2909162014925, -47, 918120963, -57}, {8535013209365, -48, 918120963, -57}}, 44}], CenteredInterval[{{{-13514896585567, -50, 865140736, -56}, {-5295826121167, -46, 865140736, -56}}, 
  44}]}}}
```

Find an eigensystem for a random representative ``mrep`` of ``m`` :

```wl
In[3]:=
ranrep[e_CenteredInterval] := e["Center"] + RandomInteger[{-1000, 1000}] / 1000 e["Radius"]
(mrep = Map[ranrep, m, {2}])//MatrixForm

Out[3]//MatrixForm=
(⁠|                                                |                                                   |                                                 |
| ---------------------------------------------- | ------------------------------------------ ... 1769003176657047/1152921504606846976000) | -(405970545404434309223/144115188075855872000)  |
| (1295482895484724051539/288230376151711744000) | -(7781348883120069132599/1152921504606846976000)  | -(2478416988811662952469/288230376151711744000) |⁠)

In[4]:= {rvals, rvecs} = Eigensystem[mrep]

Out[4]=
{{Root[369071202259841451263042514367120946204613532913603086278188381377 + 
   83547980068004018615158535034473301127724594046645097558179840000*#1 + 
   9152041454369020808951116863264015207428228599753287401472000000*#1^2 + 
   38312388521647221 ...  - 
   233632818136281030886131622819965517428901792022007433515154138437*#1 - 
   15361018696580494699116764908392230761904849747333396833129217865*#1^2 + 
   417736994936897384931041589259784242262809078914421635563696080963*#1^3 & , 2, 0], 1}}}
```

Verify that, after reordering and scaling of vectors, ``vals`` contain ``rvals`` and ``vecs`` contain ``rvecs`` :

```wl
In[5]:= MapThread[IntervalMemberQ, {vals, rvals[[{3, 1, 2}]]}]

Out[5]= {True, True, True}

In[6]:= MapThread[IntervalMemberQ, {vecs, rvecs[[{3, 1, 2}]] * (Last[#]["Center"]& /@ vecs)}, 2]

Out[6]= {{True, True, True}, {True, True, True}, {True, True, True}}
```

### Options (11)

#### Cubics (1)

A 3×3 Vandermonde matrix:

```wl
In[1]:= m = Table[{1, -1, 2} ^ j, {j, 3}]

Out[1]= {{1, -1, 2}, {1, 1, 4}, {1, -1, 8}}
```

In general, for exact 3×3 matrices the result will be given in terms of ``Root`` objects:

```wl
In[2]:= er = Eigensystem[m]

Out[2]= {{Root[-12 + 20*#1 - 10*#1^2 + #1^3 & , 1, 0], Root[-12 + 20*#1 - 10*#1^2 + #1^3 & , 3, 0], Root[-12 + 20*#1 - 10*#1^2 + #1^3 & , 2, 0]}, {{(1/2) (-18 + 10 Root[-12 + 20*#1 - 10*#1^2 + #1^3 & , 1, 0] - Root[-12 + 20*#1 - 10*#1^2 + #1^3 & , 1, 0]^2) ...  #1^3 & , 3, 0]^2), 1}, {(1/2) (-18 + 10 Root[-12 + 20*#1 - 10*#1^2 + #1^3 & , 2, 0] - Root[-12 + 20*#1 - 10*#1^2 + #1^3 & , 2, 0]^2), (1/2) (-2 + 8 Root[-12 + 20*#1 - 10*#1^2 + #1^3 & , 2, 0] - Root[-12 + 20*#1 - 10*#1^2 + #1^3 & , 2, 0]^2), 1}}}
```

To get the result in terms of radicals, use the ``Cubics`` option:

```wl
In[3]:= ec = Eigensystem[m, Cubics -> True]

Out[3]= {{-(48520/131) + (21869 (120 (9161 - 393 Sqrt[129])^1 / 3 - 2 √(3 (9161 - 393 Sqrt[129] + 400 (9161 - 393 Sqrt[129])^1 / 3 - 400 I Sqrt[3] (9161 - 393 Sqrt[129])^1 / 3 + 80 (9161 - 393 Sqrt[129])^2 / 3 + I Sqrt[3] (9161 - 393 Sqrt[129]))))) / (1572 ...  Sqrt[129])^1 / 3 - 2 √(3 (9161 - 393 Sqrt[129] + 400 (9161 - 393 Sqrt[129])^1 / 3 - 400 I Sqrt[3] (9161 - 393 Sqrt[129])^1 / 3 + 80 (9161 - 393 Sqrt[129])^2 / 3 + I Sqrt[3] (9161 - 393 Sqrt[129]))))^4 / (339552 (9161 - 393 Sqrt[129])^4 / 3), 1}}}
```

Note that the result with ``Root`` objects is better suited to subsequent numerical evaluation:

```wl
In[4]:= N[er]

Out[4]= {{7.56632  + 6.942707647121589`*^-15 I, 1.21684  - 0.324456 I, 1.21684  + 0.324456 I}, {{0.207012  - 2.2216664470789086`*^-13 I, 0.640696  - 1.1108332235394543`*^-13 I, 1.}, {-3.60351 - 1.22747 I, 3.17965  - 0.903013 I, 1.}, {-3.60351 + 1.22747 I, 3.17965  + 0.903013 I, 1.}}}

In[5]:= N[ec]

Out[5]= {{7.56632  - 2.1316282072803006`*^-14 I, 1.21684  + 0.324456 I, 1.21684  - 0.324456 I}, {{0.207012  + 5.684341886080801`*^-13 I, 0.640696  + 3.410605131648481`*^-13 I, 1.}, {-3.60351 + 1.22747 I, 3.17965  + 0.903013 I, 1.}, {-3.60351 - 1.22747 I, 3.17965  - 0.903013 I, 1.}}}
```

---

#### Method (9)

##### "Arnoldi" (6)

The Arnoldi method can be used for machine- and arbitrary-precision matrices. The implementation of the Arnoldi method is based on the "ARPACK" library. It is most useful for large sparse matrices.

The following suboptions can be specified for the method ``"Arnoldi"`` :

|                  |                                            |
| ---------------- | ------------------------------------------ |
| "BasisSize"      | the size of the Arnoldi basis              |
| "Criteria"       | which criteria to use                      |
| "MaxIterations"  | the maximum number of iterations           |
| "Shift"          | the Arnoldi shift                          |
| "StartingVector" | the initial vector to start iterations     |
| "Tolerance"      | the tolerance used to terminate iterations |

Possible settings for ``"Criteria"`` include:

|                 |                                                                        |
| --------------- | ---------------------------------------------------------------------- |
| "Magnitude"     | based on Abs                                                           |
| "RealPart"      | based on Re                                                            |
| "ImaginaryPart" | based on Im                                                            |
| "BothEnds"      | a few eigenvalues from both ends of the symmetric real matrix spectrum |

---

Compute the largest eigenpair using different ``"Criteria"`` settings. The matrix ``m`` has eigenvalues $\{1,0.8\, +0.8 i,i\}$ :

```wl
In[1]:=
d = DiagonalMatrix[{1, 0.8 + 0.8I, I}];
t = {{1, 2, 0}, {0, 3, 2}, {1, 0, 4}};
m = t.d.Inverse[t];

[image]
```

By default, ``"Criteria" -> "Magnitude"`` selects a largest-magnitude eigenpair:

```wl
In[2]:= Eigensystem[m, 1, Method -> "Arnoldi"]

Out[2]= {{0.8  + 0.8 I}, {{0.473247  - 0.289361 I, 0.70987  - 0.434042 I, 7.732758783940383`*^-16 + 2.637449508383708`*^-16 I}}}
```

Find the largest real-part eigenpair:

```wl
In[3]:= Eigensystem[m, 1, Method -> {"Arnoldi", "Criteria" -> "RealPart"}]

Out[3]= {{1.  - 8.548765893110118`*^-17 I}, {{0.603274  - 0.368864 I, 4.128257639083761`*^-16 + 2.0398200902834518`*^-16 I, 0.603274  - 0.368864 I}}}
```

Find the largest imaginary-part eigenpair:

```wl
In[4]:= Eigensystem[m, 1, Method -> {"Arnoldi", "Criteria" -> "ImaginaryPart"}]

Out[4]= {{1.1102230246251565`*^-16 + 1. I}, {{-9.49598139581231`*^-17 - 2.4343977256369846`*^-16 I, -0.381544 + 0.23329 I, -0.763088 + 0.466581 I}}}
```

---

Find two eigenpairs from both ends of the symmetric matrix spectrum:

```wl
In[1]:=
d = DiagonalMatrix[{1, 0.8, -1, 2}];
t = Orthogonalize[{{1, 2, 0, 1}, {0, 3, 2, -2}, {1, 0, 4, -1}, {0, 1, 0, 1}}];
m = t.d.Transpose[t];

In[2]:= Eigensystem[m, 2, Method -> {"Arnoldi", "Criteria" -> "BothEnds"}]

Out[2]= {{2., -1.}, {{-0.408248, 0.704361, -0.259166, -0.519656}, {2.775557561562892`*^-17, -0.528271, -0.784503, -0.324785}}}
```

---

Use ``"StartingVector"`` to avoid randomness:

```wl
In[1]:=
d = DiagonalMatrix[{1, -1, 0.}];
t = {{1, 2, 0}, {0, 3, 2}, {1, 0, 4}};
m = t.d.Inverse[t];
```

Different starting vectors may converge to different eigenpairs:

```wl
In[2]:= Eigensystem[m, 1, Method -> {"Arnoldi", "StartingVector" -> {1, 0, 0}}]

Out[2]= {{-1.}, {{-0.5547, -0.83205, -8.630319771192489`*^-17}}}

In[3]:= Eigensystem[m, 1, Method -> {"Arnoldi", "StartingVector" -> {0, 0, 1}}]

Out[3]= {{1.}, {{-0.707107, 4.738173134917611`*^-17, -0.707107}}}
```

---

Use ``"Shift" -> μ`` to shift the eigenvalues by transforming the matrix $m$ to $m-\mathbb{I} \mu$. This preserves the eigenvectors but changes the eigenvalues by ``-μ``. The method compensates for the changed eigenvalues. ``"Shift"`` is typically used to find eigenpairs when there are no criteria such as largest or smallest magnitude that can select them:

```wl
In[1]:=
d = DiagonalMatrix[{1, 2, 3}];
t = {{1, 2, 0}, {0, 3, 2}, {1, 0, 4}};
m = N@t.d.Inverse[t];

In[2]:= μ = 2.1;
```

Manually shift the matrix and adjust the resulting eigenvalue:

```wl
In[3]:= {λl1, vl1} = Eigensystem[m - μ IdentityMatrix[3], -1, Method -> "Arnoldi"]

Out[3]= {{-0.1}, {{-0.5547, -0.83205, 0.}}}

In[4]:= λl1 = λl1 + μ

Out[4]= {2.}
```

Automatically shift and adjust the eigenvalue:

```wl
In[5]:= {λl2, vl2} = Eigensystem[m, -1, Method -> {"Arnoldi", "Shift" -> μ}]

Out[5]= {{2.}, {{-0.5547, -0.83205, -1.6653345369377348`*^-16}}}
```

---

Compute two smallest generalized eigenvalues and eigenvectors with the Arnoldi method:

```wl
In[1]:=
n = 10;
s = SparseArray[{{i_, i_} -> -2., {i_, j_} /; Abs[i - j] == 1 -> 1.}, {n, n}];
m = SparseArray[{{i_, i_} -> 1}, {n, n}];
Eigensystem[{s, m}, -2, Method -> "Arnoldi"]

Out[1]= {{-0.317493, -0.0810141}, {{0.23053, 0.387868, 0.422061, 0.322253, 0.120131, -0.120131, -0.322253, -0.422061, -0.387868, -0.23053}, {-0.120131, -0.23053, -0.322253, -0.387868, -0.422061, -0.422061, -0.387868, -0.322253, -0.23053, -0.120131}}}
```

##### "Banded" (1)

The banded method can be used for real symmetric or complex Hermitian machine-precision matrices. The method is most useful for finding all eigenpairs.

Compute the two largest eigenpairs for a banded matrix:

```wl
In[1]:= m = SparseArray[{{i_, i_} -> -2., {i_, j_} /; Abs[i - j] == 1 -> 1.}, {2000, 2000}];

In[3]:= MatrixPlot[m]

Out[3]= [image]

In[3]:= {λl, vl} = Eigensystem[m, 2, Method -> "Banded"];λl

Out[3]= {-4., -3.99999}
```

Check the solution:

```wl
In[4]:= Norm[m.Transpose[vl] - Transpose[vl].DiagonalMatrix[λl]]

Out[4]= 1.89984105602231`*^-14
```

##### "FEAST" (2)

The FEAST method can be used for real symmetric or complex Hermitian machine-precision matrices. The method is most useful for finding eigenvalues in a given interval.

The following suboptions can be specified for the method ``"FEAST"`` :

|                    |                                        |
| ------------------ | -------------------------------------- |
| "ContourPoints"    | select the number of contour points    |
| "Interval"         | an interval for finding eigenvalues    |
| "MaxIterations"    | the maximum number of refinement loops |
| "NumberOfRestarts" | the maximum number of restarts         |
| "SubspaceSize"     | the initial size of subspace           |
| "Tolerance"        | the tolerance to terminate refinement  |
| "UseBandedSolver"  | whether to use a banded solver         |

---

Compute eigenpairs in the interval $(-1.,-0.9)$ :

```wl
In[1]:= s = SparseArray[{{i_, i_} -> -2., {i_, j_} /; Abs[i - j] == 1 -> 1.}, {1000, 1000}];
```

Use ``"Interval"`` to specify the interval:

```wl
In[2]:= {λl, vl} = Eigensystem[s, Method -> {"FEAST", "Interval" -> {-1.0, -0.9}}]; λl

Out[2]= {-0.996378, -0.990954, -0.985539, -0.980135, -0.97474, -0.969356, -0.963982, -0.958618, -0.953264, -0.947921, -0.942588, -0.937265, -0.931953, -0.926651, -0.92136, -0.916079, -0.91081, -0.905551, -0.900302}
```

The interval end points $(-1.,-0.9)$ are not included in the interval FEAST finds eigenvalues in.

Check the solution:

```wl
In[3]:= Norm[s.Transpose[vl] - Transpose[vl].DiagonalMatrix[λl]]

Out[3]= 9.73215450864232`*^-14
```

#### Quartics (1)

A 4×4 matrix:

```wl
In[1]:= m = Table[{1, 2, -1, -2} ^ j, {j, 4}]

Out[1]= {{1, 2, -1, -2}, {1, 4, 1, 4}, {1, 8, -1, -8}, {1, 16, 1, 16}}
```

In general, for a 4×4 matrix, the result will be given in terms of ``Root`` objects:

```wl
In[2]:= Eigensystem[m]

Out[2]= {{Root[-288 + 240*#1 - 20*#1^3 + #1^4 & , 4, 0], Root[-288 + 240*#1 - 20*#1^3 + #1^4 & , 1, 0], Root[-288 + 240*#1 - 20*#1^3 + #1^4 & , 3, 0], Root[-288 + 240*#1 - 20*#1^3 + #1^4 & , 2, 0]}, {{(6640 - 1260 Root[-288 + 240*#1 - 20*#1^3 + #1^4 & , 4, ...  2, 0]^2 + Root[-288 + 240*#1 - 20*#1^3 + #1^4 & , 2, 0]^3/2172), (-19152 - 4164 Root[-288 + 240*#1 - 20*#1^3 + #1^4 & , 2, 0] + 2338 Root[-288 + 240*#1 - 20*#1^3 + #1^4 & , 2, 0]^2 - 107 Root[-288 + 240*#1 - 20*#1^3 + #1^4 & , 2, 0]^3/4344), 1}}}
```

You can get the result in terms of radicals using the ``Cubics`` and ``Quartics`` options:

```wl
In[3]:= Short[Eigensystem[m, Cubics -> True, Quartics -> True], 10]

Out[3]//Short= {{5 + Sqrt[25 + (38 2^2 / 3/(-225 + I Sqrt[59119])^1 / 3) + (2 (-225 + I Sqrt[59119]))^1 / 3] + (1/2) √(200 - (152 2^2 / 3/(-225 + I Sqrt[59119])^1 / 3) - 4 (2 (-225 + I Sqrt[59119]))^1 / 3 + (760/Sqrt[25 + (38 2^2 / 3/(-225 + I Sqrt[59119])^1 / 3) ... ])^1 / 3) + (2 (-225 + I Sqrt[59119]))^1 / 3] + (1/2) √(200 - (152 2^2 / 3/(-225 + I Sqrt[59119])^1 / 3) - 4 (2 (-225 + I Sqrt[59119]))^1 / 3 - (760/Sqrt[25 + (38 2^2 / 3/(-225 + I Sqrt[59119])^1 / 3) + (2 (-225 + I Sqrt[59119]))^1 / 3]))}, {«1»}}
```

### Applications (16)

#### The Geometry of Eigensystem (3)

Eigenvectors with positive eigenvalues point in the same direction when acted on by the matrix:

```wl
In[1]:=
m = {{1, 2}, {2, 1}};
{{λ1, λ2}, {v1, v2}} = Eigensystem[m]

Out[1]= {{3, -1}, {{1, 1}, {-1, 1}}}

In[2]:= Graphics[{{Thick, Arrow[{{0, 0}, v1}]}, {Red, Arrow[{{0, 0}, m.v1}]}}, Axes -> True]

Out[2]= [image]
```

Eigenvectors with negative eigenvalues point in the opposite direction when acted on by the matrix:

```wl
In[3]:= Graphics[{{Thick, Arrow[{{0, 0}, v2}]}, {Red, Arrow[{{0, 0}, m.v2}]}}, Axes -> True]

Out[3]= [image]
```

---

Consider the following matrix $a$ and its associated quadratic form $q=x^T.a.x$ :

```wl
In[1]:= a = {{-2, 2}, {2, 1}};

In[2]:= q = {x, y}.a.{x, y}//Expand

Out[2]= -2 x^2 + 4 x y + y^2
```

The eigenvectors are the axes of the hyperbolas defined by $q$ :

```wl
In[3]:= {λ, v} = Eigensystem[a]

Out[3]= {{-3, 2}, {{-2, 1}, {1, 2}}}
```

The sign of the eigenvalue corresponds to the sign of the right-hand side of the hyperbola equation:

```wl
In[4]:= ContourPlot[Table[q == n, {n, {-9, -4, -1, 1, 4, 9}}]//Evaluate, {x, -3, 3}, {y, -3, 3}, Epilog -> (Arrow[{{0, 0}, #}]& /@ v), PlotLegends -> "Expressions"]

Out[4]= [image]
```

---

Here is a positive-definite quadratic form in 3 dimensions:

```wl
In[1]:= q = 20 x ^ 2 - 16 x y + 23 y ^ 2 - 12 x z - 2 y z + 17 z ^ 2;
```

Plot the surface $q=10$ :

```wl
In[2]:= cp = ContourPlot3D[q, {x, -1, 1}, {y, -1, 1}, {z, -1, 1}, Contours -> {10}, Mesh -> False, ContourStyle -> Opacity[.5]]

Out[2]= [image]
```

Get the symmetric matrix for the quadratic form, using ``CoefficientArrays`` :

```wl
In[3]:= m = Normal[CoefficientArrays[q, {x, y, z}, Symmetric -> True][[3]]]

Out[3]= {{20, -8, -6}, {-8, 23, -1}, {-6, -1, 17}}
```

Numerically compute its eigenvalues and eigenvectors:

```wl
In[4]:= {vals, vecs} = Eigensystem[N[m]]

Out[4]= {{30.4331, 20.1091, 9.45779}, {{0.675631, -0.693505, -0.250149}, {-0.301381, -0.569468, 0.764771}, {-0.672824, -0.441313, -0.593759}}}
```

Show the principal axes of the ellipsoid:

```wl
In[5]:= Show[cp, Graphics3D[{Thickness[0.015], Green, Table[Line[{{0, 0, 0}, vecs[[i]] * Sqrt[10 / vals[[i]]]}], {i, 1, 3}]}]]

Out[5]= [image]
```

#### Diagonalization (5)

Diagonalize the following matrix as $m=p.d.p^{-1}$. First, compute $m$'s eigenvalues and eigenvectors:

```wl
In[1]:= m = {{9, -7, 3}, {12, -10, 3}, {16, -16, 1}};

In[2]:= Eigensystem[m]

Out[2]= {{-3, 2, 1}, {{-1, 0, 4}, {1, 1, 0}, {-3, -3, 1}}}
```

Construct a diagonal matrix $d$ from the eigenvalues and a matrix $p$ whose columns are the eigenvectors:

```wl
In[3]:= {d = DiagonalMatrix[First[%]], p = Transpose[Last[%]]}

Out[3]= {{{-3, 0, 0}, {0, 2, 0}, {0, 0, 1}}, {{-1, 1, -3}, {0, 1, -3}, {4, 0, 1}}}
```

Confirm the identity $m=p.d.p^{-1}$ :

```wl
In[4]:= m == p.d.Inverse[p]

Out[4]= True
```

Any function of the matrix can now be computed as $f(m)=p.f(d).p^{-1}$. For example, ``MatrixPower`` :

```wl
In[5]:= MatrixPower[m, k] == p . MatrixPower[d, k].Inverse[p]

Out[5]= True
```

Similarly, ``MatrixExp`` becomes trivial, requiring only exponentiating the diagonal elements of $d$ :

```wl
In[6]:= MatrixExp[m] == p . MatrixExp[d].Inverse[p]

Out[6]= True

In[7]:= MatrixExp[d]

Out[7]= {{(1/E^3), 0, 0}, {0, E^2, 0}, {0, 0, E}}
```

---

Let $T$ be the linear transformation whose standard matrix is given by the matrix $a$. Find a basis $B$ for $\mathbb{R}^4$ with the property that the representation of $T$ in that basis $[T]_B$ is diagonal:

```wl
In[1]:=
a = (|    |    |   |   |
| -- | -- | - | - |
| -6 | 4  | 0 | 9 |
| -3 | 0  | 1 | 6 |
| -1 | -2 | 1 | 0 |
| -4 | 4  | 0 | 7 |);
```

Find the eigenvalues and eigenvectors of $a$ :

```wl
In[2]:= {λ, v} = Eigensystem[a]

Out[2]= {{5, -2, -2, 1}, {{2, 1, -1, 2}, {6, -3, 0, 4}, {1, 1, 1, 0}, {2, -1, -7, 2}}}
```

Let $B$ consist of the eigenvectors, and let $b$ be the matrix whose columns are the elements of $B$ :

```wl
In[3]:= b = Transpose[v]

Out[3]= {{2, 6, 1, 2}, {1, -3, 1, -1}, {-1, 0, 1, -7}, {2, 4, 0, 2}}
```

$b$ converts from $B$-coordinates to standard coordinates. Its inverse converts in the reverse direction:

```wl
In[4]:= bInv = Inverse[b]

Out[4]= {{-(5/14), (2/7), (1/14), (3/4)}, {(2/21), -(1/7), (1/21), 0}, {(17/21), (2/7), -(2/21), -1}, {(1/6), 0, -(1/6), -(1/4)}}
```

Thus $[T]_B$ is given by $b^{-1}.a.b$, which is diagonal:

```wl
In[5]:= bInv . a .b//MatrixForm

Out[5]//MatrixForm=
(⁠|   |    |    |   |
| - | -- | -- | - |
| 5 | 0  | 0  | 0 |
| 0 | -2 | 0  | 0 |
| 0 | 0  | -2 | 0 |
| 0 | 0  | 0  | 1 |⁠)
```

Note that this is simply the diagonal matrix whose entries are the eigenvalues:

```wl
In[6]:= % == DiagonalMatrix[λ]

Out[6]= True
```

---

A real-valued symmetric matrix is orthogonally diagonalizable as $s=o.d.o^T$, with $d$ diagonal and real valued and $o$ orthogonal. Verify that the following matrix is symmetric and then diagonalize it:

```wl
In[1]:= (s = {{1, 4, -2}, {4, 5, -3}, {-2, -3, 2}})//MatrixForm

Out[1]//MatrixForm=
(⁠|    |    |    |
| -- | -- | -- |
| 1  | 4  | -2 |
| 4  | 5  | -3 |
| -2 | -3 | 2  |⁠)

In[2]:= Transpose[s] == s

Out[2]= True
```

Compute the eigenvalues and eigenvectors:

```wl
In[3]:= {λ, v} = Eigensystem[s]

Out[3]= {{Root[3 - 12*#1 - 8*#1^2 + #1^3 & , 3, 0], Root[3 - 12*#1 - 8*#1^2 + #1^3 & , 1, 0], Root[3 - 12*#1 - 8*#1^2 + #1^3 & , 2, 0]}, {{Root[12 - 50*#1 - 53*#1^2 + 4*#1^3 & , 1, 0], Root[-13 + 19*#1 + 19*#1^2 + 2*#1^3 & , 2, 0], 1}, {Root[12 - 50*#1 - 53*#1^2 + 4*#1^3 & , 3, 0], Root[-13 + 19*#1 + 19*#1^2 + 2*#1^3 & , 1, 0], 1}, {Root[12 - 50*#1 - 53*#1^2 + 4*#1^3 & , 2, 0], Root[-13 + 19*#1 + 19*#1^2 + 2*#1^3 & , 3, 0], 1}}}
```

The matrix $d$ has the eigenvalues on the diagonal:

```wl
In[4]:= d = DiagonalMatrix[λ]

Out[4]= {{Root[3 - 12*#1 - 8*#1^2 + #1^3 & , 3, 0], 0, 0}, {0, Root[3 - 12*#1 - 8*#1^2 + #1^3 & , 1, 0], 0}, {0, 0, Root[3 - 12*#1 - 8*#1^2 + #1^3 & , 2, 0]}}
```

For an orthogonal matrix, it is necessary to normalize the eigenvectors before placing them in columns:

```wl
In[5]:= o = Transpose[FullSimplify[Normalize /@ v]]

Out[5]= {{Root[-144 + 5396*#1^2 - 27213*#1^4 + 27213*#1^6 & , 2, 0], Root[-144 + 5396*#1^2 - 27213*#1^4 + 27213*#1^6 & , 6, 0], Root[-144 + 5396*#1^2 - 27213*#1^4 + 27213*#1^6 & , 4, 0]}, {Root[-676 + 7817*#1^2 - 27213*#1^4 + 27213*#1^6 & , 1, 0], Root[-67 ...  & , 2, 0], Root[-676 + 7817*#1^2 - 27213*#1^4 + 27213*#1^6 & , 4, 0]}, {Root[-16 + 4397*#1^2 - 27213*#1^4 + 27213*#1^6 & , 5, 0], Root[-16 + 4397*#1^2 - 27213*#1^4 + 27213*#1^6 & , 4, 0], Root[-16 + 4397*#1^2 - 27213*#1^4 + 27213*#1^6 & , 6, 0]}}
```

This matrix is orthogonal:

```wl
In[6]:= OrthogonalMatrixQ[o]

Out[6]= True
```

Verify that $s=o.d.o^T$ :

```wl
In[7]:= s == o.d.Transpose[o]//FullSimplify

Out[7]= True
```

---

A matrix is called normal if $n^{\dagger }.n=n.n^{\dagger }$. Normal matrices are the most general kind of matrix that can be  diagonalized by a unitary transformation. All real symmetric matrices $s$ are normal because both sides of the equality are simply $s.s$ :

```wl
In[1]:= TensorExpand[ConjugateTranspose[s].s == s.ConjugateTranspose[s] == s.s, Assumptions -> s∈Matrices[{n, n}, Reals, Symmetric[{1, 2}]]]

Out[1]= True
```

Show that the following matrix is normal, then diagonalize it:

```wl
In[2]:= (n = {{3, -1}, {1, 3}})//MatrixForm

Out[2]//MatrixForm=
(⁠|   |    |
| - | -- |
| 3 | -1 |
| 1 | 3  |⁠)

In[3]:= n.ConjugateTranspose[n] == ConjugateTranspose[n].n

Out[3]= True
```

Confirm using ``NormalMatrixQ`` :

```wl
In[4]:= NormalMatrixQ[n]

Out[4]= True
```

Find the eigenvalues and eigenvectors:

```wl
In[5]:= {λ, v} = Eigensystem[n]

Out[5]= {{3 + I, 3 - I}, {{I, 1}, {-I, 1}}}
```

Unlike a real symmetric matrix, the diagonal matrix for this matrix is complex valued:

```wl
In[6]:= (d = DiagonalMatrix[λ])//MatrixForm

Out[6]//MatrixForm=
(⁠|       |       |
| ----- | ----- |
| 3 + I | 0     |
| 0     | 3 - I |⁠)
```

Normalizing the eigenvectors and putting them in columns gives a unitary matrix:

```wl
In[7]:=
u = Transpose[Normalize /@ v];
UnitaryMatrixQ[u]

Out[7]= True
```

Confirm the diagonalization $n=u.d.u^{\dagger }$ :

```wl
In[8]:= n == u.d.ConjugateTranspose[u]

Out[8]= True
```

---

The eigensystem of a nondiagonalizable matrix:

```wl
In[1]:=
m = {{1, 1, 0}, {0, 1, 0}, {0, 0, 2}};
{λ, vecs} = Eigensystem[m]

Out[1]= {{2, 1, 1}, {{0, 0, 1}, {1, 0, 0}, {0, 0, 0}}}
```

The dimension of the span of all the eigenvectors is less than the number of eigenvalues:

```wl
In[2]:= Count[vecs, v_ /; Norm[v] > 0] < Length[λ]

Out[2]= True
```

Estimate the probability that a random 4×4 matrix of ones and zeros is not diagonalizable:

```wl
In[3]:= Block[{nd = 0, trials = 10 ^ 4}, Do[{λ, vecs} = Eigensystem[N[RandomInteger[1, {4, 4}]]];If[Count[vecs, v_ /; Norm[v] > 0] < Length[λ], nd++], {trials}];N[nd / trials]]

Out[3]= 0.248
```

#### Differential Equations and Dynamical Systems (4)

Solve the system of ODEs $x'=y$, $y'=z$, $z'=-2 x+y+2z$. First, construct the coefficient matrix $a$ for the right-hand side:

```wl
In[1]:= a = {{0, 1, 0}, {0, 0, 1}, {-2, 1, 2}};
```

Find the eigenvalues and eigenvectors:

```wl
In[2]:= {λ, v} = Eigensystem[a]

Out[2]= {{2, -1, 1}, {{1, 2, 4}, {1, -1, 1}, {1, 1, 1}}}
```

Construct a diagonal matrix whose entries are the exponential of $t \lambda$ :

```wl
In[3]:= d = DiagonalMatrix[Exp[t λ]]

Out[3]= {{E^2 t, 0, 0}, {0, E^-t, 0}, {0, 0, E^t}}
```

Construct the matrix whose columns are the corresponding eigenvectors:

```wl
In[4]:= p = Transpose[v]

Out[4]= {{1, 1, 1}, {2, -1, 1}, {4, 1, 1}}
```

The general solution is $p.d.p^{-1}.\left\{c_1,c_2,c_3\right\}$, for three arbitrary starting values:

```wl
In[5]:= p.d.Inverse[p] . {C[1], C[2], C[3]}

Out[5]= {((E^-t/3) + E^t - (E^2 t/3)) C[1] + (-(E^-t/2) + (E^t/2)) C[2] + ((E^-t/6) - (E^t/2) + (E^2 t/3)) C[3], (-(E^-t/3) + E^t - (2 E^2 t/3)) C[1] + ((E^-t/2) + (E^t/2)) C[2] + (-(E^-t/6) - (E^t/2) + (2 E^2 t/3)) C[3], ((E^-t/3) + E^t - (4 E^2 t/3)) C[1] + (-(E^-t/2) + (E^t/2)) C[2] + ((E^-t/6) - (E^t/2) + (4 E^2 t/3)) C[3]}
```

Verify the solution using ``DSolveValue`` :

```wl
In[6]:= Simplify[% == DSolveValue[{x'[t] == y[t], y'[t] == z[t], z'[t] == -2x[t] + y[t] + 2z[t]}, {x[t], y[t], z[t]}, t]]

Out[6]= True
```

---

Suppose a particle is moving in a planar force field and its position vector $x$ satisfies $x'=A x$ and $x(0)=x_0$, where $A$ and $x_0$ are as follows. Solve this initial problem for $t\geq 0$ :

```wl
In[1]:=
a = (|    |    |
| -- | -- |
| 4  | -5 |
| -2 | 1  |);

In[2]:=
Subscript[x, 0] = (⁠|     |
| --- |
| 1.9 |
| 3.6 |⁠);
```

First, compute the eigenvalues and corresponding eigenvectors of $A$ :

```wl
In[3]:= {λ, v} = Eigensystem[a]

Out[3]= {{6, -1}, {{-5, 2}, {1, 1}}}
```

The general solution of the system is $x(t)=c_1 v_1 e^{\lambda _1 t}+c_2 v_2 e^{\lambda _2 t}$. Use ``LinearSolve`` to determine the coefficients:

```wl
In[4]:= c = LinearSolve[Transpose[v], Subscript[x, 0]]

Out[4]= {0.242857, 3.11429}
```

Construct the appropriate linear combination of the eigenvectors:

```wl
In[5]:= x[t_] = (c Exp[λ t]).v

Out[5]= {3.11429 E^-t - 1.21429 E^6 t, 3.11429 E^-t + 0.485714 E^6 t}
```

Verify the solution using ``DSolveValue`` :

```wl
In[6]:= x[t] == DSolveValue[{{u1'[t], u2'[t]} == a.{u1[t], u2[t]}, {u1[0], u2[0]} == Subscript[x, 0]}, {u1[t], u2[t]}, t] //Simplify//Chop

Out[6]= True
```

---

Produce the general solution of the dynamical system $x_{k+1}=A x_k$ when $A$ is the following stochastic matrix:

```wl
In[1]:=
a = (|     |     |     |
| --- | --- | --- |
| .90 | .01 | .09 |
| .01 | .90 | .01 |
| .09 | .09 | .90 |);
```

Find the eigenvalues and eigenvectors, using ``Chop`` to discard small numerical errors:

```wl
In[2]:= {{Subscript[λ, 1], Subscript[λ, 2], Subscript[λ, 3]}, {Subscript[v, 1], Subscript[v, 2], Subscript[v, 3]}} = Chop[Eigensystem[a]]

Out[2]= {{1., 0.89, 0.81}, {{-0.670078, -0.139906, -0.728986}, {0.707107, -0.707107, 0}, {-0.707107, 0, 0.707107}}}
```

The general solution is an arbitrary linear combination of terms of the form $\lambda _i{}^kv_i$ :

```wl
In[3]:= x[k_] = Sum[C[i] Subsuperscript[λ, i, k]Subscript[v, i], {i, 3}]

Out[3]= {-0.670078 1.^k C[1] + 0.707107 0.89^k C[2] - 0.707107 0.81^k C[3], -0.139906 1.^k C[1] - 0.707107 0.89^k C[2], -0.728986 1.^k C[1] + 0.707107 0.81^k C[3]}
```

Verify that $x$ satisfies the dynamical equation up to numerical rounding:

```wl
In[4]:= x[k] == a.x[k - 1]//Simplify//Chop

Out[4]= True
```

---

The Lorenz equations:

```wl
In[1]:= leqn = {Derivative[1][x][t] == -3(x[t] - y[t]), Derivative[1][y][t] == -x[t] z[t] + 26 x[t] - y[t], Derivative[1][z][t] == x[t] y[t] - z[t]};
```

Find the Jacobian for the right-hand side of the equations:

```wl
In[2]:= lj = D[leqn[[All, 2]], {{x[t], y[t], z[t]}}]

Out[2]= {{-3, 3, 0}, {26 - z[t], -1, -x[t]}, {y[t], x[t], -1}}
```

Find the equilibrium points:

```wl
In[3]:= eqp = Solve[leqn /. {x'[t] -> 0, y'[t] -> 0, z'[t] -> 0}, {x[t], y[t], z[t]}]

Out[3]= {{x[t] -> -5, y[t] -> -5, z[t] -> 25}, {x[t] -> 0, y[t] -> 0, z[t] -> 0}, {x[t] -> 5, y[t] -> 5, z[t] -> 25}}
```

Find the eigenvalues and eigenvectors of the Jacobian at the equilibrium point in the first octant:

```wl
In[4]:= {vals, vecs} = Eigensystem[N[lj /. eqp[[3]]]]

Out[4]= {{0.0455218  + 5.42784 I, 0.0455218  - 5.42784 I, -5.09104}, {{-0.236966 - 0.179728 I, 0.0846171  - 0.611193 I, -0.728579 + 0. I}, {-0.236966 + 0.179728 I, 0.0846171  + 0.611193 I, -0.728579 + 0. I}, {-0.784959 + 0. I, 0.547128  + 0. I, 0.290673  + 0. I}}}
```

A function that integrates backward from a small perturbation of ``pt`` in the direction ``dir`` :

```wl
In[5]:= st[pt_, dir_] := First[{x[t], y[t], z[t]} /. NDSolve[{leqn, Thread[{x[0], y[0], z[0]} == pt + 10^-6dir]}, {x, y, z}, {t, 0, -4}]];
```

Show the stable curve for the equilibrium point on the right:

```wl
In[6]:= pt = {x[t], y[t], z[t]} /. eqp[[3]];dir = vecs[[3]];sr = ParametricPlot3D[Evaluate[{st[pt, dir], st[pt, -dir]}], {t, 0, -4}, PlotStyle -> Table[{Thickness[0.015], Red}, {2}]]

Out[6]= [image]
```

Find the stable curve for the equilibrium point on the left:

```wl
In[7]:=
pt = {x[t], y[t], z[t]} /. eqp[[2]];{vals, vecs} = Eigensystem[N[lj /. eqp[[2]]]];
dir = vecs[[3]];sl = ParametricPlot3D[Evaluate[{st[pt, dir], st[pt, -dir]}], {t, 0, -4}, PlotStyle -> Table[{Thickness[0.015], Green}, {2}]];
```

Show the stable curves along with a solution of the Lorenz equations:

```wl
In[8]:= Show[{ParametricPlot3D[Evaluate[{x[t], y[t], z[t]} /. NDSolve[{leqn, x[0] == 5, y[0] == 0, z[0] == 25}, {x, y, z}, {t, 0, 100}]], {t, 0, 100}], sr, sl}, PlotRange -> {{-25, 25}, {-25, 25}, {0, 50}}]

Out[8]= [image]
```

#### Physics (4)

In quantum mechanics, states are represented by complex unit vectors and physical quantities by Hermitian linear operators. The eigenvalues represent possible observations and the squared modulus of the components with respect to eigenvectors the probabilities of those observations. For the spin operator $\sigma$ and state $\psi$ given, find the possible observations and their probabilities:

```wl
In[1]:=
σ = (ℏ/2)(|   |    |
| - | -- |
| 0 | -I |
| I | 0  |);ψ = (1/Sqrt[5])(⁠|     |
| --- |
| 1   |
| 2 I |⁠);
```

Computing the eigensystem, the possible observations are $\pm \frac{\hbar }{2}$ :

```wl
In[2]:= {λ, v} = Eigensystem[σ]

Out[2]= {{-(ℏ/2), (ℏ/2)}, {{I, 1}, {-I, 1}}}
```

Normalize the eigenvectors in order to compute proper projections:

```wl
In[3]:= {e1, e2} = Normalize /@ v

Out[3]= {{(I/Sqrt[2]), (1/Sqrt[2])}, {-(I/Sqrt[2]), (1/Sqrt[2])}}
```

The relative probabilities are $10 \%$ for $-\frac{\hbar }{2}$ and $90 \%$ for $\frac{\hbar }{2}$ :

```wl
In[4]:= {Abs[e1\[Conjugate].ψ]^2, Abs[e2\[Conjugate].ψ]^2}//Simplify

Out[4]= {(1/10), (9/10)}
```

---

In quantum mechanics, the energy operator is called the Hamiltonian $\mathcal{H}$, and a state with energy $\mathcal{E}$ evolves according to the Schrödinger equation $\psi (t)=\exp \left(-\frac{i \mathcal{E} t}{\hbar }\right)\psi (0)$. Given the Hamiltonian for a spin-1 particle in constant magnetic field in the $y$ direction, find the state at time $t$ of a particle that is initially in the state $\{1,0,0\}$ representing $S_z=-\hbar$ :

```wl
In[1]:=
ℋ = (Subscript[ω, 0]ℏ/Sqrt[2]) (|   |    |    |
| - | -- | -- |
| 0 | -I | 0  |
| I | 0  | -I |
| 0 | I  | 0  |);Subscript[ψ, 0] = {1, 0, 0};
```

Computing the eigensystem, the energy levels are $\pm \hbar \omega _0$ and $0$ :

```wl
In[2]:= {{Subscript[ℰ, 1], Subscript[ℰ, 2], Subscript[ℰ, 3]}, v} = Eigensystem[ℋ]

Out[2]= {{-ℏ Subscript[ω, 0], ℏ Subscript[ω, 0], 0}, {{-1, I Sqrt[2], 1}, {-1, -I Sqrt[2], 1}, {1, 0, 1}}}
```

Normalize the eigenvectors:

```wl
In[3]:= {Subscript[e, 1], Subscript[e, 2], Subscript[e, 3]} = Normalize /@ v

Out[3]= {{-(1/2), (I/Sqrt[2]), (1/2)}, {-(1/2), -(I/Sqrt[2]), (1/2)}, {(1/Sqrt[2]), 0, (1/Sqrt[2])}}
```

The state at time $t$ is the sum of each eigenstate evolving according to the Schrödinger equation:

```wl
In[4]:= ψ[t_] = Underoverscript[∑, i = 1, 3]Exp[-(I Subscript[ℰ, i]t/ℏ)]Subscript[e, i]\[Conjugate].Subscript[ψ, 0] Subscript[e, i] //FullSimplify

Out[4]= {Cos[(t Subscript[ω, 0]/2)]^2, (Sin[t Subscript[ω, 0]]/Sqrt[2]), Sin[(t Subscript[ω, 0]/2)]^2}
```

---

The moment of inertia is a real symmetric matrix that describes the resistance of a rigid body to rotating in different directions. The eigenvalues of this matrix are called the principal moments of inertia, and the corresponding eigenvectors (which are necessarily orthogonal) the principal axes. Find the principal moments of inertia and principal axes for the following tetrahedron:

```wl
In[1]:= tet = Tetrahedron[{{-1.25, -1, -0.75}, {3.75, -1, -0.75}, {-1.25, 3., -0.75}, {-1.25, -1, 2.25}}]

Out[1]= Tetrahedron[{{-1.25, -1, -0.75}, {3.75, -1, -0.75}, {-1.25, 3., -0.75}, {-1.25, -1, 2.25}}]
```

First compute the moment of inertia:

```wl
In[2]:= (ℐ = MomentOfInertia[tet])//MatrixForm

Out[2]//MatrixForm=
(⁠|       |       |        |
| ----- | ----- | ------ |
| 9.375 | 2.5   | 1.875  |
| 2.5   | 12.75 | 1.5    |
| 1.875 | 1.5   | 15.375 |⁠)
```

Compute the principal moments and axes:

```wl
In[3]:= {m, a} = Eigensystem[ℐ]

Out[3]= {{17.0894, 12.5, 7.91061}, {{0.350985, 0.480105, 0.803933}, {-0.26963, -0.770371, 0.577778}, {0.896721, -0.419556, -0.140938}}}
```

Verify that the axes are orthogonal:

```wl
In[4]:= a.a\[Transpose]//Chop

Out[4]= {{1., 0, 0}, {0, 1., 0}, {0, 0, 1.}}
```

The center of mass of the tetrahedron is at the origin:

```wl
In[5]:= RegionCentroid[tet]

Out[5]= {0., 0., 0.}
```

Visualize the tetrahedron and its principal axes:

```wl
In[6]:= Graphics3D[{{Opacity[.5], tet}, Arrow[{{0, 0, 0}, 3#}]& /@ a}]

Out[6]= [image]
```

---

A generalized eigensystem can be used to find normal modes of coupled oscillations that decouple the terms. Consider the system shown in the diagram:

[image]

By Hooke's law it obeys $-k_1 y_1+k_2 \left(y_2-y_1\right)=m_1 y_1''$, $k_2 \left(y_1-y_2\right)=m_2 y_2''$. Substituting in the generic solution $y(t)=x \exp (i t \omega )$ gives rise to the matrix equation $k.x=\omega ^2 m.x$, with the stiffness matrix $k$ and mass matrix $m$ as follows:

```wl
In[1]:=
k = {{k1 + k2, -k2}, {-k2, k2}};
m = {{m1, 0}, {0, m2}};
```

Find the eigenfrequencies and normal modes if $k_1=9.5$, $k_2=4.5$, $m_1=2$ and $m_2=3$ :

```wl
In[2]:=
k = k /. {k1 -> 9.5, k2 -> 4.5};
m = m /. {m1 -> 2, m2 -> 3};
```

Solve the generalized eigenvalue value problem:

```wl
In[3]:= {λ, v} = Eigensystem[{k, m}]

Out[3]= {{7.55719, 0.942811}, {{-0.970679, 0.240379}, {-0.348212, -0.937416}}}
```

The eigenfrequencies $\omega$ are the square roots of the eigenvalues:

```wl
In[4]:= {ω1, ω2} = Sqrt[λ]

Out[4]= {2.74903, 0.970984}
```

Construct the normal mode solutions as a generalized eigenvector times the corresponding exponential:

```wl
In[5]:= sol1[t_] = First[v]Exp[I ω1 t];

In[6]:= sol2[t_] = Last[v]Exp[I ω2 t];
```

Verify that both satisfy the differential equation for the system:

```wl
In[7]:= -k.sol1[t] == m.sol1''[t] && -k.sol2[t] == m.sol2''[t]//FullSimplify

Out[7]= True
```

### Properties & Relations (17)

The eigenvectors returned for a numerical matrix are unit vectors:

```wl
In[1]:= {λ, v} = Eigensystem[{{1., 2.}, {2., 1.}}]

Out[1]= {{3., -1.}, {{0.707107, 0.707107}, {-0.707107, 0.707107}}}

In[2]:= Norm /@ v

Out[2]= {1., 1.}
```

The eigenvectors returned for exact and symbolic matrices are typically not unit vectors:

```wl
In[3]:= {λ, v} = Eigensystem[{{1, 2}, {2, 1}}]

Out[3]= {{3, -1}, {{1, 1}, {-1, 1}}}

In[4]:= Norm /@ v

Out[4]= {Sqrt[2], Sqrt[2]}
```

---

``Eigensystem[m]`` is effectively equivalent to ``{Eigenvalues[m], Eigenvectors[m]}`` :

```wl
In[1]:=
m = RandomReal[1, {3, 3}];
Eigensystem[m] == {Eigenvalues[m], Eigenvectors[m]}

Out[1]= True
```

If both eigenvectors and eigenvalues are needed, it is generally more efficient to just call ``Eigensystem`` :

```wl
In[2]:= m = RandomReal[1, {1000, 1000}];

In[3]:= AbsoluteTiming[Eigensystem[m];]

Out[3]= {0.78775, Null}

In[4]:= AbsoluteTiming[{Eigenvalues[m], Eigenvectors[m]};]

Out[4]= {1.35229, Null}
```

---

Any square matrix satisfies its similarity relation:

```wl
In[1]:=
m = {{1, 2, 3}, {4, 5, 6}, {7, 8, 9}};
{vals, vecs} = Eigensystem[m];

In[2]:= Simplify[m.Transpose[vecs] - Transpose[vecs].DiagonalMatrix[vals]]

Out[2]= {{0, 0, 0}, {0, 0, 0}, {0, 0, 0}}
```

---

Any pair of square matrices satisfies the generalized similarity relation for finite eigenvalues:

```wl
In[1]:=
m = {{1, 2, 3}, {4, 5, 6}, {7, 8, 9}};
a = {{1, 1, 1}, {1, 2, 4}, {1, 3, 9}};
{gvals, gvecs} = Eigensystem[{m, a}];

In[2]:= m.Transpose[gvecs] - a.Transpose[gvecs].DiagonalMatrix[gvals]//Simplify

Out[2]= {{0, 0, 0}, {0, 0, 0}, {0, 0, 0}}
```

---

An infinite generalized eigenvalue of $\{m,a\}$ corresponds to an eigenvector of $m$ that is in the kernel of $a$ :

```wl
In[1]:=
m = {{1, 2}, {3, 2}};
a = {{1, 1}, {1, 1}};
Eigensystem[{m, a}]

Out[1]= {{∞, 2}, {{-1, 1}, {0, 1}}}
```

The vector $\{-1,1\}$ is an eigenvector of $m$ with eigenvalue $-1$ :

```wl
In[2]:= m.{-1, 1}

Out[2]= {1, -1}
```

It is also in the null space of $a$ :

```wl
In[3]:= a.{-1, 1}

Out[3]= {0, 0}
```

---

A matrix ``m`` has a complete set of eigenvectors iff ``DiagonalizableMatrixQ[m]`` is ``True`` :

```wl
In[1]:=
diag = {{17, -1 + Sqrt[3], -1 - Sqrt[3]}, {-1 + Sqrt[3], 14 + 2 Sqrt[3], 2}, {-1 - Sqrt[3], 2, 14 - 2 Sqrt[3]}};
{DiagonalizableMatrixQ[diag], Eigensystem[diag]}

Out[1]= {True, {{18, 18, 9}, {{-1 - Sqrt[3], 0, 1}, {-1 + Sqrt[3], 1, 0}, {(1/2) (-1 + Sqrt[3]), -2 + Sqrt[3], 1}}}}
```

The following matrix is missing an eigenvector for the space $\lambda =18$, symbolized by $\{0,0,0\}$ :

```wl
In[2]:=
nondiag = {{15, 3 + 3 Sqrt[3], 0}, {0, 15 + 3 Sqrt[3], 3}, {3 - 3 Sqrt[3], 0, 15 - 3 Sqrt[3]}};
{DiagonalizableMatrixQ[nondiag], Eigensystem[nondiag]}

Out[2]= {False, {{18, 18, 9}, {{-2 - Sqrt[3], (1/2) (-1 - Sqrt[3]), 1}, {0, 0, 0}, {(1/2) (-1 + Sqrt[3]), -2 + Sqrt[3], 1}}}}
```

---

For a diagonalizable matrix, ``Eigensystem`` reduces function application to application to eigenvalues:

```wl
In[1]:=
m = N[{{0, 1, 0, -1}, {-1, 0, 1, 0}, {0, -1, 0, 1}, {1, 0, -1, 0}}];
{vals, vecs} = Eigensystem[m];
```

Compute the matrix exponential using diagonalization:

```wl
In[2]:= Chop[Transpose[vecs].DiagonalMatrix[Exp[vals]].Inverse[Transpose[vecs]]]

Out[2]= {{0.291927, 0.454649, 0.708073, -0.454649}, {-0.454649, 0.291927, 0.454649, 0.708073}, {0.708073, -0.454649, 0.291927, 0.454649}, {0.454649, 0.708073, -0.454649, 0.291927}}
```

Compute the matrix exponential using ``MatrixExp`` :

```wl
In[3]:= MatrixExp[m]

Out[3]= {{0.291927, 0.454649, 0.708073, -0.454649}, {-0.454649, 0.291927, 0.454649, 0.708073}, {0.708073, -0.454649, 0.291927, 0.454649}, {0.454649, 0.708073, -0.454649, 0.291927}}
```

Note that this is not simply the exponential of each entry:

```wl
In[4]:= Exp[m]

Out[4]= {{1., 2.71828, 1., 0.367879}, {0.367879, 1., 2.71828, 1.}, {1., 0.367879, 1., 2.71828}, {2.71828, 1., 0.367879, 1.}}
```

---

For an invertible matrix $m$, $m$ and $m^{-1}$ have the same eigenvectors and reciprocal eigenvalues:

```wl
In[1]:=
m = {{0, 3, 7}, {0, -3, 1}, {-4, 2, -8}};
{λ, v} = Eigensystem[m]

Out[1]= {{Root[96 + 50*#1 + 11*#1^2 + #1^3 & , 1, 0], Root[96 + 50*#1 + 11*#1^2 + #1^3 & , 3, 0], Root[96 + 50*#1 + 11*#1^2 + #1^3 & , 2, 0]}, {{Root[97 + 185*#1 + 128*#1^2 + 36*#1^3 & , 1, 0], Root[1 + 2*#1 + 11*#1^2 + 18*#1^3 & , 1, 0], 1}, {Root[97 + 185*#1 + 128*#1^2 + 36*#1^3 & , 2, 0], Root[1 + 2*#1 + 11*#1^2 + 18*#1^3 & , 2, 0], 1}, {Root[97 + 185*#1 + 128*#1^2 + 36*#1^3 & , 3, 0], Root[1 + 2*#1 + 11*#1^2 + 18*#1^3 & , 3, 0], 1}}}
```

Because eigenvalues are sorted by absolute value, this gives the same values but in the opposite order:

```wl
In[2]:=
{λInv, vInv} = Eigensystem[Inverse[m]];
{(1/λInv), vInv}//FullSimplify

Out[2]= {{Root[96 + 50*#1 + 11*#1^2 + #1^3 & , 2, 0], Root[96 + 50*#1 + 11*#1^2 + #1^3 & , 3, 0], Root[96 + 50*#1 + 11*#1^2 + #1^3 & , 1, 0]}, {{Root[97 + 185*#1 + 128*#1^2 + 36*#1^3 & , 3, 0], Root[1 + 2*#1 + 11*#1^2 + 18*#1^3 & , 3, 0], 1}, {Root[97 + 185*#1 + 128*#1^2 + 36*#1^3 & , 2, 0], Root[1 + 2*#1 + 11*#1^2 + 18*#1^3 & , 2, 0], 1}, {Root[97 + 185*#1 + 128*#1^2 + 36*#1^3 & , 1, 0], Root[1 + 2*#1 + 11*#1^2 + 18*#1^3 & , 1, 0], 1}}}
```

---

For an analytic function $f$, the eigenvectors of $m$ are also eigenvectors of $f(m)$ with eigenvalue $f(\lambda )$ :

```wl
In[1]:=
m = {{1, 1, 0}, {-1, 1, 0}, {0, 0, 5}};
{λ, v} = Eigensystem[m]

Out[1]= {{5, 1 + I, 1 - I}, {{0, 0, 1}, {-I, 1, 0}, {I, 1, 0}}}
```

For example, the $m.m$ has the same eigenvectors with squared eigenvalues:

```wl
In[2]:= Eigensystem[m.m] == {λ^2, v}

Out[2]= True
```

Similarly, the eigenvalues of $\exp (m)$ are $\exp (\lambda )$ :

```wl
In[3]:= Eigensystem[MatrixExp[m]] == {Exp[λ], v}//FullSimplify

Out[3]= True
```

---

The eigenvalues of a real symmetric matrix are real and its eigenvectors are orthogonal:

```wl
In[1]:= (s = {{1, 4, -2}, {4, 0, -3}, {-2, -3, 2}})//MatrixForm

Out[1]//MatrixForm=
(⁠|    |    |    |
| -- | -- | -- |
| 1  | 4  | -2 |
| 4  | 0  | -3 |
| -2 | -3 | 2  |⁠)
```

This matrix is symmetric:

```wl
In[2]:= SymmetricMatrixQ[s]

Out[2]= True
```

By inspection, the eigenvalues are real:

```wl
In[3]:= {λ, v} = Eigensystem[s]

Out[3]= {{7, -2 - Sqrt[3], -2 + Sqrt[3]}, {{-1, -1, 1}, {-(-1 - 3 Sqrt[3]/-4 + Sqrt[3]), -(5 + 2 Sqrt[3]/-4 + Sqrt[3]), 1}, {-(1 - 3 Sqrt[3]/4 + Sqrt[3]), -(-5 + 2 Sqrt[3]/4 + Sqrt[3]), 1}}}
```

Confirm the eigenvectors are orthogonal to each other:

```wl
In[4]:= v.Transpose[v]//Simplify

Out[4]= {{3, 0, 0}, {0, 6 (2 + Sqrt[3]), 0}, {0, 0, -6 (-2 + Sqrt[3])}}
```

---

The eigenvalues of a real antisymmetric matrix are imaginary and its eigenvectors are orthogonal:

```wl
In[1]:= (a = {{0, 4, -2}, {-4, 0, -3}, {2, 3, 0}})//MatrixForm

Out[1]//MatrixForm=
(⁠|    |   |    |
| -- | - | -- |
| 0  | 4 | -2 |
| -4 | 0 | -3 |
| 2  | 3 | 0  |⁠)
```

This matrix is antisymmetric:

```wl
In[2]:= AntisymmetricMatrixQ[a]

Out[2]= True
```

By inspection, the eigenvalues are imaginary:

```wl
In[3]:= {λ, v} = Eigensystem[a]

Out[3]= {{I Sqrt[29], -I Sqrt[29], 0}, {{(2/13) I (-6 I + Sqrt[29]), (1/13) (-8 + 3 I Sqrt[29]), 1}, {-(2/13) I (6 I + Sqrt[29]), (1/13) (-8 - 3 I Sqrt[29]), 1}, {-3, 2, 4}}}
```

Confirm the eigenvectors are orthogonal to each other:

```wl
In[4]:= v.ConjugateTranspose[v]//Simplify

Out[4]= {{(58/13), 0, 0}, {0, (58/13), 0}, {0, 0, 29}}
```

---

The eigenvalues of a unitary matrix lie on the unit circle and its eigenvectors are orthogonal:

```wl
In[1]:= u = {{1 / Sqrt[2], I / Sqrt[2]}, {I / Sqrt[2], 1 / Sqrt[2]}};

In[2]:= UnitaryMatrixQ[u]

Out[2]= True
```

Compute the eigenvalues and eigenvectors:

```wl
In[3]:= {λ, v} = Eigensystem[u]

Out[3]= {{(1 + I/Sqrt[2]), (1 - I/Sqrt[2])}, {{1, 1}, {-1, 1}}}
```

Confirm that the eigenvalues lie on the unit circle:

```wl
In[4]:= Abs[λ]

Out[4]= {1, 1}
```

Confirm the eigenvectors are orthogonal to each other:

```wl
In[5]:= v.ConjugateTranspose[v]//Simplify

Out[5]= {{2, 0}, {0, 2}}
```

---

The eigenvectors of any normal matrix are orthogonal:

```wl
In[1]:= n = {{1, 2, -1}, {-1, 1, 2}, {2, -1, 1}};

In[2]:= NormalMatrixQ[n]

Out[2]= True
```

The eigenvalues can be arbitrary:

```wl
In[3]:= {λ, v} = Eigensystem[n]

Out[3]= {{(1/2) (1 + 3 I Sqrt[3]), (1/2) (1 - 3 I Sqrt[3]), 2}, {{(1/2) (-1 + I Sqrt[3]), (1/2) (-1 - I Sqrt[3]), 1}, {(1/2) (-1 - I Sqrt[3]), (1/2) (-1 + I Sqrt[3]), 1}, {1, 1, 1}}}
```

But the eigenvectors are orthogonal:

```wl
In[4]:= v.ConjugateTranspose[v]//Simplify

Out[4]= {{3, 0, 0}, {0, 3, 0}, {0, 0, 3}}
```

---

``SingularValueDecomposition[m]`` is built from the eigensystems of $m.m^{\dagger }$ and $m^{\dagger }.m$ :

```wl
In[1]:=
m = (⁠|      |     |      |
| ---- | --- | ---- |
| -5.2 | 0   | -5   |
| 1.7  | 5.4 | -4.3 |⁠);

In[2]:= {u, Σ, v} = SingularValueDecomposition[m];
```

Compute the eigensystem of $m.m^{\dagger }$ :

```wl
In[3]:= {Subscript[λ, L], Subscript[vec, L]} = Eigensystem[m . ConjugateTranspose[m]];
```

The columns of $u$ are the eigenvectors:

```wl
In[4]:= {u//MatrixForm, Transpose[Subscript[vec, L]]//MatrixForm}

Out[4]=
{(⁠|          |           |
| -------- | --------- |
| 0.727715 | 0.68588   |
| 0.68588  | -0.727715 |⁠), (⁠|           |           |
| --------- | --------- |
| -0.727715 | 0.68588   |
| -0.68588  | -0.727715 |⁠)}
```

Compute the eigensystem of $m^{\dagger }.m$ :

```wl
In[5]:= {Subscript[λ, R], Subscript[vec, R]} = Eigensystem[ConjugateTranspose[m].m];
```

The columns of $v$ are the eigenvectors:

```wl
In[6]:= {v//MatrixForm, Transpose[Subscript[vec, R]]//MatrixForm}

Out[6]=
{(⁠|           |            |           |
| --------- | ---------- | --------- |
| -0.327336 | -0.773103  | -0.543289 |
| 0.463069  | -0.632437  | 0.620959  |
| -0.823661 | -0.0483179 | 0.56502   |⁠), (⁠|           |           |           |
| --------- | --------- | --------- |
| 0.327336  | 0.773103  | -0.543289 |
| -0.463069 | 0.632437  | 0.620959  |
| 0.823661  | 0.0483179 | 0.56502   |⁠)}
```

Since $m$ has fewer rows than columns, the diagonal entries of $\Sigma$ are $\sqrt{\lambda _L}$ :

```wl
In[7]:= {Diagonal[Σ], Sqrt[Subscript[λ, L]]}

Out[7]= {{7.99826, 6.21352}, {7.99826, 6.21352}}
```

---

Consider a matrix $$m$$ with a complete set of eigenvectors:

```wl
In[1]:=
m = {{-(1/Sqrt[2]), -(1/Sqrt[2])}, {(1/Sqrt[2]), -(1/Sqrt[2])}};
{λ, v} = Eigensystem[m]

Out[1]= {{-(1 - I/Sqrt[2]), -(1 + I/Sqrt[2])}, {{I, 1}, {-I, 1}}}
```

``JordanDecomposition[m]`` returns matrices $$\{\text{\textit{$s$}},\text{\textit{$j$}}\}$$ built from eigenvalues and eigenvectors:

```wl
In[2]:= {s, j} = JordanDecomposition[m]

Out[2]= {{{I, -I}, {1, 1}}, {{-(1 - I/Sqrt[2]), 0}, {0, -(1 + I/Sqrt[2])}}}
```

The $$j$$ matrix is diagonal with the eigenvalues as entries:

```wl
In[3]:= j == DiagonalMatrix[λ]

Out[3]= True
```

The $$s$$ matrix has the corresponding eigenvectors as its columns:

```wl
In[4]:= s == Transpose[v]

Out[4]= True
```

---

``SchurDecomposition[n, RealBlockDiagonalForm -> False]`` for a numerical normal matrix $$n$$ :

```wl
In[1]:=
n = {{1., 3., -1.}, {-1., 1., 3.}, {3., -1., 1.}};
NormalMatrixQ[n]

Out[1]= True

In[2]:= {q, t} = SchurDecomposition[n, RealBlockDiagonalForm -> False]//Chop

Out[2]= {{{0.573031  - 0.070486 I, 0.478338  - 0.323305 I, 0.290579  - 0.498896 I}, {-0.225473 + 0.531503 I, -0.51916 - 0.2526 I, 0.290579  - 0.498896 I}, {-0.347558 - 0.461017 I, 0.0408219  + 0.575905 I, 0.290579  - 0.498896 I}}, {{0.  + 3.4641 I, 0, 0}, {0, 0.  - 3.4641 I, 0}, {0, 0, 3.}}}
```

The matrices ``{q, t}`` are built from eigenvalues and eigenvectors:

```wl
In[3]:= {λ, v} = Eigensystem[n]//Chop

Out[3]= {{0.  + 3.4641 I, 0.  - 3.4641 I, 3.}, {{0.288675  + 0.5 I, -0.57735, 0.288675  - 0.5 I}, {0.288675  - 0.5 I, -0.57735, 0.288675  + 0.5 I}, {0.57735, 0.57735, 0.57735}}}
```

The ``t`` matrix is diagonal and with eigenvalue entries, possibly in a different order from ``Eigensystem`` :

```wl
In[4]:= Diagonal[t] == λ

Out[4]= True
```

To verify that ``q`` has eigenvectors as columns, set the first entry of each vector to ``1.`` to eliminate phase differences between ``q`` and ``v`` :

```wl
In[5]:= Transpose[(#1/First[q])& /@ q] - (v / First  /@ v)//Chop

Out[5]= {{0, 0, 0}, {0, 0, 0}, {0, 0, 0}}
```

---

If matrices share a dimension‐$d$ null space, $d$ of their generalized eigenvalues will be ``Indeterminate`` and the eigenvector list will be padded with zeros:

```wl
In[1]:= a = {{1, 1, 1}, {1, 1, 1}, {1, 1, 1}};b = {{1, 2, 3}, {4, 5, 6}, {7, 8, 9}};

In[2]:= NullSpace[a]

Out[2]= {{-1, 0, 1}, {-1, 1, 0}}
```

Two generalized eigenvalues of $\{a,a\}$ are ``Indeterminate`` and produce zero vectors:

```wl
In[3]:= Eigensystem[{a, a}]
```

Eigensystem::geinsl1: Warning: a solution for the generalized eigenproblem may be incorrect.

```wl
Out[3]= {{1, Indeterminate, Indeterminate}, {{0, 0, 1}, {0, 0, 0}, {0, 0, 0}}}
```

The matrix $b$ has a one-dimensional null space:

```wl
In[4]:= NullSpace[b]

Out[4]= {{1, -2, 1}}
```

It lies in the null space of $a$ :

```wl
In[5]:= MatrixRank[Join[NullSpace[a], NullSpace[b]]]

Out[5]= 2
```

Thus, one generalized eigenvalue of $\{a,b\}$ is ``Indeterminate`` and produces a zero vector:

```wl
In[6]:= Eigensystem[{a, b}]
```

Eigensystem::geinsl1: Warning: a solution for the generalized eigenproblem may be incorrect.

```wl
Out[6]= {{0, 0, Indeterminate}, {{-1, 0, 1}, {-1, 1, 0}, {0, 0, 0}}}
```

### Possible Issues (5)

The general symbolic case very quickly gets very complicated:

```wl
In[1]:= Array[Subscript[a, ##]&, {2, 2}]

Out[1]= {{Subscript[a, 1, 1], Subscript[a, 1, 2]}, {Subscript[a, 2, 1], Subscript[a, 2, 2]}}

In[2]:= Eigensystem[%]

Out[2]= {{(1/2) (Subscript[a, 1, 1] + Subscript[a, 2, 2] - Sqrt[Subsuperscript[a, 1, 1, 2] + 4 Subscript[a, 1, 2] Subscript[a, 2, 1] - 2 Subscript[a, 1, 1] Subscript[a, 2, 2] + Subsuperscript[a, 2, 2, 2]]), (1/2) (Subscript[a, 1, 1] + Subscript[a, 2, 2] +  ... /2 Subscript[a, 2, 1]), 1}, {-(-Subscript[a, 1, 1] + Subscript[a, 2, 2] - Sqrt[Subsuperscript[a, 1, 1, 2] + 4 Subscript[a, 1, 2] Subscript[a, 2, 1] - 2 Subscript[a, 1, 1] Subscript[a, 2, 2] + Subsuperscript[a, 2, 2, 2]]/2 Subscript[a, 2, 1]), 1}}}
```

The expression sizes increase faster than exponentially:

```wl
In[3]:= Table[ByteCount[Eigensystem[Array[Subscript[a, ##]&, {n, n}]]], {n, 4}]

Out[3]= {312, 6680, 101968, 2725384}
```

---

Not all matrices have a complete set of eigenvectors:

```wl
In[1]:= Eigensystem[{{7, -1}, {9, 1}}]

Out[1]= {{4, 4}, {{1, 3}, {0, 0}}}
```

Use ``JordanDecomposition`` for exact computation:

```wl
In[2]:= JordanDecomposition[{{7, -1}, {9, 1}}]

Out[2]= {{{1, (1/3)}, {3, 0}}, {{4, 1}, {0, 4}}}
```

Use ``SchurDecomposition`` for numeric computation:

```wl
In[3]:= SchurDecomposition[N@{{7, -1}, {9, 1}}]

Out[3]= {{{0.316228, -0.948683}, {0.948683, 0.316228}}, {{4., -10.}, {0., 4.}}}
```

---

Construct a 10,000×10,000 sparse matrix:

```wl
In[1]:= n = 10^4;s = SparseArray[{Band[{1, 1}] -> 1., Band[{1, -1}, {-1, 1}, {1, -1}] -> 1., Band[{n / 2, 1}, {1, -1}, {0, 1}] -> 1., Band[{1, n / 2}, {-1, 1}, {1, 0}] -> 1.}, {n, n}];

In[2]:= ArrayPlot[s, MaxPlotPoints -> 1000, ColorFunction -> (If[# == 0, White, Black]&)]

Out[2]= [image]
```

The eigenvector matrix is a dense matrix and too large to represent:

```wl
In[3]:= MemoryConstrained[Eigensystem[s], 10 ^ 8]
```

Eigensystem::arh: Because finding 10000 out of the 10000 eigenvalues and/or eigenvectors is likely to be faster with dense matrix methods, the sparse input matrix will be converted. If fewer eigenvalues and/or eigenvectors would be sufficient, consider restricting this number using the second argument to Eigensystem.

```wl
Out[3]= $Aborted
```

Computing the few largest or smallest eigenvalues is usually possible:

```wl
In[4]:= {vals, vecs} = Eigensystem[s, 2]; vals

Out[4]= {101.496, -98.4963}
```

---

When eigenvalues are closely grouped the iterative method for sparse matrices may not converge:

```wl
In[1]:= s = SparseArray[{{x_, y_} /; Abs[x - y] ≤ 1 -> 3. Abs[x - y] - 2.}, {1000, 1000}];
```

The iteration has not converged well after 1000 iterations:

```wl
In[2]:=
{vals, vecs} = Eigensystem[s, 3];
ListPlot[vecs, PlotLabel -> vals]
```

Eigensystem::maxit2: Warning: maximum number of iterations, 1000, has been reached by the Arnoldi algorithm without convergence to the specified tolerance, but the current best computed value has been returned. You can use method options with Method -> {Arnoldi, opts} to increase the size of basis vectors, the maximum number of iterations, reduce the tolerance, or use an estimate as a shift, any of which may help.

```wl
Out[2]= [image]
```

You can give the algorithm a shift near the expected value to speed up convergence:

```wl
In[3]:= {vals, vecs} = Eigensystem[s, 3, Method -> {"Arnoldi", "Shift" -> -4}];ListPlot[vecs, PlotLabel -> vals]

Out[3]= [image]
```

---

The end points given to an interval as specified for the FEAST method are not included.

Set up a matrix with eigenvalues at 3 and 9 and find the eigenvalues and eigenvectors in the interval $(3.,9.)$ :

```wl
In[1]:=
ck = With[{n = 24}, N[SparseArray[{{i_, j_} /; Abs[i - j] == 1 :> With[{k = Min[i, j]}, Sqrt[k (n - k)]]}, {n, n}]]];
{vals, vecs} = Eigensystem[ck, Method -> {"FEAST", "Interval" -> {3., 9.}}];
vals

Out[1]= {7., 5.}
```

Enlarge the to interval $(2.9,9.1)$ such that FEAST finds the eigenvalues 3 and 9 and their corresponding eigenvectors:

```wl
In[2]:=
{vals, vecs} = Eigensystem[ck, Method -> {"FEAST", "Interval" -> {2.9, 9.1}}];
vals

Out[2]= {9., 7., 5., 3.}
```

## See Also

* [`EigenvalueDecomposition`](https://reference.wolfram.com/language/ref/EigenvalueDecomposition.en.md)
* [`Eigenvalues`](https://reference.wolfram.com/language/ref/Eigenvalues.en.md)
* [`Eigenvectors`](https://reference.wolfram.com/language/ref/Eigenvectors.en.md)
* [`NDEigensystem`](https://reference.wolfram.com/language/ref/NDEigensystem.en.md)
* [`DEigensystem`](https://reference.wolfram.com/language/ref/DEigensystem.en.md)
* [`NullSpace`](https://reference.wolfram.com/language/ref/NullSpace.en.md)
* [`JordanDecomposition`](https://reference.wolfram.com/language/ref/JordanDecomposition.en.md)
* [`SchurDecomposition`](https://reference.wolfram.com/language/ref/SchurDecomposition.en.md)
* [`SingularValueDecomposition`](https://reference.wolfram.com/language/ref/SingularValueDecomposition.en.md)
* [`QRDecomposition`](https://reference.wolfram.com/language/ref/QRDecomposition.en.md)
* [`MatrixFunction`](https://reference.wolfram.com/language/ref/MatrixFunction.en.md)

## Tech Notes

* [Eigenvalues and Eigenvectors](https://reference.wolfram.com/language/tutorial/LinearAlgebra.en.md#9501)

## Related Guides

* [Matrix Operations](https://reference.wolfram.com/language/guide/MatrixOperations.en.md)
* [Matrix Decompositions](https://reference.wolfram.com/language/guide/MatrixDecompositions.en.md)
* [Unsupervised Machine Learning](https://reference.wolfram.com/language/guide/UnsupervisedMachineLearning.en.md)

## History

* Introduced in 1988 (1.0) \| Updated in 1996 (3.0) ▪ 2003 (5.0) ▪ [2014 (10.0)](https://reference.wolfram.com/language/guide/SummaryOfNewFeaturesIn100.en.md) ▪ [2015 (10.3)](https://reference.wolfram.com/language/guide/SummaryOfNewFeaturesIn103.en.md) ▪ [2024 (14.0)](https://reference.wolfram.com/language/guide/SummaryOfNewFeaturesIn140.en.md)