Theory

Introduction

We are interested in systems of equations that can be written in the form ${\mathbf{y}}^{\prime }=A\mathbf{y}$, where $\mathbf{y}(t)$ is a vector function of $t$, and $A$ is a matrix with constant entries (independent of $t$).

Our goal is to find a complete set of solutions to ${\mathbf{y}}^{\prime }=A\mathbf{y}$. It is not enough to find just one. How many solutions constitute a complete set?

Picard-Lindelöf for first-order systems

The system: ${\mathbf{y}}^{\prime }(t)=A\mathbf{y}(t),\mathbf{y}(0)={\mathbf{y}}_0$
always has existence and uniqueness of solutions, for any initial vector ${\mathbf{y}}_0$ and constant coefficient matrix $A$.

This implies: a complete set of solutions is a vector space (by linearity) of the same dimension as ${\mathbf{y}}_0$ (by the theorem).

Eigenvectors
The eigenvector approach to solving the system starts from the observation that if $\mathbf{v}$ is an eigenvector of $A$ with eigenvalue $\lambda$, meaning that $A\mathbf{v}=\lambda \mathbf{v}$, then the vector function $\mathbf{y}(t)=\mathbf{v}{e}^{\lambda t}$ is a solution to the system.

The downside of this approach is that the solution cannot be written down until all the eigenvalues have been computed; the eigenvalues frequently involve imaginary numbers; and it is possible there are simply not enough eigenvectors to give a complete family.

Matrix exponentials
The matrix exponential approach to solving the same system starts by constructing the function $Y(t)=e^{At}$ using power series. This function is a matrix-valued function of $t$ that satisfies $Y^{\prime }=AY$, and because $e^{At}$ is always an invertible matrix, the columns of $Y(t)$ automatically give a complete set of solutions.

The downside of this approach is that it is hard to compute $e^{At}$, and the typical method requires finding the eigenvalues of $A$ anyway. On the other hand, it is not necessary to use imaginary numbers because the matrix $e^{At}$ can be put into what is called “Real Jordan Form” (with only real numbers). This form uses an extension of the idea of eigenvector that is designed specifically to handle the situation of missing eigenvectors. “Putting into” this form simply means changing variables $\mathbf{z}=P\mathbf{y}$ according to some invertible transformation $P$. The matrix $PA{P}^{-1}$ of the new equation ${\mathbf{z}}^{\prime }=(PA{P}^{-1})\mathbf{z}$ will have the Real Jordan Form.

Eigenvector approach

Suppose we have any system of the form ${\mathbf{y}}^{\prime }=A\mathbf{y}$ for a matrix with constant entries $a_{ij}$. Suppose $\mathbf{v}$ is an eigenvector of $A$ with eigenvalue $\lambda$, meaning that $A\mathbf{v}=\lambda \mathbf{v}$. (Note that $\mathbf{v}$ also has constant entries, like $A$.) Then the function $\mathbf{y}(t)=\mathbf{v}{e}^{\lambda t}$ solves the system:

$$\begin{array}{rl}{\mathbf{y}}^{\prime }(t)& =(\mathbf{v}{e}^{\lambda t}{)}^{\prime }={\mathbf{v}}^{\prime }{e}^{\lambda t}+\mathbf{v}(e^{\lambda t}{)}^{\prime }=0+\mathbf{v}\lambda e^{\lambda t}=\lambda (\mathbf{v}{e}^{\lambda t})=\lambda \mathbf{y}(t)\\ & \\ A\mathbf{y}(t)& =A(\mathbf{v}{e}^{\lambda t})=(A\mathbf{v})e^{\lambda t}=(\lambda \mathbf{v})e^{\lambda t}=\lambda (\mathbf{v}{e}^{\lambda t})=\lambda \mathbf{y}(t).\end{array}$$

Eigenvector-eigenvalue solutions

Any eigenvector-eigenvalue pair $\mathbf{v},\lambda$ determines a solution $\mathbf{y}(t)=\mathbf{v}{e}^{\lambda t}$.

Remember that if $\mathbf{v}$ is an eigenvector of $A$, then any multiple $a\mathbf{v}$ is also an eigenvector. In terms of solutions, this simply means that if $\mathbf{y}(t)$ is a solution of the form $\mathbf{v}{e}^{\lambda t}$, then $a\mathbf{y}(t)$ is also a solution. We already knew this: it is part of linearity – the fact that solutions are vectors in a vector space.

The Picard-Lindelöf theorem implies that the space of solutions has dimension $n$ equal to the dimension of the vectors $\mathbf{y}$ and ${\mathbf{y}}_0$. (The dimension of $\mathbf{y}$ is just the number of components it has.) Therefore a set of solutions ${\mathbf{y}}_1(t),\dots ,{\mathbf{y}}_k(t)$ is a complete set when (a) we have $k=n$, and (b) they are independent when considered as vectors. The Wronskian (using the determinant) can tell us when $n$ vectors having $n$ components are independent.

When there are $n$ independent eigenvectors, and all of their eigenvalues are real numbers, the eigenvector strategy works perfectly without further complication. Just find all eigenvector-eigenvalue pairs and write down all the solutions they determine, and we are done.

The main challenge of the eigenvector strategy becomes the challenge of handling situations of two types. It can happen that (a) some eigenvalues involve imaginary numbers, or (b) some systems do not have enough eigenvectors. In applications, the former tends to arise for oscillatory systems, and the latter arises for analogues of critical damping.

Complex eigenvalues

Suppose the solution to the eigenvalue equation $\det A_{\lambda }=0$ is a complex eigenvalue $\lambda =a+bi$.

Because the matrix $A_{\lambda }$ has real-valued entries, such eigenvalue roots must come in conjugate pairs, meaning we would always find both eigenvalues $\lambda =a\pm bi$. Furthermore, the corresponding eigenvectors also come in a complex conjugate pair ${\mathbf{v}}_{+}={\mathbf{v}}_a+{\mathbf{v}}_{b}i$ and ${\mathbf{v}}_{-}={\mathbf{v}}_a-{\mathbf{v}}_{b}i$:

$$\begin{array}{rl}A{\mathbf{v}}_{+}& =(a+bi){\mathbf{v}}_{+}\\ A{\mathbf{v}}_{-}& =(a-bi){\mathbf{v}}_{-}.\end{array}$$

(This follows by taking the complex conjugate of both sides of the equation, using the fact that conjugation and multiplication can be performed in either order.)

We proceed as in the $2^{\text{nd}}$ order linear theory, and combine the solutions ${\mathbf{y}}_{+}(t)={\mathbf{v}}_{+}{e}^{(a+bi)t}$ and ${\mathbf{y}}_{-}(t)={\mathbf{v}}_{-}{e}^{(a-bi)t}$ in such a way as to extract vector functions having exclusively real-valued entries. We have:

$$\begin{array}{rlrl}{\mathbf{y}}_1(t)& :=\frac{1}{2}({\mathbf{y}}_{+}+{\mathbf{y}}_{-})& & ={\mathbf{v}}_a{e}^{at}\cos (bt)-{\mathbf{v}}_b{e}^{at}\sin (bt)\\ & & & \\ {\mathbf{y}}_2(t)& :=\frac{1}{2i}({\mathbf{y}}_{+}-{\mathbf{y}}_{-})& & ={\mathbf{v}}_a{e}^{at}\sin (bt)+{\mathbf{v}}_b{e}^{at}\cos (bt).\end{array}$$

In summary: whenever $\lambda =a+bi$ is an eigenvalue for the eigenvector $\mathbf{v}={\mathbf{v}}_a+{\mathbf{v}}_{b}i$, the formulas above for ${\mathbf{y}}_1$ and ${\mathbf{y}}_2$ give two independent solutions to the system ${\mathbf{y}}^{\prime }=A\mathbf{y}$.

Insufficient eigenvectors

Suppose the solution to the eigenvalue equation $\det A_{\lambda }=0$ is a repeated eigenvalue, so for example the equation has a factor like $(\lambda -a{)}^k=0$. In this case, sometimes it happens that the eigenvalue $\lambda =a$ does not have enough eigenvectors associated to it.

Insufficient eigenvectors (part 2)

Let us consider the way around this difficulty.
The way around is not very simple! Still, we try to take the simplest path.

Suppose our situation is that $(\lambda -a{)}^2=0$ as a factor in $\det A_{\lambda }$, and that we have a single eigenvector $\mathbf{v}$ satisfying $A_a\mathbf{v}=0$. So we start with the solution $\mathbf{v}{e}^{at}$. We seek a second solution.

There is a concept of generalized eigenvector, which is any $\mathbf{x}$ that satisfies $A_{\lambda }^k\mathbf{x}=0$ for higher powers $k\ge 1$. According to linear algebra theory, whenever $(\lambda -a{)}^k=0$ from a factor in $\det A_{\lambda }$, there will be enough independent generalized eigenvectors satisfying $A_{\lambda }^k\mathbf{x}=0$. Notice that if $A_{\lambda }\mathbf{x}=0$ then $A_{\lambda }^2\mathbf{x}=0$, so regular eigenvectors always count as generalized eigenvectors.

Now let us consider the system ${\mathbf{y}}^{\prime }=A\mathbf{y}$ again to see how to extract new solutions from generalized eigenvectors with $k\ge 2$. We start by computing two linearly independent solutions to $A_{\lambda }^k\mathbf{x}=0$, called ${\mathbf{v}}_0$ and ${\mathbf{v}}_1$, which are connected according to multiplication by $A_a=(A-aI)$:

$${\mathbf{v}}_0\stackrel{(A-aI)}{⟶}{\mathbf{v}}_1\stackrel{(A-aI)}{⟶}0.$$

So ${\mathbf{v}}_1$ is a regular eigenvector, while ${\mathbf{v}}_0$ is only a generalized one with $k=2$. (In a sense, ${\mathbf{v}}_0$ is “deeper down” than ${\mathbf{v}}_1$. So regular eigenvectors do not always reach deep enough down into preimages of $A_a$.)

We already know that $({\mathbf{v}}_1{e}^{at}{)}^{\prime }=A({\mathbf{v}}_1{e}^{at})$ since ${\mathbf{v}}_1$ is a regular eigenvector with eigenvalue $a$.

Now look what happens when we plug ${\mathbf{v}}_{1}t{e}^{at}$ into the equation:

$$\begin{array}{rl}({\mathbf{v}}_{1}t{e}^{at}{)}^{\prime }& ={\mathbf{v}}_1{e}^{at}+a{\mathbf{v}}_{1}t{e}^{at},\\ & \\ A\left({\mathbf{v}}_{1}t{e}^{at}\right)& =(A-aI+aI){\mathbf{v}}_{1}t{e}^{at}\\ & =(A-aI){\mathbf{v}}_{1}t{e}^{at}+a{\mathbf{v}}_{1}t{e}^{at}\\ & =a{\mathbf{v}}_{1}t{e}^{at}.\end{array}$$

And look what happens when we plug ${\mathbf{v}}_0{e}^{at}$ into the equation:

$$\begin{array}{rl}({\mathbf{v}}_0{e}^{at}{)}^{\prime }& =a{\mathbf{v}}_0{e}^{at},\\ & \\ A\left({\mathbf{v}}_0{e}^{at}\right)& =(A-aI+aI){\mathbf{v}}_0{e}^{at}\\ & =(A-aI){\mathbf{v}}_0{e}^{at}+a{\mathbf{v}}_0{e}^{at}\\ & ={\mathbf{v}}_1{e}^{at}+a{\mathbf{v}}_0{e}^{at}.\end{array}$$

Insufficient eigenvectors

If we add the solutions to form $\mathbf{y}(t)={\mathbf{v}}_0{e}^{at}+{\mathbf{v}}_{1}t{e}^{at}$, the ODE equation holds.

Matrix exponential approach

Recall that the matrix exponential approach gives us the formula $Y(t)=e^{At}$ for the fundamental matrix solution to ${\mathbf{y}}^{\prime }=A\mathbf{y}$. Each column of $Y=\left(\begin{array}{ccc}{\mathbf{y}}_1& \cdots & {\mathbf{y}}_n\end{array}\right)$ is a solution vector function, and these columns are independent of each other. The complete family of solutions to ${\mathbf{y}}^{\prime }=A\mathbf{y}$ is:

$$Y\mathbf{c}=Y\left(\begin{array}{c}{C}_1\\ C_2\\ ⋮\\ C_n\end{array}\right)=C_1{\mathbf{y}}_1+\cdots +C_n{\mathbf{y}}_n.$$

This approach has no difficulty in finding a formula for the solution. The challenge for this approach is: How do we describe the formula $e^{At}$ in more practical terms?

A complete description of the answer to this challenge gets too complicated to cover here and requires more linear algebra than this course assumes. Instead of trying to give a complete description, our goal will be to see a glimpse of the main ideas and then let your imagination fill in the rest.

Key observation

The central observation we need is this. First, we can rewrite the series in an ‘expansion around $\lambda$’:

$$\begin{array}{rl}{e}^{At}& =e^{\lambda t}{e}^{(A-\lambda I)t}=e^{\lambda t}{e}^{A_{\lambda }t}\\ & =e^{\lambda t}\left(I+t{A}_{\lambda }+\frac{t^2}{2!}{A}_{\lambda }^2+\frac{t^3}{3!}{A}_{\lambda }^3+\cdots \right).\end{array}$$

Now suppose we have a basis of generalized eigenvectors ${\mathbf{u}}_1,\dots ,{\mathbf{u}}_k$. This means $A_{\lambda }^k{\mathbf{u}}_i=0$ for $i=1,\dots ,k$ and ${\mathbf{u}}_1,\dots ,{\mathbf{u}}_k$ are independent of each other.

Then when $e^{At}$ is applied to any of the ${\mathbf{u}}_i$, the terms of the second factor series are eventually just zeros:

$$\begin{array}{rl}{e}^{At}{\mathbf{u}}_i& =e^{\lambda t}\left(I+t{A}_{\lambda }+\frac{t^2}{2!}{A}_{\lambda }^2+\frac{t^3}{3!}{A}_{\lambda }^3+\cdots \right){\mathbf{u}}_i\\ & =e^{\lambda t}\left({\mathbf{u}}_i+t{A}_{\lambda }{\mathbf{u}}_i+\frac{t^2}{2!}{A}_{\lambda }^2{\mathbf{u}}_i+\cdots +\frac{t^{k-1}}{(k-1)!}{A}_{\lambda }^{k-1}{\mathbf{u}}_i+0+0+\cdots \right)\end{array}$$

where the zeros occur because $A_{\lambda }^j{\mathbf{u}}_i=0$ whenever $j\ge k$.

This finite series can then be calculated. No limits of matrix sums need be taken.

Matrix exponential approach (part 2)

In order to generalize the previous method to all cases you could encounter, it suffices to be able to find a set of $k$ independent generalized eigenvectors whenever $\det A_{\lambda }=(\lambda -a{)}^k$ (or has the latter as a factor). The best general method for this is to use “row reduction” from linear algebra to find a basis of the null space of $A_a$, i.e. to solve for $\mathbf{x}$ in the vector equation $(A-\lambda I)\mathbf{x}=0$ when $\lambda =a$.