Solving second order ODEs using MATLAB

Program to calculate ODE of simple pendulum:

simple pendulum:

The simple pendulum is another mechanical system that moves in an oscillatory motion. It consists of a point mass ‘m’ suspended by means of light inextensible string of length L from a fixed support

ODE:

ordinary differential equation is a differential equation containing one or more functions of one independent variable and the derivatives of those functions.ODE is used to describe the transient behaviour of a system. A simple example is a pendulum

ODE converting into matrix form

`θ = θ_ 1`

`(d θ) /(d t )= (d θ _1 )/(d t) = θ _2 ----[ 1 ]`

`(d ^2 θ)/( d t^ 2) = (d ^2 θ_1)/( d t^ 2) = d /(d t) ( (d θ _1)/( d t) ) = d ( θ 2 ) ----[2]`

therefore

`(d ^2 θ)/ (d t^ 2) + b/ m . (d θ)/ (d t )+ g/ l sin θ = 0`

put `θ = θ_ 1`

`(d ^2 θ_1)/ (d t^ 2) + b/ m . (d θ_1)/ (d t )+ g/ l sin θ_1 = 0`

put `θ_1 = θ_2`

`(d θ_2)/ (d t^ 2) + b/ m . (θ_2)+ g/ l sin θ_1 = 0` (Because equation [1])

we get ;

`(d θ _1)/( d t) = θ _2----eq^1`

`(d θ _2) /(d t )= − b/ m θ_ 2 − g/ l sin θ _1----eq^2`

matrix form

`d/ (d t) (( θ_ 1 ),(θ_ 2 )) = [ (θ_ 2),( − b /m θ_ 2 − g /l sin θ _1 ) ]`

Prrogram:

%input
b=0.02;
g=9.81;
l=1;
m=1;
%displacement and angular velocity
theta_0=[0;3];

%time
t_span = linspace(0,100,600);

%ODE
[t,result]=ode45(@(t,theta)ode_fun(t,theta,b,g,l,m),t_span,theta_0);

%plotting 
figure(1)
plot(t,result(:,1))
hold on
plot(t,result(:,2))
title('For Simple Pendulum')
xlabel('Time in sec')
ylabel('Angular Velocity(rad/s) / Displacement(rad)')
legend('Displacement','Angular velocity')

%Animation of pendulum
amplitude=result(:,1);
ct=1;
for i=1:length(amplitude)
    angle=amplitude(i);
    figure(2)
    title('Motion of Simple Pendulum')
    
    %fx point of pendulum
    x0=0;
    y0=0;
    
    x1=l*sin(angle);
    y1=-l*cos(angle);
   
    %plot
    plot([-1.1 1.1],[0 0],'linewidth',2);
    hold on 
    plot([x0 x1],[y0 y1],'color','r','linewidth',2)
    plot(x1,y1,'o','markers',20,'markerfacecolor','g')%mass ball
    title('Motion of Simple Pendulum')
    hold off 
    grid on 
    axis([-1.5 1.5 -1.5 0.5])
    pause(0.0001)
    
    M(ct) = getframe(gcf);
    ct = ct+1;
end

%animation
movie(M)
videofile = VideoWriter('simple_pen.avi','uncompressed avi');
open(videofile)
writeVideo(videofile,M)
close(videofile)

ODE function:

function [dtheta_dt]=ode_fun(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

Result:

Error faced:

Error at line number 51 that unrecognized function (gc) that can be corrected by making (gcf)