Stable Limit Cycles

We investigate the system of coupled differential equations in polar coordinates given by:
 r’ = r(1-r^2) + mu*r*cos(theta), theta’ = 1
 with initial condition r = rad and theta = 0, and mu some parameter.  
 We use NDSolve to obtain the solution and then plot the solution as a parametric plot.

mu=0;
rad = .01;

sol=NDSolve[{r'[t]==r[t]*(1-r[t]^2)+mu*r[t]*Cos[h[t]], h'[t]==1, r[0]==rad, h[0]==0},{r, h},{t,0,100}];

ParametricPlot[Evaluate[{r[t]*Cos[h[t]],r[t]*Sin[h[t]]}/.sol], {t, 0, 100}]

The applet below allows you to point and click to choose a new initial condition and enter an alternative value of mu.

Manipulate[ParametricPlot[Evaluate[{r[t]*Cos[h[t]],r[t]*Sin[h[t]]}/.Quiet@NDSolve[{r'[t]==r[t]*(1-r[t]^2)+mu*r[t]*Cos[h[t]], h'[t]-1==0, r[0]==Norm[xv0], h[0]== (1-Sign[xv0[[1]]])*Pi/2+ArcTan[xv0[[2]]/xv0[[1]]]},{r, h},{t,0,tt}]],{t,0,tt},PlotRange->{{-3,3},{-3,3}},AxesLabel->{TraditionalForm["x[t]"],TraditionalForm["v[t]"]}],{{tt,15,"time t"},.1,20,Appearance→"Labeled"},{{mu,0.75,"parameter [Mu]"},0,2,0.1,Appearance→"Labeled"},{{xv0,{1,1}},{-3,-3},{3,3},Locator}]

In the following plot, we can plot the trajectories for various values of mu together on the same graph.  
Use the parameter “startplotting” to ignore the transient behavior and locate the stable limit cycle in each case.  

Note: We need to clear the kernel here since we’ve already assigned values to these functions when we solved them previously.

Quit[]

sol[mu_?NumericQ]:=First@NDSolve[{r'[t]==r[t]*(1-r[t]^2)+mu*r[t]*Cos[h[t]], h'[t]-1==0, r[0]==.1, h[0]==0},{r, h},{t,0,100}]

ParametricPlot[Evaluate[{r[t]*Cos[h[t]],r[t]*Sin[h[t]]}/.sol[#]&/@Range[0,1,0.2]],{t,0,10}]