Week 3 - Solving second order ODEs
AIM:
To write a program that solves the following ODE. This ODE represents the equation of motion of a simple pendulum with damping.
In the above equation,
g = gravity in m/s2,
L = length of the pendulum in m,
m = mass of the ball in kg,
b=damping coefficient.
Input Values
L=1 meter,
m=0.1 kg,
b=0.05.
g=9.81 m/s2.
Convert 2nd degree ODE order into first degree ODE
`theta=theta1`
`dtheta/dt=dtheta1/dt=theta2`
`dtheta2/dt= -b/m theta2-g/2sintheta1`
Create a function of 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)*sind(theta1);
dtheta_dt=[dtheta1_dt;dtheta2_dt];
endCreate a Simple pendulum program by using ODE function
Initial conditions are provided:
-motion between 0-20 sec
- 0 is angular displacement and 0 is the angular velocity
-Solving the ODE equation by using the ode45 function
close all
clear all
clc
b=0.5;
g=9.81;
l=1;
m=0.1;
% initial conditions
theta_0=[0;5];
% time points
t_span=linspace(0,10,500);
% slove ODE
[t,results]=ode45(@(t,theta) ode_func(t,theta,b,g,l,m),t_span,theta_0);
% plotting
figure(1)
hold on
plot(t,results(:,1))
plot(t,results(:,2))
ylabel('plot')
xlabel('time')
hold off
%Animation
x0=0
y0=0
x1=-2
y1=0
x2=2
y2=0
ct=1
for i = 1:length(t)
x=l*sin(results(i));
y=-l*cos(results(i));
plot([x1,x0,x2],[y1,y0,y2],'linewidth',10)
hold on
plot([x0,x],[y0 y],'linewidth',4)
hold on
plot(x,y,'o','markersize',15,'markerfacecolor','b')
axis([-3 3 -3 0])
hold off
pause(0.1)
d(ct)=getframe(gcf);
ct=ct+1;
end
movie(d)
videofile=VideoWriter('Simple_pendulam.avi','uncompressed AVI');
open(videofile)
writeVideo(videofile,d)
close(videofile)
RESULTS:
Plotting:
Plotting of a figure has been done for plot & results
Animattion:
Conclusion:
Hence, by using MATLAB we can write a program t0 solves the ODE FUNCTION and created an animation of Simple pendulum motion.