ODE45 - How do I set up the conditions so that I'm not getting an error on 0?

1 visualización (últimos 30 días)
I have to choose an appropriate ode solver for y'(t) = te^(-1000*t) - 10*y(t) with y(0) = 0 for 0<=t<=0.5. I also have to choose appropriate values for the absolute and relative tolerance to get a solution. I am getting an error with the code I have tried to start running to plot a vector solution of it. Also, can anyone confirm that the ODE45 solver and tolerance values are the correct choice? Thanks!
a = 0;
b = .5;
y(1) = 0;
f1 = @(t,y) (t*exp(-1000*t)) - (10*y(t));
options = odeset('RelTol',10^-7,'AbsTol',10^-9);
z1 = cputime;
[T Y] = ode45(f1,[a b],0,options);
z2 = cputime;
fprintf('ODE 45 CPU time: %f seconds. Number of steps: %d\n',z2-z1,length(T))
plot(T,Y,'r')
  2 comentarios
Steven Lord
Steven Lord el 11 de Mzo. de 2016
What is the full text of the error you received?
Kaylene Widdoes
Kaylene Widdoes el 11 de Mzo. de 2016
Attempted to access y(0); index must be a positive integer or logical. Error in @(t,y)(t*exp(-1000*t))-(10*y(t))
Error in odearguments (line 87) f0 = feval(ode,t0,y0,args{:}); % ODE15I sets args{1} to yp0. Error in ode45 (line 113) [neq, tspan, ntspan, next, t0, tfinal, tdir, y0, f0, odeArgs, odeFcn, ... >>

Iniciar sesión para comentar.

Respuesta aceptada

Star Strider
Star Strider el 11 de Mzo. de 2016
Editada: Star Strider el 11 de Mzo. de 2016
This runs for me without error:
a = 0;
b = .5;
y(1) = 0;
f1 = @(t,y) (t.*exp(-1000*t)) - (10*y);
options = odeset('RelTol',10^-7,'AbsTol',10^-9);
z1 = cputime;
[T Y] = ode45(f1,[a b],0,options);
z2 = cputime;
fprintf('ODE 45 CPU time: %f seconds. Number of steps: %d\n',z2-z1,length(T))
plot(T,Y,'r')
The only changes I made were to replace ‘y(t)’ with ‘y’ in ‘f1’, and vectorise it.
  2 comentarios
Kaylene Widdoes
Kaylene Widdoes el 11 de Mzo. de 2016
Awesome! Ok - so that looks like a really rigid graph - would I be better off switching to a different solver or tolerance values to accommodate for the stiffness?
Star Strider
Star Strider el 11 de Mzo. de 2016
Thank you!
Yes. I switched it to ode15s and it gave a much better solution. (Sorry for the delay — life intrudes.)

Iniciar sesión para comentar.

Más respuestas (0)

Categorías

Más información sobre Ordinary Differential Equations en Help Center y File Exchange.

Etiquetas

Community Treasure Hunt

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

Start Hunting!

Translated by