The vanderpole equation is given as
Now I break it down into two first order differential equations.Putting x1 = x and x2= d(x1)/dt .
.
The two first order equations i get are d(x1)/dt = x2 and d(x2)/dt =myu*(1-x1^2)*x2 - x1.
I am using python to solve this differential equation.I am giving the code as below
import scipy as sp
import numpy as np
import matplotlib.pyplot as plt
from scipy.integrate import odeint
def vanderpole_oscillator(t,s):
x1,x2 = s
return [x2,myu*(1-x1**2)*x2 - x1]
t =np.linspace(-10,10,100)
myu=1
initial = (0,1)
sol = odeint(vanderpole_oscillator,y0=initial,t=t,tfirst=True)
#The solution for the two first order differential equation is given as
sol_x1=sol[:,0]
sol_x2=sol[:,1]
Now when you plot the isoclines of these two separate single order differential equation you put
x-axis as t and y-axis as x1 for one single order differential equation.
In method of isoclines what we do is the x axis is x1 and y-axis is x2.So you got to plot the isocline of d(x2)/d(x1).
I tried solving using the code below to get the quiver plot/isoclines
x1 = sol_x1
x2 = sol_x2
#d(x2)/d(x1) is given as below
u = (myu*(1-x1**2)*x2-x1)/x2
v=(myu*(1-x1**2)*x2-x1)/x2
X,Y = meshgrid(x1,x2)
qp = quiver(X,Y,u,v)
I get an error ValueError: Argument U has a size 100 which does not match 10000, the number of arrow positions . what am i doing wrong? Do i need to take x1 and x2 values arbitarily.
I then tried the code given below
x1 = np.linspace(-10,10,100)
x2 = np.linspace(-10,10,100)
#d(x2)/d(x1) is given as below
X,Y = meshgrid(x1,x2)
u = (myu*(1-X**2)*Y-X)/Y
v=(myu*(1-X**2)*Y-X)/Y
qp = quiver(X,Y,u,v)
plt.plot(t,sol_x1)
plt.plot(t,sol_x2)
I got the graph as 
How can i plot it in a better way ? Also not sure if i plotted the correct isoclines.