Any ideas why Euler method isn't running?

2 views (last 30 days)
Kaylene Widdoes
Kaylene Widdoes on 17 Feb 2016
Answered: Geoff Hayes on 18 Feb 2016
I have the following equation: (x^3)y' + 20*(x^2)y = x with y(2)=0
I have to get the Euler code to run and execute. So far this is my code, but the graph looks wrong, and the approximations aren't even showing up.
Script file:
function[T,Y] = euler(f,a,b,ya,h)
a = 0; b = 1;
T = a:h:b;
Y = zeros(1,length(T));
Y(2) = ya;
for k = 1 : length(T)-1
Y(k+1) = Y(k) + h*f(T(k),Y(k));
end
Code:
%%Part 1 - Euler
% y'(t) = 3*y(t)
% A function file euler.m is needed; you'll do this in the homework
clear all
% define the problem: function f and domain
f = @(t,y) (1/(t^2)) - ((20*y)/t);
a = 2; b = 10;
% exact solution, using a fine grid
t = a:.0001:b;
y = (1./(19.*t)) - (524288./(19.*(t.^20))); % this is a vector of values, not a function
% coarse solution
h = .01;
ya = 0;
[T1,Y1]=euler(f,a,b,ya,h);
% fine solution
h = .001;
ya = 0;
[T2,Y2]=euler(f,a,b,ya,h);
% finer solution
h = .0001;
ya = 0;
[T3,Y3]=euler(f,a,b,ya,h);
plot(t,y,'k',T1,Y1,'bo-',T2,Y2,'ro-',T3,Y3,'go-')
legend('Exact','h=0.01','h=0.001','h=0.001')
title('The Euler Method with 3 meshes')

Answers (1)

Geoff Hayes
Geoff Hayes on 18 Feb 2016
Kaylene - if I run your code, I notice that your Y array output from euler are all NaNs. Looking closely at this function
function[T,Y] = euler(f,a,b,ya,h)
a = 0; b = 1;
T = a:h:b;
Y = zeros(1,length(T));
Y(2) = ya;
for k = 1 : length(T)-1
Y(k+1) = Y(k) + h*f(T(k),Y(k));
end
the a and b inputs are overwritten with 0 and 1 respectively. Why? If I remove these initializations and change the initial value for Y (so Y(1) and not Y(2)) then the code becomes
function[T,Y] = euler(f,a,b,ya,h)
T = a:h:b;
Y = zeros(1,length(T));
Y(1) = ya;
for k = 1 : length(T)-1
Y(k+1) = Y(k) + h*f(T(k),Y(k));
end
which (at least) initializes Y to be something that can be drawn to your figure. The output using the different step sizes appears to be near identical which may not be what you are trying to show.

Categories

Find more on Linear Algebra 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!