Notes from Chapter 2 of Hirsch, Smale, Devaney (HSD).
Basics of Systems of ODEs
A system of differential equations is a collection of \(n\) interrelated differential equations: \[
\begin{aligned}
x'_1 &= f_1(t, x_1, x_2, \ldots, x_n) \\
x'_2 &= f_2(t, x_1, x_2, \ldots, x_n) \\
\vdots \\
x'_n &= f_n(t, x_1, x_2, \ldots, x_n)
\end{aligned}
\]
Each of these functions is assumed to be \(C^\infty\), which means that we can take partial derivatives of any order. This assumption is mostly for convenience. Using vector notation, we can rewrite the system compactly as \(X' = F(t, X)\), where the RHS is simply
A solution to this system would be \(X(t)\), such that \(X'(t) = F(t, X(t))\). If all the \(f_i\) inside \(F(t, X)\) do not depend on time, then the system is called autonomous and we can drop the dependence on time, such that \(X' = F(X)\). Otherwise, it is a non-autonomous system. Similar to the previous chapter, the equilibrium points for an autonomous system occur at \(X_\ast\), giving us \(F(X_\ast)=\mathbf{0}\).
A solution is a time-dependent state \(X(t)\) in the phase space. As \(t\) changes, the points \(X(t)\) trace out a curve called the trajectory.
For an autonomous system \(X' = F(X)\), an equilibrium is a state \(X_\ast\) with \(F(X_\ast)=\mathbf{0}\). It corresponds to the constant solution \(X(t)\equiv X_\ast\), so its trajectory is just a single point. The set of equilibria is the solution set of the algebraic equation \(F(X)=\mathbf{0}\), and it can be a discrete set of points or a whole curve/region.
Second Order Differential Equations
A second order ODE looks like \(x'' = f(t, x, x')\). A popular example is Newton’s equation, \(mx'' = f(x)\). Or, there is the forced harmonic oscillator, \(mx'' + bx' + kx = f(t)\). But the general, simple second order ODE is \(x'' + ax' + bx = 0\), where \(a\) and \(b\) are constants.
The upshot is that we can treat second order ODEs as systems of first order equations by defining a new variable \(y = x'\). Then, we can rewrite the second order ODE as the system: \[
\begin{aligned}
x' &= y \\
y' &= -bx -ay,
\end{aligned}
\]
which we might generalize as \[
\begin{aligned}
x' &= f(x,y) \\
y' &= g(x,y).
\end{aligned}
\]
Collectively, in a 2D system, \([f(x,y), g(x,y)]\) represents \(F(x,y)\), a vector with \(x\) and \(y\) components. As an example for an ODE, we might have the system: \[
\begin{aligned}
x' &= y, \\
y' &= -x.
\end{aligned}
\]
We can draw the direction field of this system, which produces a 2D vector \(F(x,y) = [f(x,y), g(x,y)]\) for each \((x,y)\) position on the plane. In other words, at each point \((x,y)\), we place a vector showing the direction of the system. So far, we haven’t solved the ODE; but this gives us a visual representation of how the system behaves at each point.
Code
X, Y = np.meshgrid(np.arange(-2, 2, 0.1), np.arange(-2, 2, 0.1))dx = Ydy =-Xmag = np.sqrt(dx**2+ dy**2)def sol(t, a): x = a * np.sin(t) y = a * np.cos(t)return x, yfig, ax = plt.subplots(figsize=(7,7))ax.quiver(X, Y, dx / mag, dy / mag, color="r")ax.plot(0, 0, "ro", label="Equilibrium")ax.plot(*sol(np.linspace(0, 2*np.pi, 100), a=0.5), color="r")ax.plot(*sol(np.linspace(0, 2*np.pi, 100), a=1), color="r")ax.plot(*sol(np.linspace(0, 2*np.pi, 100), a=1.5), color="r")ax.set_xlim(-2, 2)ax.set_ylim(-2, 2)ax.set_aspect("equal", adjustable="box")ax.set_xlabel("$x$")ax.set_ylabel("$y$")ax.set_title("Direction Field")fig.tight_layout()plt.show()
Direction Field
Unsurprisingly, the solution is: \[
\begin{aligned}
x(t) &= a \sin(t) \\
y(t) &= a \cos(t).
\end{aligned}
\]
We can quickly verify that this is the case \[
\begin{aligned}
x'(t) &= a\cos(t) = y(t) \\
y'(t) &= -a\sin(t) = -x(t).
\end{aligned}
\]
If we pick an \(a\), we get a circle of radius \(a\) centered at the origin. Of course, when \(a=0\), our solution is the constant function \(x(t)=0\) and \(y(t)=0\), which is the equilibrium point at the origin.
If you carry out the matrix multiplication, you’ll get back the system we started with!
In general, any 2D planar linear system can be written as \[
\begin{bmatrix}
x' \\
y'
\end{bmatrix}
=
\begin{bmatrix}
a & b \\
c & d
\end{bmatrix}
\begin{bmatrix}
x \\
y
\end{bmatrix},
\]
where we can denote the matrix of constants as \(A\). The origin is always an equilibrium (constant) solution. Other equilibria, which would also be solutions, can be found by solving the system \[
\begin{aligned}
ax + by &= 0 \\
cx + dy &= 0
\end{aligned}
\]
A way to check if there are non-zero solutions is by computing the determinant of \(A\). When \(\det(A) = 0\), the transformation “flattens” area (or volume) and the null space is nontrivial, so \(AX=\mathbf{0}\) has nonzero solutions \(X\). If \(\det(A) \ne 0\), then the only equilibrium is at the origin.
If we have a “special”, nonzero vector \(V_0\) such that \(AV_0 = \lambda_0 V_0\), then we claim that \[
X(t) = e^{\lambda_0 t} V_0
\]
is a solution (this is the solution with initial condition \(X(0)=V_0\)). These special objects have names: \(\lambda_0\) is the eigenvalue associated with the eigenvector\(V_0\). The fact that they provide a solution can be verified by computing \[
\begin{aligned}
X'(t) &= \lambda_0 e^{\lambda_0 t} V_0 \\
&= e^{\lambda_0 t} (\lambda_0 V_0)\\
&= e^{\lambda_0 t} (A V_0)\\
&= A (e^{\lambda_0 t} V_0)\\
&= A X(t)
\end{aligned}
\]
If we scale the eigenvector by some constant \(\alpha\), we get another solution \(X(t) = e^{\lambda_0 t} \alpha V_0\). For a 2D system, if we have 2 linearly independent eigenvectors, then we have two sets of solutions. Each of these eigenvector solutions is a straight-line solution.
The straight-line, here, is in the context of the phase space. The time variable \(t\) is just a parameter: as \(t\) varies, the state \(X(t)\) moves around in phase space, and the trajectory is the image set \(\{X(t): t\in\mathbb{R}\}\). For this eigenvector solution, every point has the form \(X(t)=\big(e^{\lambda_0 t}\alpha\big)V_0\), so the trajectory lies on the line through the origin in the direction of \(V_0\). In fact, since \(e^{\lambda_0 t}>0\) for all real \(t\), the trajectory stays on the same ray (half-line) from the origin passing through \(\alpha V_0\). For any \(t\), the distance from the origin is \[
\lVert X(t)\rVert = e^{\lambda_0 t}\,|\alpha|\,\lVert V_0\rVert.
\]
So, if \(\lambda_0>0\), then \(\lVert X(t)\rVert \to \infty\) as \(t\to\infty\), and \(X(t)\to (0,0)\) as \(t\to -\infty\). If \(\lambda_0<0\), the opposite happens: trajectories along this eigenvector ray decay into the origin forward in time. Finally, if \(\lambda_0=0\), then \(X(t)=\alpha V_0\) is constant, so the trajectory is just a single point.
All Solutions
Let’s suppose our 2D system has two linearly independent eigenvectors \(V_1,V_2\) with eigenvalues \(\lambda_1,\lambda_2\). The corresponding eigenvector solutions are
Since \(V_1,V_2\) form a basis for the plane, any initial value \(X(0)\) can be written as a linear combination
\[
X(0) = \alpha V_1 + \beta V_2.
\]
Now consider the function \(X(t) = \alpha X_1(t) + \beta X_2(t)\). We can verify that \(X(t)\) is also a solution by computing \[
\begin{aligned}
X'(t) &= \alpha X_1'(t) + \beta X_2'(t) \\
&= \alpha A X_1(t) + \beta A X_2(t) \\
&= A (\alpha X_1(t) + \beta X_2(t)) \\
&= A X(t).
\end{aligned}
\]
Thus, \(X(t)\) is a solution. The key property we used was just linearity: any linear combination of solutions is again a solution.
Concretely, for the homogeneous linear system \(X'=AX\), if \(Y_1(t)\) and \(Y_2(t)\) are solutions, then for any constants \(\alpha,\beta\), the function \[
X(t)=\alpha Y_1(t)+\beta Y_2(t)
\] is also a solution.
Let \(Y_1(0)\) and \(Y_2(0)\) denote the values of \(Y_1(t)\) and \(Y_2(t)\) at \(t=0\). If these initial vectors are linearly independent, then they form a basis of \(\mathbb{R}^2\). Given any initial value \(X(0)\in\mathbb{R}^2\), there are unique constants \(\alpha,\beta\) such that \[
X(0)=\alpha Y_1(0)+\beta Y_2(0),
\] and then \(X(t)=\alpha Y_1(t)+\beta Y_2(t)\) is the solution satisfying this initial condition.
