https://youtu.be/7Iso3bShAW8
%Program to define ode function
%function name ode_func
function [dtheta_dt] = ode_func(t,theta,b,g,l,m)
% defining the solution array theta contain theta1 and theta2
% from the derived ODE
% d(theta1)/dt = theta2
% d(theta2)^2/dt = -(b/m)theta2 -(g/l)*sintheta1
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
%program to simulate the pendulum motion
%damping co-efficient
b=.05;
%gravity
g=9.81;
%lenght of pendulum
l=1;
%mass of the ball
m=1;
% for stimulation the motion the
% initial angular displacement is assume as '0'
% initial angular velocity is =3 rad/sec
%initial condition is stored in array
theta_0 = [0;3]; %angular disp. 0 and angular velocity 3
%time span is taken as 0 to 20 sec.
%400 values of time taken using linspace command
t_span = linspace(0,20,400);
% solving ode by calling the function ode_func calculated in
% seperate program ode45 command
% the main variables of ode45 as ode45(@(t,theta)
% calling the function ode_func with inputs
% providing the initial condition time=t_span;
% theta_0 the initial angular displacement and velocity
% the 'result' array provide an array of two column giving angular displacement and
% velocity corresponding to the time 't'
[t,result] = ode45(@(t,theta)ode_func(t,theta,b,g,l,m),t_span,theta_0);
figure(1)
%plotting the graph of time vs angular displacement
% using first column value of result
plot(t,result(:,1));
hold on
%plotting the graph of time vs velocity
%using 2nd column value of result
plot(t,result(:,2));
xlabel("time");
ylabel("angular velocity & angular displacement");
title( "ANGULAR DISPLACEMENT,ANGULAR VELOCITY,TIME");
%Animation of motion of simple pendulum
figure(2)
ct=1;
%using for loop for plotting the image of pendulum for corresponding time
% value of i starts from 1 to no of values of displacement
for i=1:length(result(:,1))
%initialising the fixed point of pendulum as 0
x0=0;
y0=0;
%calculating the final point x1 & y1 corresponding to the time
% x1 = length of pendulum * sin value of displacement w.r.t time
x1=l*sin(result(i,1));
%(-) to get y1 value below (0,0)
y1=-l*(cos(result(i,1)));
%plotting the line through fixed point of pendulum
plot([-1 1],[0 0],'LineWidth',5,'Color','k');
hold on
%plotting the lenght of pendulum
plot([x0 x1],[y0 y1],'linewidth',2,'Color','r');
%plotting the bob of pendulum
plot(x1,y1,'-o','markersize',10,'MarkerEdgeColor','b','MarkerFaceColor','b');
hold off
axis([-1.5 1.5 -1.5 .5]);
title('MOTION OF PENDULUM');
%saving the frames for the movie
M(ct)=getframe(gcf);
ct=ct+1;
end
%creating the movie
movie(M);
%creating a avi video file
videofile=VideoWriter("simplependulum.avi","Uncompressed AVI")
%opening the video file
open(videofile)
%writing the movie to video file
writeVideo(videofile,M)
%closing the file
close(videofile)
AIM
To create a program that solves the following ODE represents the equation of motion of a simple pendulum with damping
to write a program to simulate motion of a pendulum
Equation for motion fo a pendulum
`(d^2theta)/dt^2 +b/m (d theta)/dt +g/l sintheta = 0 ->1`
In the above equation,
g = gravity in m/s2,
L = length of the pendulum in m,
m = mass of the ball in kg,
b=damping coefficient. The damping copacity of a material
The above equation is a 2nd order ODE
In MATLAB the method to solve a higher order equation is to convert to a first order equation
Assume `theta =theta1`
`(d theta)/dt =(d theta1)/dt =theta2->2`
`(d^2 theta)/dt^2 =d/dt(d theta1)/dt =(d theta2)/dt`
from eq`->1 -> d (theta2)/dt = -b/m -g/lsintheta1`
writing the equation in matrix format
`d/dt[theta1 ;theta2] = [theta2; b/m-g/l sintheta1]`
the first order ODE is solved using matlab solver ode45
the syntax for ode45 is
[t,result] = ode45(@(t,theta)ode_func(t,theta,b,g,l,m),t_span,theta_0);
the main variables of ode45 as ode45(@(t,theta)
calling the function ode_func with inputs
providing the initial condition time=t_span
theta_0 the initial angular displacement and velocity
the 'result' array provide an array of two column giving angular displacement and
velocity corresponding to the time 't'
the program of function ode_func is written seperately
function [dtheta_dt] = ode_func(t,theta,b,g,l,m)
defining the solution array theta contain theta1 and theta2
from the derived ODE
d(theta1)/dt = theta2
d(theta2)^2/dt = -(b/m)theta2 -(g/l)*sintheta1
dtheta_dt = [dtheta1_dt;dtheta2_dt]
the value for angular displacement and angular velocity is generated with respect to time
Body of content
The program will simulate the motion of a pendulum whose
Length L=1 metre,
Mass m=1 kg,
Damping coefficient b=0.05.
g=9.81 m/s2.
Simulate the motion between 0-20 sec, for angular displacement=0,
angular velocity=3 rad/sec at time t=0.
theta_0 = [0;3];
angular disp. 0 and angular velocity 3
time span is taken as 0 to 20 sec. 400 values of time taken using linspace command
t_span = linspace(0,20,400);
solving ode by calling the function ode_func using ode45 solver
plotting the graph of time vs angular displacement using first column value of result
plot(t,result(:,1));
plotting the graph of time vs velocity using 2nd column value of result
plot(t,result(:,2))
using for loop for plotting the image of pendulum for corresponding time
value of i starts from 1 to no of values of displacement
length(result(:,1))
initialising the fixed point of pendulum as x0,y0 as (0,0)
calculating the final point x1 & y1 corresponding to the time
x1 = length of pendulum * sin value of displacement w.r.t time
y2 =lenght of pendulam * cos value of displacement w.r.t time
(-) to get y1 value below (0,0)
plotting the line through fixed point of pendulum
plot([-1 1],[0 0],'LineWidth',5,'Color','k');
plotting the lenght of pendulum
plot([x0 x1],[y0 y1],'linewidth',2,'Color','r');
plotting the bob of pendulum
plot(x1,y1,'-o','markersize',10,'MarkerEdgeColor','b','MarkerFaceColor','b');
fixing the axis
axis([-1.5 1.5 -1.5 .5]);
saving the frames for the movie m using getframe command
end of for loop
creating the movie
movie(M)
creating a avi video file
videofile=VideoWriter("simplependulum.avi","Uncompressed AVI")
opening the video file
open(videofile)
writing the movie to video file
writeVideo(videofile,M)
closing the file
close(videofile)
Result
1)the graph is plotted time vs angular displacement
time vs angular velocity
2) animation is created for the motion of simple pendulum