% m34_stiffVDP.m - van der Pol oscillator for large C --
% the stiff case.  Section 10.3 of the Chebfun Guide
% describes how to invoke stiff solvers in Chebfun.

  f = @(t,u,C) [ u(2) ; -u(1)-C*(u(1)^2-1)*u(2) ] ;

  format compact
  close all
  u0 = [0.1;0];
  tspan = [0 120];
  C = input('C? (e.g., 10, 100, 1000) ');

  disp(' ')
  disp('stiff solver ode23s')
  tic
  [t,u] = ode23s(f,tspan,u0,[],C);
  toc

  plot(t,u(:,1),'linewidth',2)
  grid on
  xlabel t, ylabel('u')

  disp(' ')
  disp('non-stiff solver ode23')
  tic
  [t,u] = ode23(f,tspan,u0,[],C);
  toc
  disp(' ')