Week 3 - Solving second order ODEs

Aim:
to write a program in Matlab that solves the ODE.

Program:

% Program to solve 2nd order ODEs

  clear all
  close all
  clc

% inputs for the program
  b=0.05;
  g=9.81;
  l=1;
  m=1;

% initial condition
  theta_0=[0;6];

% time point
  t_span=linspace(0,10,500);

% solving the ODE
  [t,results]=ode45(@(t,theta)ode_func(t,theta,b,g,l,m),t_span,theta_0);

% Plotting
  plot(t,results(:,1))
  hold on
  plot(t,results(:,2))

  xlabel('time')
  ylabel('displacement')

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

% Program by simeon

Program for animation:

% Program for animation

  clear all
  close all
  clc

% inputs
  l=1;
  g=9.81;
  m=1;

% initial condition
  Theta_0=[0;6]

% time
  t_span=linspace(0,20,500);
  [t,results]=ode45('challenge_2nd ode',t_span,theta_0);

[t,results]=ode45(@(time,theta0) pendulum_damped_ode(time,theta,g,l,m),t_span,theta_0);
subplot(5,4,[13,14,17,18]);

%plotting
  plot(t,results(:,1),'color','b');
  hold on;
  plot(t,results(:,2),'color','r');
  hold off;
  subplot(5,4,[15,16,19,20]);
  plot(Results(:,1),Results(:,27));
  x label ('position_Radians');
  y label ('Angular_velocity');
  ct=1;

for i=1:length(results(:,1))
  x0=0;
  y0=0;
  x1=1*sin(results(i,1));
  y1=1*cos(results(i,1));
  subplot(5,4,[1,2,3,4,5,6,7,8,9,10,11,12]);
  plot([-1 1],[0,0],'linewidth',4,'color','k');
  hold on
  line([x0 x1],[y0 y1],'linewidth',2,'color','b');
  hold on;

% plotting circle and then mapping the center to the string
  r=0.075;
  cir_deg=linspace(0,2*pi,1000);
  xc=r*cos(cir_deg);
  yc=r*sin(cir_deg);
  %plot(x1+xc,y1+yc,'linewidth','2')
  hold off;
  
% plotting the bob using marker and face color
  plot(x1,y1,'-0','markers',20,'markerfacecolor','r');

% plot a curved rectangle and color
  rectangle('position',[x1-(r/2)),(y1-r/2),r,r],'curvature',[1 1],'faceColor','r','linewidth',3)
  hold on
  axis equal;
  grid on
  axis([-1.5 1.5 -1.5 0.5])
  pause(0.0003);
  hold off
  
  q(ct)=get frame(gcf);
  Text=springf(:pendulum%d.png",ct);
  save as(1,text)
  ct=ct+1;
  
end

  Movie(q)
  videofile=VideoWriter('pendulum_damped.avi','uncompressed.avi');
  open(videofile)
  writevideo(videofile.q)
  close(videofile)
  
% program by simeon

Explaination:

• this program is to solve the second order ode using matlab
• first clearing the screen, command line, workspace, and closing all the plots.
• then defining the inputs, lenght =1, gravity=9.81 and mass = 1kg.
• in intial condition theta_0 will vary from 0 to 3.
• and then defining time, with t_span using the linspace command.
• and then solving the ODE with the formula [t,results]=ode45(@(t,theta)ode_func(t,theta,b,g,l,m),t_span,theta_0);
• and then plotting the resulsts.
• now, declaring the function of dtheta_dt=ode_func(t,theta,b,g,l,m)
• same is applicable in the animation part
• but in the animation part we use new commands like markers, to mark the plot and markerface color to highlight the marken area.
• cir_deg and Rectangle these comands helps in the mapping the plots center to the string.
• and then with getframe, movie(q), video wrtier the matlab will present an animation.

Errors:

• Errors	Correction
  dtheta2_dt=(b/m)*theta2-(g/l)*sin(theta1)	dtheta2_dt=-(b/m)*theta2-(g/l)*sin(theta1)
  rectangle('piston',[x1-(r/2)),(y1-r/2),r,r],'curvature',[1 1],'face color','r','linewidth',3)	rectangle('position',[x1-(r/2)),(y1-r/2),r,r],'curvature',[1 1],'face color','r','linewidth',3)
  rectangle('piston',[x1-(r/2)),(y1-r/2),r,r],'curvature',[1 1],'face color','r','linewidth',3)	rectangle('piston',[x1-(r/2)),(y1-r/2),r,r],'curvature',[1 1],'faceColor','r','linewidth',3)

Result:
Thus the Matlab progrom to solve the 2nd order ODE is successfully completed and output is verified.