#!/usr/bin/octave
disp("Shooting y''=y");

f=@(x,t) [x(2);x(1)]
f([1;2]);
sol=@(t,b) (exp(t-b)+exp(-t))/(exp(-b)+1);
function s=shooting(f,b=1)
#finds s such that the solution of initial valuep orblem y''=y y(0)=1 y'(0)=s
#satisfies y(b)=1
X=lsode(f,[1;0],[0,b]);
x1=X(2,1);
X=lsode(f,[1;1],[0,b]);
x2=X(2,1);
s=(1-x1)/(x2-x1);
end

#b=1
b=1
s=shooting(f,b)
t=linspace(0,b);
X=lsode(f,[1;s],t);
plot(t,X(:,1),t,sol(t,b),";solution;");pause(3);
X(100,1)
error=abs(X(100,1)-1)

#b=14;
b=14
s=shooting(f,b)
t=linspace(0,b);
X=lsode(f,[1;s],t);
plot(t,X(:,1),t,sol(t,b),";solution;");pause(3);
X(100,1)
error=abs(X(100,1)-1)


#b=20;
b=20
s=shooting(f,b)
t=linspace(0,b);
X=lsode(f,[1;s],t);
plot(t,X(:,1),t,sol(t,b),";solution;");pause(5);
X(100,1)
error=abs(X(100,1)-1)