%%capture for i, mu in enumerate(np.exp(np.linspace(np.log(0.1), np.log(3.0), 121))):

   def system(t, xy):
       x, y = xy
       dx = y
       dy = mu * (1-x**2) * y - x
       return [dx, dy]

   vmin, vmax, wmin, wmax = -4, 4, -4, 4

   t_span = [0, 10]
   trajectory_resolution = 30
   initial_conditions = [(x, y) for x in np.linspace(vmin, vmax, trajectory_resolution) for y in np.linspace(wmin, wmax, trajectory_resolution)]
   sols = {}
   for ic in initial_conditions:
       sols[ic] = solve_ivp(system, t_span, ic, dense_output=True, max_step=0.01)

   vs = np.linspace(vmin, vmax, 200)
   v_axis = np.linspace(vmin, vmax, 20)
   w_axis = np.linspace(wmin, wmax, 20)

   v_values, w_values = np.meshgrid(v_axis, w_axis)

   dv, dw = system(0, (v_values, w_values))

   fig, ax = plt.subplots(figsize=(16,16))
   # integral curves
   for ic in initial_conditions:
       sol = sols[ic]
       ax.plot(sol.y[0], sol.y[1], color='k', alpha=0.4, linewidth=0.5)

   # vector fields
   arrow_lengths = np.sqrt(dv**2 + dw**2)
   dv = dv / arrow_lengths * np.tanh(arrow_lengths)
   dw = dw / arrow_lengths * np.tanh(arrow_lengths)
   alpha_values = 1 - (arrow_lengths / np.max(arrow_lengths))**0.2
   ax.quiver(v_values, w_values, dv, dw, color='blue', linewidth=0.5, scale=25, alpha=alpha_values)

   # nullclines
   vs = np.linspace(-0.99, 0.99, 1000)
   ax.plot(vs, vs/(1-vs**2)/mu,  color="green", alpha=0.4)
   vs = np.linspace(vmin, -1.01, 1000)
   ax.plot(vs, vs/(1-vs**2)/mu,  color="green", alpha=0.4)
   vs = np.linspace(1.01, vmax, 1000)
   ax.plot(vs, vs/(1-vs**2)/mu,  color="green", alpha=0.4)

   ax.set_title(f'van der Pol oscillator\n$\mu={mu:.2f}$')

   ax.set_xlim(vmin, vmax)
   ax.set_ylim(wmin, wmax)
   ax.set_xticks([])
   ax.set_yticks([])

   dir_path = f"./vdp"
   if not os.path.exists(dir_path):
       os.makedirs(dir_path)

   fig.savefig(f"{dir_path}/{i}.png")

import imageio.v3 as iio import os from natsort import natsorted import moviepy.editor as mp

for dir_path in ["./vdp", "./vdp2"]:

   file_names = natsorted((fn for fn in os.listdir(dir_path) if fn.endswith('.png')))

   # Create a list of image files and set the frame rate
   images = []
   fps = 12

   # Iterate over the file names and append the images to the list
   for file_name in file_names:
       file_path = os.path.join(dir_path, file_name)
       images.append(iio.imread(file_path))

   filename = dir_path[2:]
   iio.imwrite(f"{filename}.gif", images, duration=1000/fps, rewind=True)
   clip = mp.ImageSequenceClip(images, fps=fps)
   clip.write_videofile(f"{filename}.mp4")