function [t,x]=shoot_es2(L,s0,itmax,tol)
% Metodo di shooting per il problema
% Exercise 7.4, page 187 Ascher-Petzold
% Shooting NON funziona per questo esempio
% Con L=50, s0=1, itmax=1 calcola qualcosa, ma non ce la fa a fare 2 iter
% Con L=50, s0=2, opz=odeset('AbsTol',1.e-13,'RelTol',1.e-13); fa 2 iter:
% [t,x]=shoot_es2(50,2,2,1.e-5)
% Cos s0=-1 non fa neanche una iterazione
% Provare a sostituire ode45 con ode15s.

T= [1.e-8 L]; %Devo sostituire 0 con 1.e-8 (o qualcosa di piccolo, 
%perche' altrimenti ho divisione per zero nell'eq. diff. e ode45 da' NaN)

z0=[s0 0 1 0];
hold off
close all
hold on
% opz=odeset('AbsTol',1.e-13,'RelTol',1.e-13);
for k=1:itmax
%    [t,z]=ode45(@fz,T,z0,opz);
    [t,z]=ode45(@fz,T,z0);
    inc=z(end,1)/z(end,3);
    [k,inc]
    s1=s0-inc;
    if(abs(inc)<tol)break
    end
    s0=s1;
    z0(1)=s0;
    plot(t,z(:,1))
end
x=z(:,1);
plot(t,z(:,1),'--g')
plot(t,exp(-t),'-.b') % E' asintoticamente la soluzione esatta

function w=fz(t,z)
w = [ z(2); -4/t*z(2)-(t*z(1)-1)*z(1); z(4); (-2*t*z(1)+1)*z(3)-4/t*z(4)];

    