A differential equation is a mathematical equation for an unknown function of one or several variables that relates the values of the function itself and its derivatives of various orders. A matrix differential equation contains more than one function stacked into vector form with a matrix relating the functions to their derivatives.
For example, a first-order matrix ordinary differential equation is
|
=A(t)x(t)
x(t)
n x 1
t
|
A(t)
n x n
In the case where
A
x(t)=c1
λ1t | |
e |
u1+c2
λ2t | |
e |
u2+ … +cn
λnt | |
e |
un~,
More generally, if
A(t)
t | |
\int | |
a |
A(s)ds
| ||||||||||
x(t)=e |
c~,
c
n x 1
By use of the Cayley–Hamilton theorem and Vandermonde-type matrices, this formal matrix exponential solution may be reduced to a simple form.[1] Below, this solution is displayed in terms of Putzer's algorithm.[2]
The matrix equation
|
=Ax(t)+b
The steady state x* to which it converges if stable is found by setting
* | |||
|
(t)=0~,
x*=-A-1b~,
Thus, the original equation can be written in the homogeneous form in terms of deviations from the steady state,
*]~. | |||
|
An equivalent way of expressing this is that x* is a particular solution to the inhomogeneous equation, while all solutions are in the form
* | |
x | |
h+x |
~,
xh
In the n = 2 case (with two state variables), the stability conditions that the two eigenvalues of the transition matrix A each have a negative real part are equivalent to the conditions that the trace of A be negative and its determinant be positive.
The formal solution of
*] | |||
|
x(t)=x*+eAt[x(0)-x*]~,
Given a matrix A with eigenvalues
λ1,λ2,...,λn
eAt=
n-1 | |
\sum | |
j=0 |
rj+1{\left(t\right)}Pj
P0=I
Pj=
j | |
\prod | |
k=1 |
\left(A-λkI\right)=Pj-1\left(A-λjI\right), j=1,2,...,n-1
r |
1=λ1r1
r1{\left(0\right)}=1
r |
j=λjrj+rj-1, j=2,3,...,n
rj{\left(0\right)}=0, j=2,3,...,n
The equations for
ri(t)
Note the algorithm does not require that the matrix A be diagonalizable and bypasses complexities of the Jordan canonical forms normally utilized.
A first-order homogeneous matrix ordinary differential equation in two functions x(t) and y(t), when taken out of matrix form, has the following form:
dx | |
dt |
=a1x+b
|
=a2x+b2y
where
a1
a2
b1
b2
Higher order matrix ODE's may possess a much more complicated form.
The process of solving the above equations and finding the required functions of this particular order and form consists of 3 main steps. Brief descriptions of each of these steps are listed below:
The final, third, step in solving these sorts of ordinary differential equations is usually done by means of plugging in the values calculated in the two previous steps into a specialized general form equation, mentioned later in this article.
To solve a matrix ODE according to the three steps detailed above, using simple matrices in the process, let us find, say, a function and a function both in terms of the single independent variable, in the following homogeneous linear differential equation of the first order,
dx | =3x-4y, | |
dt |
dy | |
dt |
=4x-7y~.
To solve this particular ordinary differential equation system, at some point in the solution process, we shall need a set of two initial values (corresponding to the two state variables at the starting point). In this case, let us pick .
The first step, already mentioned above, is finding the eigenvalues of A in
\begin{bmatrix}x'\\y'\end{bmatrix}=\begin{bmatrix}3&-4\\4&-7\end{bmatrix}\begin{bmatrix}x\\y\end{bmatrix}~.
Once the coefficients of the two variables have been written in the matrix form A displayed above, one may evaluate the eigenvalues. To that end, one finds the determinant of the matrix that is formed when an identity matrix,
In
\det\left(\begin{bmatrix}3&-4\\4&-7\end{bmatrix}-λ\begin{bmatrix}1&0\\0&1\end{bmatrix}\right)~,
Applying further simplification and basic rules of matrix addition yields
\det\begin{bmatrix}3-λ&-4\\4&-7-λ\end{bmatrix}~.
Applying the rules of finding the determinant of a single 2×2 matrix, yields the following elementary quadratic equation,
\det\begin{bmatrix}3-λ&-4\\4&-7-λ\end{bmatrix}=0
-21-3λ+7λ+λ2+16=0
λ2+4λ-5=0~.
Now finding the two roots,
λ1
λ2
λ2+5λ-λ-5=0
λ(λ+5)-1(λ+5)=0
(λ-1)(λ+5)=0
λ=1,-5~.
λ1=1
λ2=-5
As mentioned above, this step involves finding the eigenvectors of A from the information originally provided.
For each of the eigenvalues calculated, we have an individual eigenvector. For the first eigenvalue, which is
λ1=1
\begin{bmatrix}3&-4\\4&-7\end{bmatrix}\begin{bmatrix}\alpha\\\beta\end{bmatrix}=1\begin{bmatrix}\alpha\\\beta\end{bmatrix}.
3\alpha-4\beta=\alpha
\alpha=2\beta~.
All of these calculations have been done only to obtain the last expression, which in our case is . Now taking some arbitrary value, presumably, a small insignificant value, which is much easier to work with, for either or (in most cases, it does not really matter), we substitute it into . Doing so produces a simple vector, which is the required eigenvector for this particular eigenvalue. In our case, we pick, which, in turn determines that and, using the standard vector notation, our vector looks like
\hat{v
Performing the same operation using the second eigenvalue we calculated, which is
λ=-5
\hat{v
This final step finds the required functions that are 'hidden' behind the derivatives given to us originally. There are two functions, because our differential equations deal with two variables.
The equation which involves all the pieces of information that we have previously found, has the following form:
\begin{bmatrix}x\\y\end{bmatrix}=
λ1t | |
Ae |
\hat{v
Substituting the values of eigenvalues and eigenvectors yields
\begin{bmatrix}x\\y\end{bmatrix}=Aet\begin{bmatrix}2\\1\end{bmatrix}+Be-5t\begin{bmatrix}1\\2\end{bmatrix}.
\begin{bmatrix}x\\y\end{bmatrix}=\begin{bmatrix}2&1\\1&2\end{bmatrix}\begin{bmatrix}Aet\\Be-5t\end{bmatrix}.
x=2Aet+Be-5t
y=Aet+2Be-5t.
The above equations are, in fact, the general functions sought, but they are in their general form (with unspecified values of and), whilst we want to actually find their exact forms and solutions. So now we consider the problem’s given initial conditions (the problem including given initial conditions is the so-called initial value problem). Suppose we are given
x(0)=y(0)=1
x(0)=y(0)=1
1=2A+B
1=A+2B~.
Solving these equations, we find that both constants and equal 1/3. Therefore substituting these values into the general form of these two functionsspecifies their exact forms,the two functions sought.
The above problem could have been solved with a direct application of the matrix exponential. That is, we can say that
\begin{bmatrix}x(t)\\y(t)\end{bmatrix}=\exp\left(\begin{bmatrix}3&-4\\4&-7\end{bmatrix}t\right)\begin{bmatrix}x0(t)\\y0(t)\end{bmatrix}
Given that (which can be computed using any suitable tool, such as MATLAB's expm
tool, or by performing matrix diagonalisation and leveraging the property that the matrix exponential of a diagonal matrix is the same as element-wise exponentiation of its elements)
\exp\left(\begin{bmatrix}3&-4\\4&-7\end{bmatrix}t\right)=\begin{bmatrix}4et/3-e-5t/3&2e-5t/3-2et/3\\2et/3-2e-5t/3&4e-5t/3-et/3\end{bmatrix}
the final result is
\begin{bmatrix}x(t)\\y(t)\end{bmatrix}=\begin{bmatrix}4et/3-e-5t/3&2e-5t/3-2et/3\\2et/3-2e-5t/3&4e-5t/3-et/3\end{bmatrix}\begin{bmatrix}1\\1\end{bmatrix}
\begin{bmatrix}x(t)\\y(t)\end{bmatrix}=\begin{bmatrix}e-5t/3+2et/3\ et/3+2e-5t/3\end{bmatrix}
This is the same as the eigenvector approach shown before.