Artificial Intelligence 🤖
Differential equations
Ordinary differential equations
Autonomous Equations & Stability

Autonomous equations & stability

In Chapter 11, we mention qualitative techniques which allow us to examine the behaviour of solutions. While these offer approximate solutions, they are very powerful and are especially used in solving autonomous differential equations. Their great strength comes from the structural insight they give into differential equations and their solutions so we may want to use them even when analytic methods work just fine. In this section we discuss one-dimensional (1D), autonomous ODEs as well as introduce the concept of stability of solutions.

First order autonomous differential equations take the form

dydt=f(y)\frac{d y}{d t}=f(y)

where the function ff does not explicitly depend on the independent variable. Autonomous equations of the form (12.97) are separable, yielding

1f(y)dy=t+c\int \frac{1}{f(y)} d y=t+c

where cc is a constant. While Eq. (12.98) gives an integral form of the solution, the integral itself may not be expressible as a closed-form solution resulting in the above equation being of little use. In such cases, we often resort to qualitative methods to extract information on the behaviour of solutions without actually solving the ODE. This is especially important for autonomous ODEs as they are used to model many physical processes in the real world. Let us consider a value of yy, say y=y0y=y_{0}, that renders the function f(y)f(y) zero, i.e. f(y0)=0f\left(y_{0}\right)=0. Recall from our discussion on direction fields (see Chapter 11, Subsec.11.2.1) that the arrows here have zero slope for all values of the independent variable, tt. Therefore if we start the solution at a point y=y0y=y_{0}, then, following the zero-slope direction field arrows, it is obvious that the solution will stay at that value for all tt. This means that y(t)=y0y(t)=y_{0} for all tt is a constant solution to the ODE.

For example, in Chapter 11 (Subsec.11.2.1), we considered a differential equation for the velocity v(t)v(t) of a falling object, given by dv/dt=f(v)=102vd v / d t=f(v)=10-2 v and we have shown that at v=5,f(v)=0v=5, f(v)=0 and hence v(t)=5v(t)=5 is a constant solution to the ODE. Let us revisit now the notation y(t)y(t) where yy is the dependent variable and tt is the independent variable. The points y=y0y=y_{0} that make f(y)f(y) zero are called critical points (also known as equilibrium points or fixed points). The corresponding constant solutions given by y(t)=y0y(t)=y_{0} for all values of tt are known as equilibrium solutions (or constant solutions).

Qualitative solutions & isoclines revisited

Consider the autonomous ODE,

dydt=y3\frac{d y}{d t}=y-3

Using the method of isoclines, we observe that the direction field for (12.99) shows horizontal isoclines which are independent of tt. The equation of the isoclines is y=3+cy=3+c where cc represents the slope of the arrows. A few isoclines are shown in Fig. 12.5 (coloured lines) labelled by their corresponding slope, cc. The direction field for Eq. (12.99) is also shown. At y=3y=3, the function f(y)f(y) in (12.99) is equal to zero which implies that y=3y=3 is a critical point. This corresponds to an equilibrium solution given by y(t)=3y(t)=3 which remains constant for all values of tt - this is represented by the integral curve shown in red (horizontal, dashed line). We observe from the direction field that, for y<3y<3, all solutions decrease to -\infty while for y>3y>3, all solutions increase to ++\infty. Four additional integral curves are shown in Fig. 12.5; for y>3y>3, the black, dashed curve is generated using the IC, (t,y)(t, y) as (3,6)(3,6) and the purple, dashed curve using (4,6)(4,6). The purple curve is the same curve as the black one but translated to the right. Similarly, for y<3y<3, the purple curve is the black translated to the right. This behaviour is typical for ODEs whose isoclines are horizontal.

What kind of information are we interested in?

  • Given an autonomous equation in the form of ODE (12.99), we are interested in the critical points. These are the values of yy for which
f(y)=0.f(y)=0 .

While f(y)f(y) may be a complicated function, solving the above involves solving an algebraic equation which is generally much easier than solving a differential equation. The critical points correspond to equilibrium solutions (i.e. constant solutions).

Figure 12.5: Direction field and several labelled isoclines for ODE (12.99). Five integral curves are shown and generated with the following initial conditions (t,y):(3,3)(t, y):(3,3) (red), (3,6)(3,6) (black, dashed), (4,6)(4,6) (purple, dashed), (3,1) (black, dash-dot), (4,1) (purple, dash-dot).

  • We then need to determine the behaviour of the nonequilibrium solutions. To answer this, we need to evaluate the direction field above, below, and in-between equilibria. This requires knowledge of where f(y)>0f(y)>0 and f(y)<0f(y)<0 (we already know where f(y)=0)f(y)=0).
  • if f(y)>0f(y)>0, then y>0y^{\prime}>0 and hence solutions are increasing;
  • if f(y)<0f(y)<0, then y<0y^{\prime}<0 and hence solutions are decreasing.

Going back to the example of (12.99), if the initial point in yy is any value above the equilibrium solution y=3y=3, the solutions move away from the equilibrium solution and increase indefinitely. If the initial point in yy is given at any value below the equilibrium solution, the solutions again move away from the equilibrium and decrease indefinitely. In fact, unless the initial point is given exactly at y=3y=3, all other initial conditions diverge away from the equilibrium. Further, if at any point after the initial time, yy is slightly perturbed such that y=3±εy=3 \pm \varepsilon (where ε\varepsilon is very small and positive) then, the solutions shoot away from the constant solution to ±\pm \infty. This behaviour makes the equilibrium solution unstable. We study the concept of stability further next.

Stability of solutions

First order ODEs represent 1D systems which, geometrically, define flows on a line or a circle. 21{ }^{21} We move away from the yty-t plots and focus on just the phase space. The phase

21{ }^{21} For more information, read up on the definition of the dimension of a mathematical space. space is represented by all possible real yy values 22^{22}, i.e. yRy \in \mathbb{R}. For 1D1 \mathrm{D} flows on a line, this is simply given by R\mathbb{R} (i.e. the real line). If we were to draw the phase space then for ODE (12.99), this would be the line yy, as shown next, in the schematic. The round marker at y=3y^{*}=3 indicates the critical point; note that the asterisk notation is commonly used to indicate a fixed/critical point.

At y=3y=3 we know from Fig. 12.5 that t,y\forall t, y remains constant. What happens near the critical point? Equation (12.99) has a single critical point so we only need to know the sign and magnitude of f(y)f(y) for y<3y<3 and y>3y>3. For multiple critical points, we need to determine f(y)f(y) in between equilibrium points as well. As shown before, f(y)<0f(y)<0 for y<3y<3 and f(y)>0f(y)>0 for y>3y>3. It follows that for y<3y<3, solutions decrease; this is indicated below on the phase line by an arrow pointing towards decreasing values of yy. Similarly, an arrow points towards increasing values of yy for y>3y>3.

Consider a state of the dynamical system at a particular point in time; this is represented by a point (or phase point) in the phase space. Now, as time progresses, this point (whose initial location is known to us - this is the initial condition) moves along a curve which is called the phase trajectory, beginning at the initial point. A set of phase trajectories with all possible initial conditions exhibiting different qualitative features of the flow, constitutes a phase portrait. In 1D on a line, the phase portrait looks like the figure above: the phase line marked with critical points and trajectories (shown by red arrows) which dictate the behaviour of the nonequilibruim solutions. The phase portrait is far more interesting in higher dimensions since, in 1D, trajectories only move towards or away from critical points, i.e. left or right on the phase line. We study 2D2 \mathrm{D} systems given by two first order ODEs or, equivalently, by a second order ODE later in Chapter 14. Section 14.2.

In our example above, all trajectories move away from the equilibrium y=3y=3. Such a point, where all trajectories move away from it, is called unstable. There are three types of stability:

  1. Unstable: the arrows (corresponding to the direction of trajectories) move away from the equilibrium solution, i.e. \leftarrow \circ \rightarrow.
  2. Asymptotically stable: the arrows move toward the equilibrium solution, i.e. \rightarrow \bullet \leftarrow.

22{ }^{22} Note that if yy represents a quantity like concentration for example, the phase space would be all possible, nonnegative, real values so the phase space naturally relies on the application. 3. Semi-stable: non-equilibrium solutions approach the critical point from one side but they move away from it from the other, i.e. o\rightarrow \mathbf{o} \rightarrow

As shown above, we indicate the type of stability on the phase line by a solid round marker for the stable point, a hollow round marker for the unstable point and a half-filled round marker for the semi-stable point. For ODE (12.99), we have a single, unstable point. Semi-stable points fall beyond the scope of this course. They appear in bifurcations which occur when the qualitative behaviour of the system changes. Changes in the qualitative behaviour include the creation/destruction of critical points or the exchange of stability between equilibria 23{ }^{23}.

Logistic growth example

Autonomous ODEs are often used to describe how population changes. The term 'population' may refer to people, animals, bacteria as well as the spread of a disease, a rumour etc.

A simple population growth model is given by

dydt=ky\frac{d y}{d t}=k y

where kk is a constant growth rate and y(t)y(t) denotes the population of a given species at time tt. For k>0k>0, Eq. (12.100) implies that the population keeps increasing indefinitely. This is unrealistic; at some point in time, resources like space and food run out and therefore population growth cannot be sustained. A logistic growth model gives a more realistic representation of population changes. To arrive at the logistic growth model, the growth rate, kk from Eq. (12.100), is replaced by a linear relationship, k=abyk=a-b y, where a,ba, b are positive constants. The interpretation is that, as the population yy increases, the growth rate declines. With the linear growth rate, Eq. (12.100) becomes

dydt=(aby)y\frac{d y}{d t}=(a-b y) y

In Eq. (12.101), if yy is small, i.e. y=O(ε),yayy=\mathcal{O}(\varepsilon), y^{\prime} \approx a y, since y2<<yy^{2}<<y but, as yy gets bigger, the by2-b y^{2} term dominates over the aya y term.

The objective here is to obtain qualitative information on the behaviour of solutions to Eq. (12.101) without explicitly solving the ODE. More specifically, we want to:

  • find critical points and equilibrium solutions;
  • determine the stability of critical points.

Following Subsec. 12.6.1, we first solve f(y)=0f(y)=0 for critical points, where f(y)=(aby)yf(y)=(a-b y) y from Eq. (12.101). This yields y1=0y_{1}^{*}=0 and y2=a/by_{2}^{*}=a / b as the two equilibria which represent constant solutions.

Next, we assess the stability of the two equilibria by looking at the graph of f(y)f(y). This is shown in Fig. 12.6 with the constant parameters aa and bb chosen as a=4,b=1a=4, b=1;

23{ }^{23} A truly excellent textbook on dynamical systems is Nonlinear Dynamics and Chaos: With Applications to Physics, Biology, Chemistry and Engineering (Studies in Nonlinearity) by Steven H. Strogatz. Bifurcations in 1D as well as in higher dimensions are thoroughly discussed. the slope is zero for the critical points y1=0y_{1}^{*}=0 and y2=4y_{2}^{*}=4, positive in 0<y<40<y<4 and negative elsewhere. The critical points are indicated by circle markers and the arrows on the horizontal axis (indicating yy ) show the direction nonequilibrium solutions follow. For instance, if initially, y=1y=1, the solutions approach the stable equilibrium at y2=4y_{2}^{*}=4. The physical interpretation is that the population increases and approaches the limiting value of 4 . In fact, no matter the initial population, the limiting population is always 4 . Note that while, mathematically, solutions for y<0y<0 exist, physically they do not since yy represents population. The horizontal axis in Fig. 12.6 (shown in blue) marked with the critical points and the phase trajectories is the phase portrait corresponding to ODE (12.101).

Figure 12.6: The graph of f(y)f(y) against yy where f(y)=(aby)yf(y)=(a-b y) y with a=4,b=1a=4, b=1. The critical point at y1=0y_{1}^{*}=0 is unstable (hollow marker) while the critical point at y2=4y_{2}^{*}=4 is stable (solid marker). The arrows indicate the direction nonequilibrium solutions approach.

Figure 12.7 shows the time-dependent solutions for y(t)y(t) starting at different initial conditions. The black horizontal lines indicate the constant solutions y1y_{1}^{*} and y2y_{2}^{*}. Two nonequilibrium solutions are shown in red and blue. The red curve is generated using y(0)=7y(0)=7 as the initial condition; as time evolves, the integral curve approaches the stable equilibrium at y=4y=4. The blue curve is generated using y(0)=0.1y(0)=0.1; as time evolves, population increases moving away from the unstable and approaching the stable equilibrium. The blue integral curve is an S\mathcal{S}-shape or curve (sigmoid curve), characteristic of logistic growth problems.

As mentioned earlier and as depicted in Fig. 12.7, in the example described by ODE (12.101), the population always approaches the limiting value (in this case this is equal to 4). One can argue that if y(t0)=0y\left(t_{0}\right)=0 where t0t_{0} indicates some initial time, the population remains constant at 0 . If, however, at some point in time, t>t0,y>0t>t_{0}, y>0, since y1=0y_{1}^{*}=0 is

Figure 12.7: Plots of the solutions y(t)y(t) to ODE (12.101) with a=4,b=1a=4, b=1 for 0t50 \leq t \leq 5. The black integral curves represent the constant solutions at y1=0y_{1}^{*}=0 and y2=4ty_{2}^{*}=4 \forall t while two nonequilibrium solutions are depicted in red and blue. The red integral curve is generated using y(0)=7y(0)=7 while the blue curve is generated using y(0)=0.1y(0)=0.1 as initial conditions.

unstable, the population increases to approach the limiting value according to the model given by (12.101). Another phenomenon often observed in nature is that, due to low initial species mass, the species is ultimately driven to extinction. This could be a result of environmental factors, species competition, the inability to find a mate, etc. In such a case, the y=0y^{*}=0 equilibrium, which represents extinction, is stable thus neighbouring nonequilibrium solutions decrease with time to approach extinction (see section exercises for an example of logistic growth with threshold.).

Note on stability

We end the section on dynamical systems and autonomous equations with a note on stability. We won't focus much on this but it is useful to have in mind for dynamical systems in general which you will additionally encounter in Process dynamics & control in Year II. We have already mentioned the various types of stability equilibria can take and have seen stable and unstable equilibria emerge from models on population dynamics. For a 1D system with more than one fixed point, when we talk about a stable fixed point we refer to it as being locally stable and not globally stable. To clarify this, let us consider a system with 3 critical points with y1,y3y_{1}^{*}, y_{3}^{*} as stable and y2y_{2}^{*} as unstable. The corresponding phase portrait is shown below.

Dynamical systems are often subject to perturbations or noise. It follows that when we are at an equilibrium solution, a small perturbation at some point in time, causes the solution to shift from the equilibrium before it eventually settles to equilibrium again. In population dynamics, perturbations may be a result of abrupt environmental fluctuations. This equilibrium shift may or may not influence the state of the dynamical system. Consider again the phase portrait shown above. Suppose we have an imaginary particle at a stable critical point and the system is perturbed, i.e. a disturbance is introduced into the system and we have a shift in yy from, say y1y_{1}^{*} to y2y_{2}^{*}. The particle then moves to the left approaching the point y1y_{1}^{*} since the critical point is attracting. If, however, the disturbance introduced is large and the particle shifts to y>y2y>y_{2}^{*} from y1y_{1}^{*} then, the particle will approach the next stable point at y3y_{3}^{*}. Consequently, the system is now stable at a higher equilibrium point. For this reason, we refer to y1,y3y_{1}^{*}, y_{3}^{*} as locally stable. A globally stable point in 1D1 \mathrm{D} is one for which the region of attraction is the entire real line R\mathbb{R}.

Section 12.6 exercises The population y(t)y(t) of a certain species is modelled by the following ODE:

dydt=ry(1yT)(1yK)\frac{d y}{d t}=-r y\left(1-\frac{y}{T}\right)\left(1-\frac{y}{K}\right)

where r>0,0<T<Kr>0,0<T<K. Find equilibria, assess their stability and draw the phase portrait. Sketch the equilibrium and some nonequilibrium solutions in the yy - tt plane.