15 views (last 30 days)

when i run the code below, the plots do not look correct to me. As when i make the time step smaller the eulers should get more accurate, thus closer matching the exact solution? But when i do this they eulers and exact solution seem to get further apart on the plot.

clear;

h=0.05; % time step

t=0:h:4; %time range

y=zeros(size(t)); % set y array all 0's to same size as t

y(1)=2; % set initial condition at time 0 to 2

n=numel(y);

dydt = 4*exp(0.8*t) - 0.5*y;

exact_sol=(4/1.3)*(exp(0.8*t)-exp(-0.5*t))+2*exp(-0.5*t); %This is the exact solution to dy/dt

for i=1 : n-1 %for loop to interate through y values for

y(i+1)= y(i)+ h * dydt(i); % the Euler method

end

plot(t,y) %plot Euler

hold on

plot(t,exact_sol,'red'); % plots the exact solution to this differential equation

legend('Euler','Exact'); % adds a legend

grid on

John D'Errico
on 26 Nov 2020

Edited: John D'Errico
on 26 Nov 2020

Sorry, but your code is wrong (so is your exact solution.) Why? Look at what you wrote.

h=0.05; % time step

t=0:h:4; %time range

y=zeros(size(t)); % set y array all 0's to same size as t

y(1)=2; % set initial condition at time 0 to 2

n=numel(y);

dydt = 4*exp(0.8*t) - 0.5*y;

So dydt is a VECTOR. t and y are vectors. However, y is a vector with first element 2, and EVERY other element 0.

for i=1 : n-1 %for loop to interate through y values for

y(i+1)= y(i)+ h * dydt(i); % the Euler method

end

How do you use dydt? Remember, it was fixed before the loop. Those vector elements don't change. You never update dydt when y changes. Therefore you are not solving the problem you think you are solving. You need to update dydt as a function, that changes as both y and t vary. If you don't, well then you get arbitrary garbage as a result - what you got, in fact. That is, you CANNOT create dydt in advance as avector.

So first, let me check your claim as to the solution.

syms y(t)

true_sol = dsolve(diff(y,t,1) == 4*exp(0.8*t) - 0.5*y(t),y(0) == 2)

Interestingly, that is also not the same as what you claim to be the solution. So exact_sol is also incorrect.

simplify(true_sol - (4/1.3)*(exp(0.8*t)-exp(-0.5*t))+2*exp(-0.5*t))

I'm not sure what you screwed up in the solve, but your "exact_sol" is not. Not exact, that is. And since I don't see your code, I can only guess what you did wrong on paper. I won't try to do that.

Now let me implement Euler's method. I'll still use your code as a template, but mine will be correct. :)

h = 0.05; % time step

t = 0:h:4; %time range

y = zeros(size(t)); % set y array all 0's to same size as t

y(1) = 2; % set initial condition at time 0 to 2

n = numel(y);

dydt = @(y_i,t_i) 4*exp(0.8*t_i) - 0.5*y_i;

for i=1 : n-1 % for loop to interate through y values for

y(i+1)= y(i)+ h * dydt(y(i),t(i)); % Euler update

end

plot(t,y,'-b') %plot Euler

hold on

fplot(true_sol,[0,4],'r'); % plots the exact solution to this differential equation

You can see the two curves as I plotted them are now almost on top of each other. A little larger value of h would have been good perhaps to show off the difference more. But Euler actually did pertty well here.

John D'Errico
on 26 Nov 2020

syms t

yoursol = (4/1.3)*(exp(0.8*t)-exp(-0.5*t))+2*exp(-0.5*t);

syms y(t)

mysol = dsolve(diff(y,t,1) == 4*exp(0.8*t) - 0.5*y,y(0) == 2);

expand(mysol)

simplify(mysol - yoursol)

So it looks like the two solutions are the same after all after simplification. I probably copied the wrong solution from you.

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

Start Hunting!
## 0 Comments

Sign in to comment.