  // Simple one dimension ODE
  // dy/dx=1+(y-x)^2, y(0)=0.5

  deff('[dydx]=f(x,y)','dydx=1+(y-x)^2')
  y0=0.5;x0=0;
  
  xout=linspace(0,1.9,20);
  yexact=xout+1.0./(2-xout);
 
  for j=1:4
    x=linspace(0,1.9,1+19*10^(j-1));
    y=ode('adams',y0,x0,x,f);

    xoutvec=1:10^(j-1):1+19*10^(j-1);
    out=[xout' y(xoutvec)' yexact'];
    disp('   x      Adams        exakt')
    disp(out)
  end