Blank figures appearing for plots
5 views (last 30 days)
Show older comments
I've create an adaptive time-stepping function and I'm trying to plot the graphs but it is just creating blank figures.
Here is my function:
function [tvec,yvec,dtvec]=adaptive_euler(f,ytrue,y0,t0,dt0,dtmin,dtmax,tmax,alpha,beta,gamma,TOL)
y=zeros(10^5,1);
y(1)=y0;
dt=zeros(10^5,1);
t=zeros(10^5,1);
dt(1)=dt0;
yt=y0;
yt1=y0;
t(1)=t0;
n=1;
while (t(n)<tmax && dt(n)<dtmax && dt(n)>dtmin)
yt=yt+dt(n)*f(t(n),yt);
yt1=yt1+(dt(n)/2)*f(t(n),yt1);
yt2=yt1+(dt(n)/2)*f(t(n),yt1);
Et2=(yt2-yt)/(dt(n));
c=TOL/norm(Et2);
if norm(Et2) < TOL
y(n)=yt2;
t(n+1)=t(n)+dt(n);
n=n+1;
end
dt(n) = min(max(alpha*c*dt(n-1),gamma*dt(n-1)),beta*dt(n-1));
end
tvec=zeros(n,1);
yvec=zeros(n,1);
dtvec=zeros(n,1);
for i=1:n
tvec(i)=t(i);
yvec(i)=y(i);
dtvec(i)=dt(i);
end
yreal=zeros(n,1);
error=zeros(n,1);
for i=1:n
yreal(i)=ytrue(tvec(i));
error(i)=abs(yvec(i)-yreal(i));
end
T1=table(n,error(n));
figure(1)
plot(tvec,yvec,'.',tvec,yreal,'--');
legend('Approx:','True:');
title('Graph of the Adaptive Eulers Method approximation of dy/dx=-exp(2*t)*y^2 against real values');
figure(2)
plot(tvec,dtvec)
title('Graph of the change in dt over t');
end
And here is the script:
f=@(t,y)-exp(2*t)*y^2;
t0=0;
y0=2;
dt0=0001;
dtmin=10^-10;
dtmax=0.5;
tmax=10;
alpha=0.84;
beta=4;
gamma=0.1;
TOL=0.1;
ytrue=@(t)2*exp(-2*t);
[tvec,yvec,dtvec]=adaptive_euler(f,ytrue,y0,t0,dt0,dtmin,dtmax,tmax,alpha,beta,gamma,TOL)
Thank you!
4 Comments
Adam
on 10 Mar 2020
But you set them up as vectors of length n and then you copy from 1 to n from another vector into these.
Answers (1)
Image Analyst
on 10 Mar 2020
You never enter the while loop
while (t(n)<tmax && dt(n)<dtmax && dt(n)>dtmin)
because dt(1) is 1 and dtmax is 0.5. Pass in a different dtmax.
2 Comments
Image Analyst
on 10 Mar 2020
1 is also not less than 0.001, so that won't enter the loop either. So what did you pass in?
In the meantime, see the FAQ: Subscript_indices_must_either_be_real_positive_integers_or_logicals.
See Also
Categories
Find more on Title in Help Center and File Exchange
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!