importnumpyasnp
importmatplotlib.pyplotasplt
alpha = 1.1 # prey growth rate 
beta = 0.4 # prey death rate 
gamma = 0.4 # predator death rate 
delta = 0.1 # predator growth rate 
x0 = 10 # initial prey population
y0s = [1, 2, 5, 7, 10, 12, 15] # initial predator population
fig, ax = plt.subplots(figsize=(6,6))
for y0 in y0s: 
 dt = 1e-3
 T_MAX = 20
 t = np.arange(0,T_MAX,dt)
 x = np.zeros_like(t, dtype=np.float)
 y = np.zeros_like(t, dtype=np.float)
 x[0], y[0] = x0, y0
 for i, tt in enumerate(t):
 if i == 0:
 continue
 x[i] = x[i-1] + (alpha * x[i-1] - beta * x[i-1] * y[i-1]) * dt
 y[i] = y[i-1] + (delta * x[i-1] * y[i-1] - gamma * y[i-1]) * dt
 
 ax.plot(x, y, label='$y_0=%d$' % y0)
ax.grid(True)
ax.legend()
fig.savefig('predator_prey.svg',
 bbox_inches='tight', transparent=True)