ode23s can't solve system

1 vue (au cours des 30 derniers jours)
Ryan
Ryan le 12 Avr 2014
Commenté : Star Strider le 12 Avr 2014
Hi I hope somebody can help me. I have what should be a very simple bit of code that i'm trying to use to solve an ODE. For some reason it won't run. There are no error messages but i've waited for about 3 minutes on a system that should only take a couple of seconds. Any ideas on why this is taking so long would be much appreciated. I've attached the solver and below is the code for my system of ODEs.
%antODE
function dYdt = antODE(t,Y);
gamma = 5;
mu = .6;
tau = .2;
c = .1;
dYdt = zeros(3,1);
T = Y(1);
F = Y(2);
L = Y(3);
dYdt(1) = (gamma/c)*(T - L -tau*F);
dYdt(2) = (1/(1-c))*(-F + L*T);
dYdt(3) = ((1-c)/mu)*F - (c/mu)*L;
end
  2 commentaires
Walter Roberson
Walter Roberson le 12 Avr 2014
What parameters are you calling ode23s with?
Ryan
Ryan le 12 Avr 2014
I tried to attach the file with the parameters I was using but I'm not seeing that here. The initial condition was an eigenvector tangent to the unstable manifold at a fixed point. Given Star Strider's comment on how fast Y(1) goes to negative infinity it makes sense that this wasn't working.

Connectez-vous pour commenter.

Réponses (1)

Star Strider
Star Strider le 12 Avr 2014
Modifié(e) : Star Strider le 12 Avr 2014
This works fine for me:
gamma = 5; mu = .6; tau = .2; c = .1;
antODE = @(t,y) [(gamma/c).*(y(1) - y(3) -tau.*y(2)); (1/(1-c)).*(-y(2) + y(3).*y(1)); ((1-c)/mu).*y(2) - (c/mu).*y(3)];
[t, y] = ode23s(antODE, [0 0.5], [0.1 0.1 0.1]);
figure(1)
plot(t, y)
legend('y_1', 'y_2', 'y_3', 'Location', 'SouthWest')
grid
You need to be careful with your time interval, though, because y(1) takes off to negative infinity from the outset, and reaches -1.7E+9 at t=0.5, while y(2) and y(3) behave themselves. (The output is huge, 11813 rows at t=0.5, and takes 8.102798 seconds.) I suspect the behaviour of y(1) may be the problem you are having.
  2 commentaires
Ryan
Ryan le 12 Avr 2014
This is very helpful thank you
Star Strider
Star Strider le 12 Avr 2014
My pleasure!

Connectez-vous pour commenter.

Tags

Produits

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!

Translated by