function [x,y]=aderlur(f,x0,y0,xf,h)
n=fix((xf-x0)/h);
x(1)=x0;
y(1)=y0;
for i=1:n
    x(i+1)=x0+h*i;
    yp=y(i)+h*feval(f,x(i),y(i));
    yc=y(i)+h*feval(f,x(i+1),y(i));
    y(i+1)=(yp+yc)/2;
end