Week 3 - Solving second order ODEs
AIM:- To write a program to slove a 2nd order ODEs that represent the motion of simple pendulum and to animate thats motions
Input functions:-
g= 9.81m/s Gravity
L= 1m Length of pendulum
m= 1kg mass of the ball
b= 0.05 0.05 damping coefficent
Objectives:- To simulate the motion between 0-20sec foa angular displacement = 0, angular velocity =3rad/sec at time t = 0.
function program:-
function [dtheta_dt] = ode_func(t ,theta, b, g, l, m)
theta1 = theta(1)
theta2 = theta(2)
dtheta1_dt = theta2;
dtheta2_dt = -(b/m)*theta2 - (g/l)*sin(theta1);
dtheta_dt = [dtheta1_dt; dtheta2_dt];
end
main program:-
clear all
close all
clc
b=0.05; %damping coefficient
g=9.81; %gravity
l=1; %length of pendulum
m=1; %mass of the bob
%intial condition
theta_0 = [0; 5]; % theta= angular displlacement of the pendulum
% time points
t_span = linspace(0,10,500); %t = time
% solve ODE
[t, results] = ode45(@(t,theta) ode_func(t, theta, b,g,l,m),t_span, theta_0);
figure(1)
plot(t,results(:,1))
hold on
plot(t,results(:,2))
ylabel('angular velocity and angular displacement')
xlabel('time')
grid on
z = results(:,1);
ct = 1;
for i = 1:length(z)
THETA1 = z(i);
%joint on base
x0 = 0;
y0 = 0;
x1 = l.*sin(THETA1);
y1 = -l.*cos(THETA1);
figure(2)
plot(x0, y0, 'o', 'MarkerSize', 10, 'MarkerFaceColor', [0 0 0], 'MarkerEdgeColor', [1 1 0])
hold on
plot([x0 x1], [y0 y1], '-bo', 'LineWidth', 3, 'MarkerSize', 20, 'MarkerFaceColor', 'r', 'MarkerEdgeColor', 'r', 'MarkerIndices', 2)
axis([-1.5 1.5 -1.5 0.5]);
xlabel('Horizental Displacement in m')
ylabel('Vertical Displacement in m')
pause(0.03)
title('Pendulum Motion with matlab ode45')
M(ct)= getframe(gcf);
ct = ct+1;
hold off
end
%crate the movie of pendulum motion
movie(M)
videofile = VideoWriter('pendulum.avi', 'Uncompressed AVI');
open(videofile)
writeVideo(videofile, M)
close(videofile)
plot:-
Pendulum:-
Conclusion:- 2nd order ODEs slove by ode45 command 2nd order odes can be sloved and prectical use of matlab can be seen and achived.
link simple pendulum video
youtube link :- https://youtu.be/yMbCz3ulgtY