Week 3 - Solving second order ODEs

dtheta2_dt = -(b/k)*theta2-(g/l)*sin(theta1)

we have to all the function during main code run. so separate script is created for above 2nd order equation

by using a inbuild command called 'ODE45'

ODE solver can solve systems of equations y' = f(t,y)

%function
function [dtheta_dt] = ode_fun2(t, theta, b,g,l,k)
theta1 = theta(1);
theta2 = theta(2);
dtheta1_dt = theta2;
dtheta2_dt = -(b/k)*theta2-(g/l)*sin(theta1);
dtheta_dt = [dtheta1_dt; dtheta2_dt];
end

%main script
% 2ODE
b = 0.05;
g =9.81;
l = 1;
k = 1;

%initial_condition i.e solution matrix (ang dis. = 0, ang vel.=5)
theta_0 = [0; 3];

%time_points

t_span = linspace(0,10,500);

%solve ODE

[t, results] = ode45(@(t,theta) ode_fun2(t, theta, b,g,l,k),t_span, theta_0);
%results = [angular displacement, angular velocity]
plot(t,results(:,1))
hold on
plot(t,results(:,2))
xlabel('plot')
ylabel('time')
%%Animation
x0 = 0;
y0 = 0;

ct = 1; % frame for every second

for i = 1:length(t)
     x1 = l*sin(results(i,1));
     y1 = -l*cos(results(i,1));
     figure(2)
     plot([-1 1],[0 0],'linewidth',6,'color','g');
     axis([-2 2 -2 2]);
     line([x0 x1],[y0 y1],'linewidth',2,'color','b');
     hold on;
     plot(x0,y0,'o','markersize',3,'markerfacecolor','y');
     plot(x1,y1,'o','markersize',5,'markerfacecolor','r');
     title('damped simple pendulum');
     legend('support rod','length of pendulum','fixed point','bob of pendulum');
     grid on;
     hold off;
     M(ct)=getframe(gcf);
     ct=ct+1;
end
%creating a movie for pendulum oscillations
movie(M);
videofile=videowriter('pendulum_avi','uncompressed AVI');
open(videofile);
writeVideo(videofile,M);
close(videofile);





   

RESULT:- Animation for damped pendulum is obtained