Blank figures appearing for plots
3 visualizaciones (últimos 30 días)
Mostrar comentarios más antiguos
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 comentarios
Adam
el 10 de Mzo. de 2020
But you set them up as vectors of length n and then you copy from 1 to n from another vector into these.
Respuestas (1)
Image Analyst
el 10 de Mzo. de 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 comentarios
Image Analyst
el 10 de Mzo. de 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.
Ver también
Categorías
Más información sobre Title en Help Center y File Exchange.
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!