from numpy import *
from scipy import *
from scipy.integrate import odeint
from matplotlib.pyplot import *
from mpl_toolkits.axes_grid.axislines import SubplotZero

def myFun(u,t=0.,mu=.5):
	x = u[0]
	v = u[1]
	dx = v
	dv = mu*(1.-x**2)*v-x
	return (dx,dv)

t = linspace(-5.5,5.5,300)
u0 = array([1.,1.])
mu = [.1, .2, .5, 1., 1.5, 2.0]
fig = figure(figsize=(5.5,7))
ax = SubplotZero(fig,111)
fig.add_subplot(ax)
ax.grid(True)
for direction in ["right","top"]:
        ax.axis[direction].set_visible(False)
for m in mu:
	u = odeint(myFun,u0,t,args=(m,))
	ax.plot(u[:,0],u[:,1],lw=1.5,label=r'$\mu=%.1f$'%m)
x = linspace(-3,3,15)
y = linspace(-4,4,15)
x,y = meshgrid(x,y)
X,Y = myFun([x,y])
M = (hypot(X,Y))
M[M==0]=1.
X,Y = X/M, Y/M
ax.quiver(x,y,X,Y,M,pivot='mid',cmap=cm.jet)
ax.minorticks_on()
ax.legend(handletextpad=0,loc='upper left')
setp(ax.get_legend().get_texts(),fontsize=12)
fig.savefig("Van_der_pols_equation_phase_portrait.svg",bbox_inches="tight",\
        pad_inches=.15)