Finding average step size and slope with ode45

I want to be able to calculate the average step size and slope for each value in the loop.
I first want to plot the log values of dT and TOL.
I believe my average step size (dT) should be the average of 1/T.
t = [0, 5]; % time span
y0 = 1; % initial condition
dy0 = 2; % initial condition first derivative
dydt = @(t,y) [y(2); ((1 - y(1)^2))*y(2) - y(1)];
for j = 1:10
TOL(j) = 1*exp(-1-j);
options = odeset('AbsTol',TOL(j),'RelTol',TOL(j));
[T{j},Y{j}] = ode45(dydt, t, [y0, dy0], options);
dT = mean(1/length(T)); % average step size
plot(log(dT),log(TOL));
hold on;
grid
j = j+1;
end
hold off;
Errors I'm facing:
  • my plot doesn't show up
  • my slope receives the error 'Undefined function 'log' for input arguments of type 'cell'.'
  • not sure if I am thinking of average step size correctly
Side note: Would I be able to generate multiple slopes using polyfit in one graph?

 Accepted Answer

Subscript ‘dT’ as well. The code will be more efficient if you move the plot outside the loop:
t = [0, 5]; % time span
y0 = 1; % initial condition
dy0 = 2; % initial condition first derivative
dydt = @(t,y) [y(2); ((1 - y(1)^2))*y(2) - y(1)];
for j = 1:10
TOL(j) = 1*exp(-1-j);
options = odeset('AbsTol',TOL(j),'RelTol',TOL(j));
[T{j},Y{j}] = ode45(dydt, t, [y0, dy0], options);
dT(j) = mean(1/length(T)); % average step size
end
figure
plot(log(dT),log(TOL));
grid

3 Comments

Changing the name of ‘dT’ isn’t necessary. You just need to subscript it as ‘dT(j)’ to create a vector from it. The plot then shows up for me with the code I posted.
Re-posting:
t = [0, 5]; % time span
y0 = 1; % initial condition
dy0 = 2; % initial condition first derivative
dydt = @(t,y) [y(2); ((1 - y(1)^2))*y(2) - y(1)];
for j = 1:10
TOL(j) = 1*exp(-1-j);
options = odeset('AbsTol',TOL(j),'RelTol',TOL(j));
[T{j},Y{j}] = ode45(dydt, t, [y0, dy0], options);
dT(j) = mean(1/length(T)); % average step size
end
figure
plot(log(dT),log(TOL));
grid
EDIT —
The plot shows up when I run my code, although the plot image doesn’t seem to be appearing when I post it here. I am attaching it as well to work around that problem.
I realized my error.
Thanks again!!
As always, my pleasure!

Sign in to comment.

More Answers (0)

Categories

Find more on Programming in Help Center and File Exchange

Products

Release

R2016a

Tags

Asked:

EB
on 9 Oct 2019

Commented:

on 12 Oct 2019

Community Treasure Hunt

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

Start Hunting!