Pitchfork bifurcations

In this lecture, we're going to take a look at pitchfork bifurcations. In this kind of bifurcation, changing a parameter of the problem leads to the formation or destruction of fixed points. We will begin with the equation $\dot{x} = rx-x^3$.

In [1]:
import numpy as np;
import matplotlib.pyplot as plt;
plt.rcParams.update({"text.usetex":True});
%config InlineBackend.figure_format = "svg"
from ipywidgets import interactive
from scipy.optimize import fsolve
In [2]:
def f(x, r):
    return r*x- x**3;

def dfdx(x, r):
    return r-3*x**2;

x = np.linspace(-2, 2); r_a = np.linspace(-2, 2);

for r in r_a:
    sol = fsolve(f, [-2,0, 2], args=(r), full_output=True)
    if sol[2] == 1: # Solution has convered
        for i in np.arange(0, np.size(sol[0])):
            root = sol[0][i];
            slope_at_root = dfdx(root, r)
            if slope_at_root > 0:
                plt.plot(r, root, 'bx')
            else:
                plt.plot(r, root, 'ro')
ax = plt.gca(); ax.set_aspect(1);
plt.xlim(np.min(r_a), np.max(r_a))
plt.axhline(0); plt.axvline(0)
plt.xlabel("r");
plt.ylabel("x*");
plt.title("Supercritical pitchfork bifurcation");

As seen from the figure, as the parameter $r$ is changed we go from having one fixed point to three fixed points. The stability is asertained from the slope at the fixed point. The behaviour can be ascertained from plotting the RHS of the function, as shown below.

In [3]:
def show_fx(r):

    def f(x, r):
        return r*x- x**3;
    
    x = np.linspace(-2, 2);
    y = f(x, r);
    
    plt.plot(x, y);
    plt.axhline(0,color='k');
    
w = interactive(show_fx, r = (-2, 2, 0.1));
w
In [4]:
def show_fx(r):

    def f(x, r):
        return -x + r*np.tanh(x);
    
    x = np.linspace(-2, 2);
    y = f(x, r);
    
    plt.plot(x, y);
    plt.axhline(0,color='k');
    
w = interactive(show_fx, r = (-2, 2, 0.1));
w

In the snippet that follows, we aim at looking at the bifurcation diagram $\dot{x} = -x + r\tanh(x)$

In [5]:
def f(x, r):
    return -x + r*np.tanh(x);

def dfdx(x, r):
    return -1 + r*1/(np.cosh(x))**2;

x = np.linspace(-2, 2); r_a = np.linspace(-2, 2);

for r in r_a:
    sol = fsolve(f, [-2,0, 2], args=(r), full_output=True)
    if sol[2] == 1: # Solution has convered
        for i in np.arange(0, np.size(sol[0])):
            root = sol[0][i];
            slope_at_root = dfdx(root, r)
            if slope_at_root > 0:
                plt.plot(r, root, 'bx')
            else:
                plt.plot(r, root, 'ro')
ax = plt.gca(); ax.set_aspect(1);
plt.xlim(np.min(r_a), np.max(r_a))
plt.axhline(0); plt.axvline(0)
plt.xlabel("r");
plt.ylabel("x*");
plt.title("Supercritical pitchfork bifurcation");

The equation above has the same normal form as the equation we had worked with in the first example, albeit shifted. This results in the same kind of bifurcation behaviour. In what follows, we can explain with the help of an interactive plot, what the reason behind the observations of the bifurcation diagram are.

In [6]:
def show_fx(r):
    
    x = np.linspace(-2, 2);
    y = f(x, r);
    
    plt.plot(x, x);
    plt.plot(x, np.tanh(x)*r)
    plt.axhline(0,color='k');
    
w = interactive(show_fx, r = (-2, 2, 0.1));
w

Let's shift our attention to a different prototypical problem $\dot{x} = rx+x^3$.

In [7]:
def show_fx(r):

    def f(x, r):
        return r*x + x**3;
    
    x = np.linspace(-2, 2);
    y = f(x, r);
    
    plt.plot(x, y);
    plt.axhline(0,color='k');
    plt.title("$\\dot{x} = rx + x^3$")
    
w = interactive(show_fx, r = (-2, 2, 0.1));
w

What about the case where we can have many more roots, all with different stability characteristics? We see that they lead to a much more involved bifurcation diagram - one which shows hysterisis as we will show by solving the initial value problem and look at the equilibrium point (large time value of x)

In [8]:
def show_fx(r):

    def f(x, r):
        return r*x + x**3 - x**5;
    
    x = np.linspace(-2, 2, 200);
    y = f(x, r);
    
    plt.plot(x, y);
    plt.axhline(0,color='k');
    plt.title("$\\dot{x} = rx + x^3 - x^5$");
    plt.ylim(-0.25, 0.25)
    
w = interactive(show_fx, r = (-2, 2, 0.005));
w
In [9]:
def f(x, r):
    return r*x + x**3 - x**5;

def dfdx(x, r):
    return r + 3*x**2 - 5*x**4;

x = np.linspace(-2, 2); r_a = np.linspace(-0.5, 1,200);
xg = np.linspace(-2,2,20);

for r in r_a:
    sol = fsolve(f, xg, args=(r), full_output=True)
    if sol[2] == 1: # Solution has convered
        for i in np.arange(0, np.size(sol[0])):
            root = sol[0][i];
            slope_at_root = dfdx(root, r)
            if slope_at_root > 0:
                plt.plot(r, root, 'bx')
            else:
                plt.plot(r, root, 'ro')
ax = plt.gca(); ax.set_aspect(0.4);
plt.xlim(np.min(r_a), np.max(r_a))
plt.axhline(0); plt.axvline(0)
plt.xlabel("r");
plt.ylabel("x*");
plt.title("Subcritical pitchfork bifurcation");
In [10]:
from scipy.integrate import solve_ivp
In [11]:
def plot_subcritical(r=0.3, x0 = 0.2, t_end =30):
    
    def xd(t, x, r):
        return r*x + x**3 - x**5
    
    sol = solve_ivp(xd, [0, t_end], [x0], args=[r], rtol=1e-6);
    plt.plot(sol.t, sol.y.T);
    plt.xlabel("$t$"); plt.ylabel("x");
    plt.title("$x = %f$"%(sol.y.T[-1]))
    plt.ylim([-1.2, 1.2])
    
w = interactive(plot_subcritical,  r = (-0.5, 0.4, 0.01), x0 = (-0.1,0.5,0.1), t_end=(10, 50, 10));
w
In [ ]: