In :

%matplotlib inline
import numpy as np
from scipy.integrate import odeint
from matplotlib import pyplot as plt
plt.style.use('bmh')




In :

# 4th-order Runge-Kutta
def rk4(x, t, h, f):
# x is coordinates (as a vector)
# h is timestep
# f(x) is a function that returns the derivative
# "Slopes"
k1  = f(x,          t)
k2  = f(x + k1*h/2, t + h/2)
k3  = f(x + k2*h/2, t + h/2)
k4  = f(x + k3*h,   t + h)
# Update time and position
x_  = x + h*(k1 + 2*k2 + 2*k3 + k4)/6
return x_

def streamline(x, T, f, h = 0.05):
# Return a trajectory of length T (in time units)
# starting at position given by the vector x
# Calculated by Runge-Kutta, with timestep h

# Number of timesteps
Nt = int(T / h)
traj = np.zeros((2, Nt))
traj[:,0] = x
for i in range(Nt-1):
t = i*h
traj[:,i+1] = rk4(traj[:,i], h, t, f)
return traj

def get_streamlines(f, T, h, N, xlims, ylims):
Nt = int(T/h)
trajectories = np.zeros((2, Nt, N))
trajectories[0,0,:] = np.random.uniform(low = xlims, high = xlims, size = N)
trajectories[1,0,:] = np.random.uniform(low = ylims, high = ylims, size = N)
for i in range(1, Nt):
t = i*h
trajectories[:,i,:] = rk4(trajectories[:,i-1,:], t, h, f)
return trajectories

def plot_streamlines(f, T, h, N, xlims, ylims, figsize=(9,9)):
# Calculate N streamlines from random initial positions
lines = get_streamlines(f, T, h, N, xlims, ylims)
# Plot all the lines
fig = plt.figure(figsize = figsize)
plt.plot(lines[0,:,:], lines[1,:,:], lw = 1, c = '#A60628', alpha = 0.25)
# Set plot limits
plt.xlim(*xlims)
plt.ylim(*ylims)




In :

def f(x, t):
dx = x[0,:] + np.exp(-x)
dy = -x
return np.array([dx, dy])

T = 0.5
h = 0.005
N = 1000

xlims = (-3, 3)
ylims = (-2, 2)

plot_streamlines(f, T, h, N, xlims, ylims)







In :

def f(x, t):
a = 0
dx = -x + a*x*(x**2 + x**2)
dy =  x + a*x*(x**2 + x**2)
return np.array([dx, dy])

T = 0.5
h = 0.005
N = 1000

xlims = (-3, 3)
ylims = (-3, 3)

plot_streamlines(f, T, h, N, xlims, ylims)







In :

def f(X, t):
x = X
y = X
dx = x*(3 - x - 2*y)
dy = y*(2 - x - y)
return np.array([dx, dy])

T = 2
h = 0.002
N = 2000

xlims = (0, 4)
ylims = (0, 3)

plot_streamlines(f, T, h, N, xlims, ylims)