function y = bob_vdpol(T,y)
%  Van der Pol equation in relaxation oscillation.
%  The parameters here were specified by Shampine.  The limit cycle has
%  portions where the solution components change slowly and the problem
%  is quite stiff alternating with regions of very sharp change,
%  quasi-discontinuities, where it is not stiff.  The initial
%  conditions are close to a portion of slow change so as to test
%  schemes for the selection of the initial step size.
%  
%   Evaluation of a test set for stiff ODE solvers, 
 
  

  x1 = y(1);                                  % x value
  x2 = y(2);                                  % y value

  x1 = x1 + T*x2;                             % Compute x_hat(k_1)
  x2 = x2 + x2*(1-x1^2)*T - x1*T;             %                                              
  y = [x1 x2]';                               %
  
