Estimation of the eigenvalues of a matrix

Question

Given the following matrix \begin{pmatrix}0 & 1 & 0 & 0\\ -k & -\lambda_1 & -k & 0\\ 0 & 0 & 0 & 1\\ -k & 0 & -k & -\lambda_2\end{pmatrix} where $k,\lambda_1$ and $\lambda_2$ are strictly positive, is there a way to estimate the eigenvalues, or at least say if they are positive or negative, complex or real? One eigenvalue is zero, but to find the others you need to solve a parametric cubic equation, which is not easy. Note that this matrix is the representation of a mechanical system, so a physical context may come in handy.

Answer 1

The eigenvalues of the matrix are the roots of the characteristic polynomial.

1. Compute the characteristic polynomial. Mathematica tells me that it is: $$p(s) = s^4+s^3(\lambda_1+\lambda_2)+s^2(\lambda_1\lambda_2+2k)+s k(\lambda_1+\lambda_2)=s q(s)$$ with $$q(s) = s^3 +s^2(\lambda_1+\lambda_2)+s(\lambda_1\lambda_2+2k)+ k(\lambda_1+\lambda_2)$$
2. Apply Routh's Stability Criterion to $q(s)$.

Routh's array should be as follows: $$\begin{matrix} s^3 & 1 & \lambda_1\lambda_2+2k\\ s^2 & \lambda_1+\lambda_2 & k(\lambda_1+\lambda_2)\\ s & \lambda_1\lambda_2+k\\ s^0 & k(\lambda_1+\lambda_2) \end{matrix}.$$ All the eigenvalues have negative real part if and only if the first column has only positive numbers (which is the case, given your assumptions).

Answer 2

Here is a treatment that takes advantage of the observation that this matrix with parameters can be considered as a simpler one with a rather simple "rank-one perturbation". Here is how :

$$B=\underbrace{\begin{pmatrix}0 & 1 & 0 & 0\\ 0 & -\lambda_1 & 0 & 0\\ 0 & 0 & 0 & 1\\ 0 & 0 & 0 & -\lambda_2\end{pmatrix}}_A+\underbrace{ \begin{pmatrix}0 & 0 & 0 & 0\\ -k & 0 & -k & 0\\ 0 & 0 & 0 & 0\\ -k & 0 & -k & 0\end{pmatrix}}_P$$

Observe that $A$ is a triangular matrix and $$P=cd^T \ \text{with} \ c=\begin{pmatrix} \ \ 0\\-k\\ \ \ 0\\-k\end{pmatrix} \ \ \text{and} \ \ d=\begin{pmatrix}1\\0\\1\\0\end{pmatrix} $$

It is matrix $P$ which is the "rank-one perturbation matrix" (column vector times row vector).

Why is this interesting ? Because we have a nice formula :

$$\det(A + \mathbf{c}\mathbf{d}^T) = \det(A)(1 + \mathbf{d}^T A^{-1}\mathbf{c})$$

(see here for a proof) Let us apply it, not to $A$ itself, but to $A-\lambda I$, giving

$$\det((A-\lambda I) + \mathbf{c}\mathbf{d}^T) = \det(A-\lambda I)(1 + \mathbf{d}^T (A-\lambda I)^{-1}\mathbf{c})$$

I have used a CAS (see computational details at the bottom of this anwer), giving the following final expression to (1) :

$$\det(B-\lambda I)=\underbrace{\lambda^2(\lambda+\lambda_1)(\lambda+\lambda_2)}_{\det(A-\lambda I)}\underbrace{\left(1+\frac{k}{\lambda(\lambda+\lambda_1)}+\frac{k}{\lambda(\lambda+\lambda_2)}\right)}_{f(\lambda)}$$

Instead of expanding the products in the above formula, we will leave it like that. Indeed, the second factor $f(\lambda)$ holds all the information we need about the perturbation on the characteristic equation, under the form of a factor to the original characteristic equation having the form $(1+\text{something})$.

Due to the positivity of parameters $k, \lambda_1, \lambda_2$, if we take any $\lambda > 0$ into $f(\lambda)$, it is impossible that the result is zero : therefore if $f(\lambda)=0$, it can only be with negative value(s) of $\lambda$. QED.

A different understanding of this phenomena can be obtained through graphical representations as those given below : the roots of $f(\lambda)$ are in between its poles ( materialized by vertical asymptotes $\lambda = 0, \lambda=-\lambda_1, \lambda=-\lambda_2$) constraining the real roots to be negative.

Fig. 1: A case with a single real eigenvalue $\lambda \approx -4.66$. Please note the horizontal asymptate $y=1$. Besides, a conjugate pair of eigenvalues can be "guessed" with negative real part $\approx -1.8$ in the vicinity of the summit of the third branch.

Fig. 2: A case with three real eigenvalues, all of them negative.

---

Edit (2 days after) : I agree that, as remarked by Dan, Descartes rule of signs gives the result with less calculations. The main interest of the method I have developed is that it allows to bracket the solutions, for example in Fig. 2, with 3 solutions $\lambda, \lambda', \lambda''$, we have :

$$-\lambda_1 < \lambda < -\lambda_2 < \lambda' \le \lambda'' < 0 $$

In particular, all solutions are $ > - \max(\lambda_1,\lambda_2)$

---

Computational detail :

$$(A-\lambda I)^{-1}=\begin{pmatrix}-\tfrac{1}{\lambda} & -\tfrac{1}{\lambda(\lambda+\lambda_1)} & 0 & 0\\ 0 & -\tfrac{1}{(\lambda+\lambda_1)} & 0 & 0\\ 0 & 0 & -\tfrac{1}{\lambda} & -\tfrac{1}{\lambda(\lambda+\lambda_2)}\\ 0 & 0 & 0 & -\tfrac{1}{(\lambda+\lambda_2)}\end{pmatrix}.$$

(In fact, it is not so difficult to invert $A-\lambda I$ even "by hand" because it is formed by two $2 \times 2$ triangular blocks).