ODE simple pendulum

• Program to simulate the transient behaviour of a simple pendulum and to create an animation of it\'s motion.

clear all
close all
clc
%damping coefficient
b=0.05
%gravity
g=9.81
%length of pendulum
l=1
%mass of ball
m=1
%initial condition
theta_0=[0;3];

%time points
t_span=linspace(0,10,500)
  
%solving ODE
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)*sin(theta1);
  dtheta_dt=[dtheta1_dt;dtheta2_dt];
end

[t,results]=ode45(@(t,theta) ode_func(t,theta,b,g,l,m),t_span,theta_0);
figure(1);
plot(t,results(:,1),\'color\',\'r\');
hold on
plot(t,results(:,2),\'color\',\'b\');
hold off
figure(2);
plot(results(:,1),results(:,2));
xlabel(\'position-radian\');
ylabel(\'angular velocity\');


figure(3);

ct=1
%loop to create different position of pendulum
for i=1:length(results(:,1));
  x0=0
  y0=0
  x1=l*sin(results(i,1));
  y1=l*cos(results(i,1));
  
  %support for pendulum
  plot([-1 1],[0 0],\'linewidth\',10,\'color\', \'black\')
  hold on
  %plotting pendulum rop by using coordinate
  plot([x0 -x1],[y0 -y1],\'linewidth\',5,\'color\',\'b\')
   %creatting the pendulum ball using markersize
  plot(-x1,-y1,\'-o\',\'markers\',20,\'markerfacecolor\',\'r\');

  %defining axis
  axis([-1.5 1.5 -1.5 1]);
  pause(0.001);
  hold off
  ct=ct+1;
text=sprintf(\"pendulum_animation_ode%d.png\",ct);
saveas(gca,text);
end

• Programing for function used.

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

• Link for the animation.

 https://youtu.be/-STzVjzxrso

• Explanation
• Starting the programe by using the command (clear all,close all,clc) for getting A fresh start.
• Giving the initial input that are already given to us such as damping coefficient,mass,length,gravity.
• We have taken the value of theta as of choice shown in the programe.
• Defining the function as given the ODE equation which is in the quadratic form so we break the qudratic equations in different form such as theta1 and theta2 as shown in the formula.

(d^2 θ)/(dt^2 )+b/m dθ/dt+g/l sin⁡θ=0

• The equation is been resolved to theta1 and theta2 shown in the function.
• As for solving this we use ODE45 which is inbuilt function of octave and has its own syntex that we have followed by typing in command window help(\'ode45\'),we got our syntex that we have tu put the values.
• As we have 2degree equation we have to plot the two results that are shown in the programe.
• output of that two plot is shown bellow.

• Now we plot the relation between angular velocity and the position radian.

• For the pendulum section we just created the for loop
• Whear we just given the values to our variable such as (x0,x1,y0,y1).
• Plot the given variables for our animation of pendulum.
• Give support to the pendulum.
• PLoting pendulum rop by using coordinates.
• Creating the pendulum ball using marker size
• We have ploted these things and defined the axis which is most important thing.
• Give the pause command for the better animation.
• We get our desired animation of pendulum while runing this.
• Small picture of this is shown bellow.

• As we have 500 plots so to save those plots into animation we uses the following commands (text,sprintf,saveas)for saving each plot at every angle for getting the desired animation.
• ERRORS.

• In octave we used the ode45 or lsode for the programing of code so we used the ode45 instead of ode.

• In this we have not competed the sytex as we miss the bracket of results so put the bracket.