# BifurcationsΒΆ

The behavior of the solution of a dynamical system may change abruptly as a function of some control parameter. The most commonly observed transitions in dynamical states are bifurcations. This page discusses some of the most important classes of bifurcations on the example of a handful of simple equations, which can be viewed as archetypical, and to which many more complex dynamical systems reduce close to the transition point.

The following examples are taken - and slightly varied - from the excellent textbook Introduction to the Modeling and Analysis of Complex Systems of Hiroki Sayama, section 8.1., pp.131.

A bifurcation is a qualitative, topological change of a systemβs phase space that occurs when some parameters are slightly varied across a critical thresholds and the stability of an equilibrium point changes between stable and unstable. This is the case when the eigenvalues $\lambda_i$ of the Jacobian matrix at an equilibrium point satisfy the following:

• For discrete-time models: $\lambda_i = 1$ for some $i$, while $\lambda_i < 1$ for the rest.
• For continuous-time models: $Re(\lambda_i) = 0$ for some $i$, while $Re(\lambda_i) < 0$ for the rest.

Consider the example of a:

defined by $\frac{dx}{dt} = r-x^2$

First find the equilibrium points by letting $\frac{dx}{dt}=0$. This gives $x_{eq} = \pm \sqrt r$ which means that equilibrium points exist only for non-negative $r$. The critical condition for the occurence of a bifurcation is found by taking the derivative of $\frac{dx}{dt}=F(x)$:

$\frac{dF}{dx}=-2x$

and looking at it at the equilibrium point

$\frac{dF}{dx}\lvert_{x=x_{eq}}=\pm 2\sqrt r = 0$

Thus, $r = 0$ is when the bifurcation occurs.

Plugging each solution for $x_{eq} = \pm \sqrt r$ into $\frac{dF}{dx}=-2x$ shows that $x_{eq} = \sqrt r$ is stable while $x_{eq} = -\sqrt r$ is unstable.

This can be visually demonstated by drawing a bifurcation diagram in which the parameter used - in this case $r$ - is varied along the horizontal axis, while the location(s) of the equilibrium point(s) of the system are shown on the vertical axis. Bifurcation diagrams work only for systems with one variable and one parameter, but are still conceptually helpful in understanding the nature of bifurcations.

In [1]:
from numpy import *
import matplotlib.pyplot as plt
% matplotlib inline

# prepare plot
fig = plt.figure(figsize=(9,6))

# stable equilibrium
def xeq1(r):
return sqrt(r)

# unstable equilibrium
def xeq2(r):
return -sqrt(r)

domain = linspace(0, 10)
ax1.plot(domain, xeq1(domain), 'b-', label = 'stable equilibrium', linewidth = 3)
ax1.plot(domain, xeq2(domain), 'r--', label = 'unstable equilibrium', linewidth = 3)
ax1.legend(loc='upper left')
#neutral equilibrium point
ax1.plot([0], [0], 'go')
ax1.axis([-10, 10, -5, 5])
ax1.set_xlabel('r')
ax1.set_ylabel('x_eq')

#black arrows indicating the attracting dynamics of the stable and the repelling dynamics of the unstable equilibrium point:
ax1.annotate('', xy=(-7, -4), xytext=(-7, 4), arrowprops=dict(arrowstyle="->",connectionstyle="arc3",lw=1),)
ax1.annotate('', xy=(-5, -4), xytext=(-5, 4), arrowprops=dict(arrowstyle="->",connectionstyle="arc3",lw=1),)
ax1.annotate('', xy=(-3, -4), xytext=(-3, 4), arrowprops=dict(arrowstyle="->",connectionstyle="arc3",lw=1),)
ax1.annotate('', xy=(-1, -4), xytext=(-1, 4), arrowprops=dict(arrowstyle="->",connectionstyle="arc3",lw=1),)
ax1.annotate('', xy=(1, -4), xytext=(1, -1.5), arrowprops=dict(arrowstyle="->",connectionstyle="arc3",lw=1),)
ax1.annotate('', xy=(1, 0.7), xytext=(1, -0.7), arrowprops=dict(arrowstyle="->",connectionstyle="arc3",lw=1),)
ax1.annotate('', xy=(1, 1.5), xytext=(1, 4), arrowprops=dict(arrowstyle="->",connectionstyle="arc3",lw=1),)
ax1.annotate('', xy=(3, -4), xytext=(3, -2), arrowprops=dict(arrowstyle="->",connectionstyle="arc3",lw=1),)
ax1.annotate('', xy=(3, 1.5), xytext=(3, -1.5), arrowprops=dict(arrowstyle="->",connectionstyle="arc3",lw=1),)
ax1.annotate('', xy=(3, 2), xytext=(3, 4), arrowprops=dict(arrowstyle="->",connectionstyle="arc3",lw=1),)
ax1.annotate('', xy=(5, -4), xytext=(5, -2.5), arrowprops=dict(arrowstyle="->",connectionstyle="arc3",lw=1),)
ax1.annotate('', xy=(5, 2), xytext=(5, -2), arrowprops=dict(arrowstyle="->",connectionstyle="arc3",lw=1),)
ax1.annotate('', xy=(5, 2.5), xytext=(5, 4), arrowprops=dict(arrowstyle="->",connectionstyle="arc3",lw=1),)
ax1.annotate('', xy=(7, -4), xytext=(7, -3), arrowprops=dict(arrowstyle="->",connectionstyle="arc3",lw=1),)
ax1.annotate('', xy=(7, 2.3), xytext=(7, -2.3), arrowprops=dict(arrowstyle="->",connectionstyle="arc3",lw=1),)
ax1.annotate('', xy=(7, 3), xytext=(7, 4), arrowprops=dict(arrowstyle="->",connectionstyle="arc3",lw=1),)

Out[1]:
<matplotlib.text.Annotation at 0x7f04d10>

Each vertical slice of the bifurcation diagram (left plot) for a particular parameter value depicts a phase space of the dynamical system. For example, for $r = 5$ in the diagram above, there are two equilibrium points, one stable (blue/solid) and the other unstable (red/dashed). Flows of the systemβs state are visualized by downward and upward arrows.

### A transcritical bifurcationΒΆ

is a bifurcation where one equilibrium point βpasses throughβ another, exchanging their stabilities

The system $\frac{dx}{dt} = rx-x^2$ has the equilibrium points $x_{eq} = 0$ and $x_{eq} = r$, with the exception of $r = 0$ where equilibrium points collide and swap their stabilities.

In [2]:
# stable equilibrium > 0
def xeq1(r):
return r

# unstable equilibrium > 0
def xeq2(r):
return r - r  # = 0

# prepare plot
fig = plt.figure(figsize=(9,6))

domain1 = linspace(-10, 0)
domain2 = linspace(0, 10)
plt.plot(domain1, xeq1(domain1), 'r--', linewidth = 3)
plt.plot(domain1, xeq2(domain1), 'b-', linewidth = 3)
plt.plot(domain2, xeq1(domain2), 'b-', linewidth = 3)
plt.plot(domain2, xeq2(domain2), 'r--', linewidth = 3)
#neutral equilibrium point
plt.plot([0], [0], 'go')
plt.axis([-10, 10, -10, 10])
plt.xlabel('r')
plt.ylabel('x_eq')
plt.title('Transcritical bifurcation')

#black arrows:
plt.annotate('', xy=(0, -1), xytext=(0, -9), arrowprops=dict(arrowstyle="->",connectionstyle="arc3",lw=1),)
plt.annotate('', xy=(0, 1), xytext=(0, 9), arrowprops=dict(arrowstyle="->",connectionstyle="arc3",lw=1),)
plt.annotate('', xy=(-5, -1), xytext=(-5, -4), arrowprops=dict(arrowstyle="->",connectionstyle="arc3",lw=1),)
plt.annotate('', xy=(-5, 1), xytext=(-5, 9), arrowprops=dict(arrowstyle="->",connectionstyle="arc3",lw=1),)
plt.annotate('', xy=(-5, -9), xytext=(-5, -6), arrowprops=dict(arrowstyle="->",connectionstyle="arc3",lw=1),)
plt.annotate('', xy=(5, 4), xytext=(5, 1), arrowprops=dict(arrowstyle="->",connectionstyle="arc3",lw=1),)
plt.annotate('', xy=(5, -9), xytext=(5, -1), arrowprops=dict(arrowstyle="->",connectionstyle="arc3",lw=1),)
plt.annotate('', xy=(5, 6), xytext=(5, 9), arrowprops=dict(arrowstyle="->",connectionstyle="arc3",lw=1),)

Out[2]:
<matplotlib.text.Annotation at 0x80155d0>

## A pitchfork bifurcationΒΆ

is a bifurcation where an equilibrium point splits into three. Two of these (the outermost two) have the same stability as the original equilibrium point, while the one between them has a stability opposite to the original stability.

There are two types of pitchfork bifurcations.

### A supercritical pitchfork bifurcationΒΆ

makes a stable equilibrium point split into two stable and one unstable equilibrium.

For example, the system $\frac{dx}{dt} = rx-x^3$ has the three equilibrium points $x_{eq} = 0$ and $x_{eq} = \pm \sqrt r$, with the last two existing only for $r \ge 0$.

In [3]:
# 1 equilibrium
def xeq1(r):
return r - r      # = 0

# 2 equilibrium
def xeq2(r):
return sqrt(r)

# 3 equilibrium
def xeq3(r):
return -sqrt(r)

domain1 = linspace(-10, 0)
domain2 = linspace(0, 10)
plt.plot(domain1, xeq1(domain1), 'b-', linewidth = 3)
plt.plot(domain2, xeq1(domain2), 'r--', linewidth = 3)
plt.plot(domain2, xeq2(domain2), 'b-', linewidth = 3)
plt.plot(domain2, xeq3(domain2), 'b-', linewidth = 3)
#neutral equilibrium point
plt.plot([0], [0], 'go')
plt.axis([-10, 10, -5, 5])
plt.xlabel('r')
plt.ylabel('x_eq')
plt.title('Supercritical pitchfork bifurcation')

#black arrows:
plt.annotate('', xy=(0, -1), xytext=(0, -4), arrowprops=dict(arrowstyle="->",connectionstyle="arc3",lw=1),)
plt.annotate('', xy=(0, 1), xytext=(0, 4), arrowprops=dict(arrowstyle="->",connectionstyle="arc3",lw=1),)
plt.annotate('', xy=(-5, -0.5), xytext=(-5, -4), arrowprops=dict(arrowstyle="->",connectionstyle="arc3",lw=1),)
plt.annotate('', xy=(-5, 0.5), xytext=(-5, 4), arrowprops=dict(arrowstyle="->",connectionstyle="arc3",lw=1),)
plt.annotate('', xy=(3, 1.5), xytext=(3, 0.5), arrowprops=dict(arrowstyle="->",connectionstyle="arc3",lw=1),)
plt.annotate('', xy=(3, -1.5), xytext=(3, -0.5), arrowprops=dict(arrowstyle="->",connectionstyle="arc3",lw=1),)
plt.annotate('', xy=(3, 2.2), xytext=(3, 4), arrowprops=dict(arrowstyle="->",connectionstyle="arc3",lw=1),)
plt.annotate('', xy=(3, -2.2), xytext=(3, -4), arrowprops=dict(arrowstyle="->",connectionstyle="arc3",lw=1),)
plt.annotate('', xy=(7, 2), xytext=(7, 0.5), arrowprops=dict(arrowstyle="->",connectionstyle="arc3",lw=1),)
plt.annotate('', xy=(7, -2), xytext=(7, -0.5), arrowprops=dict(arrowstyle="->",connectionstyle="arc3",lw=1),)
plt.annotate('', xy=(7, 3), xytext=(7, 4), arrowprops=dict(arrowstyle="->",connectionstyle="arc3",lw=1),)
plt.annotate('', xy=(7, -3), xytext=(7, -4), arrowprops=dict(arrowstyle="->",connectionstyle="arc3",lw=1),)

Out[3]:
<matplotlib.text.Annotation at 0x92bc7f0>

### A subcritical pitchfork bifurcationΒΆ

makes an unstable equilibrium point split into two unstable and one stable equilibrium.

The system $\frac{dx}{dt} = rx+x^3$ has the three equilibrium points $x_{eq} = 0$ and $x_{eq} = \pm \sqrt -r$, with the last two existing only for $r \le 0$.

In [4]:
# 1 equilibrium
def xeq1(r):
return r - r      # = 0

# 2 equilibrium
def xeq2(r):
return sqrt(-r)

# 3 equilibrium
def xeq3(r):
return -sqrt(-r)

domain1 = linspace(-10, 0)
domain2 = linspace(0, 10)
plt.plot(domain1, xeq1(domain1), 'b-', linewidth = 3)
plt.plot(domain1, xeq2(domain1), 'r--', linewidth = 3)
plt.plot(domain1, xeq3(domain1), 'r--', linewidth = 3)
plt.plot(domain2, xeq1(domain2), 'r--', linewidth = 3)
#neutral equilibrium point
plt.plot([0], [0], 'go')
plt.axis([-10, 10, -5, 5])
plt.xlabel('r')
plt.ylabel('x_eq')
plt.title('Subcritical pitchfork bifurcation')

#black arrows:
plt.annotate('', xy=(1, -4), xytext=(1, -1), arrowprops=dict(arrowstyle="->",connectionstyle="arc3",lw=1),)
plt.annotate('', xy=(1, 4), xytext=(1, 1), arrowprops=dict(arrowstyle="->",connectionstyle="arc3",lw=1),)
plt.annotate('', xy=(5, -4), xytext=(5, -0.5), arrowprops=dict(arrowstyle="->",connectionstyle="arc3",lw=1),)
plt.annotate('', xy=(5, 4), xytext=(5, 0.5), arrowprops=dict(arrowstyle="->",connectionstyle="arc3",lw=1),)
plt.annotate('', xy=(-3, 0.5), xytext=(-3, 1.5), arrowprops=dict(arrowstyle="->",connectionstyle="arc3",lw=1),)
plt.annotate('', xy=(-3, -0.5), xytext=(-3, -1.5), arrowprops=dict(arrowstyle="->",connectionstyle="arc3",lw=1),)
plt.annotate('', xy=(-3, 4), xytext=(-3, 2.2), arrowprops=dict(arrowstyle="->",connectionstyle="arc3",lw=1),)
plt.annotate('', xy=(-3, -4), xytext=(-3, -2.2), arrowprops=dict(arrowstyle="->",connectionstyle="arc3",lw=1),)
plt.annotate('', xy=(-7, 0.5), xytext=(-7, 2), arrowprops=dict(arrowstyle="->",connectionstyle="arc3",lw=1),)
plt.annotate('', xy=(-7, -0.5), xytext=(-7, -2), arrowprops=dict(arrowstyle="->",connectionstyle="arc3",lw=1),)
plt.annotate('', xy=(-7, 4), xytext=(-7, 3), arrowprops=dict(arrowstyle="->",connectionstyle="arc3",lw=1),)
plt.annotate('', xy=(-7, -4), xytext=(-7, -3), arrowprops=dict(arrowstyle="->",connectionstyle="arc3",lw=1),)

Out[4]:
<matplotlib.text.Annotation at 0x941cef0>

## Combined bifurcationsΒΆ

The system $\frac{dx}{dt} = r+x-x^3$ has three equilibrium points, which are rather complicated to calculate. Solving $\frac{dx}{dt} = 0$ in terms of $r$ however, gives $r=-x+x^3$, which is enough for drawing a bifurcation diagram.

Information on stability can also be otained from the Jacobian matrix, which in this case is a $1 x 1$ matrix whose eigenvalue is its content itself (because it is a scalar). The critical condition at which a bifurcation occurs is given by $RE(\frac{dF}{dx})\mid_{x=x_{eq}} = 1-3x^2$, which indicates that equilibrium points are stable when $x^2>\frac{1}{3}$, otherwise they are unstable.

In [5]:
# 1 equilibrium
def xeq1(r):
return -r + r**3

domain1 = linspace(-1.3, -sqrt(1/3.))
domain2 = linspace(-sqrt(1/3.), sqrt(1/3.))
domain3 = linspace(sqrt(1/3.), 1.3)
plt.plot(xeq1(domain1), domain1, 'b-', linewidth = 3)
plt.plot(xeq1(domain2), domain2, 'r--', linewidth = 3)
plt.plot(xeq1(domain3), domain3, 'b-', linewidth = 3)
plt.axis([-1, 1, -1.5, 1.5])
plt.xlabel('r')
plt.ylabel('x_eq')

plt.annotate('', xy=(0.75, 1.2), xytext=(0.75, -1.4), arrowprops=dict(arrowstyle="->",connectionstyle="arc3",lw=1),)
plt.annotate('', xy=(0.5, 1.1), xytext=(0.5, -1.4), arrowprops=dict(arrowstyle="->",connectionstyle="arc3",lw=1),)
plt.annotate('', xy=(0.5, 1.25), xytext=(0.5, 1.4), arrowprops=dict(arrowstyle="->",connectionstyle="arc3",lw=1),)
plt.annotate('', xy=(0.25, -0.9), xytext=(0.25, -1.4), arrowprops=dict(arrowstyle="->",connectionstyle="arc3",lw=1),)
plt.annotate('', xy=(0.25, -0.8), xytext=(0.25, -0.3), arrowprops=dict(arrowstyle="->",connectionstyle="arc3",lw=1),)
plt.annotate('', xy=(0.25, 1), xytext=(0.25, -0.1), arrowprops=dict(arrowstyle="->",connectionstyle="arc3",lw=1),)
plt.annotate('', xy=(0.25, 1.15), xytext=(0.25, 1.4), arrowprops=dict(arrowstyle="->",connectionstyle="arc3",lw=1),)
plt.annotate('', xy=(0, -1.05), xytext=(0, -1.4), arrowprops=dict(arrowstyle="->",connectionstyle="arc3",lw=1),)
plt.annotate('', xy=(0, -0.9), xytext=(0, -0.1), arrowprops=dict(arrowstyle="->",connectionstyle="arc3",lw=1),)
plt.annotate('', xy=(0, 0.9), xytext=(0, 0.1), arrowprops=dict(arrowstyle="->",connectionstyle="arc3",lw=1),)
plt.annotate('', xy=(0, 1.05), xytext=(0, 1.4), arrowprops=dict(arrowstyle="->",connectionstyle="arc3",lw=1),)
plt.annotate('', xy=(-0.75, -1.2), xytext=(-0.75, 1.4), arrowprops=dict(arrowstyle="->",connectionstyle="arc3",lw=1),)
plt.annotate('', xy=(-0.5, -1.1), xytext=(-0.5, 1.4), arrowprops=dict(arrowstyle="->",connectionstyle="arc3",lw=1),)
plt.annotate('', xy=(-0.5, -1.25), xytext=(-0.5, -1.4), arrowprops=dict(arrowstyle="->",connectionstyle="arc3",lw=1),)
plt.annotate('', xy=(-0.25, 0.9), xytext=(-0.25, 1.4), arrowprops=dict(arrowstyle="->",connectionstyle="arc3",lw=1),)
plt.annotate('', xy=(-0.25, 0.8), xytext=(-0.25, 0.3), arrowprops=dict(arrowstyle="->",connectionstyle="arc3",lw=1),)
plt.annotate('', xy=(-0.25, -1), xytext=(-0.25, 0.1), arrowprops=dict(arrowstyle="->",connectionstyle="arc3",lw=1),)
plt.annotate('', xy=(-0.25, -1.15), xytext=(-0.25, -1.4), arrowprops=dict(arrowstyle="->",connectionstyle="arc3",lw=1),)

Out[5]:
<matplotlib.text.Annotation at 0x948c170>

Obviously, this system shows hysteresis, that is, a dependence not only on the input parameter $r$ in this case, but also on the system's history.

## Hopf bifurcationΒΆ

A Hopf bifurcation - named after Mathematician Eberhard Hopf - causes the appearance (or disappearance) of a limit cycle around the equilibrium point, a cyclic, closed trajectory in the phase space. An example for a Hopf bifurcation gives the dynamical model of a nonlinear oscillator called the van der Pol oscillator.

This model is defined by the second-order differential equation $\frac{d^2 x}{dt^2} + r(x^2-1)\frac{dx}{dt}+x=0$.

Introducing the additional variable $y = \frac{dx}{dt}$ allows to make the system first-order

$\frac{dx}{dt} = y$

$\frac{dy}{dt} = -r(x^2-1)y-x$

$(0, 0)$ is the only equilibrium point of this system, with the Jacobian matrix given as follows:

$\left(\begin{matrix} 0 & 1 \ -2rxy-1 & -r(x^2-1)\end{matrix}\right)\mid_{(x,y)=(0,0)} = \left(\begin{matrix} 0 & 1 \ -1 & r \end{matrix}\right)$

Eigenvalues of this matrix can be calculated as follows:

$\left\lvert\begin{matrix} 0-\lambda & 1 \ -1 & r-\lambda \end{matrix}\right\lvert = 0$

$\ -\lambda (r-\lambda )+1=\lambda^2-r\lambda + 1=0$

$\lambda = \frac{r \pm\sqrt{r^2 - 4}}{2}$

The critical condition for a bifurcation to occur thus is $Re(\lambda) = 0$

You can check whether the bifurcation is Hopf or not by looking at the imaginary components of the dominant eigenvalues of which the real parts are at a critical value (zero); if there are non-zero imaginary components, it must be a Hopf bifurcation.

The following code iterates through five values of $r$ and shows the appearance of the limit cycle at $r = 0$.

In [6]:
dt = 0.01

# prepare plots
fig = plt.figure(figsize=(18,6))

def plot_phase_space():
x = y = 0.1
xresult = [x]
yresult = [y]
for t in range(10000):
nextx = x + y * dt
nexty = y + (-r * (x**2 - 1) * y - x) * dt
x, y = nextx, nexty
xresult.append(x)
yresult.append(y)
plt.plot(xresult, yresult)
plt.axis('image')
plt.axis([-3, 3, -3, 3])
plt.title('r = ' + str(r))

rs = [-1, -0.1, 0, .1, 1]
for i in range(len(rs)):