Compare results of different step size using Euler's method

15 views (last 30 days)
x1=zeros(10000,1);
x2=zeros(10000,1);
x3=zeros(10000,1);
Time=zeros(10000,1);
for i=1:10001
t=(i-1)*0.0005;
Time(i)=t;
if t==0
x1(i)=0;
x2(i)=0;
x3(i)=2;
else
x1(i)=x1(i-1)+x2(i-1)*0.0005;
x3(i-1)=2-9*sin(x1(i-1))-x2(i-1);
x2(i)=x2(i-1)+x3(i-1)*0.0005;
end
Xa=x1*180/pi;
plot(Time,Xa,'red')
I can get the plot for one time interval but I want to get it for different time intervals such as at 0.001, 0.01,0.1 and also compare the results on the same plot. how do i do this? please help??

Accepted Answer

ME
ME on 30 Oct 2019
I have made a few tweaks to your code and it will now allow you to alter one parameter value (the time step size) and produce a new approximation that will be added to your previous plot.
clear all
h = 0.00001;
x1=zeros(round(5/h),1);
x2=zeros(round(5/h),1);
x3=zeros(round(5/h),1);
Time=zeros(round(5/h),1);
for i=1:(5/h)+1
t=(i-1)*h;
Time(i)=t;
if t==0
x1(i)=0;
x2(i)=0;
x3(i)=2;
else
x1(i)=x1(i-1)+x2(i-1)*h;
x3(i-1)=2-9*sin(x1(i-1))-x2(i-1);
x2(i)=x2(i-1)+x3(i-1)*h;
end
end
Xa=x1*180/pi;
hold on
plot(Time,Xa)
Here h is your time step and I have just modified everything else to produce the correct length arrays for storing the results. I also took out the 'red' from your plotting command so that each approximation will automatically plot in a different colour.
  2 Comments
ME
ME on 30 Oct 2019
Edited: ME on 30 Oct 2019
I guess you'd have to adjust this to something like:
h = sort([0.00001 0.005 0.01],'descend');
for step = 1:numel(h)
x1=zeros(round(5/h(step)),1);
x2=zeros(round(5/h(step)),1);
x3=zeros(round(5/h(step)),1);
Time=zeros(round(5/h(step)),1);
for i=1:(5/h(step))+1
t=(i-1)*h(step);
Time(i)=t;
if t==0
x1(i)=0;
x2(i)=0;
x3(i)=2;
else
x1(i)=x1(i-1)+x2(i-1)*h(step);
x3(i-1)=2-9*sin(x1(i-1))-x2(i-1);
x2(i)=x2(i-1)+x3(i-1)*h(step);
end
end
Xa=x1*180/pi;
hold on
plot(Time,Xa)
end
Here the sort is to make sure that you start with the largest step size, otherwise you'll run into issues when you come to plotting.

Sign in to comment.

More Answers (1)

Lucky
Lucky on 30 Oct 2019
Thanks a lot again, but it still doesn't solve my problem.
I need all the plots on a single figure.
clear all
h = sort([0.00001 0.005 0.01 1],'descend');
for step = 1:numel(h)
x1=zeros(round(5/h(step)),1);
x2=zeros(round(5/h(step)),1);
x3=zeros(round(5/h(step)),1);
Time=zeros(round(5/h(step)),1);
for i=1:(5/h(step))+1
t=(i-1)*h(step);
Time(i)=t;
if t==0
x1(i)=0;
x2(i)=0;
x3(i)=2;
else
x1(i)=x1(i-1)+x2(i-1)*h(step);
x3(i-1)=2-9*sin(x1(i-1))-x2(i-1);
x2(i)=x2(i-1)+x3(i-1)*h(step);
end
end
end
Xa=x1*180/pi;
figure(1)
subplot(2,1,1)
plot(Time,Xa)
hold on
plot(Time,Xa)
hold off
  3 Comments
Jan Smid
Jan Smid on 15 Mar 2021
Hello, may I ask what is the reasoning behind the "5/h(step)" while predefining the correct length arrays? It's not clear to me what exactly does this part of the code do.

Sign in to comment.

Categories

Find more on Matrices and Arrays 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!

Translated by