Question

t0 = 0; tf = 20; y0 = [10;60];

a = .8; b = .01; c = .6; d = .1;

[t,y] = ode45(@f,[t0,tf],y0,[],a,b,c,d);

u1 = y(:,1); u2 = y(:,2); % y in output has 2 columns corresponding to u1 and u2

figure(1);

subplot(2,1,1); plot(t,u1,'b-+'); ylabel('u1'); subplot(2,1,2); plot(t,u2,'ro-'); ylabel('u2');

figure(2) plot(u1,u2); axis square; xlabel('u_1'); ylabel('u_2'); % plot the phase plot

%----------------------------------------------------------------------

function dydt = f(t,y,a,b,c,d)

u1 = y(1); u2 = y(2);

dydt = [ a*u1-b*u1*u2 ; -c*u2+d*u1*u2 ];

end

Only Part B please.

B .)By reading the matrix Y and the vector t, find (approximately) the last three values of t in the interval 0 ≤ t ≤ 45 at which y reaches a local maximum. Note that, because the M-file LAB04ex1.m is a function file, all the variables are local and thus not available in the Command Window. To read the matrix Y and the vector t, you need to modify the M-file by adding the line [t, Y(:,1), Y(:,2)]. Do not include the whole output in your lab write-up. Include only the values necessary to answer the question, i.e. just the rows of [t, y, v] with local y-maxima and the adjacent rows. To quickly locate the desired rows, recall that the local maxima of a differentiable function appear where its derivative changes sign from positive to negative. (Note: Due to numerical approximations and the fact that the numerical solution is not necessarily computed at the exact t-values where the maxima occur, you should not expect v(= y 0 ) to be exactly 0 at local maxima, but only close to 0)

Instructions: Just like for all the other lab reports, unless otherwise specified, include in your lab report all M-files, figures, MATLAB input commands, the corresponding output, and hte answers to the questions. 1. (a) Modify the function ex_with_2eqs to solve the IVP (4) for 0 <t < 45 using the MATLAB routine ode45. Call the new function LABO4ex1. Let [t,Y] (note the upper case Y) be the output of ode45 and y and v the unknown functions. Use the following commands to define the ODE: function dYdt= f(t,Y) y=Y(1); v=Y(2); dydt = [v; -4*sin(t)-4*v-3*y]; Plot y(t) and v(t) in the same window (do not use subplot), and the phase plot showing v vs y in a separate window. Add a legend to the first plot. (Note: to display v(t) = y(t), use 'v(t)=y'' (t)'). Add a grid. Use the command ylim([-2.9,2.9]) to adjust the y-limits for both plots. Adjust the c-limits in the phase plot so as to reproduce the pictures in Figure 7. 0 10 20 30 40 -2 -1 0 1 2 Figure 7: Time series y = y(t) and v = v(t) = y(t) (left), and phase plot v = y vs. y for (4). (b) By reading the matrix Y and the vector t, find (approximately) the last three values of t in the interval 0 <t < 45 at which y reaches a local maximum. Note that, because the M-file LAB04ex1.m is a function file, all the variables are local and thus not available in the Command Window. To read the matrix Y and the vector t, you need to modify the M-file by adding the line [t, Y(:,1), Y(:,2)]. Do not include the whole output in your lab write-up. Include only the values necessary to answer the question, i.e. just the rows of [t, y, v] with local y-maxima and the adjacent rows. To quickly locate the desired rows, recall that the local maxima of a differentiable function appear where its derivative changes sign from positive to negative. (Note: Due to numerical approximations and the fact that the numerical solution is not necessarily computed at the exact t-values where the maxima occur, you should not expect v=y) to be exactly 0 at local maxima, but only close to 0). (c) What seems to be the long term behavior of y? (d) Modify the initial conditions to y(0) = 1.6, v(0) = 1.2 and run the file LABO4ex1.m with the modified initial conditions. Based on the new graphs, determine whether the long term behavior of the solution changes. Explain. Include the pictures with the modified initial conditions to support your answer.

Answer 1

8&Matlab code for solving ode clear all close all #Answering question b. Initial conditions for ode yO=(1.6;1.2]; {minimum and maximum x tspan=[0 45]; Solution of ODEs using ode45 matlab function soll= ode45(e(t,y) odefcni(t,y), tspan, y0); Equally splitting x into small intervals tl = linspace(tspan(1),tspan (end),501); tyy is the corresponding x y v and z Yy1 = deval(soli,t1); figure (1) hold on plot (t1, yyl(1,:), 'linewidth', 2) plot(ti,yy1(2,:), 'linewidth', 2) ylim( [-2.9 2.9]) title('Plot for y(t) and v(t) vs. t solution') xlabel('t') ylabel('y(t),v=y'') legend ('y(t)', 'v(t)) box on grid on figure (2) plot(yyl(1,:),yy1(2, :), 'linewidth', 2) xlabel('y') ylabel('v=y''') title('Phase potrait plot') grid on; box on; ylim([-2.9 2.9]); xlim([-2.9 2.9]) Function for evaluating the ODE function dydt = odefcnl(t,y) eql= y(2); eq2= -4.*sin(t)-4.*Y(2)-3. *Y(1); Evaluate the ODE for our present problem dydt = [eql; eq2]; end ----End of Code----

Plot for y(y and v(t) vs. t solution FLA A yt),v=y" 1 VA HA 5 10 15 20 25 30 35 40 45 Phase potrait plot V=y' -2.5 -2 -1.5 -1 -0. 5 0 0.5 1 1.5 2 25 у

%%Matlab code for solving ode
clear all
close all
%Answering question b.
%Initial conditions for ode
y0=[1.6;1.2];
      
        %minimum and maximum x
        tspan=[0 45];
        %Solution of ODEs using ode45 matlab function
        sol1= ode45(@(t,y) odefcn1(t,y), tspan, y0);
        %Equally splitting x into small intervals
        t1 = linspace(tspan(1),tspan(end),501);
        %yy is the corresponding x y v and z
        yy1 = deval(sol1,t1);
        figure(1)
        hold on
        plot(t1,yy1(1,:),'linewidth',2)
        plot(t1,yy1(2,:),'linewidth',2)
        ylim([-2.9 2.9])
        title('Plot for y(t) and v(t) vs. t solution')
        xlabel('t')
        ylabel('y(t),v=y''')
        legend('y(t)','v(t)')
        box on
        grid on
        figure(2)
        plot(yy1(1,:),yy1(2,:),'linewidth',2)
        xlabel('y')
        ylabel('v=y''')
        title('Phase potrait plot')
        grid on; box on;
        ylim([-2.9 2.9]); xlim([-2.9 2.9])
%Function for evaluating the ODE
function dydt = odefcn1(t,y)

    eq1= y(2);
    eq2= -4.*sin(t)-4.*y(2)-3.*y(1);
  

    %Evaluate the ODE for our present problem
    dydt = [eq1;eq2];
  
end
%----------------------------End of Code-----------------------------%