[x,y]=meshgrid([-5:0.5:5],[-5:0.5:5]);
u=ones(size(x));
v=2*x*y./(x^2 + y^2);
quiver(x,y,u,v);
r=sqrt(u.^2+v.^2);
quiver(x,y,u./r,v./r,'r');
hold on
plot([-5,5],[0,0])
plot([0,0],[-5,5])
axis([-5,5,-5,5])
sol=dsolve('(t^2 + y^2)*Dy=2*t*y');
for k=-5:5
sol1=subs(sol,'C2',k);
xx=-5:.1:5;
yy=subs(sol1,'t',xx);
 figure(3);
    axis([-5, 5, -5, 5]);
    hold on
    plot(xx, yy, 'k');  hold on
    hold off;
end
