% phase portrait for u'' + f(u)=0 where f(u) = u(1-u);
clear all
set(0,                           ...
   'defaultaxesfontsize', 20,   ...
   'defaultaxeslinewidth', 1.0, ...
   'defaultlinelinewidth', 1.2, ...
   'defaultpatchlinewidth', 0.7);
% nullclines  
x=[-1:.0001:1.2];  % specify the x- coordinate values
N = length(x);
 
f = x.*(1-x); %the rhs function
figure(1)
plot(x,zeros(N,1),'--', [0,0],[-1,1],'--', [1,1],[-1,1],'--','linewidth',2)
hold on
plot([0, 1],[0, 0],'*','linewidth',2)
axis([-1 1.5 -1 1])

hold off
 

%integral curves
F =  x.^2 - 2*x.^3/3;
% create 4 curves

xj = [0,0.25,0.5,1,-0.6]; % specify  initial points
for j = 1:5
    xs = xj(j);
    Fj = xs.^2 - 2*xs.^3/3;;
    ndx = find(Fj-F>=0);
    
    
 
plot(x(ndx),sqrt(Fj-F(ndx)),x(ndx),-sqrt(Fj-F(ndx)),'linewidth',2)
hold on
end
plot(x,zeros(N,1),'--', [0,0],[-1,1],'--', [1,1],[-1,1],'--','linewidth',2)
 
plot([0,1],[0,0],'*','linewidth',2)
annotation('textarrow',[0.35, 0.44],[0.85 0.85],'linewidth',2)
annotation('textarrow',[0.44 0.35],[ 0.2  0.2],'linewidth',2)
axis([-0.6 1.2 -1 1])
set(gca, 'XTick', -0.6 : 0.3 : 1.2);
set(gca, 'YTick', -1 : 0.5 : 1);
hold off
xlabel('u','fontsize',20)
ylabel('v','fontsize',20)
hold off
