# What is the fastest way to find the characteristic polynomial of a matrix?

Finding the characteristic polynomial of a matrix of order $n$ is a tedious and boring task for $n > 2$.

I know that:

• the coefficient of $\lambda^n$ is $(-1)^n$,
• the coefficient of $\lambda^{n-1}$ is $(-1)^{n-1}(a_{11} + a_{22} + \dots + a_{nn})$,
• the constant term is $\det{A}$.

When finding the coefficient of the linear term $\lambda$ of the characteristic polynomial of a $3\times 3$ matrix, one has to calculate the determinant of the matrix $A - \lambda I_n$ anyway. (But you don't have to sum all the terms, only the linear terms.)

Does anybody know a faster way?

• You could expand the determinant in terms of the variables, and you get it as a function of permutations. It would be similar to what you calculated, i.e that $M(\lambda) = (-1)^n \lambda^n +(-1)^n \sum_{|\sigma|=n-1} \prod a_{i\sigma(i)} + \ldots$ – Calvin Lin May 29 '13 at 15:33
• Do you mean to perform row operations on the matrix $A$ or the matrix $A - \lambda I_n$? If you're talking about the later, I don't see any way how to zero out the first column under the main diagonal. Could you please develop that point? – David May 29 '13 at 15:48
• @CalvinLin I answered to Julien who had withdrawn his comment. I will try that nonetheless, thanks. – David May 29 '13 at 15:55

Once upon a less enlightened time, when people were less knowledgeable in the intricacies of algorithmically computing eigenvalues, methods for generating the coefficients of a matrix's eigenpolynomial were quite widespread. One of the more prominent methods for computing the coefficients was a method ascribed to both the Frenchman Leverrier, and the Russian Faddeev (who was an (co-)author of one of the oldest references on the practice of numerical linear algebra).

The (Faddeev-)Leverrier method is a method that will require you to do a number of matrix multiplications to generate the coefficients of the characteristic polynomial. Letting the $n\times n$ matrix $\mathbf A$ have the monic characteristic polynomial $(-1)^n \det(\mathbf A-\lambda\mathbf I)=\lambda^n+c_{n-1}\lambda^{n-1}+\cdots+c_0$, the algorithm proceeds like so:

$\mathbf C=\mathbf A;$
$\text{for }k=1,\dots,n$

$\text{if }k>1$
$\qquad \mathbf C=\mathbf A\cdot(\mathbf C+c_{n-k+1}\mathbf I);$

$c_{n-k}=-\dfrac{\mathrm{tr}(\mathbf C)}{k};$

$\text{end for}$

If your computing environment can multiply matrices, or take their trace (sum of the diagonal elements, $\mathrm{tr}(\cdot)$), then you can easily program (Faddeev-)Leverrier. The method works nicely in exact arithmetic, or in hand calculation (assuming you have the stamina to repeatedly multiply matrices), but is piss-poor in inexact arithmetic, as the method tends to greatly magnify rounding errors in the matrix, ever yielding coefficients that become increasingly inaccurate as the iteration proceeds. But, for the simple $3\times 3$ case envisioned by the OP, this should work nicely.

People interested in this old, retired method might want to see this paper.

• For searching purposes: the method is also sometimes named after Frame and Souriau. – J. M. isn't a mathematician May 30 '13 at 14:48
• The computational complexity is suboptimal, but the method is nevertheless useful for symbolic arithmetic because there are no data-dependent branches. It is also useful for computing the adjugate matrix. – ccorn Oct 25 '16 at 1:23
• Oh, and the method gives useful formulae for the invariants, which is a nice thing to have e.g. when you need derivatives of those with respect to the matrix object as a whole. Needless to say, I like it very much. Glad to find it mentioned here. – ccorn Oct 25 '16 at 1:33

By hand, I think Ted Shifrin's method is fastest.

By computer, evaluate $\det (\lambda \cdot \mathrm{Id}- A)$ for $n$ values of $\lambda$, and find the characteristic polynomial by Lagrange interpolation.

• Yes, but then you are in front of solving a not so nice equation. – AlainD May 30 '13 at 10:01
• @user1938185: What do you mean? Can you explain? – Giuseppe Negro May 30 '13 at 10:14
• @user, No need to solve polynomials in this method; you only need to choose $n$ good points to evaluate your determinant on, and then use polynomial interpolation methods to retrieve the coefficients. In particular, if one is able to evaluate the determinant at the roots of unity, then you can use FFT to retrieve the coefficients. – J. M. isn't a mathematician May 30 '13 at 14:52
• This requires $n$ determinant evaluations, for a total cost of $O(n^4)$ operations. I believe this can be done in $O(n^3)$; for instance, if the coefficients are floating point numbers, you can compute a Schur factorization $A=QTQ^*$ of the matrix in $O(n^3)$ and then compute $\det (\lambda I - A ) = \det(\lambda I - T)$ in $O(n^2)$. – Federico Poloni Feb 16 '16 at 13:32
• @Federico, you don't really need to go up to Schur; a similarity transformation to Hessenberg form ought to suffice. – J. M. isn't a mathematician Dec 10 '16 at 11:52

This may or may not make you happy. The coefficient of $\lambda^{n-k}$ is $(-1)^{n-k}$ times the sum of all the principal $k\times k$ minors, i.e., the determinants of all $k\times k$ submatrices built using the same rows and columns. Note that this generalizes your observations for $k=n$ and $k=1$.

• A bit unrelated, but thanks for your text on Diff. Geometry. – user12802 May 29 '13 at 17:55

You can use the Cayley-Hamilton theorem, which says that the matrix $A$ is a root of the minimal polynomial, which divides the characteristic polynomial. In facts, the minimal polynomial has the same roots than the characteristic, except those roots may have lower multiplicity. In particular, if the eigenvalues of $A$ are distinct, then the two polynomials are equals.

You compute the $n=\dim$ powers: $I, A, A^2, A^3,...,A^n$, and then look for coefficients $c_i$ such that $A^n = c_{n-1}A^{n-1}+...+c_1A^1+c_0I$, then the minimal polynomial is $\det(\lambda I - A) = P(X) = \lambda^n-c_{n-1}\lambda^{n-1}-...-c_1\lambda-c_0$.

For example, if $A=\begin{pmatrix} 0&1&0\\1&0&0\\0&0&2\end{pmatrix}$, then $A^2=\begin{pmatrix} 1&0&1\\0&1&0\\0&0&4\end{pmatrix}$, and $A^3=\begin{pmatrix} 0&1&0\\1&0&0\\0&0&8\end{pmatrix}$,

so that $A^2-I=\begin{pmatrix} 0&0&0\\0&0&0\\0&0&3\end{pmatrix}$ and $A^3-A=\begin{pmatrix} 0&0&0\\0&0&0\\0&0&6\end{pmatrix}$, thus $A^3-A=2(A^2-I)$, and

the minimal polynomial is $\lambda^3-\lambda-2(\lambda^2-1)=(\lambda^2-1)(\lambda-2)=(\lambda+1)(\lambda-1)(\lambda-2)$. All the roots being simple, the minimal polynomial is also the characteristic polynomial, and $\chi_A(\lambda)=\det(\lambda I-A)=\lambda^3-2\lambda^2-\lambda+2$.

The search of the coefficients $c_i$ may be done in a systematic way, similar to the Gauss method. However the system is overdetermined, with $n \times n$ equations, and is usually quicker to solve than you think, especially for matrices found in homework.

The eigenvalues will almost always be distinct (so the characteristic and minimal polynomials will be the same), except maybe for matrices found in homework. For these, is is better to look to the power matrix $A^k, k \le n$ as you compute them, and look if they are linear combinations of the previous ones.

For example if $A=\begin{pmatrix} 0&1&1\\1&0&1\\1&1&0\end{pmatrix}$, then $A^2=\begin{pmatrix} 2&1&1\\1&2&1\\1&1&2\end{pmatrix}$, and $A^2=A+2I$. So the minimal polynomial is $\lambda^2-\lambda-2 = (\lambda+1)(\lambda-2)$. The characteristic polynomial being a polynomial of degree 3 with the same roots, it can either be $(\lambda+1)^2(\lambda-2)$ or $(\lambda+1)(\lambda-2)^2$. The multiplicity $\nu_i$ of $(x-\lambda_i)$ in $\chi_A(x) = \prod (x-\lambda_i)^{\nu_i}$, is the dimension of the associated eigenspace $E_{\lambda_i}=\ker (A-\lambda_i I)=\{x\mid Ax=\lambda_i x\}$. In our case, $E_{-1} = \ker(A+I) = \ker \begin{pmatrix} 1&1&1\\1&1&1\\1&1&1\end{pmatrix}=2$, so $\chi_A(\lambda)=(\lambda+1)^2(\lambda-2)$.

I don't know if the method is less tedious, but I find it less boring.

• Less systematic than the Faddeev-Leverrier method, though. – AlainD May 30 '13 at 10:33
• @J-M thanks for the edit. – AlainD May 30 '13 at 19:12
• No problem; nice, but requires some thinking from the user... :) – J. M. isn't a mathematician May 30 '13 at 19:17
• The original question was to find a method that was neither tedious nor boring. All of the proposed method, except the one of @David Speyer, are tedious. And all of the proposed methods, except mine are boring, precisely because they don't require thinking. – AlainD May 30 '13 at 19:23

I believe the fastest way known is the Pernet-Storjohann randomized algorithm with complexity $O(n^\theta)$, where $\theta$ is an admissible exponent for the complexity of matrix multiplication in terms of field operations. On the other hand, it is known that the complexity of the product of two matrices is bounded above by the complexity of the the characteristic polynomial, so the Pernet-Storjohann complexity is "sharp".