Mini-Lecture 5-3 - Exponential growth in continuous time and discrete time.ipynb

Kernel: SageMath 9.3

Recall exponential growth/decay in continuous time:
Remember, continuous time model means differential equation, the things we've been studying since the beginning of LS 30A.

Differential equation: $X'(t) = r X(t)$
In other words, a constant percentage (per-capita, per-mass, per-molecule) growth rate

Solution of that differential equation: $X(t) = X(0) e^{rt}$
where $X(0) =$ the initial value of $X$ (at time $t = 0$).

Behavior:

• If $r > 0$: Exponential growth

• If $r < 0$: Exponential decay

(If $r = 0$, then $X(t) =$ constant for all time.)

choices = {"Actual solution": 1, "Simulated using desolve_odeint": 2, "Both": 3}

@interact(r=slider(-0.3, 0.31, 0.1, default=0.1, label="$r$"), 
          which=selector(choices, label="Which:"))
def continuous_exponential(r, which):
    p = Graphics()
    tmax = 20
    X0 = 5
    if which // 2:
        X = var("X")
        vectorfield(X) = r*X
        t_range = srange(0, tmax + 0.01, 0.1)
        solution = desolve_odeint(vectorfield, X0, t_range, X)
        p += list_plot(list(zip(t_range, solution)), plotjoined=True, color="blue", thickness=3)
    if which % 2:
        X(t) = X0 * e^(r*t)
        p += plot(X(t), (t, 0, tmax), color="red")
    p.show(ymin=0, ymax=100, axes_labels=("$t$", "$X$"))

Exponential growth/decay in discrete time:

Discrete time model: $X_{t+1} = r X_t$
(a.k.a. “difference equation”)

Start with $X_0$.

The difference equation then says $X_1 = r X_0$.

For the next step, it says $X_2 = r X_1$...
                             ... $= r \cdot (r X_0) = r^2 X_0$.

For the next step, it says $X_3 = r X_2$...
                             ... $= r \cdot (r^2 X_0) = r^3 X_0$.

For the next step, it says $X_4 = r X_3$...
                             ... $= r \cdot (r^3 X_0) = r^4 X_0$.

Conclusion: $X_t = r^t X_0$.

Exponential growth/decay in discrete time:

Discrete time model: $X_{t+1} = r X_t$
(a.k.a. “difference equation”)

Solution of that difference equation: $X_t = X_0 r^t$
where $X_0 =$ the initial value of $X$ (at time $t = 0$).

Behavior:

• If $r > 1$: Exponential growth

• If $r < 1$: Exponential decay

(If $r = 1$, then $X(t) =$ constant for all time.)

choices = {"Actual solution": 1, "Simulated using iteration": 2, "Both": 3}

@interact(r=slider(0.7, 1.3, 0.1, default=1.1, label="$r$"), 
          which=selector(choices, label="Which:"))
def discrete_exponential(r, which):
    p = Graphics()
    tmax = 20
    X0 = 5
    if which // 2:
        f(X) = r*X
        solution = [X0]
        for t in srange(20):
            solution.append(f(solution[-1]))
        p += list_plot(solution, color="blue", size=30)
    if which % 2:
        solution = [X0 * r^t for t in srange(0, 21)]
        p += list_plot(solution, color="red", size=20)
    p.show(ymin=0, ymax=100, axes_labels=("$t$", "$X$"))

What if $X(0)$ is negative?

choices = {"Actual solution": 1, "Simulated using desolve_odeint": 2, "Both": 3}

@interact(r=slider(-0.3, 0.31, 0.1, default=0.1, label="$r$"), 
          which=selector(choices, label="Which:"))
def continuous_exponential(r, which):
    p = Graphics()
    tmax = 20
    X0 = -5
    if which // 2:
        X = var("X")
        vectorfield(X) = r*X
        t_range = srange(0, tmax + 0.01, 0.1)
        solution = desolve_odeint(vectorfield, X0, t_range, X)
        p += list_plot(list(zip(t_range, solution)), plotjoined=True, color="blue", thickness=3)
    if which % 2:
        X(t) = X0 * e^(r*t)
        p += plot(X(t), (t, 0, tmax), color="red")
    p.show(ymin=-100, ymax=0, axes_labels=("$t$", "$X$"))

The same thing, in discrete time:

choices = {"Actual solution": 1, "Simulated using iteration": 2, "Both": 3}

@interact(r=slider(0.7, 1.3, 0.1, default=1.1, label="$r$"), 
          which=selector(choices, label="Which:"))
def discrete_exponential(r, which):
    p = Graphics()
    tmax = 20
    X0 = -5
    if which // 2:
        f(X) = r*X
        solution = [X0]
        for t in srange(20):
            solution.append(f(solution[-1]))
        p += list_plot(solution, color="blue", size=30)
    if which % 2:
        solution = [X0 * r^t for t in srange(0, 21)]
        p += list_plot(solution, color="red", size=20)
    p.show(ymin=-100, ymax=0, axes_labels=("$t$", "$X$"))

Question: What are the equilibrium points of the exponential growth/decay model $X' = rX$?

Answer: The only equilibrium point is at $X = 0$.

Now, is it stable or unstable?

@interact(r=slider(0.2, 1, 0.2, default=0.2, label="$r$"))
def continuous_growth(r):
    initial_states = [0, 1, 5, 10, -1, -5, -10]
    p = Graphics()
    for X0 in initial_states:
        f(t) = X0*e^(r*t)
        color = "purple" if X0 == 0 else "red"
        p += plot(f(t), t, 0, 5, thickness=2, color=color)
    p.show(axes_labels=("$t$", "$X$"), ymin=-100, ymax=100)

@interact(r=slider(-1, -0.199, 0.2, default=-0.2, label="$r$"))
def continuous_decay(r):
    initial_states = [0, 10, 50, 100, -10, -50, -100]
    p = Graphics()
    for X0 in initial_states:
        f(t) = X0*e^(r*t)
        color = "purple" if X0 == 0 else "green"
        p += plot(f(t), t, 0, 5, thickness=2, color=color)
    p.show(axes_labels=("$t$", "$X$"), ymin=-100, ymax=100)

@interact(r=slider(1.2, 2, 0.2, default=1.2, label="$r$"))
def discrete_growth(r):
    initial_states = [0, 1, 5, 10, -1, -5, -10]
    p = Graphics()
    for X0 in initial_states:
        X_list = [X0*r^t for t in srange(0, 6)]
        color = "purple" if X0 == 0 else "red"
        p += list_plot(X_list, size=30, color=color)
        p += list_plot(X_list, plotjoined=True, color=colors[color].lighter(0.8))
    p.show(axes_labels=("$t$", "$X$"), ymin=-100, ymax=100)

@interact(r=slider(0.2, 0.8, 0.2, default=0.8, label="$r$"))
def discrete_decay(r):
    initial_states = [0, 10, 50, 100, -10, -50, -100]
    p = Graphics()
    for X0 in initial_states:
        X_list = [X0*r^t for t in srange(0, 6)]
        color = "purple" if X0 == 0 else "green"
        p += list_plot(X_list, size=30, color=color)
        p += list_plot(X_list, plotjoined=True, color=colors[color].lighter(0.8))
    p.show(axes_labels=("$t$", "$X$"), ymin=-100, ymax=100)