y=dsolve('Dy=y/x-y^2/x','y(1)=0.5','x')
y =
1/(1/x + 1)
梯形公式：
y(1)=0.5;x(1)=1;h=0.1;
for n=1:99
x(n+1)=x(n)+h;
y(n+1)=y(n)+h*( y(n)/x(n)-y(n)^2/x(n));
y(n+1)=y(n)+h/2*( y(n)/x(n)-y(n)^2/x(n) +y(n+1)/x(n+1)-y(n+1)^2/x(n+1) );
end
plot(x,y,'y')