clc; clear;
syms y(t) c
sol = diff(y, t) == -t*y^3
cond=y(0) == 1/sqrt(c);
sol2=(dsolve(sol,cond));
cond2=y(0) == -1/sqrt(c);
sol3=(dsolve(sol,cond2));
for c=linspace(0.25,10,8);
t=linspace(-3,3,100);
plot(t,subs(sol2),'b','LineWidth',1.5);
hold on
plot(t,subs(sol3),'b','LineWidth',1.5);grid on;
end
[t, y] = meshgrid(-3:0.2:3, -2:0.2:2);
S = -t.*y.^3;
L = sqrt(1 + S.^2)
quiver(t, y, 1./L, S./L, 0.5)
set(gca, 'XLim', [-3 3], 'YLim', [-2 2]);