3.4: Linear Systems of Differential Equations

Reading: Notes on Diffy Q’s Section 3.2, Section 3.3, and Section 3.4

An Initial Model

Question 1:

Let \(x(t)\) denote the number of bacteria in a colony (labeled Colony 1) at time \(t\) hours since noon. Let \(y(t)\) denote the number of bacteria in a separate colony (labeled Colony 2) at time \(t\) hours since noon. The size of populations \(x\) and \(y\) can be modeled by system of differential equations

\[\begin{align} \frac{dx}{dt} &= 3x\\ \frac{dy}{dt} &= -2y \end{align}\]

Question 1a:

Explain in practical terms how each of the populations is changing over time. Are the two populations interacting? What happens to the size of each colony in the long run?

Solution to Question 1a:

Question 1b:

Solve each of the differential equations and give general solutions for \(x\) and \(y\).

Solution to Question 1b:

Modeling Interacting Populations

Question 2:

Now imagine Colony 1 and Colony 2 are arranged such that the size of populations \(x\) and \(y\) can be modeled by system of differential equations

\[\begin{align} \frac{dx}{dt} &= 3x+10y\\ \frac{dy}{dt} &= -2y \end{align}\]

Question 2a:

Explain in practical terms how each of the populations is changing over time. Are the two populations interacting? What happens to the size of each colony in the long run?

Solution to Question 2a:

Question 2b:

Notice the solution for \(y\) is the same for this modified setup. Give an educated guess for the general form of the expression for the new solution for \(x(t)\).

Solution to Question 2b:

Question 2c:

Substitute your answer in Question 2b as well as the original solution for \(y(t)\) in Question 1b into the differential equation for \(x\),

\[\frac{dx}{dt} = 3x+10y,\]

and find the general solution for \(x\).

Solution to Question 2c:

Linear Systems of Differential Equations

Question 3:

If we consider a general system of differential equations of the form

\[\begin{align} x' &= ax+by \\ y' &= cx+dy \end{align}\]

then we can generalize our approach in the previous example by guessing solutions of the form

\[x(t) = C_1 e^{r_1t} + C_2 e^{r_2t} \qquad \mbox{and} \qquad y(t) = C_3e^{r_1t} + C_4 e^{r_2t},\]

where we note there is some dependence between constant \(C_1\) and \(C_3\) and constants \(C_2\) and \(C_4\). Explain how we might interpret the meaning of the constants \(r_1\) and \(r_2\) in practical terms.

Solution to Question 3:

Question 4:

By plugging the guess for \(x\) and \(y\) into the differential equation \(x'= ax+by\) you can derive two equations that the values \(r_1\) and \(r_2\) must satisfy. Replace each ?? with an appropriate expression to give the two equations.

Tip

Run the code cell below to help save some time with the calculus and algebraic manipulations.

Solution to Question 4:

\(r_1 C_1 = ??\)

\(r_2 C_2 = ??\)

import sympy as sym

# Creating symbols
r1, r2, C1, C2, C3, C4, a, b, t, x, y = sym.symbols('r1, r2, C1, C2, C3, C4, a, b, t, x, y')  

x = C1 * sym.exp(r1 * t) + C2 * sym.exp(r2 * t)  # General solution for x
y = C3 * sym.exp(r1 * t) + C4 * sym.exp(r2 * t)  # General solution for y


left = x.diff(t,1)  # find formula for x' on left side
right = sym.simplify(a*x + b*y)  # simplify right side

print("The left side is ", left,
     "\n \n The right side side is ", right)

Question 5:

Similar to Question 4, by plugging the guess for \(x\) and \(y\) into the differential equation \(y'= cx+dy\) you can derive two equations that the values \(r_1\) and \(r_2\) must satisfy. Replace each ?? with an appropriate expression to give the two equations.

Tip

Run the code cell below to help save some time with the calculus and algebraic manipulations.

Solution to Question 5:

\(r_1 C_3 = ??\)

\(r_2 C_4 = ??\)

#import sympy as sym  # only need to import once

# Creating symbols
r1, r2, C1, C2, C3, C4, c, d, t, x, y = sym.symbols('r1, r2, C1, C2, C3, C4, c, d, t, x, y')  

x = C1 * sym.exp(r1 * t) + C2 * sym.exp(r2 * t)  # General solution for x
y = C3 * sym.exp(r1 * t) + C4 * sym.exp(r2 * t)  # General solution for y


left = y.diff(t,1)  # find formula for x' on left side
right = sym.simplify(c*x + d*y)  # simplify right side

print("The left side is ", left,
     "\n \n The right side side is ", right)

Question 6:

From Question 4 and Question 5, we get a system of two equations that \(r_1\) must simultaneously satisfy. Write the resulting system of two equations below.

Solution to Question 6:

\(r_1 C_1 = ??\)

\(r_1 C_3 = ??\)

Essentials of Linear Algebra

An \(m \times n\) matrix is a rectangular array of numbers with \(m\) rows and \(n\) columns. For example

\[\begin{bmatrix} 3 & -2 & 7 \\ -12 & 0 & 5 \end{bmatrix}\]

is a 2 by 3 matrix since it has 2 rows and 3 columns.

Multiplication of a Matrix by a Scalar

We can multiple a scalar (a regular number) and a matrix by simply multiplying each value in the matrix by the scalar. For example:

\[2 \begin{bmatrix} 3 & -2 & 7 \\ -12 & 0 & 5 \end{bmatrix} = \begin{bmatrix} 6 & -4 & 14 \\ -24 & 0 & 10 \end{bmatrix} \qquad \mbox{and} \qquad \lambda \begin{bmatrix} C_1 \\ C_3 \end{bmatrix} = \begin{bmatrix} \lambda C_1 \\ \lambda C_3 \end{bmatrix}\]

Adding Matrices

We can add two \(m\) by \(n\) matrices \(A\) and \(B\) by adding values in the same row and column, for example

\[\begin{bmatrix} 3 & -2 & 7 \\ -12 & 0 & 5 \end{bmatrix} + \begin{bmatrix} 0 & 2 & -4 \\ 10 & -3 & 4 \end{bmatrix} = \begin{bmatrix} 3 & 0 & 3 \\ -2 & -3 & 9 \end{bmatrix} .\]

Note

As a result of this definition we cannot add two matrices if they have different dimensions.

Multiplication of Matrices

Multiplication of two matrices is a little more complicated. Let \(A\) denote an \(m\) by \(n\) matrix and \(B\) denote an \(n\) by \(p\) matrix, then we have

\[AB = \begin{bmatrix} \color{dodgerblue}{\mathbf{a_{11}}} & \color{dodgerblue}{\mathbf{a_{12}}} & \color{dodgerblue}{\mathbf{\ldots}} & \color{dodgerblue}{\mathbf{a_{1n}}} \\ a_{21} & a_{22} & \ldots & a_{2n} \\ \vdots & \vdots & & \vdots \\ a_{m1} & a_{m2} & \mathbf{\ldots} & a_{mn} \end{bmatrix} \begin{bmatrix} b_{11} & \color{dodgerblue}{\mathbf{b_{12}}} & \ldots & b_{1p} \\ b_{21} & \color{dodgerblue}{\mathbf{b_{22}}} & \ldots & b_{2p} \\ \vdots & \color{dodgerblue}{\mathbf{ \vdots }} & & \vdots \\ b_{n1} & \color{dodgerblue}{\mathbf{b_{n2}}} & \ldots & b_{np} \end{bmatrix} = \begin{bmatrix} c_{11} & \color{dodgerblue}{\mathbf{c_{12}}} & \ldots & c_{1p} \\ c_{21} & c_{22} & \ldots & c_{2p} \\ \vdots & \vdots & & \vdots \\ c_{m1} & c_{m2} & \ldots & b_{mp} \end{bmatrix}\]

where entry \(c_{ij}\) in the \(i^{\mbox{th}}\) row and \(j^{\mbox{th}}\) column of the product is found by multiplying and adding entries of the \(i^{\mbox{th}}\) row of \(A\) and the \(j^{\mbox{th}}\) column of \(B\) as follows:

\[c_{ij} = a_{i1}b_{1j} + a_{i2}b_{2j} + \ldots + a_{in}b_{nj} = \sum_{k = 1}^n a_{ik}b_{kj}.\]

Note

The product \(\mathbf{AB}\) is only defined if the number of columns of \(\mathbf{A}\) matches the number of rows of \(\mathbf{B}\).

Matrix Arithmetic in Python Using Sympy

The sympy package can be used to perform matrix operations such as multiplication and addition.

We enter matrices using sym.Matrix([ [row1], [row2], ... , [last_row] ]).
We can find the size of a matrix using sym.shape(M).
We can add two matrices \(M\) and \(N\) (assuming they are the same size) using M + N.
We can multiply scalar \(c\) and matrix \(M\) with c * M.
We can multiply two matrices \(M\) and \(N\) (assuming the have the appropriate sizes) with M * N.

The code cells below verify the calculations in the examples above as well as an example of matrix multiplication.

Creating Matrices

#import sympy as sym  # we have already imported sympy

# Create a 2 by 3 matrix M
M = sym.Matrix([[3, -2, 7],
          [-12, 0, 5]])

# Create a 2 by 3 matrix N
N = sym.Matrix([[0, 2, -4],
          [10, -3, 4]])

# Create a 3 by 2 matrix P
P = sym.Matrix([[5, 2], 
          [-1, 8],
          [6, -3]])

Scalar Multiplication

# Scalar multiplication example

2 * M

Adding Matrices

# Matrix Addition

M + N

# See what happens if you try the command below
# M + P

Matrix Multiplication

# Matrix Multiplication

M * P

# See what happens if you try the command below
# M * N

Question 7:

Compute the product \(AB\) for the matrices:

\[A = \begin{bmatrix} 3 & 2 & -1 \\ -7 & 0 & 5 \end{bmatrix} \qquad \mbox{and} \qquad B = \begin{bmatrix} 1 \\ 2 \\ -1 \end{bmatrix}\]

Feel free to compute by hand or using Python by inserting a code cell to do the calculations.

Solution to Question 7:

\[AB = \begin{bmatrix} ?? \\ ?? \end{bmatrix}\]

Question 8:

Compute the product \(AB\) for the matrices:

\[A = \begin{bmatrix} a & b \\ c & d \end{bmatrix} \qquad \mbox{and} \qquad B = \begin{bmatrix} C_1 \\ C_2 \end{bmatrix}\]

Feel free to compute by hand or using Python by inserting a code cell to do the calculations.

Solution to Question 8:

\[AB = \begin{bmatrix} ?? \\ ?? \end{bmatrix}\]

Matrix Terminology

Some important terminology when working with matrices:

A matrix that has only one column is often called a column vector.
A column vector that is a column of all zeros is called the zero vector and denoted \(\mathbf{0}\) or \(\vec{0}\) in order to distinguish itself from the scalar value 0.
An \(n \times n\) matrix (same number of rows and columns) is called a square matrix.

Question 9:

Give the matrix \(A\) such that the system of equations for \(r_1\) in Question 6 can be written in matrix form as

\[A \begin{bmatrix} C_1 \\ C_3 \end{bmatrix} = r_1 \begin{bmatrix} C_1 \\ C_3 \end{bmatrix}.\]

Solution to Question 9:

\[A = \begin{bmatrix} ?? & ??\\ ?? & ?? \end{bmatrix}\]

Expressing Systems in Matrix Form

Let \(\mathbf{x'}\) denote a column vector of derivatives such as

\[\mathbf{x'} = \begin{bmatrix} \dfrac{dx_1}{dt} \\ \dfrac{dx_2}{dt} \\ \dfrac{dx_3}{dt} \end{bmatrix}.\]

Question 10:

Give the matrix \(A\) such that the system of differential equations below can be written in the form \(\mathbf{x'} = A \mathbf{x}\).

\[\begin{array}{rl} \dfrac{dx_1}{dt} &= 4 x_1 + 7 x_2 - x_3 \\ \dfrac{dx_2}{dt} &= -2x_1-11x_3\\ \dfrac{dx_3}{dt} &= 8x_3-x_2 \end{array}\]

Solution to Question 10:

\[\begin{bmatrix} \dfrac{dx_1}{dt} \\ \dfrac{dx_2}{dt} \\ \dfrac{dx_3}{dt} \end{bmatrix} = \begin{bmatrix} ?? & ?? & ?? \\ ?? & ?? & ?? \\ ?? & ?? & ?? \end{bmatrix} \begin{bmatrix} x_1 \\ x_2 \\ x_3 \end{bmatrix}.\]

Question 11:

Express the system below in matrix form.

\[\begin{array}{rl} x' &= \cos{(2t)}x + \sin{(2t)}z\\ y' &= e^ty \\ z' &= \sin{(2t)}x + \cos{(2t)}z \end{array}\]

Solution to Question 11:

\[\begin{bmatrix} x' \\ y' \\ z' \end{bmatrix} = \begin{bmatrix} ?? & ?? & ?? \\ ?? & ?? & ?? \\ ?? & ?? & ?? \end{bmatrix} \begin{bmatrix} x \\ y \\ z \end{bmatrix}.\]

Eigenvalues of a \(2 \times 2\) Matrix

We say a scalar \(\lambda\) is an eigenvalue of a square matrix \(A\) if there exists a nonzero vector \(\mathbf{v}\) such that \(A\mathbf{v} = \lambda \mathbf{v}\). The eigenvalues of matrix have many useful applications and interpretations. We will see in the context of differential equations, they can tell us very important information about how the solutions behave.

We will mostly be working with a system of two differential equations, so it will suffice to restrict our attention to the case where \(A\) is a \(2 \times 2\) matrix, but the discussion below can be generalized to deal with much larger systems. In the \(2 \times 2\) case, \(\lambda\) is an eigenvalue of \(A\) if there exists a nonzero vector \(\mathbf{v}\) such that

\[\begin{align} \lambda \mathbf{v} &= \begin{bmatrix} a & b \\ c & d \end{bmatrix} \mathbf{v} \\ \textbf{0} &= \begin{bmatrix} a & b \\ c & d \end{bmatrix} \mathbf{v} - \lambda \mathbf{v} \\ \textbf{0} &= \begin{bmatrix} a & b \\ c & d \end{bmatrix} \mathbf{v} - \begin{bmatrix} \lambda & 0 \\ 0 & \lambda \end{bmatrix} \mathbf{v} \\ \textbf{0} &= \begin{bmatrix} a-\lambda & b \\ c & d-\lambda \end{bmatrix} \mathbf{v} \end{align}\]

Using linear algebra, we can show that \(\lambda\) is an eigenvalue of \(A\) if and only if

\[ \mbox{det} \left( \begin{bmatrix} a-\lambda & b \\ c & d-\lambda \end{bmatrix} \right) = \color{dodgerblue}{(a-\lambda)(d-\lambda) - bc = 0}.\]

The resulting quadratic equation (in \(\lambda\)) is called the characteristic equation for \(A\).

Note

We can only find eigenvalues of square matrices. We cannot find eigenvalues for matrices whose number of rows and columns are not equal.

Finding Eigenvalues with Sympy

Let \(M\) denote a square matrix.

To find the characteristic polynomial corresponding \(M\), use M.charpoly().
To find the eigenvalues of \(M\), use the command M.eigenvals().

Example: Finding Eigenvalues of \(3 \times 3\) Matrix

Consider the \(3 \times 3\) matrix

\[M = \begin{bmatrix} -13 & -8 & -4 \\ 12 & 7 & 4\\ 24 & 16 & 7 \end{bmatrix}\]

In the first code cell below, we find the eigenvalues by finding the zeros of the characteristic equation.
In the second code cell below, we check our eigenvalues using sym.eigenvals(M).

# Finding Eigenvalues with Characteristic Polynomials

# import sympy as sym  # we have already imported

# Note the lambda is reserved, so we misspell on purpose
lamda = sym.symbols('lamda')  

# Define matrix M
M = sym.Matrix([[-13, -8, -4],
                [12, 7, 4], 
                [24, 16, 7]])

# Find characteristic polynomial of M
p = M.charpoly(lamda)

# Factor the characteristic equation
sym.factor(p.as_expr())

# Checking our previous result

M.eigenvals()

A Step by Step Process For Solving Linear Systems of ODES

Recall the model for bacteria populations \(x\) and \(y\) in Colony 1 and Colony 2 in Question 2. Follow the steps below to find general solutions to the system

\[\begin{align} \frac{dx}{dt} &= 3x+10y \\ \frac{dy}{dt} &= -2y \end{align}.\]

Step 1:

Write the system in matrix form and identify the matrix \(A\).

Solution to Step 1:

\[A = \begin{bmatrix} ?? & ?? \\ ?? & ?? \end{bmatrix}\]

Step 2:

Determine the characteristic equation.

Tip

You may find the equation by hand or by inserting a code cell to help with the calculation.

Solution to Step 2:

Step 3:

Solve the characteristic equation.

Tip

You may solve the equation by hand or by inserting a code cell to help with the calculation.

Solution to Step 3:

Step 4:

Give a general solution for \(x(t)\).

Solution to Step 4:

Step 5:

Plug your solution for \(x(t)\) into the differential equation \(x'\) and solve for \(y(t)\).