Matlab ode45

darii · #1 $Old$ October 24th, 2009, 07:23 AM

Hello!
i want to resolve this equation in matlab:
y'(1)=-y(2)
y'(2)=y(1)

I have do it with:
%%%%%%%%
function dy=test1(t,y)

dy1=-y(2);
dy2=y(1);
dy=[dy1;dy2];
%%%%%%%%%%%
t0=0:0.0001:2*pi;
t0=double(t0);
for x1=1:5
x0=[0,x1];
x0=double(x0);

[t,y]=ode45('test1',t0,x0)
t=double(t);
y=double(y);

hold on
plot(y(:,1), y(:,2))
title('Gráfica general')
grid on
end

But now i have two problems:
First i want stop the plot when i give a time, for example t=pi, but if i write t0=[0,pi] the plot is wrong.
After, i want that when y(2)=0 the plot stop, and it will appear a semicircle, but i can't do it :S.

Thanks!

elbarto · #2 $Old$ October 24th, 2009, 08:24 AM

Try use logical indexing. I am not familiar with ODE45 so there might be a nicer way to do this but in terms of plotting you could try this.

Code:

clear;clc;cla;
t0=0:0.01:2*pi;
t0=double(t0);
hold on
for x1=1:5
    x0=[0,x1];
    x0=double(x0);
    
    [t,y]=ode45('test1',t0,x0);
    t=double(t);
    y=double(y);
    
    Lind = (y(:,2) < 0 & t < pi);%use logical indexing
    plot(y(Lind,1), y(Lind,2))
end
title('Gráfica general')
grid on
hold off
axis equal

The plot I get is only a quadrant tho so you might need to play with the code a bit. I did not fully understand your question so I hope that helps a bit.

Elbarto

elbarto · #3 $Old$ October 24th, 2009, 08:28 AM

also there is no need to use "double" when defining variables. All variables typed into the workspace are doubles by default. You can check this with the "whos" command.

Eg:

Code:

EDU>> a = 1.06

a =

    1.0600

EDU>> b = 1

b =

     1

EDU>> whos
  Name      Size            Bytes  Class     Attributes

  a         1x1                 8  double              
  b         1x1                 8  double              

EDU>>

Elbarto

darii · #4 $Old$ October 24th, 2009, 08:52 AM

But...this is a "small trap", no?
I mean, i want to plot the graphic until y(2)=0...but here you say that i plot all the graphic and only i show one quadrant, no?

darii · #5 $Old$ October 24th, 2009, 11:30 AM

I have other problem...(i do it how you say me, thanks!).
I have other differential equations:
dx = y(4)+y(2);
dy = y(5)-y(1);
dz = y(6);
dpx = y(5) - ((1-mu)/(r1.^3))*(y(1)-mu) - (mu/(r2.^3))*(y(1)-mu+1);
dpy = -y(4) -((1-mu)/(r1.^3) + mu/(r2.^3))*y(2);
dpz = -((1-mu)/(r1.^3) + mu/(r2.^3))*y(2);

and i want to resolve the system, i do:
[t,y]= ode45('test2',tf, cond_initials);

but, i want plot this..."y" have 6 columns...i want to do something like that:
plot(y(:,1),y(:,2),y(:,3),y(:,4),y(:,5),y(:,6))

but it is wrong...what can i do?

thank you very much!

Mush · #6 $Old$ October 24th, 2009, 11:38 AM

Quote:

Originally Posted by darii $View Post$

Hello!
i want to resolve this equation in matlab:
y'(1)=-y(2)
y'(2)=y(1)

I have do it with:
%%%%%%%%
function dy=test1(t,y)

dy1=-y(2);
dy2=y(1);
dy=[dy1;dy2];
%%%%%%%%%%%
t0=0:0.0001:2*pi;
t0=double(t0);
for x1=1:5
x0=[0,x1];
x0=double(x0);

[t,y]=ode45('test1',t0,x0)
t=double(t);
y=double(y);

hold on
plot(y(:,1), y(:,2))
title('Gráfica general')
grid on
end

But now i have two problems:
First i want stop the plot when i give a time, for example t=pi, but if i write t0=[0,pi] the plot is wrong.
After, i want that when y(2)=0 the plot stop, and it will appear a semicircle, but i can't do it :S.

Thanks!

What is it you want to plot in the first instance? Is it the function against time? Because if so, that's not what you're doing. You're plotting the results of the two equations against each other when you should be plotting the result of the appropriate equation against time.

Secondly, I'm not sure that you can stop ODE45 using a particular condition. ODE45 runs through the time interval that you give it, and I'm not sure you can stop it given a particular condition without editing the code for ODE45.

Mush · #7 $Old$ October 24th, 2009, 11:41 AM

Quote:

Originally Posted by darii $View Post$

I have other problem...(i do it how you say me, thanks!).
I have other differential equations:
dx = y(4)+y(2);
dy = y(5)-y(1);
dz = y(6);
dpx = y(5) - ((1-mu)/(r1.^3))*(y(1)-mu) - (mu/(r2.^3))*(y(1)-mu+1);
dpy = -y(4) -((1-mu)/(r1.^3) + mu/(r2.^3))*y(2);
dpz = -((1-mu)/(r1.^3) + mu/(r2.^3))*y(2);

and i want to resolve the system, i do:
[t,y]= ode45('test2',tf, cond_initials);

but, i want plot this..."y" have 6 columns...i want to do something like that:
plot(y(:,1),y(:,2),y(:,3),y(:,4),y(:,5),y(:,6))

but it is wrong...what can i do?

thank you very much!

I wouldn't recommending attempting to plot in 6 dimensions lol :S. The plot command is for plotting in 2 dimensions. You must use mesh, contour, surf or plot3 to plot in three dimensions... and I don't reckon you're able to plot in anything more than that!

darii · #8 $Old$ October 24th, 2009, 12:34 PM

I have to plot the function projected in (x,y), for that i are plotting the results.
After, i have to integrate until y=0, but i don't know how do it.

And...for the 6 equations, it is the same, but i can't do it :S

darii · #9 $Old$ October 24th, 2009, 01:22 PM

I just saw that with "events" in ode45 could be posible integrate until y=0.
Somebody has used events? because i never did it