how to code simple pendulum motion using ode45

how to code simple pendulum motion (displacement vs time) using ode45

Answer 1

[Ameer Hamza]

The equation of simple pendulum is , which can be converted to two first order ODEs

and then using ode45 like this

theta_ic = [0.5; 0]; % initial conditions: theta(t=0)=0.5; dtheta(t=0)=0.
tspan = [0 10];
[t, theta] = ode45(@odeFun, tspan, theta_ic);
plot(t, theta);
legend({'$\theta$', '$\dot{\theta}$'}, ...
    'Location', 'best', ...
    'Interpreter', 'latex', ...
    'FontSize', 16)

function dtheta = odeFun(t, theta)
    g = 9.8;
    l = 1;
    % theta(1) = theta, theta(2) = dtheta
    
    dtheta = zeros(2, 1);
    dtheta(1) = theta(2);
    dtheta(2) = -g/l*theta(1);
end

[Ashan Cooray]

thanks.

[Ameer Hamza]

I am glad to be of help!

[Raphaël Candelier]

Isn't there a sin() function missing in the solution ?

[Sam Chak]

Yes, @Raphaël Candelier, you're right. However, at 0.5 rad or 28.65°, the difference is generally not noticeable at the beginning.

0.5 - sin(0.5)

ans =

0.020574461395797

fv1 = @(t, x, y) y;
fv2 = @(t, x, y) - (9.8/1)*sin(x);
opt = odeset('RelTol', 1e-4, 'AbsTol', 1e-6);
[t, v] = ode45(@(t, x) ([fv1(t, x(1), x(2)); fv2(t, x(1), x(2))]), [0 10], [0.5 0], opt);
plot(t, v, 'linewidth', 1.5)

At 10 seconds, you should be able notice the slight difference in phase.