%PENDBVP1: MATLAB script M-file that solves a 
%boundary value problem for a simple pendulum.
g=9.81; l=1;
rhs = @(t,y) [y(2);-(g/l)*sin(y(1))];
bc = @(yleft,yright) [yleft(1)-pi/4;yright(1)];
%Here, 'yleft' and 'yright' are variable names, and 
%can be changed, but order matters. I.e., in the 
%syntax @(yleft, yright) the first variable will 
%refer to y at the left endpoint. The designation
%yleft(1) refers to the first component of y 
%(which is theta in this case). If we had a condition 
%on the derivative of theta at the left endpoint, we 
%would use yleft(2). Finally, MATLAB sets each 
%component of the vector defined in bc to 0.
initial_partition = linspace(0,.2,25);
%MATLAB's collocation methods require an initial 
%partition of the independent variable, which is 
%specifed here as 'initial_partition.' (The variable
%name can be changed.) If MATLAB requires a more
%refined partition, it will update this choice. 
initial_solution = [pi/4;-1];
%This is an approximation of the solution along 
%the entire interval specified in 'initial_partition.'
%Here, the first component is theta (t) and the 
%second is theta'(t). For this example, this approximation
%is very crude, but we'll see that it's sufficient.
%For examples with more precise approximations see
%pendbvp1a and pendbvp3.
solinit = bvpinit(initial_partition,initial_solution);
sol = bvp4c(rhs,bc,solinit)
pause
[sol.x' sol.y(1,:)' sol.y(2,:)'] %Note the use of sol.x, even though 
                                 %our independent variable is t
pause
plot(sol.x, sol.y(1,:))
title('Pendulum angle as a function of time t','interpreter','latex','FontSize',16)
xlabel('Time $t$','interpreter','latex','FontSize',14)
ylabel('Pendulum angle $\theta$','interpreter','latex','FontSize',14)