Question

MATLAB question!

4. (a) Modify the code, to compute and plot the quantity E = 1/2mv^2 + 1/2ky^2 as a function of time. What do you observe? Is the energy conserved in this case?

(b) Show analytically that dE/dt < 0 for c > 0 while dE/dt > 0 for c < 0.

(c) Modify the code to plot v vs y (phase plot). Comment on the behavior of the curve in the context of the motion of the spring. Does the graph ever get close to the origin? Why or why not?

given code
---------------------------------------------------------------

clear all;   
m = 4; % mass [kg]
k = 9; % spring constant [N/m]
c = 4; % friction coefficient [Ns/m]
omega0 = sqrt(k/m); p = c/(2*m);
y0 =-0.8; v0 = 0.3; % initial conditions
[t,Y] = ode45(@f,[0,15],[y0,v0],[],omega0, p); % solve for 0<t<15
y = Y(:,1); v = Y(:,2); % retrieve y, v from Y
figure(1); plot(t,y,'ro-',t,v,'b+-');% time series for y and v
grid on; axis tight;
%---------------------------------------------------
function dYdt = f(t,Y,omega0,p); % function defining the DE
y = Y(1); v = Y(2);
dYdt=[ v ; -omega0^2*y-2*p*v]; % fill-in dv/dt
end

-----------------------------------------------------------------------------------------

Answer 1

clear all;
m = 4; % mass [kg]
k = 9; % spring constant [N/m]
c = 4; % friction coefficient [Ns/m]
omega0 = sqrt(k/m); p = c/(2*m);
y0 =-0.8; v0 = 0.3; % initial conditions
[t,Y] = ode45(@f,[0,15],[y0,v0],[],omega0, p); % solve for 0<t<15
y = Y(:,1); v = Y(:,2); % retrieve y, v from Y
figure(1); plot(t,y,'ro-',t,v,'b+-');% time series for y and v
title('y(t) and v(t) plot')
xlabel('time (t)')
ylabel('y(t) and v(t) plot')
grid on; axis tight;

%Answering question 4a.
E=((1/2).*m.*v.^2)+((1/2).*k.*y.^2);
figure(2); plot(t,E,'ro-');% time series for E
title(sprintf('E(t) vs. time plot for c=%d',c))
xlabel('time (t)')
ylabel('E(t) plot')
grid on; axis tight;

%Answering question 4a.

figure(5); plot(y,v,'ro-');% time series for E
title(sprintf('Phase space plot for c=%d',c))
xlabel('y(t)')
ylabel('v(t)')
grid on; axis tight;

%Answering question 4b.
E=((1/2).*m.*v.^2)+((1/2).*k.*y.^2);
figure(3); plot(t,E,'ro-');% time series for E
title('E(t) vs. time plot for c>0')
xlabel('time (t)')
ylabel('E(t) plot')

grid on; axis tight;
m = 4; % mass [kg]
k = 9; % spring constant [N/m]
c = -4; % friction coefficient [Ns/m]
omega0 = sqrt(k/m); p = c/(2*m);
y0 =-0.8; v0 = 0.3; % initial conditions
[t,Y] = ode45(@f,[0,15],[y0,v0],[],omega0, p); % solve for 0<t<15
y = Y(:,1); v = Y(:,2); % retrieve y, v from Y

E=((1/2).*m.*v.^2)+((1/2).*k.*y.^2);
figure(4); plot(t,E,'bo-');% time series for E
title('E(t) vs. time plot for c<0')
xlabel('time (t)')
ylabel('E(t) plot')
grid on; axis tight;

fprintf('\nFor question 4a)\n')
fprintf('\tFrom the energy plot it can be shown that the energy is decreasing gradually and E conserved.\n')
fprintf('Due to friction energy is dissipated.\n')

fprintf('\nFor question 4b)\n')
fprintf('\tFor C>0 Energy gradually decreases so dE/dt<0.\n')
fprintf('\tFor C<0 Energy gradually decreases so dE/dt>0.\n')

fprintf('\nFor question 4b)\n')
fprintf('\tFrom phase plot, graph get close to origin\n')
fprintf('As velocity decreases and energy dissipates it gets close to origin.\n\n')
%---------------------------------------------------
function dYdt = f(t,Y,omega0,p) % function defining the DE
y = Y(1); v = Y(2);
dYdt=[ v ; -omega0^2*y-2*p*v]; % fill-in dv/dt
end

%%%%%%%%%%%%%%%%%% End of Code %%%%%%%%%%%%%%%%%%%%