close all
clear all
clc
% Objective- Program to solve ODE which represents the equation of motion of simple pendulum with damping ,and simulate the motion of simple pendulum
% The way the pendulum moves depends on the Newtons second law. When this law is written down, we get a second order Ordinary Differential Equation that describes the position of the pendulum as well as the position of ball w.r.t time.
% equation for motion of simple pendulum
% d^2theta/dt^2+(b/m)dtheta/dt+(g/l)sin(theta)=0
% Inputs
% 'b' is the damping co-efficient
b=0.05;
% 'g' is the acceleration due to gravity in m/s^2
g=9.81;
% 'l' is the length of pendulum in m
l=1;
% 'm' is the mass of ball in kg
m=1;
% This is a second order differential equation.
% For solving this second order differential equation we split it into a set of two first order diffrential equation
% After solving the equation we get the values angular displacement and angular velocity
% initial conditions of pendulum
% initial displacement is zero and initial velocity is 3 rad/sec
theta_0=[0;3];
% time interval -we simulate motion between 0-20 sec
t_span=linspace(0,20,500);
% Solving ODE
% For solving the equation we have to call the ODE function
% calling
[t,results]= ode45(@(t,theta) ode_func(t,theta,b,g,l,m),t_span,theta_0);
% ode45 is an in built solver in matlab and here its main inputs are time 't' and theta
% Here we are saying that ODe is a function of t and theta ,
% then we are defining the ODE function by calling it and passing inputs t, theta,b,g,l and m
% By solving the equation we get the integration time 't' and results
% Since we are solving two ODE the results is a column vector having two columns
% First column gives the value of displacement and second column gives the value of velocity
% Animation of motion of simple pendulum
% We have to animate the movement of pendulum with a ball attached to its end
% coordinates of pendulum
% x0=0
% y0=0
% The values of x1 and y1 are finds by using pythagorean theorem
% x1=l*sin(theta)
% y1=-l*cos(theta)
% theta is the angular displacement and the values of dispalcement is the first column of results
% val=theta
val=results(:,1);
% Animation done by using for loop
% Counter ct
ct=1;
for i=1:length(val);
x0=0;
y0=0;
x1=l*sin(val(i));
y1=-l*cos(val(i));
% ploting
% ploting the line which to assumed be the line from which the pendulum is suspended
plot([-0.8 0.8],[0 0],'linewidth',5,'color','b')
hold on
% ploting pendulum
plot([x0 x1],[y0 y1],'LineWidth',3,'color','g')
hold on
% ploting ball attached to the end of pendulum
plot(x1,y1,'o','markersize',15,'MarkerFaceColor',[0.5 0.5 0.5],'MarkerEdgeColor','r')
title('motion of simple pendulum')
axis([-1,1,-1.5,0.2])
hold off
% getting each frame and storing it to a variable M
M(ct)=getframe(gcf);
ct=ct+1;
pause(0.05)
end
% video creation ,writing and saving
% 'movie(M)' for creating movie with all frames
movie(M)
% video file creation and video writing
videofile =VideoWriter('simple_pendulum.avi','Uncompressed AVI');
open(videofile)
writeVideo(videofile,M)
close(videofile)
% youtube link of the video https://youtu.be/OoU29fH4RWU
% Codes for ODE function
% equation for motion of simple pendulum
% d^2theta/dt^2+(b/m)dtheta/dt+(g/l)sin(theta)=0
% Here this second order differential equation we split it into a set of first order diffrential equations
% such that,
% theta1=theta(1)
%theta2=theta(2)
%dtheta1_dt=theta2
%dtheta2_dt=-(b/m)*theta2-(g/l)*sin(theta1);
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