theta_final = -0.1;
t_final = 3.5;

sol = shooting_pendulum(θfinal,tfinal,b,l,m);

p1 = plot(sol,xlabel="Time",ylabel="Angle",label="θ",vars=(0,1))
scatter!([tfinal],[θfinal],label="θ required")
p2 = plot(sol,xlabel="Time",ylabel="Angular rate",label="ω",vars=(0,2))
plot(p1,p2,layout=(2,1))


x0 = 1;
tspan = [0,1];

[t,x] = ode45(@(t,x)cos(x),tspan,x0);

figure(1)
plot(t,x)
xlabel("t")
ylabel("x(t)")


function dxdt=f(t,x)
    g = 9.81;
    l = 1.0;
    m = 1.0;
    b = 0.1;
    dxdt = zeros(2,1);
    theta = x(1);
    omega = x(2);
    dxdt(1) = theta;
    dxdt(2) = -g/l*sin(theta) - b/(m*l^2)*theta;
end