0% found this document useful (0 votes)
8 views103 pages

Dynamical Systems Notes

The document contains lecture notes for the Dynamical Systems A course at the University of Johannesburg, covering topics such as the introduction to dynamical systems, first-order linear discrete systems, and the classification of difference equations. It includes definitions, theorems, examples, and exercises to illustrate key concepts. The notes emphasize the qualitative behavior of systems and provide various applications in fields like finance and biology.

Uploaded by

7868bm24tr
Copyright
© All Rights Reserved
We take content rights seriously. If you suspect this is your content, claim it here.
Available Formats
Download as PDF, TXT or read online on Scribd
0% found this document useful (0 votes)
8 views103 pages

Dynamical Systems Notes

The document contains lecture notes for the Dynamical Systems A course at the University of Johannesburg, covering topics such as the introduction to dynamical systems, first-order linear discrete systems, and the classification of difference equations. It includes definitions, theorems, examples, and exercises to illustrate key concepts. The notes emphasize the qualitative behavior of systems and provide various applications in fields like finance and biology.

Uploaded by

7868bm24tr
Copyright
© All Rights Reserved
We take content rights seriously. If you suspect this is your content, claim it here.
Available Formats
Download as PDF, TXT or read online on Scribd

Dynamical Systems A

APM08X1 Honours Course Notes–WEEK 1

Prof F Nyabadza

University of Johannesburg
Department of Applied Mathematics

These notes were prepared by Prof Nyabadza for the Dynamical Systems A course for
honours students at UJ.
Dynamical Systems A Lecture Notes

Contents
1 Introduction to Dynamical Systems 3

2 First-Order Linear Discrete Systems 4

3 Classification of Difference Equations 6

4 Exercises and Problem Sets 8

2
Dynamical Systems A Lecture Notes

1 Introduction to Dynamical Systems


What is a Dynamical System?
Definition 1.
A dynamical system is a mathematical formalism for describing how the state of a
system evolves over time according to deterministic rules. It typically takes the form:
dx
ẋ = = f (x, t, µ)
dt
where:
• x ∈ Rn is the state vector
• t ∈ R is the independent variable (time)
• µ ∈ Rm is a vector of parameters
• f : Rn × R × Rm → Rn is the vector field
Key Concept
Classification of Dynamical Systems:

• Continuous vs Discrete: Time can be continuous (t ∈ R) or discrete


(t ∈ Z)

• Autonomous vs Non-autonomous: Does the system depend explicitly on


time?

• Linear vs Nonlinear: Is the vector field a linear function of x?

• Dimensionality: Number of state variables (n)

Continuous vs Discrete: A Philosophical Distinction


Important
ODE vs Dynamical Systems:

• ODE Approach: Focus on finding explicit solutions (analytical methods,


integration techniques)

• Dynamical Systems Approach: Focus on qualitative behavior without


necessarily finding explicit solutions (stability analysis, bifurcations, chaos
theory)

Discrete Dynamical Systems: The Focus of This Course


For a discrete-time, one-dimensional, autonomous system:
xt+1 = f (xt , µ), t = 0, 1, 2, . . .
where xt ∈ R and µ ∈ R is a parameter.

3
Dynamical Systems A Lecture Notes

2 First-Order Linear Discrete Systems


The General Linear Equation
Consider the system:
xt+1 = αxt + β, t = 0, 1, 2, . . . (1)
with initial condition x0 .

Theorem 1.
The general solution to (1) is given by:
  t

t 1 − α
α x0 + β if α ̸= 1

xt = 1−α

x + βt
0 if α = 1

Proof. Proceed by iteration:

x1 = αx0 + β
x2 = αx1 + β = α2 x0 + αβ + β
x3 = αx2 + β = α3 x0 + α2 β + αβ + β
..
.
t−1
X
xt = α t x0 + β αk
k=0

1−αt
Pt−1
The geometric series k=0 αk evaluates to 1−α
when α ̸= 1 and t when α = 1.

Fixed Points and Stability


Definition 2.
A fixed point (or equilibrium point) of xt+1 = f (xt ) is a point x∗ such that:

f (x∗ ) = x∗

For equation (1), the fixed point is:

β
x∗ = (for α ̸= 1)
1−α
Definition 3.
Let x∗ be a fixed point of xt+1 = f (xt ).

• x∗ is stable if for every ϵ > 0, there exists δ > 0 such that |x0 − x∗ | < δ implies
|xt − x∗ | < ϵ for all t ≥ 0.

• x∗ is asymptotically stable if it is stable and there exists δ > 0 such that |x0 −
x∗ | < δ implies limt→∞ xt = x∗ .

• x∗ is unstable if it is not stable.

4
Dynamical Systems A Lecture Notes

Theorem 2.
For xt+1 = αxt + β:

• If |α| < 1, the fixed point x∗ is asymptotically stable.

• If |α| > 1, the fixed point is unstable.

• If |α| = 1, the fixed point is stable but not asymptotically stable.

Applications and Examples


Financial Mathematics: Compound Interest
Example 1. ‘
Let An be the amount in a bank account after n years with annual interest rate r and
initial deposit A0 .
An+1 = (1 + r)An
The solution is An = (1 + r)n A0 .

Solution
Parameter Analysis:

• Growth factor: α = 1 + r

• Fixed point: A∗ = 0 (only fixed point)

• Stability: |α| = |1 + r| > 1 for r > 0 (the fixed point is unstable), it is < 1
for r < 0 (the fixed point is stable)

Biological Systems: Population Growth


Example 2. ‘
Let Mn be the number of cells after n generations, with each cell producing h daughter
cells:
Mn+1 = hMn
Solution: Mn = M0 hn .

Important
Dynamical Interpretation:

• |h| > 1: Exponential growth (unstable fixed point at 0)

• |h| < 1: Exponential decay (asymptotically stable fixed point at 0)

• h = 1: Constant population (neutrally stable)

5
Dynamical Systems A Lecture Notes

Numerical Methods: Finding Square Roots


Example√ 3.
To find a, iterate:  
1 a
xn+1 = xn +
2 xn

√ a nonlinear map that is quadratic with fixed points at x = ± a.
This is
Find 3.

Solution

Analysis for 3 with x0 = 1:

x0 = 1.000000
 
1 3
x1 = 1+ = 2.000000
2 1
 
1 3
x2 = 2+ = 1.750000
2 2
 
1 3
x3 = 1.75 + ≈ 1.732143
2 1.75
x4 ≈ 1.732051
x5 ≈ 1.732051

Rapid convergence to 3 ≈ 1.732051.

3 Classification of Difference Equations


Definition 4.

(i) Homogeneous linear first-order: xn+1 = αxn

(ii) Non-homogeneous linear first-order: xn+1 = αxn + β

(iii) Homogeneous linear second-order: axn+2 + bxn+1 + cxn = 0

(iv) Nonlinear: xn+1 = f (xn ) where f is nonlinear

Higher-Order Linear Equations


Example 4. ‘
For the second-order equation:

xn+2 + bxn+1 + cxn = 0

Find the general solution by assuming a solution of the form xn = λn .

6
Dynamical Systems A Lecture Notes

Solution
We begin by forming the Characteristic Equation Assume a solution of the form
xn = λn where λ ̸= 0.
Substituting into the recurrence relation:

λn+2 + bλn+1 + cλn = 0

Factor out λn (since λ ̸= 0 gives

λn (λ2 + bλ + c) = 0

This gives us the characteristic equation

λ2 + bλ + c = 0

Use the quadratic formula we obtain the following



−b ± b2 − 4c
λ=
2
Let λ1 and λ2 denote the two roots of the characteristic equation. The form of the
general solution depends on the discriminant ∆ = b2 − 4c. We have the following
cases;

• Case 1: If ∆ > 0 we have two distinct real roots λ1 ̸= λ2 . So we have

xn = Aλn1 + Bλn2 .

• Case 2: If ∆ = 0 the we have one repeated real root λ1 = λ2 , and thus

xn = (A + Bn)λn .

• Case 3: If ∆ < 0 we have complex conjugate roots givenpby λ1,2 = α ± βi.



So, expressing in polar form gives λ = ρe where ρ = α2 + β 2 and θ =
arctan(β/α). So,
xn = ρn (A cos(nθ) + B sin(nθ))

Example 5. ‘
The Fibonacci sequence is generated by:

Fn+2 = Fn+1 + Fn , F0 = 0, F1 = 1

Characteristic equation:
√ λ2 − λ − 1 = 0
1± 5
Roots: λ1,2 = 2 h √ n  √ n i
√1 1+ 5 1− 5
General solution: Fn = 5 2
− 2

7
Dynamical Systems A Lecture Notes

Solution
Binet’s Formula:

ϕn − (−ϕ)−n 1+ 5
Fn = √ , ϕ= ≈ 1.618034
5 2

This shows exponential growth with rate ϕ (the golden ratio).

4 Exercises and Problem Sets


Tutorial Problems
1. Transformation to Homogeneous Form:
β
(a) For xn+1 = αxn + β with fixed point x∗ = 1−α
, show that the transformation
yn = xn − x∗ yields yn+1 = αyn .
(b) Find the general solution for yn .
(c) Transform back to obtain the general solution for xn .

2. Graphical Analysis: For xn+1 = 12 xn + 1 with x0 = 2.5:

(a) Plot xn vs n for n = 0, 1, . . . , 10.


(b) Identify the fixed point.
(c) Comment on the stability and long-term behavior.

3. Verification of Solutions:

(a) Verify that Pn = −(n + 1) is a particular solution of xn+1 = 3xn + 2n + 1.


(b) Show that the general solution with x0 = 1 is xn = 2(3n ) − n − 1.

4. Financial Modeling: Consider a savings account with a monthly compounding


interest rate r%. Let an be the amount at the end of month n.
r

(a) Derive the recurrence relation: an+1 = 1 + 1200 an + p n .
(b) For constant monthly deposit pn = p, find the closed-form solution.
(c) Analyze the long-term behaviour for different values of r and p.

5. Second-Order Equations: Find the general solution of:

(a) xn+2 + 2xn+1 + xn = 0


(b) xn+2 + 2xn+1 + 5xn = 0

6. Classification Exercise: Classify the following as linear/nonlinear, homogeneous/non-


homogeneous, and determine their order:

(a) xn+1 = sin(xn−1 )


(b) xn+2 = 3xn+1 − xn
(c) xn+1 = xn−1 + 2

8
Dynamical Systems A Lecture Notes

(d) xn+3 = xn+1 + xn−3 + n2

7. Particular Solutions: Find particular solutions for:

(a) xn+1 = −xn + cos(2n) (Hint: Try pn = A cos n + B sin n)


n
(b) xn+1 = x2n + n 21
(c) xn+1 = xn + 1

8. Computational Exploration: Investigate the logistic map:

Pn+1 = Pn + αPn (1 − Pn )

for α = 0.1, 0.7, 1.2 with initial conditions P0 = 0, 0.001, 2.

(a) Find all equilibrium solutions.


(b) Plot time series for each combination of parameters and initial conditions.
(c) Analyze stability and bifurcation behavior.

Summary

Equation Type Solution


Homogeneous 1st order: xn+1 = αxn xn = α n x0
n
Non-homogeneous 1st order: xn+1 = αxn + β xn = αn x0 + β 1−α
1−α

Characteristic equation for 2nd order: λ2 + bλ + c = 0

Distinct real roots: λ1 , λ2 xn = Aλn1 + Bλn2

Repeated root: λ xn = (A + Bn)λn

Complex roots: r(cos θ ± i sin θ) xn = rn (A cos nθ + B sin nθ)

Stability Criteria
• For xn+1 = αxn : Fixed point at 0 is stable if |α| < 1

• For xn+1 = f (x): Fixed point x∗ is stable if |f ′ (x∗ )| < 1

• For second-order systems: Jury conditions determine stability from the character-
istic polynomial

9
APM08X1 – Dynamical Systems A
Week 2 Notes: Orbits in Discrete Dynamical Systems

Prepared for Students

1 Introduction to Orbits
In discrete dynamical systems, we study how a quantity evolves when a rule is applied
repeatedly. Rather than solving equations explicitly, we often investigate the long-term
behaviour of repeated updates. This naturally leads to the concept of an orbit.

Orbit of a Point
Let F : R → R be a function and let x0 ∈ R be an initial value. The orbit of x0
under F is the sequence

x0 , x1 = F (x0 ), x2 = F (x1 ) = F 2 (x0 ), ..., xn = F n (x0 ).

The set
O(x0 ) = {x0 , x1 , x2 , . . . }

is called the orbit (or trajectory) starting at x0 . The initial value x0 is called the
seed.

Interpretation. Each step applies the same rule again:

xn+1 = F (xn ).

Thus an orbit describes how the system evolves in time.


Remark.

• In discrete-time systems (maps), an orbit is a sequence.

• In continuous-time systems (flows), an orbit is a curve.

2 Computing Orbits
2.1 Basic Iteration

1
Example 1: Square-Root Map

Let F (x) = x and x0 = 144. Compute the first five points of the orbit.
Solution.

x1 = 144 = 12,

x2 = 12 ≈ 3.4641,

x3 = 3.4641 ≈ 1.8612,

x4 = 1.8612 ≈ 1.3643.

The orbit decreases toward 1. Since 1 = 1, the point x = 1 is a fixed point.

2.2 Iterating Functions Algebraically

Example 2: Iterating a Polynomial


Let F (x) = x2 − 2. Compute F 2 (x) and F 3 (x).
Solution.
First iterate:
F (x) = x2 − 2.

Second iterate:

F 2 (x) = F (F (x)) = (x2 − 2)2 − 2 = x4 − 4x2 + 2.

Third iterate:

F 3 (x) = F (F 2 (x)) = (x4 − 4x2 + 2)2 − 2.

Expanding:
F 3 (x) = x8 − 8x6 + 20x4 − 16x2 + 2.

This shows how algebraic complexity grows rapidly under iteration.

2.3 Closed Forms of Iterates


Example 3: Exponential Growth of Iterates
Let F (x) = x2 .
Then
F (x) = x2 , F 2 (x) = x4 , F 3 (x) = x8 .

Pattern:
n
F n (x) = x2 .

2
Thus iteration doubles the exponent each time.

2.4 Bounded Orbits


Example 4: Iterating Cosine
Let F (x) = cos x and x0 = 1.

x1 = cos(1) ≈ 0.5403,

x2 = cos(0.5403) ≈ 0.8576,

x3 = cos(0.8576) ≈ 0.6543,

x4 = cos(0.6543) ≈ 0.7935.

The orbit oscillates but approaches a fixed point near

x∗ ≈ 0.739.

Indeed, solving x = cos x gives this value.

3 Types of Orbits
3.1 Fixed Points
Fixed Point
A point x∗ is a fixed point of F if

F (x∗ ) = x∗ .

The orbit is constant:


x∗ , x∗ , x∗ , . . .

Example 5: Fixed Points of x3


Solve
x3 = x.

x3 − x = 0 ⇒ x(x2 − 1) = 0.

Hence
x = 0, x = 1, x = −1.

3
Example 6: Fixed Points of x sin x
Solve
x sin x = x.

x(sin x − 1) = 0.

Thus either
x=0

or
π
sin x = 1 ⇒ x = + 2kπ.
2

3.2 Geometric Interpretation


Fixed points occur where the graphs intersect:

y = F (x) and y = x.

Example 7: F (x) = x2
Solve
x2 = x ⇒ x(x − 1) = 0.

Hence fixed points are 0 and 1.

3.3 Periodic Orbits


Periodic Point
A point x is periodic of period k if

F k (x) = x.

If k is the smallest such integer, it has prime period k.

Example 8: Period-2 Map


Let F (x) = −x.

F (1) = −1, F (−1) = 1.

Thus {1, −1} is a period-2 orbit.

4
Example 9: Period-2 Cycle
Let F (x) = x2 − 1.
Orbit of 0:
0, −1, 0, −1, . . .

Thus 0 and −1 form a 2-cycle.

Example 10: Period-3 Cycle


Let
3 5
F (x) = − x2 + x + 1.
2 2
Compute:
F (0) = 1, F (1) = 2, F (2) = 0.

Hence
0→1→2→0

is a period-3 orbit.

4 Existence of Fixed Points


Fixed Point Theorem
Let J = [a, b] and let f : J → J be continuous. Then f has at least one fixed
point in J.

Idea. If the graph of f lies inside the square [a, b] × [a, b], it must cross the diagonal
y = x.
Define
g(x) = f (x) − x.

Since f maps [a, b] into itself:

f (a) ≥ a, f (b) ≤ b.

Thus
g(a) ≥ 0, g(b) ≤ 0.

If either equals zero we have a fixed point. Otherwise signs differ, so by the Interme-
diate Value Theorem there exists c ∈ (a, b) such that

g(c) = 0 ⇒ f (c) = c.

5
Example 11
Show f (x) = 1 − x2 has a fixed point in [0, 1].
Solve:
1 − x2 = x.

x2 + x − 1 = 0.


−1 + 5
x= ≈ 0.618.
2
Since this lies in [0, 1], a fixed point exists.

5 Tutorial 2
Basic Iteration

1. Compute the first five iterates of x0 = 0.5 for f (x) = x2 .

2. Compute four iterates of 0 under f (x) = x2 + 1.



3. Compute the orbit of 2 under f (x) = x + 2.

Algebraic Iterates

4. Compute g 2 (x) and g 3 (x) for g(x) = sin x.

5. Find a general expression for F n (x) if F (x) = x3 .

Fixed Points

6. Find fixed points of f (x) = x2 − 2.

7. Find fixed points of f (x) = x3 − 3x.

8. Show f (x) = cos x has a fixed point in [0, 1].

Periodic Orbits

9. Show f (x) = 1 − x2 has a period-2 orbit.

10. Determine whether x = 1 is periodic for f (x) = x2 − 1.

11. Construct a function f : [0, 1] → [0, 1] with exactly three fixed points.

6
Conceptual

12. Explain the difference between fixed and periodic points.

13. Explain why iterates of x2 grow rapidly.

14. Give an example of a bounded orbit and an unbounded orbit.

7
Week 3 Notes:
Graphical Analysis of Discrete Dynamical Systems
APM08X1 – Dynamical Systems A

1 Graphical Analysis of Discrete Dynamical Systems


In this section we study one-dimensional discrete dynamical systems of the form

xn+1 = F (xn ),

where F : R → R is a nonlinear function. Our main goal is to understand the long-term


behaviour of solutions without explicitly computing formulas for xn .
A powerful graphical technique for this purpose is called cobwebbing. Cobweb dia-
grams allow us to visualise how iterates of a function evolve and whether they converge
to steady states, oscillate, or diverge.

Cobwebbing Process
Given a function F (x) and an initial value (seed) x0 , the cobwebbing procedure
is as follows:

1. Plot the graphs of y = F (x) and the diagonal line y = x on the same axes.
Fixed points occur at intersections where F (x) = x.

2. Start at the point (x0 , x0 ) on the line y = x.

3. Draw a vertical line from (x0 , x0 ) to the curve y = F (x), reaching (x0 , F (x0 )).

4. Draw a horizontal line from (x0 , F (x0 )) to the line y = x, reaching (F (x0 ), F (x0 )).
This point corresponds to x1 = F (x0 ).

5. Repeat this vertical–horizontal process to generate

x0 , x1 = F (x0 ), x2 = F (x1 ), x3 = F (x2 ), . . .

The resulting “cobweb” visually represents the orbit of the system.

1
Illustrative Example
Consider the function
F (x) = sin x,

on the interval 0 ≤ x ≤ π. The graphs of y = sin x and y = x intersect at x = 0,


indicating a fixed point. Starting from a seed such as x0 = 0.5, the cobweb diagram
shows that successive iterates approach the fixed point at the origin.
This illustrates how graphical methods can immediately reveal stability properties.

Stability Criterion for Fixed Points


Let x∗ be a fixed point of F , i.e., F (x∗ ) = x∗ . Suppose F is differentiable at x∗ .
Then:

• If |F ′ (x∗ )| < 1, then x∗ is asymptotically stable (attracting).

• If |F ′ (x∗ )| > 1, then x∗ is unstable (repelling).

• If |F ′ (x∗ )| = 1, the test is inconclusive; further analysis is required.

When F ′ (x∗ ) < 0 and |F ′ (x∗ )| < 1, orbits converge to x∗ via oscillatory (alter-
nating) approach. When F ′ (x∗ ) > 0 and |F ′ (x∗ )| < 1, convergence is monotone.

2 Assignment: Plotting in Mathematica


The following Mathematica commands may be used to plot y = sin x and y = x on the
same axes.

h = Plot[Sin[x], {x, 0, Pi}, Axes -> True, AxesOrigin -> {0, 0}]

g = Plot[x, {x, 0, Pi}, Axes -> True, AxesOrigin -> {0, 0}]

Show[h, g]

You should observe that the two curves intersect at the origin, confirming the fixed
point.

Exercises
Use Mathematica to plot the following pairs of functions and interpret their intersections
and stability:

1. y = rx(1 − x) and y = x, for r > 1, 0 ≤ x ≤ 1

2
2. y = rx(1 − x) and y = x, for r < 1, 0 ≤ x ≤ 1

3. y = 2x and y = x, 0 ≤ x < 5

4. y = 0.5x and y = x, 0 ≤ x ≤ 2
x
5. y = − and y = x, 0 ≤ x ≤ 2
1+x

Example 1: Stability of Fixed Points


Consider F (x) = x3 .
The fixed points satisfy x3 = x, giving

x∗ = −1, 0, 1.

Cobweb diagrams show that:

• If |x0 | < 1, then xn → 0 as n → ∞

• If |x0 | > 1, the orbit diverges

Thus, x = 0 is attracting, while x = ±1 are unstable.

Example 2: Linear Maps


Consider F (x) = 12 x+ 23 . There is a unique fixed point, and since the slope satisfies
|F ′ (x)| < 1, all orbits converge to this fixed point.
Similarly, for F (x) = − 12 x + 92 , oscillatory convergence occurs due to the negative
slope.

Example 3 (Harder ): The Logistic Map and Period Doubling


Consider the logistic map

F (x) = rx(1 − x), x ∈ [0, 1],

with parameter r > 0.


Step 1: Find fixed points. Setting F (x∗ ) = x∗ gives

rx∗ (1 − x∗ ) = x∗ =⇒ x∗ r(1 − x∗ ) − 1 = 0.


Thus x∗0 = 0 and x∗1 = 1 − 1r (valid for r > 1).


Step 2: Stability of fixed points. We have F ′ (x) = r(1 − 2x).

• At x∗0 = 0: F ′ (0) = r. So x∗0 is stable iff r < 1.

3
• At x∗1 = 1 − 1r : F ′ 1 − 1r 1
= 2 − r. So x∗1 is stable iff
 
= r 1−2 1− r
|2 − r| < 1, i.e., 1 < r < 3.

Step 3: Period-2 orbits. For r > 3, the fixed point x∗1 loses stability. A period-
2 orbit {p, q} bifurcates, satisfying F (F (p)) = p. Solving F 2 (x) = x and discarding
the fixed points yields:
p
(r + 1) ± (r + 1)(r − 3)
p, q = .
2r

This period-2 orbit is stable for 3 < r < 1 + 6 ≈ 3.449.
Conclusion: As r increases from 1 to 4, the logistic map undergoes a period-
doubling cascade: fixed point → period-2 → period-4 → · · · → chaos. This is a
hallmark of complex dynamical behaviour arising from a simple quadratic rule.

Example 4 (Harder ): The Tent Map and Sensitive Dependence


Define the tent map on [0, 1] by

2x, 0 ≤ x ≤ 12 ,
T (x) =
1
2(1 − x), < x ≤ 1.

2

Fixed points. Solving T (x) = x: on [0, 12 ], 2x = x ⇒ x∗ = 0; on ( 21 , 1],


2 − 2x = x ⇒ x∗ = 23 .
Stability. Since |T ′ (x)| = 2 > 1 everywhere (where the derivative exists), both
fixed points are unstable.
Sensitive dependence. Suppose two initial conditions differ by δ0 ≪ 1. After
n iterations, the separation grows as |δn | ≈ 2n |δ0 |, as long as orbits remain in [0, 1].
This exponential divergence of nearby orbits is the defining hallmark of chaos.
Conjugacy with logistic map. The tent map is topologically conjugate to the
logistic map at r = 4 via the change of variables x = sin2 (πθ/2). This means their
dynamical properties — including chaos — are qualitatively identical.
Density of periodic orbits. One can show that the periodic points of T are
exactly the rationals in [0, 1], which are dense. Yet the set of non-periodic orbits has
full measure. This coexistence of dense periodic orbits with generic chaos is a subtle
and important feature.

Example 5 (Harder ): Schwarzian Derivative and Global Stability

4
For maps on an interval, the Schwarzian derivative of F is defined as
2
F ′′′ (x) 3 F ′′ (x)

SF (x) = ′ − ,
F (x) 2 F ′ (x)

wherever F ′ (x) ̸= 0.

Singer’s Theorem (1978). Let F : [a, b] → [a, b] be a C 3 map with SF (x) < 0
for all x where F ′ (x) ̸= 0. Then each stable periodic orbit of F attracts at least
one critical point of F .

Consequence for the logistic map. For F (x) = rx(1 − x), we compute:

F ′ (x) = r(1 − 2x), F ′′ (x) = −2r, F ′′′ (x) = 0.

Then
3 4r2 −6
SF (x) = 0 − · 2 2
= 2
< 0 for all x ̸= 21 .
2 r (1 − 2x) (1 − 2x)
Since F has a unique critical point at xc = 12 , Singer’s theorem implies that the
logistic map has at most one stable periodic orbit for any given r. This pow-
erful constraint, invisible from cobweb diagrams alone, severely restricts the possible
attracting behaviour.

3 Phase Portraits
A phase portrait is a global graphical representation of a discrete dynamical system on
the real line. It shows:

• Fixed points

• Direction of motion

• Stability information

Phase portraits provide a “zoomed-out” view compared to cobweb diagrams.

Example: F (x) = x3
The fixed points are x∗ = 0, ±1.

• Arrows point towards 0 for −1 < x < 1

• Arrows point away from ±1

5
This confirms that 0 is attracting, while ±1 are repelling.

Example: F (x) = x2
Fixed points are x∗ = 0 and x∗ = 1.

• If 0 < x < 1, then xn → 0

• If x > 1, then xn → ∞

• If −1 < x < 0, the orbit enters (0, 1) and then converges to 0

• If x < −1, the orbit diverges

Harder Example: Saddle-Node Bifurcation in F (x) = x2 + c


Consider Fc (x) = x2 + c with parameter c ∈ R.
Fixed points. Solving x2 + c = x gives x2 − x + c = 0, so

1± 1 − 4c
x∗ = .
2

Fixed points exist (as real numbers) iff c ≤ 14 .


Stability. Fc′ (x) = 2x. At each fixed point x∗ :

1±1 − 4c √
Fc′ (x∗ ) =2· = 1 ± 1 − 4c.
2

• Lower fixed point (− branch): F ′ = 1 − 1 − 4c. For c < 1
4
, this satisfies
|F ′ | < 1 when 0 < c < 14 ; the fixed point is stable.

• Upper fixed point (+ branch): F ′ = 1 + 1 − 4c > 1 always. This fixed point
is always unstable.

Bifurcation. At c = 14 , the two fixed points collide and annihilate at x∗ = 12 .


For c > 14 , no fixed points exist and all orbits escape to −∞. This is a saddle-node
bifurcation (also called a fold bifurcation) — a fundamental mechanism by which
fixed points are created or destroyed as a parameter varies.
Phase portrait interpretation. Draw Fc (x) = x2 + c and y = x for three cases:

• c < 14 : two intersections, one stable and one unstable

• c = 14 : tangency (one semi-stable fixed point)

• c > 14 : no intersection, all orbits diverge

6
1
The cobweb at c = 4
shows a characteristic “staircase” slowdown near the tangency.

4 Tutorial 3
1. Use colour coding to describe orbit behaviour for:

F (x) = 3x, F (x) = 12 x, F (x) = 3x + 1, F (x) = x − x2 , F (x) = sin x

2. Sketch phase portraits for each function.

3. Perform a full orbit analysis for:

F (x) = 12 x − 3, F (x) = |x|, F (x) = −x5 , F (x) = 1


x

4. For w(x) = rx(1 − x):

• Compute iterates for r = 2, 3, 4


• Compare short-term and long-term behaviour

5. Draw phase portraits for:

f (x) = x3 , f (x) = x1 , f (x) = cos x

6. For the logistic map F (x) = rx(1 − x), find all period-2 orbits analytically by
solving F (F (x)) = x. For which values of r are these period-2 orbits stable?
Verify your answer using the derivative criterion applied to G(x) = F (F (x)) at
the period-2 points.

7. Consider F (x) = λ sin(πx) on [0, 1] for λ ∈ (0, 1].

• Find all fixed points.


• Use the stability criterion to determine for which λ the nonzero fixed point
is stable.
• Compute SF (x) (the Schwarzian derivative) and confirm its sign.
• Use Singer’s theorem to deduce an upper bound on the number of stable
periodic orbits.

8. For the map Fc (x) = x2 + c:

7
• Identify all values of c at which a bifurcation occurs among fixed points
and period-2 orbits.
• Sketch a bifurcation diagram: plot the location and stability of fixed points
and period-2 orbits as a function of c ∈ [−2, 14 ].
• Determine the value of c at which the period-2 orbit itself loses stability.

8
Week 4 Notes: Fixed Points and Their
Stability
Dynamical Systems A - APM08X1

What to expect: This week we explore the stability analysis of fixed points in discrete
dynamical systems, with applications to population models and the logistic map. We
derive stability conditions, explore bifurcations, and introduce cobweb diagrams.

1 Introduction to Fixed Points and Stability


Definition

Fixed Point (Equilibrium Point): A point x∗ is called a fixed point or steady


state of the map f (x) if it satisfies:

f (x∗ ) = x∗

At a fixed point, the system remains unchanged under iteration.

1.1 Intuitive Understanding of Stability


Example

Physical Analogy - Ball on a Landscape:


Consider a ball placed on three different positions on a hilly landscape:

1. Point A (Valley bottom): If you slightly displace the ball, it returns to point
A. This represents a stable fixed point (attractor).

2. Point B (Hilltop): If you slightly displace the ball, it rolls away from point
B. This represents an unstable fixed point (repeller).

3. Point C (Flat plateau): The ball stays wherever placed. This represents a
neutral or indifferent fixed point.
Week 4: Fixed Points and Stability Dynamical Systems A

B (Unstable)

A (Stable) C (Neutral)

1.2 Mathematical Formulation of Stability


Given a discrete dynamical system:

xn+1 = f (xn ), n = 0, 1, 2, . . .

with a fixed point x∗ satisfying f (x∗ ) = x∗ .


The key question: If we start near x∗ (i.e., x0 = x∗ + ϵ with small ϵ), will the orbit
converge to x∗ or diverge away from it?

Example

Example: f (x) = x2
Fixed points: x∗ = 0 and x∗ = 1 (since 02 = 0 and 12 = 1).
Behavior analysis:

• For 0 < x0 < 1: xn → 0 as n → ∞

• For x0 > 1: xn → ∞ as n → ∞

• For x0 = 1: xn = 1 for all n (but this is unstable to perturbations)

Conclusion: x∗ = 0 is stable (attracting), while x∗ = 1 is unstable (repelling).

2 Linear Stability Analysis


2.1 Derivation of Stability Condition
Let xn = x∗ + ηn , where ηn is a small perturbation from the fixed point. Then:

xn+1 = f (xn ) = f (x∗ + ηn )


Using Taylor expansion around x∗ :

f ′′ (x∗ ) 2
f (x∗ + ηn ) = f (x∗ ) + f ′ (x∗ )ηn + ηn + · · ·
2!
Since f (x∗ ) = x∗ and neglecting higher-order terms for small ηn :

x∗ + ηn+1 ≈ x∗ + f ′ (x∗ )ηn


Thus:

ηn+1 ≈ f ′ (x∗ )ηn = ληn


where λ = f ′ (x∗ ) is called the multiplier or eigenvalue.
The solution to this linearized equation is:

2
Week 4: Fixed Points and Stability Dynamical Systems A

ηn = λn η0
Therefore:

• If |λ| < 1, then ηn → 0: perturbations decay ⇒ stable

• If |λ| > 1, then |ηn | → ∞: perturbations grow ⇒ unstable

• If |λ| = 1, linear analysis is inconclusive: neutral/marginally stable

Theorem
Linear Stability Theorem for 1D Maps:
Let x∗ be a fixed point of f (x) (i.e., f (x∗ ) = x∗ ). Then:

1. x∗ is locally asymptotically stable if |f ′ (x∗ )| < 1

2. x∗ is unstable if |f ′ (x∗ )| > 1

3. If |f ′ (x∗ )| = 1, the fixed point is non-hyperbolic and linear analysis is incon-


clusive; higher-order terms must be considered.

Preliminary Concepts
Part 1: Setup and Notation
Consider the discrete dynamical system:

xn+1 = f (xn ), n = 0, 1, 2, . . .

where f : R → R is a differentiable function, and x∗ is a fixed point satisfying f (x∗ ) = x∗ .


Step 1: Introduce Perturbation Variable
Define the deviation or perturbation from the fixed point:

η n = xn − x∗

This represents how far the system is from equilibrium at time step n.
Our goal is to understand how ηn evolves over time. Specifically:

• If |ηn | → 0 as n → ∞, the fixed point is stable

• If |ηn | grows without bound, the fixed point is unstable

Step 2: Express the Dynamics in Terms of ηn


Starting from xn+1 = f (xn ) and substituting xn = x∗ + ηn :

x∗ + ηn+1 = f (x∗ + ηn )

Since f (x∗ ) = x∗ , we can write:

ηn+1 = f (x∗ + ηn ) − f (x∗ )

3
Week 4: Fixed Points and Stability Dynamical Systems A

Step 3: Linearization Using Taylor Expansion


Assume ηn is small (we’re considering local stability near x∗ ). Expand f (x∗ + ηn ) in
a Taylor series around x∗ :

∗ ∗ ′ ∗ f ′′ (x∗ ) 2 f ′′′ (x∗ ) 3


f (x + ηn ) = f (x ) + f (x )ηn + ηn + ηn + · · ·
2! 3!
Substituting f (x∗ ) = x∗ :
f ′′ (x∗ ) 2 f ′′′ (x∗ ) 3
f (x∗ + ηn ) = x∗ + f ′ (x∗ )ηn + ηn + ηn + · · ·
2 6
Therefore:
f ′′ (x∗ ) 2
 
∗ ′ ∗
ηn+1 = x + f (x )ηn + ηn + · · · − x∗
2
f ′′ (x∗ ) 2 f ′′′ (x∗ ) 3
ηn+1 = f ′ (x∗ )ηn + ηn + ηn + · · ·
2 6
Step 4: Linear Approximation (Key Simplification)
For sufficiently small ηn , the higher-order terms (ηn2 , ηn3 , . . . ) are much smaller than
the linear term f ′ (x∗ )ηn .
More precisely, if ηn is small enough that:
f ′′ (x∗ ) f ′′′ (x∗ ) 2
ηn + ηn + · · · ≪ |f ′ (x∗ )|
2 6
then we can make the linear approximation:

ηn+1 ≈ f ′ (x∗ )ηn

Let λ = f ′ (x∗ ) (the multiplier or eigenvalue at the fixed point). Then:

ηn+1 ≈ ληn

Part 2: Analysis of the Linearized System


Step 5: Iteration of the Linearized Equation
The linearized equation ηn+1 = ληn is a simple geometric progression:

η1 = λη0
η2 = λη1 = λ2 η0
η3 = λη2 = λ3 η0
..
.
η n = λn η 0

Thus, the behavior of ηn depends entirely on λ = f ′ (x∗ ).


Step 6: Stability Analysis Based on |λ|
The solution ηn = λn η0 leads to three distinct cases:

Case 1: |λ| < 1 (Stable Case):

• Since |λ| < 1, we have |λ|n → 0 as n → ∞

4
Week 4: Fixed Points and Stability Dynamical Systems A

• Therefore, |ηn | = |λ|n |η0 | → 0 as n → ∞


• This means the perturbation decays exponentially: |ηn | ≤ C|λ|n for some
constant C
• The fixed point is locally asymptotically stable

Case 2: |λ| > 1 (Unstable Case):

• Since |λ| > 1, we have |λ|n → ∞ as n → ∞


• Therefore, |ηn | = |λ|n |η0 | → ∞ as n → ∞ (unless η0 = 0)
• The perturbation grows exponentially, no matter how small the initial η0
• The fixed point is unstable

Case 3: |λ| = 1 (Borderline Case):

• When |λ| = 1, the linear term alone predicts |ηn | = |η0 | for all n
• However, this is misleading because higher-order terms (ηn2 , ηn3 , . . . ) become
important
• These higher-order terms can determine stability or instability
• The fixed point is non-hyperbolic or critically stable
• Linear analysis is inconclusive; we need nonlinear analysis

Part 3: Justifying the Linear Approximation


Step 7: Formal Justification (for |λ| = ̸ 1)
For a rigorous proof, we need to show that the linear approximation accurately predicts
stability/instability when |λ| = ̸ 1.
Lemma (Local Linearization Lemma): If f is continuously differentiable at x∗
with |f ′ (x∗ )| =
̸ 1, then there exists a neighborhood U of x∗ such that:

|f (x) − x∗ | ≤ ρ|x − x∗ | for all x ∈ U

where ρ is a constant satisfying:

• If |f ′ (x∗ )| < 1, then ρ < 1 (contraction)

• If |f ′ (x∗ )| > 1, then ρ > 1 (expansion)

Proof of Lemma: Since f is continuously differentiable, for any ϵ > 0, there exists
δ > 0 such that for |x − x∗ | < δ:

f (x) − f (x∗ )
− f ′ (x∗ ) < ϵ
x − x∗

This implies:
|f (x) − x∗ | ≤ (|f ′ (x∗ )| + ϵ)|x − x∗ |
Choose ϵ small enough so that |f ′ (x∗ )| + ϵ has the desired property.
Step 8: Contraction Mapping Argument (for |λ| < 1)

5
Week 4: Fixed Points and Stability Dynamical Systems A

When |f ′ (x∗ )| < 1, choose ϵ > 0 such that ρ = |f ′ (x∗ )| + ϵ < 1. Then for x ∈
(x∗ − δ, x∗ + δ):
|f (x) − x∗ | ≤ ρ|x − x∗ |
This shows f is a contraction mapping near x∗ . By the Contraction Mapping
Theorem:
1. Starting from any x0 in this neighborhood, the sequence xn converges to x∗
2. The convergence is at least linear: |xn − x∗ | ≤ ρn |x0 − x∗ |
This proves local asymptotic stability.
Step 9: Instability Proof (for |λ| > 1)
When |f ′ (x∗ )| > 1, consider the inverse map (if it exists locally). Alternatively, use a
direct argument:
Choose ϵ > 0 such that ρ = |f ′ (x∗ )| − ϵ > 1. For sufficiently small ηn :

|ηn+1 | = |f (x∗ + ηn ) − x∗ | ≥ ρ|ηn |

By induction:
|ηn | ≥ ρn |η0 |
Since ρ > 1, |ηn | → ∞ as n → ∞ for any η0 ̸= 0.
This proves instability.

Part 4: The Borderline Case |λ| = 1


Step 10: Why Linear Analysis Fails When |λ| = 1
When |f ′ (x∗ )| = 1, the linear term gives:

ηn+1 ≈ ±ηn

This suggests the perturbation neither grows nor decays at linear order. However, the
higher-order terms can dominate:

• If f (x) = x∗ + (x − x∗ ) + a(x − x∗ )2 + · · · , then:

ηn+1 = ηn + aηn2 + · · ·

The quadratic term can cause growth or decay depending on the sign of a.
• If f (x) = x∗ + (x − x∗ ) + b(x − x∗ )3 + · · · , then:

ηn+1 = ηn + bηn3 + · · ·

The cubic term determines stability (superstable if b < 0, unstable if b > 0).

Examples of borderline cases:


1. f (x) = x + x2 at x∗ = 0: f ′ (0) = 1, but unstable (positive quadratic term)
2. f (x) = x − x3 at x∗ = 0: f ′ (0) = 1, but stable (negative cubic term)
3. f (x) = x + x3 at x∗ = 0: f ′ (0) = 1, but unstable (positive cubic term)
This demonstrates why we must consider nonlinear terms when |f ′ (x∗ )| = 1.

6
Week 4: Fixed Points and Stability Dynamical Systems A

Part 5: Geometric Interpretation


Step 11: Visual Understanding
Consider the graph of y = f (x) and the line y = x:
y
Geometric Interpretation:
Slope at fixed point determines stability y=x
Slope = 1<(borderline)
Slope 1 (stable) x∗
Slope > 1 (unstable)

The slope f ′ (x∗ ) determines how f (x) crosses the line y = x:


• If 0 < f ′ (x∗ ) < 1: f (x) crosses y = x from above to below (stable)
• If −1 < f ′ (x∗ ) < 0: Oscillatory convergence (stable)
• If f ′ (x∗ ) > 1: f (x) crosses y = x from below to above (unstable)
• If f ′ (x∗ ) < −1: Oscillatory divergence (unstable)
• If f ′ (x∗ ) = ±1: Tangent crossing (borderline case)

Part 6: Summary and Conclusion


Step 12: Theorem Summary
We have proven the Linear Stability Theorem through the following logical steps:

1. Defined perturbation ηn = xn − x∗ from fixed point x∗


2. Used Taylor expansion to linearize: ηn+1 ≈ f ′ (x∗ )ηn
3. Solved linearized equation: ηn = [f ′ (x∗ )]n η0
4. Analyzed three cases:
• |f ′ (x∗ )| < 1: ηn → 0 (stable)
• |f ′ (x∗ )| > 1: |ηn | → ∞ (unstable)
• |f ′ (x∗ )| = 1: Linear term inconclusive
5. Provided rigorous justification using contraction mapping
6. Explained why nonlinear terms matter in borderline case
7. Gave geometric interpretation

The theorem is fundamental because it reduces stability analysis to computing a single


derivative at the fixed point, making it computationally tractable for many applications.

7
Week 4: Fixed Points and Stability Dynamical Systems A

Applications and Examples


Example

Example 1: Logistic Map f (x) = rx(1 − x) at x∗ = 0

• f ′ (x) = r(1 − 2x)

• f ′ (0) = r

• Stability condition: |r| < 1 (since r > 0, this means 0 < r < 1)

• For r > 1, x∗ = 0 is unstable

• At r = 1, we need nonlinear analysis (transcritical bifurcation)

Example

Example 2: Sine Map f (x) = sin(x) at x∗ = 0

• f ′ (x) = cos(x)

• f ′ (0) = 1 (borderline case!)

• Linear analysis inconclusive


x3
• Need nonlinear terms: sin(x) = x − 6
+ ···

• Since the cubic term is negative, x∗ = 0 is actually stable (though not linearly
determined)

Example

Example 3: Superstable Fixed Point f (x) = x2 at x∗ = 0

• f ′ (x) = 2x

• f ′ (0) = 0 (less than 1 in absolute value)

• This is not just stable, but superstable (fastest possible convergence)

• Nearby points converge to 0 very quickly

Important Notes and Limitations


• Local Nature: The theorem only guarantees local stability. A fixed point can be
locally stable but globally unstable.

• Continuous Differentiability: The proof assumes f is continuously differentiable


near x∗ .

8
Week 4: Fixed Points and Stability Dynamical Systems A

• Multi-dimensional Analogue: In higher dimensions, replace |f ′ (x∗ )| with the


spectral radius of the Jacobian matrix.
• Practical Application: Always check |f ′ (x∗ )| first. If it equals 1, proceed to
nonlinear analysis (center manifold theory, normal forms).
• Historical Context: This theorem is a discrete analogue of the Hartman-Grobman
theorem for continuous systems.

Key Takeaway: The derivative at a fixed point serves as a linear stability indicator
|f ′ (x∗ )| < 1: Green light (stable) |f ′ (x∗ )| > 1: Red light (unstable) |f ′ (x∗ )| = 1:
Yellow light (proceed with caution)

2.2 Classification of Fixed Points


Based on the value of λ = f ′ (x∗ ):
• Stable node: 0 < λ < 1
• Stable alternating: −1 < λ < 0 (oscillatory convergence)
• Unstable node: λ > 1
• Unstable alternating: λ < −1 (oscillatory divergence)
• Neutral: λ = ±1 (requires nonlinear analysis)
• Superstable: λ = 0 (fastest convergence)

3 Examples and Applications


3.1 Example 1: Michaelis-Menten Type Population Model
Example
kx
Model: f (x) = b+x , where b, k > 0 (population growth with saturation)
Step 1: Find fixed points
kx
x= ⇒ x(b + x) = kx ⇒ x(b + x − k) = 0
b+x
So: x∗1 = 0 or x∗2 = k − b
Note: x∗2 > 0 only if k > b.
Step 2: Compute derivative
bk
f ′ (x) = (using quotient rule)
(b + x)2

Step 3: Evaluate stability

9
Week 4: Fixed Points and Stability Dynamical Systems A

1. At x∗1 = 0: f ′ (0) = k
b

• Stable if k
b
<1⇒k<b
• Unstable if k > b

2. At x∗2 = k − b (exists only if k > b):

bk bk b
f ′ (k − b) = 2
= 2 =
(b + (k − b)) k k

• Stable if b
k
< 1 ⇒ k > b (always true when it exists)
• So x∗2 is stable whenever it exists

Biological Interpretation:

• When growth rate k is less than saturation constant b, population goes extinct
(x∗ = 0 is stable)

• When k > b, population stabilizes at x∗ = k − b

• This is a classic transcritical bifurcation at k = b

3.2 Example 2: Logistic Map


Example

Model: f (x) = rx(1 − x), with r > 0 (classic logistic map)


Step 1: Find fixed points

x = rx(1 − x) ⇒ x[r(1 − x) − 1] = 0

So: x∗1 = 0 or x∗2 = 1 − 1r


Note: x∗2 ∈ (0, 1) only if r > 1.
Step 2: Compute derivative

f ′ (x) = r(1 − 2x)

Step 3: Evaluate stability

1. At x∗1 = 0: f ′ (0) = r

• Stable if |r| < 1 ⇒ 0 < r < 1


• Unstable if r > 1

2. At x∗2 = 1 − 1
r
(exists only if r > 1):
      
′ 1 1 2
f 1− =r 1−2 1− =r 1−2+ =2−r
r r r

• Stable if |2 − r| < 1 ⇒ 1 < r < 3

10
Week 4: Fixed Points and Stability Dynamical Systems A

• Unstable if r > 3 or r < 1 (but r < 1 means x∗2 doesn’t exist)

Summary of stability for logistic map:

• 0 < r < 1: Only x∗1 = 0 exists and is stable

• 1 < r < 3: x∗1 = 0 is unstable, x∗2 = 1 − 1


r
is stable

• r > 3: Both fixed points are unstable

Bifurcation points:

• r = 1: Transcritical bifurcation (exchange of stability)

• r = 3: Period-doubling bifurcation (stable fixed point becomes unstable, giving


birth to a period-2 cycle)

4 Period-2 Cycles and Higher Periods


When r > 3 in the logistic map, neither fixed point is stable. What happens?
Definition

Period-k Cycle: A set of points {x1 , x2 , . . . , xk } forms a period-k cycle if:

f (x1 ) = x2 , f (x2 ) = x3 , ..., f (xk−1 ) = xk , f (xk ) = x1

Each xi is a fixed point of the k-th iterate f k (x) = f (f (. . . f (x) . . . )) (k times).

4.1 Finding Period-2 Cycles for Logistic Map


For f (x) = rx(1 − x), we find period-2 points by solving:
f 2 (x) = f (f (x)) = x
but excluding fixed points (which satisfy f (x) = x).
Example

Derivation of period-2 points:

f 2 (x) = f (f (x)) = f (rx(1 − x)) = r[rx(1 − x)][1 − rx(1 − x)]

= r2 x(1 − x)[1 − rx(1 − x)]


Solve f 2 (x) = x:
x = r2 x(1 − x)[1 − rx(1 − x)]
One solution is x = 0 (fixed point). Factor out x:

0 = x r2 (1 − x)[1 − rx(1 − x)] − 1




1
The bracketed term gives a cubic equation. Since we know x = 1 − r
is also a

11
Week 4: Fixed Points and Stability Dynamical Systems A

fixed point, we can factor (x − (1 − 1r )) out. After factorization, we get a quadratic:


 
2 1 1+r
x − 1+ x+ 2 =0
r r

The solutions are:


 
  s 2
1 1 1 1+r
x1,2 =  1+ ± 1+ −4 2 
2 r r r

p
r+1± (r + 1)(r − 3)
=
2r
These are real when (r + 1)(r − 3) ≥ 0. Since r > 0, this means r ≥ 3.
At r = 3, the period-2 cycle is born (bifurcation point), with x1,2 = 23 (initially
equal to the fixed point).
For r > 3, we have two distinct points forming a period-2 cycle.

4.2 Stability of Period-2 Cycles


The stability of a period-k cycle {x1 , x2 , . . . , xk } is determined by the product of deriva-
tives:

λ = f ′ (x1 ) · f ′ (x2 ) · · · · · f ′ (xk )


The cycle is stable if |λ| < 1. √
For the logistic map’s period-2 cycle when 3 < r < 1 + 6 ≈ 3.449, the cycle is stable.
At r ≈ 3.449, another period-doubling occurs to a period-4 cycle, beginning the famous
period-doubling route to chaos.

5 Cobweb Diagrams: Graphical Analysis


Cobweb diagrams provide a geometric way to visualize iterations of a map.

Definition
Constructing a cobweb diagram:

1. Plot y = f (x) and y = x on the same axes

2. Start at x0 on the x-axis

3. Draw vertical line to y = f (x) to get (x0 , f (x0 ))

4. Draw horizontal line to y = x to get (f (x0 ), f (x0 ))

5. Draw vertical line to y = f (x) to get (f (x0 ), f (f (x0 )))

6. Repeat

12
Week 4: Fixed Points and Stability Dynamical Systems A

Example

Cobweb examples for logistic map with different r:

1. r = 2.8 (stable fixed point): Cobweb converges to x∗ = 1 − 1


2.8
≈ 0.643

2. r = 3.3 (stable period-2): Cobweb oscillates between two values

3. r = 3.55 (period-4 or chaotic): More complex behavior

4. r = 3.829 (period-3 window in chaos): Has a stable period-3 cycle

6 Tutorial Problems
Tutorial Problems
1. Classification and Fixed Points: For each equation, indicate whether it is
linear or nonlinear, and find all fixed points.

(a) f (x) = (1 − x)(1 − βx), β constant


x
(b) g(x) = 1+x
αx
(c) h(x) = xe , α constant
(d) f (x) = x[2x2 − 2x + 1]1/2 + α, α constant
k
(e) F (x) = k1 +k2 /x
, k, k1 , k2 constants

2. Stability Analysis: Determine when the following steady states are stable:

(a) f (x) = mx(1 − x), x∗ = 0



(b) g(x) = −x2 (1 − x), x∗ = (1 + 5)/2

1
(c) h(x) = 2+x , x∗ = 2 − 1
(d) j(x) = x ln x2 , x∗ = ek/2

For each, sketch the function and use cobwebbing with appropriate x0 to show
solution behavior.

3. Population Model (Ricker Model): A population model is given by:

f (N ) = αN e−βN , α, β > 0

(a) Sketch f (N ) vs N
(b) Interpret parameters α and β biologically
(c) Show that N ∗ = ln α
β
is a fixed point
(d) Show it is stable provided |1 − ln α| < 1
(e) Find the range of α for stability

4. Generalized Growth Model: Plot f (x) = λx(1 + ax)−b for λ, a, b > 0

13
Week 4: Fixed Points and Stability Dynamical Systems A

(a) How do a and b affect the graph shape?


(b) Find all fixed points
(c) Analyze their stability

5. Exploring the Logistic Map: Using a computer or graphing calculator:

(a) For r = 2.5, start with x0 = 0.1 and iterate 20 times. What do you
observe?
(b) Repeat for r = 3.2, r = 3.5, and r = 3.8
(c) For r = 3.83, find the period-3 cycle
(d) Create a bifurcation diagram for 2.5 ≤ r ≤ 4

7 Additional Examples and Applications


7.1 Example 3: A Model with Multiple Fixed Points
Example

Model: f (x) = x3 − 3x2 + 3x


Fixed points: Solve x = x3 − 3x2 + 3x

0 = x3 − 3x2 + 2x = x(x − 1)(x − 2)

So: x∗1 = 0, x∗2 = 1, x∗3 = 2


Stability: f ′ (x) = 3x2 − 6x + 3 = 3(x − 1)2

• At x∗1 = 0: f ′ (0) = 3 (unstable)

• At x∗2 = 1: f ′ (1) = 0 (superstable)

• At x∗3 = 2: f ′ (2) = 3 (unstable)

Despite x∗1 and x∗3 having the same derivative magnitude, their stability is the
same (unstable).

7.2 Example 4: Sine Map


Example

Model: f (x) = r sin(πx), with 0 ≤ x ≤ 1, r > 0


Fixed points: x = r sin(πx)

• Always x∗ = 0 is a fixed point

• For r > 1/π ≈ 0.318, additional fixed points appear

14
Week 4: Fixed Points and Stability Dynamical Systems A

Stability at x∗ = 0: f ′ (x) = rπ cos(πx)

f ′ (0) = rπ

So x∗ = 0 is stable if r < 1/π, unstable if r > 1/π.


This is another example of a bifurcation (pitchfork bifurcation).

8 Summary and Key Takeaways


• Fixed points satisfy f (x∗ ) = x∗

• Stability is determined by |f ′ (x∗ )|:

– |f ′ (x∗ )| < 1: Stable (attracting)


– |f ′ (x∗ )| > 1: Unstable (repelling)
– |f ′ (x∗ )| = 1: Neutral (requires nonlinear analysis)

• Period-k cycles are fixed points of f k (x)

• Cobweb diagrams provide visual iteration insight

• Bifurcations occur when stability changes as parameters vary

• The logistic map exhibits rich dynamics: fixed points → period-2 → period-4 →
chaos

Next Week: Bifurcation theory and introduction to chaos in 1D maps

15
The Logistic Map Dynamical Systems — Week 5

Dynamical Systems
Week 5 — Lecture Notes & Tutorial
The Logistic Map: Cobwebbing, Bifurcations and Chaos

Part I Lecture Notes

1 Derivation of the Logistic Map


Recall the logistic growth equation for a population N (t):
 
dN N
= rN 1 − , (1)
dt K
where r > 0 is the intrinsic growth rate and K > 0 is the carrying capacity.

1.1 Non-dimensionalisation
Introduce the dimensionless variable
N
x= =⇒ N = Kx.
K
dx 1 dN dx dN
Differentiating: = , so K = .
dt K dt dt dt
Substituting into (1):
dx dx
K = r(Kx)(1 − x) =⇒ = rx(1 − x) =: f (x).
dt dt
Since x = N/K is a ratio (population per carrying capacity), it is dimensionless. The
logistic growth equation becomes:
dx
= rx(1 − x). (2)
dt

1.2 The Logistic Map


The orbits of f (x) = rx(1 − x) form the logistic map. For an initial value x0 we iterate:
xn+1 = f (xn ) = rxn (1 − xn ), n = 0, 1, 2, . . . (3)
with x1 = f (x0 ), x2 = f (x1 ) = f 2 (x0 ), x3 = f 3 (x0 ), and so on.
We restrict to 0 < x < 1 and 1 < r < 4 to keep the population bounded; outside this
range the population goes to extinction.

1
The Logistic Map Dynamical Systems — Week 5

2 Fixed Points and Stability


2.1 Finding the Fixed Points
Fixed points x∗ satisfy f (x∗ ) = x∗ :
rx∗ (1 − x∗ ) = x∗ x∗ r(1 − x∗ ) − 1 = 0.
 
=⇒
This yields two fixed points:
1
x∗1 = 0, x∗2 = 1 − (exists only for r > 1).
r

2.2 Stability Criterion


Stability of fixed points
A fixed point x∗ of f is
• stable (attracting) if |f ′ (x∗ )| < 1,
• unstable (repelling) if |f ′ (x∗ )| > 1,
• marginal if |f ′ (x∗ )| = 1 (higher-order analysis needed).

For the logistic map f (x) = rx(1 − x):


f ′ (x) = r(1 − 2x).

2.2.1 Stability of x∗1 = 0


f ′ (0) = r.
Hence x∗1 = 0 is stable for 0 < r < 1 and unstable for r > 1.

2.2.2 Stability of x∗2 = 1 − 1/r


f ′ (x∗2 ) = r 1 − 2 1 − 1 2
 
r
=r r
− 1 = 2 − r.
Therefore x∗2 is stable when |2 − r| < 1, i.e. 1 < r < 3.
Fixed point stability summary
• x∗1 = 0: stable for 0 < r < 1; unstable for r > 1.
• x∗2 = 1 − 1/r: exists for r > 1; stable for 1 < r < 3.
• At r = 3: f ′ (x∗2 ) = −1 (tangent to graph has slope −1, perpendicular to diagonal
y = x). Period-doubling onset.
• For r > 3: |f ′ (x∗2 )| > 1, so x∗2 is unstable.

Remark. At r = 1 the two fixed points x∗1 and x∗2 coincide. The moment r passes 1, x∗2
comes into existence. If r increases to 3 then f ′ (x∗2 ) = −1, which means the gradient of
the tangent line is −1 and it is at right angles to the line y = x.

3 Cobweb Diagrams
A cobweb (staircase) diagram simultaneously displays the graph of y = f (x), the
diagonal y = x, and the orbit x0 , x1 , x2 , . . . encoded as an alternating vertical/horizontal
path.

2
The Logistic Map Dynamical Systems — Week 5

3.1 Algorithm
Starting at the point (x0 , x0 ) on the diagonal, each iteration consists of:
1. Vertical move: draw a vertical segment from (xn , xn ) to (xn , f (xn )), reading off the
new iterate xn+1 = f (xn ) on the curve.
2. Horizontal move: draw a horizontal segment from (xn , xn+1 ) to (xn+1 , xn+1 ), trans-
ferring the new value to the diagonal ready for the next step.
Reading a cobweb diagram
• Converging spiral/staircase inward ⇒ attracting fixed point.
• Diverging path outward ⇒ repelling fixed point.
• Closed rectangle with 2 corners on curve ⇒ stable 2-cycle.
• Closed path with 4 corners on curve ⇒ stable 4-cycle.
• Irregular, non-repeating dense path ⇒ chaotic orbit.

3.2 Computed Examples for the Logistic Map


The following cases were run with x0 = 0.1 (unless noted). They can be reproduced with
cobweb diagram.m.

Case r x0 N Behaviour

(a) 0.5 0.5 10 Rapid decay to x∗1 = 0 (stable fixed point)


(b) 1.5 0.1 10 Monotone convergence to x∗2 = 1/3
(c) 2.5 0.1 10 Oscillating convergence to x∗2 = 0.6
(d) 3.0 0.05 2000 Marginal case; very slow convergence
(e) 3.2 0.1 1000 Period-2 cycle
(f) 3.3 0.1 2000 Period-2 cycle (similar behaviour)
(g) 3.5 0.05 100 Period-4 cycle
(h) 3.5699 0.05 100 Period-8 cycle (near Feigenbaum point)
(i) 3.6 0.2 1000 Chaotic behaviour begins
(j) 3.8 0.05 100 Fully chaotic orbit

4 Bifurcation Diagrams
A bifurcation diagram plots the long-run attractor xn (after discarding transients)
against the parameter r. Successive period-doublings are visible in the following windows:

3
The Logistic Map Dynamical Systems — Week 5

Plot r range Feature

(a) 0.5 ≤ r ≤ 2.5 Forward transcritical bifurcation at r = 1


(b) 0.8 < r < 3.1 Supercritical pitchfork bifurcation near r = 3
(c) 2.0 < r < 3.35 Stable x∗2 branch; onset of 2-cycle
(d) 2.9 < r < 3.65 Period-doubling: 2-cycle and 4-cycle
(e) 2.92 < r < 3.56995 Full cascade: 2, 4, 8 cycles
(f) 2.9 < r < 3.6 Period-8 cycle clearly visible
(g) 2.9 < r < 3.9 Full cascade into chaos

Period-doubling cascade — Feigenbaum universality


The period-doubling bifurcation values r1 = 3, r2 ≈ 3.449, r3 ≈ 3.544, . . . accumulate
at
r∞ ≈ 3.56995 (the Feigenbaum point).
The ratio of successive intervals converges to the Feigenbaum constant:
rn − rn−1
−→ δ ≈ 4.6692.
rn+1 − rn
This ratio is universal: it appears in every unimodal map undergoing period-
doubling, regardless of the specific functional form.

4
The Logistic Map Dynamical Systems — Week 5

Part II Tutorial Problems

5 Problem 1 — Tumour Growth: the Gompertz Model


The growth of cancerous tumours can be modelled by the equation
dN
= −aN ln(bN ), (4)
dt
where N (t) is proportional to the number of cells in the tumour, and a, b > 0 are param-
eters.

5.1 Part (a) — Direction Field


Sketch the right-hand side of (4) and show the direction flows for varying initial condi-
tions.

Let g(N ) = −aN ln(bN ). Key observations:


• g(N ) = 0 when ln(bN ) = 0, i.e. bN = 1, giving the fixed point N ∗ = 1/b.
• g(N ) = 0 also at N = 0 (since N ln(bN ) → 0 as N → 0+ ).
• g(N ) > 0 (growth) when bN < 1, i.e. N < 1/b.
• g(N ) < 0 (shrinkage) when bN > 1, i.e. N > 1/b.
Direction flows: for N < 1/b the vector field points right (N increases); for N > 1/b it
points left (N decreases); both flows drive N toward N ∗ = 1/b.

5.2 Part (b) — Fixed Points and Stability


Find the fixed points and determine their stability.

dN
Fixed points satisfy = 0:
dt
1
−aN ln(bN ) = 0 =⇒ N = 0 or N = .
b

Stability: Let g(N ) = −aN ln(bN ). Then:

g ′ (N ) = −a ln(bN ) + 1 .
 

• At N ∗ = 0: g ′ (0+ ) → +∞ (unstable; any positive perturbation grows away from


zero).
• At N ∗ = 1/b: g ′ (1/b) = −a[ln(1) + 1] = −a < 0 (stable; perturbations decay back
to 1/b).

5
The Logistic Map Dynamical Systems — Week 5

Conclusion
N ∗ = 0 is an unstable fixed point. N ∗ = 1/b is a stable fixed point. All positive
initial conditions converge to the equilibrium 1/b.

5.3 Part (c) — Biological Interpretation of a and b


Interpret the parameters a and b biologically.
• b determines the tumour’s carrying capacity: the equilibrium cell count is N ∗ =
1/b. A larger b corresponds to a smaller maximum tumour size.
• a controls the rate of convergence toward equilibrium. Larger a gives faster
growth (or regression) toward N ∗ = 1/b.
• Together, a and b set the time-scale and the equilibrium size of the tumour.

6 Problem 2 — The Allee Effect


dN
For some species the effective per-capita growth rate /N is highest at intermediate
dt
N . This is called the Allee effect (e.g. species needing a minimum population to locate
mates or resist predation).

6.1 Part (a) — Demonstrating the Allee Effect


Show that the model
1 dN
= r − a(N − b)2 (5)
N dt
provides an example of the Allee effect if r, a, b satisfy certain conditions.

The per-capita growth rate h(N ) = r − a(N − b)2 is a downward-opening parabola with:
• Maximum value r attained at N = b.
• h(N ) > 0 (population grows) for |N − b| <
p p p 
r/a, i.e. N ∈ b − r/a, b + r/a .
• h(N ) < 0 (population declines) outside this interval.
The Allee effect requires a positive maximum
p at intermediate N , achieved when b > 0,
r > 0, a > 0. If additionally b > r/a, then h(N ) < 0 near N = 0: very small
populations decline — the strong Allee effect.

6.2 Part (b) — Fixed Points and Classification


Find all fixed points and classify them.

dN
= N r − a(N − b)2 = 0:
 
Fixed points of
dt
r
r
N1∗ = 0, ∗
N2,3 =b± .
a

6
The Logistic Map Dynamical Systems — Week 5

Let F (N ) = N [r − a(N − b)2 ]. Then:

F ′ (N ) = r − a(N − b)2 − 2aN (N − b).

• At N1∗ = 0: F ′ (0) = r − ab2 . Stable if r < ab2 (strong Allee effect); unstable if
r > ab2 .
• At N2∗ = b − r/a: unstable — this is the Allee threshold; populations below it
p

decline to extinction.
• At N3∗ = b + r/a: stable — this is the carrying capacity K = b + r/a.
p p

6.3 Part (c) — Number of Fixed Points


Show that depending on the values of r, a and b there are either 2 or 3 fixed points.
• 3 fixed points: N1∗ = 0, b − r/a, b + r/a when b > r/a (so that N2∗ > 0).
p p p

• 2 fixed points: when p


p
b = r/a the twopnon-trivial roots merge in a saddle-node
bifurcation. When b < r/a the root b− r/a < 0 is biologically irrelevant, leaving
N ∗ = 0 and the carrying capacity K.

6.4 Part (d) — Discussion of the Full Model


Discuss the growth model
dN
= N r − a(N − b)2 .
 
(6)
dt
Key features of model (6):
1. N = 0 is always a fixed point. Its stability depends on whether r < ab2 (stable:
small populations go extinct — strong Allee effect) or r > ab2 (unstable: all positive
populations grow initially).
2. When three
p fixed points exist, the interior fixed points give the Allee
p threshold
A = b − r/a (lower, unstable) and the carrying capacity K = b + r/a (upper,
stable).
3. Biological implication: if the initial population N0 lies below the Allee threshold A,
the population declines to extinction. If N0 > A, the population grows toward K.
4. Bifurcation structure: as r increases through ab2 , a saddle-node bifurcation cre-
ates two new fixed points (the threshold and the carrying capacity) from a single
tangency.
5. Conservation biology: a minimum viable population must exceed the Allee threshold
for a species to persist.

7
The Logistic Map Dynamical Systems — Week 5

Summary: Allee effect conditions

• Strong Allee effect (decline near N = 0): requires r < ab2 . p


• Three biologically meaningfulrfixed points exist when: b > r/a.
r
• Carrying capacity: K = b + .
r a
r
• Allee threshold: A = b −
p
(positive only when b > r/a).
a

End of Week 5 Notes & Tutorial

8
Theorem (Dynamics of fc (x) = x2 + c)

(a) For c > 1/4: all orbits tend to +∞.


(b) For c = 1/4: the unique fixed point x∗1 = x∗2 = 1/2 is neutral (saddle-node
bifurcation).
(c) For −3/4 < c < 1/4: x∗2 is attracting, x∗1 is repelling; no 2-cycles exist.
(d) For c = −3/4: x∗2 is neutral (|fc′ (x∗2 )| = 1); a period-doubling bifurcation
occurs.
(e) For c < −3/4: x∗2 is repelling; an attracting 2-cycle {x∗3 , x∗4 } exists.

7.2 Period-doubling bifurcation


When c < −3/4, the fixed point x∗2 loses stability and a 2-cycle (period-2 orbit) is born.
Points of period 2 satisfy fc2 (x) = x but fc (x) ̸= x. Since fc (x∗1 ) = x∗1 and fc (x∗2 ) = x∗2 ,
the factor (x − x∗1 )(x − x∗2 ) = x2 − x + c divides fc2 (x) − x:

fc2 (x) − x = (x2 + c)2 + c − x = x4 + 2cx2 − x + c2 + c.

Polynomial long division by (x2 − x + c) gives quotient x2 + x + c + 1, so the period-2


points satisfy √
2 ∗ −1 ± −3 − 4c
x + x + c + 1 = 0 =⇒ x3,4 = .
2
These are real iff −3 − 4c ≥ 0, i.e. c ≤ −3/4, confirming the bifurcation at c = −3/4.
Stability of the 2-cycle. By the chain rule:

(fc2 )′ (x∗3 ) = fc′ (x∗4 ) · fc′ (x∗3 ) = 2x∗4 · 2x∗3 = 4x∗3 x∗4 .

Since x∗3 x∗4 = c + 1 (product of roots of x2 + x + c + 1 = 0):

(fc2 )′ (x∗3 ) = 4(c + 1).

The 2-cycle is attracting iff |4(c + 1)| < 1, i.e. −5/4 < c < −3/4. At c = −5/4:
(fc2 )′ (x∗3 ) = −1 and another period-doubling occurs, generating a 4-cycle.
Exercise 7.1 (Period-doubling cascade). Investigate the dynamics of fc (x) = x2 + c
for the following ranges. In each case, state whether there are fixed points, 2-cycles, or
neither, and whether they are attracting or repelling.
(a) −3/4 < c < 1/4
(b) c = −3/4
(c) −5/4 < c < −3/4
(d) c = −5/4
(e) c < −5/4
Use the cobwebbing method to illustrate your conclusions for representative values of c
in each range.

10
Applied Mathematics – Honours

n
Applied Mathematics – Honours


n ∈

n − n

n n

n n n n ∈

n
Applied Mathematics – Honours

n − n

n n n − n

− n n − n − n

n n − n − n

 
n n − n − n n − n ✓

n · n · n − n · − n

n n − n · − n

n n − n − n

n n − n
 
n n − n

n n ≥
Applied Mathematics – Honours



≤ ≤
Applied Mathematics – Honours

n n − n − n

n · n − n − n

n n − n − n

n · n − n n n − n
n n
n n − n − n n n − n − n

n
 
n n − n n
n
n n − n − n n n − n − n

n
n n ≜
− n − − n

n n

× ×
× ×
Applied Mathematics – Honours

×

∗ ∗ ∗
n



− ∗ − − ∗
h i
∗ ∗ ∗ ∗
− − −

̸
h i
∗ ∗ ∗
− − − −


∗ ∗

− − − − − − − − −

− − − − −

̸
√ p
∗ − ± − − − ± − − − −
,

∗ − −
̸
− −
p∗ / p∗3 ∈ ,
Applied Mathematics – Honours

β < x, δ

• ∗

• ∗

• ∗ ∗

′ ∗ ∗
| | ⇒
′ ∗ ∗
| | ⇒



′ ′
′ − −

≤ ≤
• ∗

• ∗

• ∗
̸
∗ ∗

n
n
Applied Mathematics – Honours



n n − n
n n
n n − n − n

∗ ∗ ∗ ∗ ∗

∗ ∗

p
∗ − ± − − −
,

− − − −

∗ − n n
− n − n
Applied Mathematics – Honours

• ∗

• ∗

• ∗

∗ ′ ∗
| |

• ∗

• ∗

• ∗
inconclusive ′

• ∗
| ′ ∗
|

heterozygote advantage


Applied Mathematics – Honours

n →
n →

n → ∈

n n



′ ∗
Applied Mathematics – Honours

∗ ∗


∗ ∗
′ ∗
| |


∗ ′ ∗

n − n

n → →∞


2D Dynamical Systems

These are of the form

Xn+1 = f (Xn ) where Xn = (x1n , x2n ), Xn ∈ R2 .

Consider a continuous system of interacting species in which species x is the prey and
species y is the predator. The interaction is a predator–prey interaction.
So the model is given by
dx
= rx − a12 xy ,
dt |{z} | {z }
growth of prey decline of prey

dy
= a21 xy − δy .
dt | {z } |{z}
growth of predator decline of predator

Here:

• r — growth rate of the prey,

• a12 — predation rate,

• a21 — growth rate of the predator as a result of feeding on the prey,

• δ — death rate of the predator.

Consistency of Units
Consider
dx
= rx − a12 xy,
dt
 
prey 1 1
= ·prey − ·prey · pred,
|time
{z } |time{z } prey · time
| {z }
dx r a12

dt
1 1
so r = and a12 = .
time pred · time

1
Converting to a Discrete System
Converting the continuous system into a discrete system results in the following:

xn+1 = xn + r xn − a12 xn yn = f (xn , yn ),


yn+1 = yn + a21 xn yn − δ yn = g(xn , yn ).

In general, a 2-D nonlinear discrete system has the form


)
xn+1 = f (xn , yn )
(∗)
yn+1 = g(xn , yn )

The steady states of (∗) are obtained from the solutions of

x∗ = f (x∗ , y ∗ ),
y ∗ = g(x∗ , y ∗ ).

To determine the stability of the steady state we set

xn = x∗ + un =⇒ xn+1 = x∗ + un+1 ,

yn = y ∗ + vn =⇒ yn+1 = y ∗ + vn+1 ,
where un , vn are small perturbations around the steady state.

So

xn+1 = x∗ + un+1 = f (x∗ + un , y ∗ + vn ) = f (xn , yn ),


yn+1 = y ∗ + vn+1 = g(x∗ + un , y ∗ + vn ) = g(xn , yn ).

We can use Taylor series expansion around (x∗ , y ∗ ) so that


∂f (x∗ , y ∗ ) ∂f (x∗ , y ∗ )
x∗ + un+1 = f (x∗ , y ∗ ) + un + vn + O(un , vn ),
∂x ∂y
∂g(x∗ , y ∗ ) ∂g(x∗ , y ∗ )
y ∗ + vn+1 = g(x∗ , y ∗ ) + un + vn + O(un , vn ).
∂x ∂y
Neglecting higher order terms we have:
∂f ∂f
un+1 = un + vn ,
∂x ∂y
∂g ∂g
vn+1 = un + vn .
∂x ∂y
In matrix form:
∂f ∂f
 
   
un+1  ∂x
 un .
∂y 
=
vn+1  ∂g ∂g  vn
∂x ∂y
 
un
Hence Wn+1 = J Wn , where J is the Jacobian matrix and Wn = .
vn

2
Note that
∂f (x∗ , y ∗ ) ∂f (x∗ , y ∗ )
 
 ∂x ∂y 
J = ∂g(x∗ , y ∗ ) ∂g(x∗ , y ∗ )  .

∂x ∂y
The stability of (∗) is determined by the eigenvalues of J. If the eigenvalues of J, say λi ,
satisfy:
(a) |λi | < 1 for i = 1, 2, then (∗) is stable,
(b) |λi | > 1, then (∗) is unstable,
(c) |λi | = 1, we cannot conclude.

The Stability Theorem of a 2-D Discrete System


Theorem 1. Let x∗ be a fixed point of (∗) where x ∈ R2 , f : R2 → R2 , and (∗) is such
that Xn+1 = f (Xn ) (where Xn is a vector and f is vector-valued, with components x11 , x12
and f1 , f2 ). Let J be the Jacobian matrix evaluated at x∗ . Then:
(i) x∗ is stable if all eigenvalues of J have magnitude less than 1.
(ii) x∗ is unstable if at least one of the eigenvalues has magnitude greater than 1.
(iii) Otherwise, if the magnitude is equal to one, we cannot conclude.

For 2 × 2 systems (2-D), the eigenvalues are obtained from:


|J − λI| = 0 (characteristic polynomial).
 
a b
Let J = . Then
c d
a−λ b
|J − λI| = = 0,
c d−λ
λ2 − (a + d) λ + ad − bc = 0.
Note that a + d is the trace of J, written tr J, and ad − bc is the determinant of J, written
det J.
So we have
λ2 − tr J λ + det J = 0,
p
tr J ± (tr J)2 − 4 det J
λ1,2 = .
2

For stability we require that |λ1,2 | < 1:


p
tr J + (tr J)2 − 4 det J
< 1, (∗∗)
2
p
tr J − (tr J)2 − 4 det J
< 1. (∗ ∗ ∗)
2

3
Consider (∗∗):
p p
tr J + (tr J)2 − 4 det J tr J + (tr J)2 − 4 det J
< 1 =⇒ −1 < < 1,
2 2
p
=⇒ −2 < tr J + (tr J)2 − 4 det J < 2,
p
=⇒ −2 − tr J < (tr J)2 − 4 det J < |2 −{ztr J} .
| {z }
(a)
p
For (a): (tr J)2 − 4 det J < 2 − tr J. Squaring both sides:

tr J < 1 + det J.

Considering (b) and squaring:


p
−2 − tr J < (tr J)2 − 4 det J,

=⇒ tr J < −(1 + det J),


=⇒ |tr J| < 1 + det J.

Exercise: Use (∗ ∗ ∗) to obtain further relationships between the trace and determinant.

However, all four cases together give:

|tr J| < 1 + det J < 2. (A)

(A) provides the necessary and sufficient condition for the eigenvalues of J to have mag-
nitude less than 1.
The conditions in (A) are called the Jury conditions.

So a visualisation of (A) is best done on a graph of det J vs tr J:

det J

(−3, 1) (2, 1)
1 det J = 1

Stable
tr J
−3 −2 −1 1 2 3

−1
(0, −1)

4
Love Affairs
Consider two lovers, Romeo and Juliet, where the love of Romeo at time step n is Rn and
that of Juliet is Jn .
The dynamics of their love affair is given by
)
Rn+1 = ar Rn + pr Jn
(1)
Jn+1 = aj Jn + pj Rn

At equilibrium:

R∗ = ar R∗ + pr J ∗ =⇒ (ar − 1) R∗ + pr J ∗ = 0,
J ∗ = aj J ∗ + pj R∗ =⇒ pj R∗ + (aj − 1) J ∗ = 0.

In matrix form:    ∗  
ar − 1 pr R 0
∗ = .
pj aj − 1 J 0
At steady state (R∗ , J ∗ ) = (0, 0).
System (1) has a unique solution (0, 0) provided det A ̸= 0, where
 
ar − 1 pr
A= .
pj aj − 1

The Jacobian matrix of (1) is  


ar p r
J= .
p j aj
So
tr J = ar + aj , det J = ar aj − pr pj .
Using the Jury conditions, (0, 0) is stable if:

|ar + aj | < 1 + (ar aj − pr pj ) < 2.

A special case is when Romeo and Juliet preserve their love/hate from day to day, i.e.

ar + pj = 1 and pr + aj = 1.

Exercise: Show that det A = 1.

5
Mathematical Modelling

Introduction
Modelling involves translating real-life scenarios into models.
Definition:
A model is a representation of a real-life scenario.
Models are classified by their representational form:
ˆ If represented physically, we have a physical model.
ˆ If represented mathematically, we have a mathematical model.
ˆ If represented statistically, we have a statistical model, and so on.

Mathematical modelling is therefore the process of expressing real-life scenarios using


mathematical formulations.

The Modelling Process


The diagram below illustrates the cyclic nature of the mathematical modelling process.
When the model’s solutions do not adequately describe the real-life scenario, the cycle
returns to the conceptualisation phase for refinement.

Results do not
match the problem
Conceptualisation
Interpretation &Application
Real-lifeScenario Conceptualisation:Identify Problem

Interpretation &Application

Analyse /Link to Real-lifeScenario FormulateEquations(M

Steps of the Modelling Process


Step 1 — Problem Identification
Have a real-life scenario/problem. Identify the problem or research question.
Problem identification is a skill. Ask yourself:
“What is the problem?” “What is the research question?”

1
Step 2 — Conceptualise the Model
ˆ Identify the state variables.
ˆ Understand the dynamics/interplay and linkages between the state variables.
ˆ Set parameters — rates of flow, growth, decline, etc.

Step 3 — Write the Model


Write the model using model assumptions.

Step 4 — Solve the Model


Solve the model. The solutions must say something meaningful about the problem. If
they do not, return to the conceptualisation phase and re-examine the dynamics, variables,
and assumptions.

Step 5 — Interpret and Link to Reality


Link the solution to the real-life scenario to answer the research question.

Example 1: Dynamics of Red Blood Cells (RBCs)


Background
RBCs carry oxygen via blood flow. Abnormalities can cause irregularities in blood flow
and the supply of oxygen. Common abnormalities include anaemias (e.g. sickle-cell
anaemia, haemolysis).
ˆ RBCs are produced in bone marrow at a rate of approximately 200 billion per
day.
ˆ Each cell has a half-life of 55–60 days; cells are lost after about 120 days.
Goal: Model the production dynamics of red blood cells.

State Variables

Rn = amount of RBCs at stage n,


Bn = amount of bone marrow at stage n.

Model

Rn+1 = Rn − aRn + βBn


|{z} |{z}
parameter measuring rate of
loss of RBCs production
(1)
Bn+1 = Bn − εBn + δ(Re − Rn )
|{z} | {z }
rate of loss of rate of growth of BM;
bone marrow baseline level of RBCs

2
System (1) constitutes a mathematical model of RBC production dynamics.

Parameter summary:

Parameter Description
a Rate of loss of RBCs
β Rate of RBC production from bone marrow
ε Rate of loss of bone marrow
δ Rate of growth of bone marrow
Re Baseline (equilibrium) level of RBCs

Example 2: Drug Uptake by the Body


Model Setup
Let Bn and Tn be the amounts of drug in the blood vessels and tissues on day n,
respectively. The rate of change of the drug in the blood system depends on:
ˆ the ingestion rate Dn ,
ˆ the excretion rates µ (blood) and s (tissue), and
ˆ the exchange rates with body tissues (k21 , k12 ).

Model Equations

Bn+1 = Bn + Dn − µBn − k21 Bn + k12 Tn (2)


|{z} | {z } | {z }
excretion transfer rate transfer rate
blood→tissue tissue→blood

Tn+1 = Tn − sTn + k21 Bn − k12 Tn (3)


|{z} | {z } | {z }
excretion rate transfer transfer
from tissue blood→tissue tissue→blood

Example 3: Spread of Influenza (Flu) — SIR Model


Population Classes
The population is divided into three classes:
S Susceptible — not currently diseased, but can be infected.
I Infected — infectious; spreads infection through contact.
R Recovered/Immune — has recovered and is now immune.

3
SIR Model Equations

Sn+1 = Sn − βSn In (4)


| {z }
rate of infection;
interaction of S and I

In+1 = In + βSn In − γIn (5)


| {z } |{z}
new infections recovery rate

Rn+1 = Rn + γIn (6)

where β is the rate of infection and γ is the recovery rate.

Tutorial Problems
Problem 1 — Red Blood Cell Precursor Model
Let Mt be the concentration of red blood cell precursors in the bone marrow, and Rt the
concentration of mature red blood cells in circulation. Consider the discrete model:

Mt+1 = aMt (1 − bRt ),


Rt+1 = cMt + (1 − d)Rt ,

where a, b, c, d are positive constants.


(a) Interpret the parameters a, b, c, and d.
(b) Find the steady state(s).
(c) Perform the stability analysis of the non-trivial equilibrium.
(d) What happens biologically if b = 0?

Problem 2 — Drug Concentration Model


Let Dt be the drug concentration in the blood and Tt the concentration in the tissues.
The dynamical system

Dt+1 = Dt − αDt + βTt ,


Tt+1 = Tt + αDt − βTt ,

models these concentrations.


(a) Interpret each term and each parameter in this model.
(b) Show that the total drug concentration is conserved.
(c) Find the steady state and its stability.

4
Problem 3 — Influenza Spread (SIR with Population)
The spread of influenza between susceptibles S, infectives I, and the recovered R is
modelled by:
St It
St+1 = St − β ,
N
St It
It+1 = It + β − δIt ,
N
Rt+1 = Rt + γIt ,

where all parameters are assumed positive and N is the total population.
(a) Interpret the parameters β and δ.
(b) Show that the total population remains constant over time.
(c) If the disease-free equilibrium is (S ∗ , 0, 0), derive the conditions for the stability of
this steady state.

Problem 4 — Analysis of Examples 1 and 2


For Examples 1 and 2:
(a) Find the steady states.
(b) Determine the stability of each steady state.
(c) Discuss the biological implications of the conditions for stability of each steady state.

5
Population Dynamics Weeks 1011

Population Dynamics
Weeks 1011: HostParasitoid Systems and Interacting Species

1. HostParasitoid Systems

Denition 1.1. A parasitoid is an insect that lays its eggs inside a host. The host is
then responsible for generating the next generation of parasitoids.

Let

ˆ Hn = the host population at generation n,

ˆ Pn = the parasitoid population at generation n.

1.1 The NicholsonBailey Model

A classical example of a hostparasitoid model is the NicholsonBailey model, given


by

NicholsonBailey Model (∗)



Hn+1 = kHn e−aPn ,
(∗)
= cHn 1 − e−aPn ,
P 
n+1

where

ˆ k = growth rate of the host,

ˆ c = growth rate of the parasitoid,

ˆ a = searching eciency of the parasitoid.

The model can equivalently be written in the form


 
kHn 1
f (Hn , Pn ) = aPn , g(Hn , Pn ) = cHn 1 − aPn .
e e

1.2 Steady States of the NicholsonBailey Model

Setting Hn+1 = Hn = H ∗ and Pn+1 = Pn = P ∗ , the steady-state equations are:

1
Population Dynamics Weeks 1011


H ∗ = kH ∗ e−aP , (1)
−aP ∗
P ∗ = cH ∗ 1 − e (2)

.

From (1):
∗ ∗
H ∗ 1 − ke−aP = 0 =⇒ H ∗ = 0 or 1 = ke−aP .


Case 1: H ∗ = 0. Substituting into (2) gives P ∗ = 0. This yields the trivial steady
state
E0 = (0, 0).

E0 represents the extinction of both populations and is of no biological interest.

Case 2: 1 = ke−aP . Solving for P ∗ :


1 ∗ ∗ ln k
= e−aP =⇒ k = eaP =⇒ ln k = aP ∗ =⇒ P∗ = .
k a

ln k
Substituting P ∗ = into (2):
a
 
ln k ∗

−a· lnak

∗ − ln k ∗ 1 k−1
= cH ∗ ·

= cH 1 − e = cH 1 − e = cH 1 − .
a k k

Therefore
k ln k
H∗ = .
ac(k − 1)

This gives the non-trivial (coexistence) steady state


 
k ln k ln k
E1 = , .
ac(k − 1) a

Existence Condition for E1


The coexistence steady state E1 exists if and only if k > 1.

1.3 Stability Analysis of E1

To study the stability of the non-trivial steady state, we compute the Jacobian matrix of
system (∗) at E1 .

2
Population Dynamics Weeks 1011

The Jacobian of (f, g) with respect to (H, P ) is:


 
∂f ∂f
ke−aP −akH e−aP
!
 ∂H ∂P 
J =
 ∂g
= .
∂g  c 1 − e−aP ace−aP H

∂H ∂P

∗ ∗ k ln k
Evaluating at E1 using ke−aP = 1 (i.e. e−aP = 1/k ) and H ∗ = :
ac(k − 1)


J11 = ke−aP = 1,
∗ k ln k 1 k ln k
J12 = −akH ∗ e−aP = −ak · · =− ,
ac(k − 1) k c(k − 1)
 
−aP ∗
 1 c(k − 1)
J21 =c 1−e =c 1− = ,
k k
∗ 1 k ln k ln k
J22 = ace−aP H ∗ = ac · · = .
k ac(k − 1) k−1

Therefore the Jacobian evaluated at E1 is:

k ln k
 
1 −
c(k − 1) 
(3)

JE1 =
 .
c(k − 1) ln k 
k k−1
The trace and determinant of JE1 are:

ln k
tr JE1 = 1 + , (4)
k−1
ln k ln k
det JE1 = + ln k = ln k + . (5)
k−1 k−1

1.4 Jury Stability Conditions

For a 2 × 2 discrete-time map, the non-trivial steady state is locally asymptotically stable
if and only if the Jury conditions hold:

Jury Stability Conditions

|tr JE1 | < 1 + det JE1 < 2.

3
Population Dynamics Weeks 1011

Condition 1: | tr JE1 | < 1 + det JE1 . Substituting (4) and (5):


ln k ln k
1+ < 1 + ln k + .
k−1 k−1

This simplies to 0 < ln k , which holds for all k > 1. Condition 1 is satised.

Condition 2: 1 + det JE1 < 2. This requires det JE1 < 1, i.e.
ln k
ln k + < 1.
k−1

We ask: is the opposite true? Suppose

ln k
ln k + > 1.
k−1

Multiplying through by (k − 1):

(k − 1) ln k + ln k > k − 1 =⇒ k ln k − ln k + ln k > k − 1 =⇒ k ln k − k + 1 > 0.

Dene f (k) = k ln k − k + 1. We need to show f (k) > 0 for k > 1.

ˆ f (1) = 1 · 0 − 1 + 1 = 0.
1
ˆ f ′ (k) = ln k + k · − 1 = ln k > 0 for k > 1.
k

So f is strictly increasing for k > 1, hence f (k) > f (1) = 0 for all k > 1. Therefore

ln k
ln k + > 1 for all k > 1,
k−1

which means det JE1 > 1 and the Jury condition 1 + det JE1 < 2 is not satised.

Instability of E1
The Jury condition 1 + det JE1 < 2 fails for all k > 1. Consequently, the coexistence
equilibrium E1 of the NicholsonBailey model is unstable.

1.5 A Modied HostParasitoid Model

To correct the instability of the coexistence equilibrium, the NicholsonBailey model can
be modied by introducing logistic self-regulation in the host. The adjusted model is:

4
Population Dynamics Weeks 1011

Modied HostParasitoid Model (∗∗)



Hn+1 = Hn e r(1−Hn /k) · e−aPn ,
(∗∗)
−aPn

n+1 = cHn 1 − e
P ,

Exercise
Carry out the full stability analysis of model (∗∗).

2. Interacting Species Models

2.1 Competing Species

Example 2.1. Consider a model of interacting species given by



An+1 = µ1 An − µ3 An Bn ,
(6)
B
n+1 = µ2 Bn − µ4 An Bn ,

where µi > 0, i = 1, . . . , 4, are positive constants.

(a) Type of Interaction

This is a competing species model. The interaction of the two species leads to a
decline in both populations: the terms −µ3 An Bn and −µ4 An Bn represent the negative
eect each species has on the other's growth.

(b) Fixed Points

At a steady state, An+1 = An = A∗ and Bn+1 = Bn = B ∗ :

A∗ = µ1 A∗ − µ3 A∗ B ∗ = A∗ [µ1 − 1 − µ3 B ∗ ] ,
B ∗ = µ2 B ∗ − µ4 A∗ B ∗ = B ∗ [µ2 − 1 − µ4 A∗ ] .

µ1 − 1
From the rst equation: A∗ = 0 or B ∗ = .
µ3
µ2 − 1
From the second equation: B ∗ = 0 or A∗ = .
µ4

5
Population Dynamics Weeks 1011

Combining the cases gives three xed points:


 
µ2 − 1 µ1 − 1
(0, 0), , .
µ4 µ3

The coexistence xed point exists subject to µ1 > 1 and µ2 > 1.

(c) Stability Analysis

The Jacobian of system (6) is


!
µ1 − µ3 B −µ3 A
J= .
−µ4 B µ 2 − µ4 A

Stability at (0, 0). !


µ1 0
J(0,0) = .
0 µ2
The eigenvalues are λ1 = µ1 and λ2 = µ2 . Given µ1 = 1.2 > 1 and µ2 = 1.3 > 1, we have
|λ1 | > 1 and |λ2 | > 1, so (0, 0) is unstable.

Stability at the coexistence point. With µ1 = 1.2, µ2 = 1.3, µ3 = 0.001, µ4 = 0.002,


the coexistence point is
   
µ2 − 1 µ1 − 1 0.3 0.2
, = , = (150, 200).
µ4 µ3 0.002 0.001

Evaluating the Jacobian at this point:


! !
µ1 − µ3 (µ1 − 1)/µ3 −µ3 · (µ2 − 1)/µ4 1 −0.15
J= = .
−µ4 · (µ1 − 1)/µ3 µ2 − µ4 · (µ2 − 1)/µ4 −0.4 1

The characteristic equation is (λ − 1)2 − (−0.15)(−0.4) = 0, giving



(λ − 1)2 = 0.06 =⇒ λ = 1 ± 0.06.

Thus
λ1 ≈ 1.245, |λ1 | > 1; λ2 ≈ 0.755, |λ2 | < 1.

Since one eigenvalue lies outside the unit circle, the coexistence point is a saddle point
and hence unstable.

6
Population Dynamics Weeks 1011

2.2 A Symbiotic Interaction Model

Example 2.2. Consider the following population model:


  
xn
xn+1 = xn + r1 xn 1 − ,


k1 + αyn
  (7)
yn
yn+1 = yn + r2 yn 1 −

 .
k2 + βxn

The presence of each species increases the eective carrying capacity of the other
(through the terms αyn and βxn in the denominators). This is a model of symbiosis
(mutualism).

3. Tutorial Problems

Exercise
Problem 1. Consider a hostparasitoid model given by
  −k
aPn
Hn+1 = f Hn 1 + ,


k (8)
H
Pn+1 = Hn − n+1 ,


f

where f , k , a are positive constants. Find the steady states and determine their
stabilities.

Exercise
Problem 2. Consider a model for leaf-eating herbivores whose population size at
stage n is Hn , living on a tree with leaf mass Vn :

k(1−Hn )
Vn+1 = Vn e
 ,

1 Hn
 k, b > 0. (9)
Hn+1 = bHn 1 + −
 ,
b Vn

(a) What are the possible interpretations of the parameters k and b?

(b) Find the steady states.

(c) Determine the stability of the steady states.

7
Population Dynamics Weeks 1011

Exercise
Problem 3. Consider the model

xn+1 = α xn ynβ ,

σyn α, β, σ, δ > 0. (10)
yn+1 =
 ,
δ + xn

(a) What kind of interactions are depicted by this model?

(b) Find the steady states and the conditions under which they exist.
xn 1
(c) By dening Xn = and Yn = α β yn , show that model (10) transforms to
σ−δ

β
Xn+1 = Xn Yn ,

Yn (∗)
Yn+1 =
 ,
ε + Xn (1 − ε)

where ε is to be determined.

(d) Find the steady states of (∗).

Exercise
Problem 4. Consider the model
   
Nt

Nt+1 = Nt exp r 1 −
 − aPt ,
k r, k, a > 0. (11)

Pt+1 = Nt 1 − e−aPt ,
 

(a) Describe the type of interaction modelled, explaining your reasoning.

(b) Find the trivial steady state and determine its stability.

8
Practical 1 - Discrete model
Design and simulating a model for Flu
Objective: The main objective is to develop a flu model taking into account the epidemiology
of the disease and fit the designed model to data in Excel.

Model
Influenza or flu is a viral infection that is readily transmitted through the air and
causes respiratory problems in humans and other animals. It occurs seasonally
and may result in deaths among the young and old.

Important features of Flu that we would want to capture in a model are:


n Influenza (flu)
n Is a contagious respiratory illness.
n Caused by influenza viruses (3 types).
n The illness can be mild, severe and sometimes fatal.
n Not really related to common cold.
n Transmission process
n Infectious individuals cough and sneeze viruses in the
air.
n Viruses enter the human body through the nose and
mouth.
n If infected you will get the symptoms 1-4 days later.
n One can spread the illness before the onset of symptoms
and 3-4 days after the onset of symptoms.
n Recovery begin 2-3 days after becoming infectious.
n Overall recovery takes at least 7 days.
n Recovery is with immunity.
n Can be treated with antiviral drugs
n Prevention is usually by flue vaccines administered on a
yearly basis.
Symptoms: include headaches, coughing, fever, body and muscle aches, soar throat
and chilles.
n The transmission process is key
n We have people at risk.
n If we introduce a single infective, there is a chance of
new infections being generated.
n Once infection takes place, symptoms manifest1-4 days
later, meaning there is a short period of latency.
n Once symptoms appear, individuals become highly
infectious.
n After a few days individuals recover.
n Recovery is with immunity for that particular strain.
THE MODEL
n Consider a closed system, boarding school, prison, refugee camps, Isolation camps.
n The time frame is usually a few weeks to a year at national level for an influenza
epidemic.
n No births, deaths, immigration and emigrations are assumed.
n Individuals are grouped according to their state of infection, creating classes or
compartments.
n Each compartment contains individuals of the same characteristic, in this case, the
same stage of infection.

Draw the flow diagram of the model that takes into account the described transmission
process. SUGGEST YOUR OWN MODEL. I suggest you draw the model in power point
and paste it here.
Attach parameters to the flow diagram and then work out the DISCRETE algebraic
expressions for the word equations above. Rewrite the equations in algebraic form in the
table below:

List all the major assumptions made in this model (At least 5)

1……………………………………………………………………………………………..

2………………………………………………………………………………………………

3……………………………………………………………………………………………….

4……………………………………………………………………………………………….

5………………………………………………………………………………………………
Model Simulation

Suppose the proposed model is of the `SEIR’ type. The equation in R is often called a
redundant equation, since those in R are not involved in the transmission of the disease.

Note that 𝑁 = 𝑆 + 𝐸 + 𝐼 + 𝑅 and the equation in R can always be obtained from

𝑅 = 𝑁 − 𝑆 − 𝐸 − 𝐼. In the simulation below we shall consider only three state variables 𝑆, 𝐸


and 𝐼. In the box below write down the set of discrete equations, using standard incidence
function.

We want to simulate the SEIR model and try to fit it to data, the 1978 English boarding school
influenza outbreak. The data is for the daily number of children confined to bed in a boarding school
in North England during the outbreak in 1978 of the reemerging A/H1N1 strain. The school had 763
boys of which 512 boys were confined to bed sometime during the outbreak. The data for the
confirmed cases is shown as follows:
Assume that the control population consists of N = 763, which is based on the number of students.
Taking I0 = 3, (E0 = 4 assumed) and S0 = 756 = N -E0 - [Link] would like to simulate the model and
determine the best values of the parameters, that match the data from the model. These parameters are
key to understanding the nature of a disease.

Epidemiologists often examine what is called the basic reproduction ratio given by

R0 = (duration spent in E)x(infection rate)x(duration of infectivity)

which provides a measure of how rapidly a disease will spread and how much of the population will
be affected by a particular disease.

We simulate the model, finding the Least Sum of Square Errors (LSSE) by changing the
parameters, infection rate, progression rate and recovery rate. The SSE is found from the difference
between the data on the infecteds and the model simulation of In.

Instructions

1. In column A1 type Day and B1 type Exposed and C1 type Infected


2. Fill the columns A and C with the data in the table above
3. In column H1 type N (the total number of people)
4. In column H2 type S0 (the initial number of susceptible people)
5. In column H3 type E0 (the initial number of susceptible people)
6. In column H4 type I0 (the initial number of infected people)
7. In column H5 type (name of parameter)(the infection rate)
8. In column H6 type (name of parameter)(the progression rate)
9. In column H7 type (name of parameter)( (the recovery rate)
10. In column I1 type 763 (size of N)
11. In column I2 type =H1-H3-H4 (size of S0)
12. In column I3 type 4 (the initial number of exposed people)
13. In column I4 type 3 (the initial number of infected people)
14. In column I5 type 3.99 (the infection rate)
15. In column I6 type 0,72 (the progression rate)
16. In column I7 type 3.5 (the recovery rate)
17. In column K4 type R0 (the reproduction number)
18. In column L4 type =I5/(I6*I7) (the reproduction number)
19. In column D1, type S(n)
20. In column E1, type I(n)
21. In column F1, type Error
22. In column B2, type =I$3
23. In column D2, type =I$2
24. In column E2, type =I$4
25. In column F2, type =(C2-E2)^2
26. In column D3, type =D2-(I$5/I$1)*D2*E2
27. In column E3, type =(1-I$7)*E2+I$6 *B2
28. In column B3, type =(1-I$5)*B2+(I$5/I$1)*D2*E2
29. In column F3, type = (C3-E3)^2
30. Run down the excel sheet values for columns B, D, E and F to the last day.
31. In column F16 provide the sum of all the errors (=SUM(F2:F15))
32. In Column G, provide the number of the recovered
33. Plot columns A, C and E, on a scatter plot.
34. Play around with the parameters for beta and gamma and so that the sum of the errors
is as small as possible (Hint: you may need to increase the decimal points on beta and
gamma).
Dynamical Systems – Discrete Maps in
1D
Cobweb Diagrams: Tutorial & MATLAB Laboratory
Ten Worked Exercises with Complete Code

1 Background: Discrete Dynamical Systems


A one-dimensional discrete dynamical system is defined by an iteration rule of the form
x_{n+1} = f(x_n), n = 0, 1, 2, ...
where f : I → I is a continuous map on an interval I ⊆ ℝ. Starting from an initial condition x₀ ∈ I, we
obtain the orbit (trajectory):
x₀, f(x₀), f²(x₀), f³(x₀), ...
The long-run behaviour of this orbit is the central question of the theory: does it settle to a fixed point, a
periodic orbit, or exhibit sensitive dependence on initial conditions (chaos)?

1.1 Key Definitions


Fixed point: A point x* satisfying f(x*) = x*. Stability is determined by |f'(x*)| < 1 (stable) or |f'(x*)| > 1
(unstable).
Period-p orbit: A set {x₀, x₁, ..., x_{p-1}} with f^p(x₀) = x₀ and p the smallest such integer.
Chaos: Sensitive dependence on initial conditions together with topological mixing. For the logistic
map this occurs (roughly) for r > 3.57.

Stability criterion for fixed points


If x* is a fixed point of f, then:
• |fʹ(x*)| < 1 ⇒ x* is attracting (stable)
• |fʹ(x*)| > 1 ⇒ x* is repelling (unstable)
• |fʹ(x*)| = 1 ⇒ marginal case; higher derivatives needed

2 The Cobweb Construction


The cobweb (or staircase) diagram is a graphical tool that simultaneously displays:
1. The graph of y = f(x)
2. The diagonal y = x
3. The orbit x₀, x₁, x₂, ... encoded as an alternating vertical/horizontal path

Dynamical Systems – Cobwebbing Tutorial | Page 1


2.1 Algorithm
Begin at the point (x₀, 0) on the x-axis. Each iteration consists of two moves:
4. Vertical move: draw a vertical segment from (x_n, x_n) up (or down) to (x_n, f(x_n)) — this
reads off the new iterate x_{n+1} = f(x_n) on the curve.
5. Horizontal move: draw a horizontal segment from (x_n, x_{n+1}) to (x_{n+1}, x_{n+1}) — this
transfers the new value to the diagonal, ready for the next step.

Fixed points appear as intersections of y = f(x) with y = x. An attracting fixed point "pulls" the cobweb
inward (a spiral or staircase that tightens); a repelling one pushes it outward. Periodic orbits produce
closed polygonal paths that cycle repeatedly through p segments.

Reading a cobweb diagram


Converging spiral/staircase inward → attracting fixed point
Diverging spiral/staircase outward → repelling fixed point (orbit escapes)
Closed box with 2 corners on curve → stable 2-cycle
Closed box with 4 corners on curve → stable 4-cycle
Irregular, non-repeating path → chaotic orbit

Dynamical Systems – Cobwebbing Tutorial | Page 2


3 MATLAB Implementation
3.1 Main function: logistic_cobweb_and_bifurcation.m
Save the following code as cobweb_diagram.m in your working directory (or on the MATLAB path).

function logistic_cobweb_and_bifurcation()
% LOGISTIC_COBWEB_AND_BIFURCATION
% Runs a cobweb diagram (for a chosen r and x0) AND a bifurcation diagram
% for the logistic map x_{n+1} = r x_n (1 - x_n), showing chaos.
%
% Just run:
% logistic_cobweb_and_bifurcation

clc; close all;

% ---------------- USER SETTINGS ----------------


r_cobweb = 3.2; % try 3.2 (periodic), 3.5, 3.7, 3.9 (chaotic)
x0 = 0.2; % initial condition for cobweb
N_cobweb = 80; % iterations for cobweb
x_range = [0, 1]; % plotting range

r_min = 2.5; % bifurcation range start


r_max = 3.5; % bifurcation range end
nr = 4000; % number of r samples (increase for sharper plot)
burn_in = 1000; % transient iterations to discard
keep = 200; % points to keep per r (increase to see structure)

% Logistic map (vector-safe)


f = @(x, r) r .* x .* (1 - x);

% ---------------- 1) COBWEB DIAGRAM ----------------


cobweb_diagram_internal(f, r_cobweb, x0, N_cobweb, x_range);

% ---------------- 2) BIFURCATION DIAGRAM ----------------


bifurcation_diagram_logistic(f, r_min, r_max, nr, burn_in, keep);
end

% =====================================================================
% COBWEB (internal)
% =====================================================================
function cobweb_diagram_internal(f, r, x0, N, x_range)

% Orbit
x = zeros(1, N+1);
x(1) = x0;
for k = 1:N
x(k+1) = f(x(k), r);
if ~isfinite(x(k+1))
error('Non-finite orbit at k=%d.', k);
end
end

% Cobweb path
cx = zeros(1, 2*N + 1);
cy = zeros(1, 2*N + 1);

Dynamical Systems – Cobwebbing Tutorial | Page 3


cx(1) = x0; cy(1) = x0;

for k = 1:N
cx(2*k) = x(k); cy(2*k) = x(k+1);
cx(2*k+1) = x(k+1); cy(2*k+1) = x(k+1);
end

% Plot
figure('Color','w','Name','Cobweb (Logistic Map)');
xp = linspace(x_range(1), x_range(2), 1000);
yp = f(xp, r);

plot(xp, yp, 'LineWidth', 2); hold on;


plot(x_range, x_range, '--', 'LineWidth', 1);
plot(cx, cy, 'LineWidth', 1.2);
plot(x0, x0, 'o', 'MarkerSize', 7, 'MarkerFaceColor', 'auto');

xlabel('x_n', 'FontSize', 13);


ylabel('x_{n+1}', 'FontSize', 13);
title(sprintf('Cobweb: Logistic map x_{n+1}=r x_n(1-x_n), r=%.4f, x_0=%.4f',
r, x0), ...
'FontSize', 12);

legend({'f(x)', 'y=x', 'cobweb', 'x_0'}, 'Location', 'best');


axis([x_range(1) x_range(2) x_range(1) x_range(2)]);
grid on; box on; hold off;
end

% =====================================================================
% BIFURCATION DIAGRAM (Logistic Map)
% =====================================================================
function bifurcation_diagram_logistic(f, r_min, r_max, nr, burn_in, keep)

r_vals = linspace(r_min, r_max, nr);

% Start all r's with same x0 (can be random too)


x = 0.5 * ones(1, nr);

% Burn-in (discard transient)


for i = 1:burn_in
x = f(x, r_vals);
end

% Collect last "keep" iterates for each r


R_plot = zeros(keep, nr);
X_plot = zeros(keep, nr);

for i = 1:keep
x = f(x, r_vals);
R_plot(i, :) = r_vals;
X_plot(i, :) = x;
end

% Plot (scatter-like using dots)


figure('Color','w','Name','Bifurcation Diagram (Logistic Map)');
plot(R_plot(:), X_plot(:), '.', 'MarkerSize', 1);

xlabel('r', 'FontSize', 13);

Dynamical Systems – Cobwebbing Tutorial | Page 4


ylabel('Asymptotic x_n', 'FontSize', 13);
title(sprintf('Bifurcation diagram: Logistic map (burn-in=%d, keep=%d, nr=%d)',
...
burn_in, keep, nr), 'FontSize', 12);

xlim([r_min r_max]);
ylim([0 1]);
grid on; box on;
end

3.2 Code disection


The function is structured in four clearly commented stages:
6. Build the orbit: a simple for-loop applies f repeatedly to produce the sequence x₀, x₁, ..., x_N.
7. Assemble the cobweb path: the alternating vertical/horizontal segments are packed into two
coordinate vectors cx, cy for a single plot call.
8. Plot: the curve of f, the diagonal, the cobweb path, and the initial condition are drawn with
distinct colours.
9. Labels & formatting: axis limits, grid, legend, and a descriptive title are applied.

Usage tip
The function accepts any anonymous function handle for f.
Example: cobweb_diagram(@(x) 3.5*x.*(1-x), 0.1, 100, [0 1], 'Logistic r=3.5')
Use vectorised operations (dot operators: .*, ./, .^) so that f can be evaluated on the plotting grid xp as
well as on scalars.

Dynamical Systems – Cobwebbing Tutorial | Page 5


4 Ten Exercises
The table below summarises all ten exercises. Detailed instructions follow. The demo script (Section
5) runs all exercises automatically.

# Map f(x) Parameter Initial x₀ Task


1 f(x) = rx(1−x) r = 2.5 x₀ = 0.2 Show convergence to fixed
point
2 f(x) = rx(1−x) r = 3.2 x₀ = 0.5 Identify period-2 orbit
3 f(x) = rx(1−x) r = 3.5 x₀ = 0.1 Identify period-4 orbit
4 f(x) = rx(1−x) r = 3.9 x₀ = 0.4 Observe chaotic behaviour
5 f(x) = x²+ c c = −0.5 x₀ = 0.3 Find stable 2-cycle
6 f(x) = μ − x² μ = 1.0 x₀ = 0.0 Stable fixed point analysis
7 f(x) = x·e^(r(1−x)) r = 1.5 x₀ = 0.5 Ricker map — fixed point
8 f(x) = x·e^(r(1−x)) r = 2.5 x₀ = 0.5 Ricker map — 2-cycle
9 f(x) = ax/(1+x²) a = 2.5 x₀ = 0.5 Convergence to fixed point
10 f(x) = sin(πx) — x₀ = 0.7 Stable orbit of sine map

Exercise 1 – Logistic Map: Stable Fixed Point (r = 2.5)


The logistic map f(x) = rx(1−x) models bounded population growth. For 1 < r < 3 every orbit converges
to the unique positive fixed point x* = 1 − 1/r.

cobweb_diagram(@(x) 2.5*x.*(1-x), 0.2, 60, [0 1], 'Logistic r=2.5');

• Verify analytically: fixed points satisfy rx*(1−x*) = x*. Find x* and confirm |fʹ(x*)| < 1.
• Describe the shape of the cobweb (monotone staircase or oscillating spiral).
• Repeat with x₀ = 0.9. Does the orbit converge to the same fixed point?

Exercise 2 – Logistic Map: Period-2 Orbit (r = 3.2)


For 3 < r < 1+√6 ≈ 3.449 the fixed point loses stability via a period-doubling bifurcation and the orbit
converges to a 2-cycle.

cobweb_diagram(@(x) 3.2*x.*(1-x), 0.5, 80, [0 1], 'Logistic r=3.2');

• Identify the two values p, q that the orbit alternates between.


• Observe the closed rectangular cobweb box. Sketch it and label the four corners.
• Compute f²(x) = f(f(x)) symbolically and verify that p and q are fixed points of f².

Dynamical Systems – Cobwebbing Tutorial | Page 6


Exercise 3 – Logistic Map: Period-4 Orbit (r = 3.5)
A second period-doubling occurs near r ≈ 3.449, giving a stable 4-cycle for r ∈ (3.449, 3.544).

cobweb_diagram(@(x) 3.5*x.*(1-x), 0.1, 80, [0 1], 'Logistic r=3.5');

• Count the distinct values visited. Confirm there are exactly four.
• Print x(end-8:end) in MATLAB to see the last nine iterates.
• How does the cobweb box differ visually from the period-2 case?

Exercise 4 – Logistic Map: Chaos (r = 3.9)


For r > r_∞ ≈ 3.5699 (the Feigenbaum accumulation point) the system exhibits chaotic behaviour for
most initial conditions.

cobweb_diagram(@(x) 3.9*x.*(1-x), 0.4, 100, [0 1], 'Logistic r=3.9');

• Describe the qualitative difference from Exercises 1–3.


• Run with x₀ = 0.40001 (a tiny perturbation). After 50 iterates, how different are the two orbits?
This illustrates sensitive dependence.
• Increase N to 200 and observe that the cobweb densely fills the curve.

Exercise 5 – Quadratic Map f(x) = x² − 0.5 (Stable 2-Cycle)


The family f_c(x) = x² + c (the real slice of the Mandelbrot family) exhibits period-doubling as c
decreases. At c = −0.5 there is a stable 2-cycle.

cobweb_diagram(@(x) x.^2 - 0.5, 0.3, 80, [-1.5 1.5], 'Quadratic c=-0.5');

• Solve x² − 0.5 = x for the fixed points. Determine their stability.


• The 2-cycle satisfies f(f(x)) = x but f(x) ≠ x. Find the two cycle points by solving f(f(x)) = x.
• Change the domain to [−2, 2] and experiment with c = −2. What happens to the orbit?

Exercise 6 – Map f(x) = μ − x² with μ = 1.0


This variant emphasises asymmetric cobweb geometry. Find the fixed points algebraically and
compare with the graphical result.

cobweb_diagram(@(x) 1.0 - x.^2, 0.0, 60, [-1.5 1.5], 'mu-x^2, mu=1.0');

Dynamical Systems – Cobwebbing Tutorial | Page 7


• Find both fixed points analytically (roots of μ − x² = x).
• Compute the derivative at each fixed point and classify as stable or unstable.
• Try μ = 1.9. Does the qualitative behaviour change? What type of orbit do you observe?

Exercise 7 – Ricker Map: Stable Fixed Point (r = 1.5)


The Ricker map f(x) = x·exp(r(1−x)) is widely used in ecology (salmon population dynamics). It has a
non-polynomial nonlinearity.

cobweb_diagram(@(x) x.*exp(1.5*(1-x)), 0.5, 80, [0 3], 'Ricker r=1.5');

• Show that x* = 1 is always a fixed point of the Ricker map for any r.
• Compute fʹ(1) = 1−r. For r = 1.5, is the fixed point stable?
• What is the critical value of r at which x* = 1 loses stability?

Exercise 8 – Ricker Map: Period-2 Orbit (r = 2.5)


Above r = 2 the fixed point x* = 1 of the Ricker map undergoes a period-doubling bifurcation.

cobweb_diagram(@(x) x.*exp(2.5*(1-x)), 0.5, 80, [0 3.5], 'Ricker r=2.5');

• Since fʹ(1) = 1−r = −1.5, confirm the instability of x* = 1.


• Use the cobweb to estimate the two cycle values numerically.
• Compare the Ricker and logistic cobwebs for a 2-cycle — note the different curve shapes.

Exercise 9 – Rational Map f(x) = ax/(1+x²) with a = 2.5


This map, related to models of neural firing and optics, has a more complex graph shape with a
maximum interior to [0, ∞).

cobweb_diagram(@(x) 2.5*x./(1+x.^2), 0.5, 80, [0 2.5], 'ax/(1+x^2) a=2.5');

• Find all fixed points by solving ax/(1+x²) = x. One solution is x = 0 — find the others.
• Classify each fixed point by computing fʹ(x*).
• Try x₀ = 2.0 (close to the unstable fixed point). How does the orbit behave?

Exercise 10 – Sine Map f(x) = sin(πx)


The sine map on [0,1] shares topological properties with the logistic map. It is conjugate to the logistic
map with r = 4 via a change of variables.

Dynamical Systems – Cobwebbing Tutorial | Page 8


cobweb_diagram(@(x) sin(pi*x), 0.7, 80, [0 1], 'Sine map sin(pi x)');

• Find the nonzero fixed point x* of sin(πx) = x numerically (use MATLAB’s fzero or Newton’s
method).
• Compute |fʹ(x*)| = |π cos(πx*)| and determine its stability.
• The topological conjugacy h(x) = sin²(πx/2) maps the logistic r=4 orbit to the sine map orbit.
Verify this numerically for x₀ = 0.3.

End of Tutorial

Dynamical Systems – Cobwebbing Tutorial | Page 9

You might also like