Week 3 - Solving second order ODEs

simulate pendulam motion:

GIVEN PROBLEM:

simulate animation of pendulam motion.

solution:

we have to solve the following ode.this ode represents the equation of motion of a simple pendulam

d^(theta)/d(t^2) + (b/m)*(d(theta)/d(t)) +(g/l)*sin(theta)

b= damping ,m=mass,g=gravity, l=length

A pendulam is a weight suspended from a fixed pivot. when a pendulam is displaced from its equilibrium position its try to restore its eqilibrium due to action of gravity and duo to inertia it goes beyond equilibrium position and same restoration of equilibrium position happen . Thus it tends to swing like this from one point to another and comes to that point after some time ,this time is called 'Time period' .

MATHEMATICAL EQUATION

(d^(theta)/d(t^2) )+ ( (b/m)*(d(theta)/d(t))) +((g/l)*sin(theta)) .......(1)

let , theta= theta1; d(theta)/dt = d(theta1)/dt=theta2

so from eqn (1)

d(theta2)/dt =(-b/m)*theta2-(g/l)*sin(theta1)

(d/dt)*[theta1 theta2] =[ theta2 -(b/m)*theta2-(g/l)*sin(theta1) ]

CODE:

First of all to simulate the projet we have to simulate a functioin

FUNCTION CODE

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)*(theta1);
dtheta_dt=[dtheta1_dt;dtheta2_dt];
end

PROGRAM CODE:

clear all
close all
clc
 %inputs
 b=.05;
 g= 9.81;
 l=1;
 m=.1;
 %initial condition
 theta_0 = [0;3];
 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);
%pendulam
ct=1;
for i=1:length(results(:,1));
    x0=0;
    y0=0;
    x1= l1*sin(results(i,1));
    y1=-l*cos(results(i,1));
    plot([-.5,.5],[0,0],'linewidth',3)
    line([x0 x1],[y0 y1],'linewidth',3);
    hold on
    plot(x1,y1,'o','markers',15,'markerfacecolor','r');
    grid on
    axis([-1.5 1.5 -1.5 .5]);
    pause(.01)
    h(ct)=getframe(gcf);
    
ct=ct+1;
end
movie(h)
videofile=VideoWriter('pendulam.avi','Uncompressed AVI');
open(videofile)
writeVideo(videofile,h)
close(videofile)

error occured at line 19 . Modifying error x1= l*sin(results(i,1));

Modified code:

clear all
close all
clc
 %inputs
 b=.05;
 g= 9.81;
 l=1;
 m=.1;
 %initial condition
 theta_0 = [0;3];
 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);
%pendulam
ct=1;
for i=1:length(results(:,1));
    x0=0;
    y0=0;
    x1= l*sin(results(i,1));
    y1=-l*cos(results(i,1));
    plot([-.5,.5],[0,0],'linewidth',3)
    line([x0 x1],[y0 y1],'linewidth',3);
    hold on
    plot(x1,y1,'o','markers',15,'markerfacecolor','r');
    grid on
    axis([-1.5 1.5 -1.5 .5]);
    pause(.01)
    h(ct)=getframe(gcf);
    
ct=ct+1;
end
movie(h)
videofile=VideoWriter('pendulam.avi','Uncompressed AVI');
open(videofile)
writeVideo(videofile,h)
close(videofile)

RESULT:

still animation is not proper so we have to adopt a hold off comand

FINAL CODE:

clear all
close all
clc
 %inputs
 b=.05;
 g= 9.81;
 l=1;
 m=.1;
 %initial condition
 theta_0 = [0;3];
 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);
%pendulam
ct=1;
for i=1:length(results(:,1));
    x0=0;
    y0=0;
    x1= l*sin(results(i,1));
    y1=-l*cos(results(i,1));
    plot([-.5,.5],[0,0],'linewidth',3)
    line([x0 x1],[y0 y1],'linewidth',3);
    hold on
    plot(x1,y1,'o','markers',15,'markerfacecolor','r');
    grid on
    axis([-1.5 1.5 -1.5 .5]);
    pause(.01)
    h(ct)=getframe(gcf);
    hold off
ct=ct+1;
end
movie(h)
videofile=VideoWriter('pendulam.avi','Uncompressed AVI');
open(videofile)
writeVideo(videofile,h)
close(videofile)

FINAL RESULT

link of the final animation given below

https://youtu.be/HThKq3LjvVA

conclusion:

Here we analyse a differential equation of simple pendulam. By altering variables different results are obtained . By changing damping factor time taken to come to equilibrium position is obtained.