How to obtain the period of a function created by ODE45

I am trying to create a code to measure the evolution of the period of a simple pendulum. Using ODE45 I was able to find the numeric solution to this problem. However, I am still not sure on how to find how the period changes. I was hoping on finding the 0 values with ODE45 but for some reason I get an error of "varargout{4} not assigned during ODE45".

function [] = pendulum()
tspan=[0 70]; 
h0=[pi/18, 0]; 
[t,z,te,ye,ie]=ode45(@pendulum,tspan,h0, @Events);
a = pi/18*cos(sqrt(9.81/0.8)*t);
y = plot(t,h(:,1), 'b'); hold on;
w = plot(t, a, 'r');
legend('Numerically', 'Theory')
function dhdt = pendulum(t,h)
G=9.81; L=0.8; b=6.8809*10^(-3); 
h1=h(1); 
h2=h(2); 
dhdt = [h2 ; -G/L*h1+b*h2^2;];
end
end

I would very much appreciate any help in this project

Answer 1

[Star Strider]

Try this:

function [] = pendulum()
    tspan=[0 70];
    h0=[pi/18, 0];
    opts = odeset('Events', @Events);
    [t,z,te,ye,ie]=ode45(@pendulum,tspan,h0, opts);
    a = pi/18*cos(sqrt(9.81/0.8)*t);
    y = plot(t,z(:,1), 'b'); hold on;
    w = plot(t, a, 'r');
    e = plot(te, ye(:,1), '+g')
    legend('Numerically', 'Theory', 'Zero Crossings')
    function dhdt = pendulum(t,h)
        G=9.81; L=0.8; b=6.8809E-3;
        h1=h(1);
        h2=h(2);
        dhdt = [h2 ; -G/L*h1+b*h2^2;];
    end
    function [value,isterminal,direction] = Events(t,y)
        value = y(1);
        isterminal = 0;
        direction = 0;
    end
end

Note that ‘te’ are the times of the zero-crossings, and the plot call in the ‘e’ assignment plots them. The rest I leave to you.