Pendulum oscillation simulation using Octave

The following program solves the ODE which represents the equation of motion of a simple pendulum with damping and simulates the pendulum motion using Octave. The initial conditions are length of string =1m,mass of the bob = 1kg,damping coefficient=0.05, and acceleration of gravity=9.81m/s2. The motion of pendulum is simulated for a time period of 20 seconds with initial angular displacement taken as 0 and initial angular velocity taken as 3 rad/sec2.

The main body of the program is given below, for the sake of modularity the ODE is taken as seperate function and is solved using ode45:

clear all
close all
clc

l=1;
m=1;
b=0.05;
g=9.81;
t=linspace(0,20,750);
%y0(1) is initial angular displacement and y0(2) is initial angular velocity
y0=[0 3];


[t,result]=ode45(@(t,teta) second_ode(t,teta,b,g,l,m),t,y0);
result
figure(1)
plot(t,result(:,1))
hold on
plot(t,result(:,2))
xlabel('Time-->')
ylabel('Angular velocity and Angular acceleration-->')
legend('Angular veocity','Angular acceleration')

figure(2)
plot(result(:,1),result(:,2))
xlabel('Angular velocity-->')
ylabel('Angular acceleration-->')

ang_vel=result(:,1);
counter=1;
for i=1:length(ang_vel)
  ang_vel1=ang_vel(i);
  x0=0;
  y0=0;
  
  x1=l*sin(ang_vel1);
  y1=-l*cos(ang_vel1);
  
  figure(3)
  plot([-1 1],[0 0],'linewidth',2,'color','k')
  hold on
  plot([0 0],[1 -1],'linewidth',1,'linestyle','--')
  hold on
  plot([x0 x1],[y0 y1],'linewidth',1.5,'color','r')
  hold on
  plot(x1,y1,'o','markers',20,'markerfacecolor','y')
  hold off
  axis([-1.5 1.5 -1.5 1.5])
  pause(0.0001)
  Text=sprintf("pendulum%d.png",counter);
  saveas(1,Text);
  counter=counter+1;
  
end

The ODE is called as a function into main body of the program, the code of the function of ODE is given below:

function [tetadot]=second_ode(t,y,b,g,l,m)
 tetadot(1)=y(2);
 tetadot(2)=-((b/m)*y(2))-((g/l)*sin(y(1))); 
end

Using ode45 the ODE is solved to obtain the angular velocity and angular displacement which is then plotted against time for a time period of 20 seconds, the resulting plot is given below:

Then the angular velocity is plotted against angular acceleration ,the resulting plot is given below:

The pendulum is simulated using Octave to obtain image file of pendulum at different position of its oscillation as the time advances. The image files are then converted into .avi, the resulting video which shows the oscillations of pendulum is attached below: