function doall=van_der_pol
clear all; close all;
% Eli Tziperman: example for course in nonlinear dynamics and
% chaos:  Solve van der pol equation:
% \ddot{x}+mu*(1-x^2)\dot{x}+x=0

% define parameters:
% ------------------
mu=10.2;
% Initial values:
x(1:2,1)=[0.1,0.1]';
% time step: needs to be small enough for a given problem:
delta_t=0.01;
% Print (and save) every nprint steps:
nprint=1;
% Integration time:
T=50.0;

% Number of steps:
N=floor(T/delta_t);
NP=floor(N/nprint);


% Integrate:
for n=1:1:NP
    xp=x(:,n)';
    for j=1:1:nprint
      x_tilde=xp+delta_t*rhs(xp,mu);
      xp=xp+0.5*( rhs(xp,mu) + rhs(x_tilde,mu) )*delta_t;
    end
    x(:,n+1)=xp';
    time(n)=n*delta_t;
end

% Plot simple time series:
plot(time(1:n),x(1,1:n))
xlabel('time')
ylabel('x')


%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
function ff=rhs(x,mu)
% return the rhs of the dynamical system:

ff(1)=x(2);
ff(2)=mu*(1-x(1)^2)*x(2)-x(1);
