AIM
To write a programme that solves the following ODE.This ODE represents the equation of motion of a simple pendulam with damping.
(ⅆ^2 θ)/(ⅆt^2 )+b/m ⅆθ/ⅆt+g/L sinθ =0
Where g=gravity in m/s^2
L=length of pendulam=1meter
m=mass of pendulam in kg=1
b=damping coefficient=.05
solving the ODE
let θ=θ_1
so that , ⅆθ/ⅆt=(ⅆθ_1)/ⅆt=θ_2------------(1)
therefore (ⅆ^2 θ)/(ⅆt^2 )=(ⅆ^2 θ_1)/(ⅆt^2 )=ⅆ/ⅆt ((dθ_1)/dt)=ⅆ/ⅆt (θ_2 )--------(2)
substituting values of equation 1 and 2 in the main ODE
(ⅆ^2 θ1)/(ⅆt^2 )+b/m ⅆθ1/ⅆt+g/L sinθ1 =0
Equal to, (ⅆ^ θ2)/ⅆt+b/m θ2/+g/L sinθ1 =0
So that, (ⅆ^ θ2)/ⅆt=-b/m θ2/-g/L sinθ1
So in matrix form we can represent the solution of above equations as below,
ⅆ/ⅆt [█(θ_1@θ2)]=[█(θ_2@-b/m θ_2-(g/l 〖)sin〗θ )]
So making the above solving to a matlab function
function [dtheta_dt]=ode_func(t,theta,b,g,l,m)
theta1=theta(1);
theta2=theta(2);
dtheta1_dt=theta(2);
dtheta2_dt=-(b/m)*theta2-(g/l)*sin(theta1);
dtheta_dt=[dtheta1_dt;dtheta2_dt];
end
Then the programme for simple pendulam animation is as follows using the above function
clear all
close all
clc
b=0.05;
g=9.81;
l=1;
m=1;
%initial condition
theta_0=[0;3];
%time points
t_span=linspace(0,20,500);
%solve ODE
[t,results]=ode45(@(t,theta)ode_func(t,theta,b,g,l,m),t_span,theta_0);
ang_disp=results(:,1);
ct=1;
for i=1:length(ang_disp)
x0=0;
y0=0;
x(i)=l*sin(ang_disp(i));
y(i)=-l*cos(ang_disp(i));
% plotting
plot([x0,x(i)],[y0,y(i)],'linewidth',3);
axis([-1.5 1.5 -1.2 .2]);
ylabel('plot');
xlabel('time');
hold on
plot(x(i),y(i),'o','MarkerFaceColor','red');
hold off
M(ct)=getframe(gcf);
ct=ct+1;
pause(0.3);
end
movie(M)
videofile=VideoWriter('simple_pendulam01a.avi','uncompressed AVI');
open(videofile)
writeVideo(videofile,M)
close(videofile)