Question

1. Consider the equation of motion governing large deflections of a simple pendulum do + mulsin Mu. ED) where m-mass of the bob - Ikg. I-length, c-damping constant, acceleration due to gravity, MO external torque, e - angular deflection, and t-time. MO) (a) Linearize the equation for small angular displacements using the approximation, sin = 0. Find the solution of the linearized equation for the following conditions: 1=1 4,3 =981./s.c=M,0)=0,00 = 0) =10", and 4 = 0) = 0 (b) Find the solution of the linearised equation for the following conditions: de 1-1 m,8-981m/s",e-M/(-0,01-0)-48 (0) 0 (c) Find the solution of the nonlinear equation, Eq. (EI), for the conditions indicated in part (a). (d) Find the solution of the nonlinear equation, E. (El) for the conditions indicated in part (b) (e plot the degree () vs. time (0) figures by using MATLAB for all the above problems

IMPORTANT NOTES:

Using the Classical Fouth-order Runge-Kutta method to solve all the following problems, with step size h = 0.01, and t = [0:1]

Please use MATLAB to solve the problem. Thanks!

Answer 1

The equation is first written in a way required for RK-4 method. The equations (2) and (3) below are finally used in matlab code:

Matlab code below:

clear;
clc;

%independent variable
tfinal = 1;
h = 0.01;
nsteps = tfinal/h;

%Initial conditions
%question (a)
y0 = 10*pi/180;
z0 = 0;
x0 = 0;
yrk_sola = rk4(x0,y0,z0,nsteps,h,1);

%question (b)
y0 = 45*pi/180;
z0 = 0;
x0 = 0;
yrk_solb = rk4(x0,y0,z0,nsteps,h,1);

%question (c)
y0 = 10*pi/180;
z0 = 0;
x0 = 0;
yrk_solc = rk4(x0,y0,z0,nsteps,h,0);

%question (d)
y0 = 45*pi/180;
z0 = 0;
x0 = 0;
yrk_sold = rk4(x0,y0,z0,nsteps,h,0);

%Plotting
tspan = x0:h:tfinal;
plot(tspan,yrk_sola);
hold on
plot(tspan,yrk_solb);
hold on
plot(tspan,yrk_solc, 'LineStyle',"--");
hold on
plot(tspan,yrk_sold,'LineStyle',"--");
hold off
legend('question-(a)','question-(b)','question-(c)','question-(d)','Location','best')
title('Numerical Methods')
xlabel('time, sec')
ylabel('theta, deg')
grid on

%Function f(x,y,z)
function f = function_f(x1,y1,z1)
    f = z1;
end

%Function phi(x,y,z)
function phi = function_phi(x1,y1,z1,order)
    Mt = 0;
    m = 1;
    l = 1;
    c = 0;
    g = 9.81;
    if order==1
        phi = (Mt - c*z1 - m*g*l*y1)/(m*l*l);
    else
        phi = (Mt - c*z1 - m*g*l*sin(y1))/(m*l*l);
    end
end

function yrk_deg = rk4(x0,y0,z0,nsteps,h,order)
    %RK-4 Method
    x = x0;
    yoldrk = y0;
    zoldrk = z0;
    
    yrk = [y0];
    zrk = [z0];

    for n=1:1:nsteps
        
        k1 = h*function_f(x,yoldrk,zoldrk);
        l1 = h*function_phi(x,yoldrk,zoldrk,order);
        k2 = h*function_f(x+0.5*h,yoldrk+0.5*k1, zoldrk+0.5*l1);
        l2 = h*function_phi(x+0.5*h,yoldrk+0.5*k1, zoldrk+0.5*l1,order);
        k3 = h*function_f(x+0.5*h,yoldrk+0.5*k2, zoldrk+0.5*l2);
        l3 = h*function_phi(x+0.5*h,yoldrk+0.5*k2, zoldrk+0.5*l2,order);
        k4 = h*function_f(x+h , yoldrk+k3, zoldrk+l3);
        l4 = h*function_phi(x+h , yoldrk+k3, zoldrk+l3,order);
        yrk(n+1) = yoldrk + (k1+2*k2+2*k3+k4)/6;
        zrk(n+1) = zoldrk + (l1+2*l2+2*l3+l4)/6;
        yoldrk = yrk(n+1);
        zoldrk = zrk(n+1);
        
        x= n*h;
    end
    
    yrk_deg = yrk*180/pi;
end

Output:

Numerical Methods 50 40 question-a) question-(b) question-c) - question-(d) 30 20 10 theta, deg 0 -10 -20 -30 -40 -50 0 0.2 0.8 1 0.4 0.6 time, sec