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