Week 3 - Solving second order ODEs

https://drive.google.com/file/d/1De3PEgajegZQYvHxh__Jjih9lCM_TY55/view?usp=sharing

%clearing line screen,command line,work space and closing all plots
clear all
close all
clc

%inputs
b=0.05;
g=9.81; % specific gravity in m/s^2
l=1; %length of string in m
m=1; % mass in kg

% initial condition
% for solving first order ode it is required value of displacement and velocity
% initial condition of pendulum is 0
% initial velocity is 3m/s
theta_0=[0;3];

%time intervals
t_span=linspace(0,20,300);

% solving ODE
%[time,results]=ode45('challenge_2nd_ode',t_spam,theta_0);
[t,results]=ode45(@(time,theta)ode_function(time,theta,b,g,l,m),t_span,theta_0);
% plotting
subplot(5,4,[13,14,17,18]);
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(:,2));
ylabel('position-radian');
xlabel('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');
  
  %  r=0.075; % plotting circle  and then mapping the center of string
  %  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;
  
  plot (x1,y1,'-o','markers',20,'markerfacecolor','r'); %plotting the bob using marker and facecolor
  
%  r=0.25; % plot a curved rectangle snd color its inner
%  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 1]);
pause(0.0003);
hold off
##  q(ct)=getframe(gcf);
Text=sprintf("pendulum%d.png",ct);
saveas(1,Text)
ct=ct+1;

end 
% movie(q)
% video file=videowriter('pendulum_damped.avi','uncompressed avi');
% open giffile
% writegif(giffile,q)
%close(giffile)


% function of pendulum oscillation with for loop
function[dtheta_dt]=pendulum_damped_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

%clearing line screen,command line,work space and closing all plots
clear all
close all
clc

%inputs of pendulum geometric
b=0.05;
g=9.81; % specific gravity in m/s^2
l=1; %length of string in m
m=1; % mass in kg

% initial condition
% for solving first order ode it is required value of displacement and velocity
% initial condition of pendulum is 0
% initial velocity is 3m/s
theta_0=[0;3];

%time intervals
t_span=linspace(0,20,300);

% solving ODE
%[time,results]=ode45('challenge_2nd_ode',t_spam,theta_0);
[t,results]=ode45(@(time,theta)ode_function(time,theta,b,g,l,m),t_span,theta_0);
% plotting
subplot(5,4,[13,14,17,18]);
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(:,2));
ylabel('position-radian');
xlabel('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');

% r=0.075; % plotting circle and then mapping the center of string
% 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;

plot (x1,y1,'-o','markers',20,'markerfacecolor','r'); %plotting the bob using marker and facecolor

% r=0.25; % plot a curved rectangle snd color its inner
% rectangle('position',[(x1-(r/2)),(y1-(r/2)), r,r],'curvature',[1,1],'facecolor','r',linewidth,',3)
% hold on
% axis equal

To save output of pendulum images as png

grid on;
axis([-1.5 1.5 -1.5 1]);
pause(0.0003);
hold off
## q(ct)=getframe(gcf);
Text=sprintf("pendulum%d.png",ct);
saveas(1,Text)
ct=ct+1;

end

To make gif file the command as shown
% movie(q)
% video file=videowriter('pendulum_damped.avi','uncompressed avi');
% open giffile
% writegif(giffile,q)
%close(giffile)