Week 3 - Solving second order ODEs

SOLVING SECOND ORDER DIFFERENTIAL EQUATION, SIMULATING PENDULUM MOTION

AIM: - To solve a second order ODE for depicting motion of pendulum

ODE45 syntax

[t,y] = ode45(odefun,tspan,y0)

The command ode45 is used for solving second order differential. Similarly ode23 is used for single order equations

1. odefun - The function handle where the second order differential is loaded.
2. tspan – timespan of the program.
3. y0- The initial conditions provided.

PROCEDURE:-

1. Saving the input parameters like damping coefficient, mass of the bob, length of the pendulum and gravity.
2. The second order differential formula used here is d2θdt2+(b/m)(dθ/dt)+gLsinθ=0.
3. Using this formula create a function handle with name “my_ode”.
4. Provide the initial conditions (displacement, velocity) and the time span that will be used.
5. Use the function handle according to the syntax given to generate results.
6. To check if the results produced are correct plot t vs results graph for confirmation.
7. Create the animation for oscillation by forming a forloop
8. The forloop should be the length of the number of displacement
9. Since we get 2 input values (Displacement, Velocity) one needs to mention that its the first column which needs to be taken into consideration.
10. Fix a point and name them as (x0,y0) and with the help of trignometric functions allot (x1,y1) for generating the movement
11. Plot the animation by using suitable commands (Plot, line).
12. The forloop enables to generate every movement of pendulum with the help of ODE45 and displacement values.
13. Allot a particular axis to make a central appearance, give a pause count as low as possible to make the animation smooth.
14. Create a movie out of it
15. Using the videoWriter save the video.

Animation Link:-https://www.youtube.com/watch?v=cOjO05yRJws

function [dtheta_dt] = my_ode(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 
-----------------------------------------------------------------------------------------
clear all
close all
clc

%Inputs
%damping coefficient
b = 0.05;
%gravitational force
g = 9.81;
%length
l = 1;
%mass of the bob
m = 0.1;

%Equation of motion of a single pendulam
%(d2Theta/dt2)+b/m(dtheta/dt)+g/l(sintheta)=0

%time_span
time_space = linspace(0,20,500);

%initial_condition
theta_0 = [0;5];

%ode_operation
[t,results] = ode45(@(t,theta) my_ode(t,theta,b,g,l,m),time_space,theta_0);

%Animation of pendulum
ct=1;
for i=1:length(results)
    THETA = results(i,1);
    
    x0 = 0;
    y0 = 0;
    x1 = -(l*sin(THETA));
    y1 = -(l*cos(THETA));
    
    
   plot(x1,y1,'O','markersize',15,'markerfacecolor','r')
   line([-1 1],[0 0],'linewidth',5,'color','g');
   line([x0,x1],[y0,y1],'linewidth',2);
   axis([-3 3 -3 3])
   legend('BOB')
   title('OSCILLATION OF PENDULUM')
   
   
  
  
   pause(0.0001);
   M(ct) = getframe(gcf);
   ct=ct+1;
end

movie(M);%creating the movie
videofile = VideoWriter('Oscilation of pendulam2.avi','uncompressed AVI');
open(videofile);%opening the file 
writeVideo(videofile,M);
close(videofile);