Academic Integrity: tutoring, explanations, and feedback — we don’t complete graded work or submit on a student’s behalf.

I am having trouble answering these two parts on Matlab Diff eq: In MATLAB, clic

ID: 3605389 • Letter: I

Question

I am having trouble answering these two parts on Matlab Diff eq:

In MATLAB, click ‘New’ and then ‘Function’. Create the function:
function dydt = pendSystem(t, y, g, L)
% t is time, g is gravitational constant
% L is pendulum length
% y is the column vector (theta, thetaPrime)’
% the output (dydt) is the column
% vector (dydt(1), dydt(2))’
theta = y(1);
thetaPrime = y(2);
dydt = zeros(2,1); %create a 2x1 zero vector
dydt(1) = thetaPrime;
dydt(2) = -g / L * sin(theta);
end

In MATLAB, click ‘New’ and then ‘Script’. Create the script:
figure, hold on % creates figure window and holds open
g = 9.8; % gravitational acceleration constant (m/sec^2)
L = 9.8; % length of pendulum in meters
thetaVelocity = 0; % initial theta' value (rad/sec)
tspan = [0 11*pi]; % time span
% ode45 numerically solves the system y' = pendSystem(t, y)
% the output tfor is an nx1 column vector of times
% the output yfor is a nx2 matrix whose rows correspond to times
% and whose columns correspond to theta and theta'
thetaPosition = 0.5; % initial theta position in radians
[tfor, yfor] = ode45(@(t,y) pendSystem(t, y, g, L), tspan,
[thetaPosition, thetaVelocity]);
plot(tfor, yfor(:,1))
% plot points using first column of yfor (theta-values)
thetaPosition = 1; % initial theta position in radians
[tfor, yfor] = ode45(@(t,y) pendSystem(t, y, g, L), tspan,
[thetaPosition, thetaVelocity]);
plot(tfor, yfor(:,1))
thetaPosition = 3; % initial theta position in radians
[tfor, yfor] = ode45(@(t,y) pendSystem(t, y, g, L), tspan,
[thetaPosition, thetaVelocity]);
plot(tfor, yfor(:,1))
title('Period of a Pendulum verses Amplitude')
ylabel('Amplitude')
xlabel('time')
grid on
grid minor
ax = gca;
ax.GridAlpha = 0.8;
ax.MinorGridAlpha = 0.6;
legend('A = 0.5','A = 1.0', 'A = 3.0');
hold off

Part 5: Answer the following questions using the MATLAB solution to equation (1)
Using the graph from Part 4, what are the approximate periods corresponding to each
amplitude: A 0 0.5, 1.0, 3.0 ? Does the period depend on the amplitude?
How do these periods compare with the periods from Part 3?
Give the percent relative error for each. (|true – approx.|/|true| * 100).
Does the period depend on mass?
Change the program to see what happens to the period if g is 4 times larger?
What happens to the period if L is 4 times larger?

Theory Assume a pendulum consists of a mass m is suspended from a ceiling by a rigid massless rod of constant length L rotating about a pivot point (where is the angular displacement in radians measured from vertica) Assume there is no damping or other external forces except for gravity (where g is acceleration due to gravity) Recall that the length s of an arc of a circle cut by an angle (in radians) is given by s-r6(here it would be sL). After differentiating twice with respect to time we get the tangential acceleration sde So by Newton's second law,the tangential force is F But we also dt, 2. . Newton dt dt see that the tangential force is F=-m sin (the negative sign is because the direction of the force is opposite of , i.e. when > 0 the force acts to decrease and vis-versa). So,-mg sin = mL which after dividing by m and rearranging yields dt The equation (1) is not linear and cannot be solved using elementary functions. However for small e-values, sin ~ and we can approximate the differential equation (1) by the linear differential equation de g

Explanation / Answer

function dydt = pendSystem(t, y, g, L)

% t is time, g is gravitational constant

% L is pendulum length

% y is the column vector (theta, thetaPrime)’

% the output (dydt) is the column

% vector (dydt(1), dydt(2))’

theta = y(1);

thetaPrime = y(2);

dydt = zeros(2,1); %create a 2x1 zero vector

dydt(1) = thetaPrime;

dydt(2) = -g / L * sin(theta);

end

In MATLAB, click ‘New’ and then ‘Script’. Create the script:

figure, hold on % creates figure window and holds open

g = 9.8; % gravitational acceleration constant (m/sec^2)

L = 9.8; % length of pendulum in meters

thetaVelocity = 0; % initial theta' value (rad/sec)

tspan = [0 11*pi]; % time span

% ode45 numerically solves the system y' = pendSystem(t, y)

% the output tfor is an nx1 column vector of times

% the output yfor is a nx2 matrix whose rows correspond to times

% and whose columns correspond to theta and theta'

thetaPosition = 0.5; % initial theta position in radians

[tfor, yfor] = ode45(@(t,y) pendSystem(t, y, g, L), tspan,

[thetaPosition, thetaVelocity]);

plot(tfor, yfor(:,1))

% plot points using first column of yfor (theta-values)

thetaPosition = 1; % initial theta position in radians

[tfor, yfor] = ode45(@(t,y) pendSystem(t, y, g, L), tspan,

[thetaPosition, thetaVelocity]);

plot(tfor, yfor(:,1))

thetaPosition = 3; % initial theta position in radians

[tfor, yfor] = ode45(@(t,y) pendSystem(t, y, g, L), tspan,

[thetaPosition, thetaVelocity]);

plot(tfor, yfor(:,1))

title('Period of a Pendulum verses Amplitude')

ylabel('Amplitude')

xlabel('time')

grid on

grid minor

ax = gca;

ax.GridAlpha = 0.8;

ax.MinorGridAlpha = 0.6;

legend('A = 0.5','A = 1.0', 'A = 3.0');

hold off