0
$\begingroup$

I am given an example file for solving a second order ODE using ODE 45

function ex_with_2eqs

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 ];

I am supposed to alter this function to solve another equation.

$y''+ 4y'+ 3y = \cos(t)$ with $y(0) = −1; y'(0) = 0$

They gave me the following to define the function

function dYdt = f(t,Y)

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

dYdt = [ v ; cos(t)-4*v-3*y ];

end

So, what I have now is

function LAB04ex1

t0 = 0; tf = 40; y0 = [-1;0];

[t,Y] = ode45(@f,[t0,tf],y0,[]);

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
end

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

function dYdt = f(t,Y)
u1 = Y(1); u2 = Y(2);
dYdt = [ v ; cos(t)-4*v-3*y ];
end

But when I run the file absolutely nothing happens. I think I am way off on my part. I figured the parameters a,b,c,d were no longer needed. I think I changed the right y's to Y's. But overall, I am completely lost.

  • 0
    As a programming question this is probably more appropriate on stack overflow.2012-11-08

1 Answers 1

1

It's because your v and y are undefined in the scope of the function dYdy.

Change that function to this:

function dYdt = f(t,Y)
u1 = Y(1); u2 = Y(2);
dYdt = [ u2 ; cos(t)-4*u2-3*u1 ];
end

and you should be fine.

(v changed to u2, y changed to u1).

Here's figure 1:

enter image description here

Very nice stable LCO behavior there.