Chapter 4: The Logistic ODE for Urban Growth

4.1 From Malthus to Verhulst

The simplest model of urban population growth begins with exponential (Malthusian) growth$dN/dt = rN$, where$r$ is the intrinsic growth rate. But no city grows exponentially forever—resources, land, infrastructure, and social carrying capacity impose limits. The logistic equation, introduced by Pierre-François Verhulst in 1838, captures this saturation:

$$\frac{dN}{dt} = r \, N \left(1 - \frac{N}{K}\right)$$

Here $K$ is the carrying capacity—the maximum sustainable population. The term$(1 - N/K)$ represents density-dependent feedback: as $N \to K$, growth slows to zero.

4.2 Analytical Solution

The logistic equation is separable. We rearrange and integrate:

$$\int \frac{dN}{N(1 - N/K)} = \int r \, dt$$

Using partial fractions on the left side:

$$\frac{1}{N(1 - N/K)} = \frac{1}{N} + \frac{1/K}{1 - N/K}$$

Integrating both sides:

$$\ln N - \ln\left(1 - \frac{N}{K}\right) = rt + C$$

Applying the initial condition $N(0) = N_0$ and solving for $N(t)$:

$$\boxed{N(t) = \frac{K}{1 + \left(\frac{K - N_0}{N_0}\right) e^{-rt}}}$$

This is the logistic sigmoid. Key properties:

• $N(0) = N_0$
• $N(t) \to K$ as $t \to \infty$
• Inflection point at $N = K/2$, where growth is fastest
• Maximum growth rate: $dN/dt|_{\max} = rK/4$

4.3 Stability Analysis

Define $f(N) = rN(1 - N/K)$. The fixed points satisfy $f(N^*) = 0$:

$$N_1^* = 0, \qquad N_2^* = K$$

Linear stability is determined by the sign of$f'(N^*)$:

$$f'(N) = r\left(1 - \frac{2N}{K}\right)$$

• At $N_1^* = 0$: $f'(0) = r > 0$ — unstable (any small population grows)
• At $N_2^* = K$: $f'(K) = -r < 0$ — stable (perturbations decay)

Phase Portrait

The phase portrait plots $dN/dt$ vs$N$. It is a downward parabola with roots at 0 and$K$, maximum at$N = K/2$. Arrows point right for$0 < N < K$ (growth) and left for$N > K$ (decline).

Bifurcation

If we treat $r$ as a bifurcation parameter and allow it to become negative (e.g., during urban decline), we get a transcritical bifurcation at $r = 0$: the stability of $N^* = 0$ and$N^* = K$ exchange.

$$\text{For } r < 0: \quad N^* = 0 \text{ is stable}, \quad N^* = K \text{ is unstable (population collapses)}$$

4.4 Python: Logistic Growth Simulation

We solve the logistic ODE for multiple cities with different growth rates and carrying capacities, and display the phase portrait alongside the time series.

Logistic ODE: Multi-City Growth & Phase Portrait

Python

script.py67 lines

import numpy as np

# ── RK4 solver for dN/dt = r*N*(1 - N/K) ──
def rk4_logistic(N0, r, K, dt, n_steps):
    """Fourth-order Runge-Kutta for logistic ODE."""
    t = np.zeros(n_steps + 1)
    N = np.zeros(n_steps + 1)
    N[0] = N0
    for i in range(n_steps):
        f = lambda n: r * n * (1 - n / K)
        k1 = dt * f(N[i])
        k2 = dt * f(N[i] + k1/2)
        k3 = dt * f(N[i] + k2/2)
        k4 = dt * f(N[i] + k3)
        N[i+1] = N[i] + (k1 + 2*k2 + 2*k3 + k4) / 6
        t[i+1] = t[i] + dt
    return t, N

# Analytical solution for comparison
def logistic_exact(t, N0, r, K):
    return K / (1 + ((K - N0) / N0) * np.exp(-r * t))

print("=" * 65)
print("  Logistic Growth: Multiple City Scenarios")
print("=" * 65)
print()

# Three cities with different parameters
cities = [
    {"name": "Fast-growing metropolis", "N0": 0.5e6, "r": 0.04, "K": 12e6},
    {"name": "Medium city",             "N0": 0.2e6, "r": 0.03, "K": 5e6},
    {"name": "Slow-growing town",       "N0": 0.1e6, "r": 0.02, "K": 1.5e6},
]

dt = 0.5  # years
n_steps = 400  # 200 years

for city in cities:
    t, N_num = rk4_logistic(city["N0"], city["r"], city["K"], dt, n_steps)
    N_exact = logistic_exact(t, city["N0"], city["r"], city["K"])
    max_error = np.max(np.abs(N_num - N_exact))

# Find inflection point
    dNdt = city["r"] * N_num * (1 - N_num / city["K"])
    t_inflect = t[np.argmax(dNdt)]
    N_inflect = N_num[np.argmax(dNdt)]

print(f"  {city['name']}:")
    print(f"    N0 = {city['N0']/1e6:.1f}M, r = {city['r']}, K = {city['K']/1e6:.1f}M")
    print(f"    Inflection at t = {t_inflect:.1f} yr, N = {N_inflect/1e6:.2f}M (K/2 = {city['K']/2e6:.2f}M)")
    print(f"    RK4 max error vs analytical: {max_error:.2f}")
    print()

# Phase portrait data
print("  Phase Portrait: dN/dt vs N  (K = 10M, r = 0.03)")
K, r = 10e6, 0.03
N_vals = np.linspace(0, 12e6, 25)
print(f"  {'N (millions)':>14} {'dN/dt (millions/yr)':>22} {'Direction':>12}")
print("-" * 52)
for Nv in N_vals:
    dNdt_v = r * Nv * (1 - Nv / K)
    direction = ">>>" if dNdt_v > 0 else ("<<<" if dNdt_v < 0 else "---")
    print(f"  {Nv/1e6:14.2f} {dNdt_v/1e6:22.4f} {direction:>12}")

print()
print("  Fixed points: N* = 0 (unstable), N* = K (stable)")
print(f"  Max growth rate = rK/4 = {r*K/4/1e6:.4f} M/yr at N = K/2 = {K/2/1e6:.1f}M")

Click Run to execute the Python code

Code will be executed with Python 3 on the server

4.5 Fortran: RK4 Logistic Solver

A Fortran implementation of the fourth-order Runge-Kutta scheme for the logistic ODE.

Fortran: RK4 Logistic ODE Solver

Fortran

program.f9056 lines

program logistic_rk4
  implicit none
  integer, parameter :: n_steps = 500
  real(8) :: N_val, t, dt, r, K, N0
  real(8) :: k1, k2, k3, k4
  real(8) :: N_exact, err_max, err
  integer :: i

! Parameters
  r  = 0.03d0       ! growth rate (1/year)
  K  = 10.0d6       ! carrying capacity
  N0 = 0.5d6        ! initial population
  dt = 0.4d0        ! time step (years)

write(*,'(A)')      '================================================='
  write(*,'(A)')      '  Fortran RK4 Logistic Growth Solver'
  write(*,'(A)')      '================================================='
  write(*,'(A,F10.0)') '  N0 = ', N0
  write(*,'(A,F8.4)')  '  r  = ', r
  write(*,'(A,F12.0)') '  K  = ', K
  write(*,'(A,F8.4)')  '  dt = ', dt
  write(*,'(A)')      '-------------------------------------------------'
  write(*,'(A)') '    t(yr)    N(millions)   N_exact(M)    Error'
  write(*,'(A)') '-------------------------------------------------'

N_val = N0
  t = 0.0d0
  err_max = 0.0d0

do i = 0, n_steps
    ! Analytical solution
    N_exact = K / (1.0d0 + ((K - N0)/N0) * exp(-r * t))
    err = abs(N_val - N_exact)
    if (err > err_max) err_max = err

! Print every 50 steps
    if (mod(i, 50) == 0) then
      write(*,'(F10.1, F14.4, F14.4, ES14.3)') t, N_val/1.0d6, N_exact/1.0d6, err
    end if

! RK4 step
    k1 = dt * r * N_val * (1.0d0 - N_val / K)
    k2 = dt * r * (N_val + k1/2) * (1.0d0 - (N_val + k1/2) / K)
    k3 = dt * r * (N_val + k2/2) * (1.0d0 - (N_val + k2/2) / K)
    k4 = dt * r * (N_val + k3) * (1.0d0 - (N_val + k3) / K)
    N_val = N_val + (k1 + 2*k2 + 2*k3 + k4) / 6.0d0
    t = t + dt
  end do

write(*,'(A)')      '-------------------------------------------------'
  write(*,'(A,ES12.3)') '  Maximum RK4 error: ', err_max
  write(*,'(A,F12.4)')  '  Final N (M):       ', N_val / 1.0d6
  write(*,'(A,F12.4)')  '  K (M):             ', K / 1.0d6
  write(*,'(A)')      '================================================='

end program logistic_rk4

Click Run to execute the Fortran code

Code will be compiled with gfortran and executed on the server

4.6 Extensions and Urban Applications

Time-Varying Carrying Capacity

In reality, carrying capacity changes over time as technology, policy, and infrastructure evolve. A common extension is:

$$\frac{dN}{dt} = rN\left(1 - \frac{N}{K(t)}\right), \qquad K(t) = K_0(1 + \alpha t)$$

Allee Effect

Very small cities may struggle to attract resources. The Allee effect introduces a minimum viable population $A$:

$$\frac{dN}{dt} = rN\left(\frac{N}{A} - 1\right)\left(1 - \frac{N}{K}\right)$$

This creates three fixed points: $N^* = 0$ (stable),$N^* = A$ (unstable threshold), and$N^* = K$ (stable). Cities below the Allee threshold collapse; above it, they grow to capacity. This is a bistable system.

4.7 Summary & Key Takeaways

• The logistic ODE $dN/dt = rN(1 - N/K)$ captures saturation-limited urban growth
• Exact solution: $N(t) = K/[1 + ((K-N_0)/N_0)e^{-rt}]$ — the logistic sigmoid
• Two fixed points: $N^* = 0$ (unstable) and $N^* = K$ (stable)
• Inflection at $N = K/2$ marks the transition from accelerating to decelerating growth
• The Allee effect extension creates bistability, relevant for understanding city collapse thresholds

← Entropy Production Next: Fisher-KPP Waves →