Trying to make euler's method program work

2 views (last 30 days)
ce
ce on 5 Dec 2022
Edited: Torsten on 6 Dec 2022
function [x,y] = euler_for(x0,y0,x_end,h,fcn)
%Solve the ivp
% y' = f(x,y), x0 <= x <= b, y(x0) = y0
%Use euler's method with a stepsize h. The
%user must supply a program with some name, say, deriv,
%and a first line of the form
% function ans = deriv(x,y)
%A sample call would be
% [t,z] = euler_for(t0,z0,b,delta,'deriv')
%
%output:
%the routine eulercls will return two vectors, x and y.
%The vector x will contain the node points
% x(1) = x0, x(j) = x0 + (j-1)*h, j=1,2,...,N
%with
% x(N) <= x_end - h, x(N) + h > x_end - h
%The vector y will contain the estimates of the solution Y
%at the node points in x.
%
n = fix((x_end-x0)/h)+1;
x = linspace(x0,x_end,n);
y = zeros(n,1);
y(1) = y0;
for i = 2:n
y(i) = y(i-1)+h*feval(fcn,x(i-1),y(i-1));
end
I'm not sure what is is meant by the function that I need to supply. I know what function I'm supposed to use, but I don't know how to add it to the program. I tried making a function called fcn but when I run the program it says there is an error in fcn, although it doesn't say what it is. This is my fcn code:
function f = fcn(x,y)
f = (cos(atan(x)))^2;
end

Answers (1)

Alan Stevens
Alan Stevens on 5 Dec 2022
Like this:
% Arbitrary values here - insert your own
x0 = 0;
x_end = 10;
h = 0.1;
y0 = 1;
n = fix((x_end-x0)/h)+1;
x = linspace(x0,x_end,n);
y = zeros(n,1);
y(1) = y0;
for i = 2:n
y(i) = y(i-1)+h*fcn(x(i-1));
end
plot(x,y),grid
xlabel('x'),ylabel('y')
function f = fcn(x) % Your function doesn't depend on y
f = (cos(atan(x)))^2;
end
  2 Comments
ce
ce on 5 Dec 2022
I think I might have misunderstood the way this program works. The goal of it is to solve differential equations. For example, the function I've been trying to solve in this code is Y'(x) = cos(Y(x))^2, but I'm not sure how to input this into the program. I know the the solution is Y(x) = arctan(x). The code I posted is from my textbook but it doesn't expalain how to use it and it's very confusing to me. I tried the code you gave me, but now I don't know how to enter the function since its a DE.
Torsten
Torsten on 6 Dec 2022
Edited: Torsten on 6 Dec 2022
x0 = 0;
x_end = 10; % Integration interval [x0,x_end]
h = 0.1; % Step size of integration
fcn = @(x,y) (cos(y))^2; % Y' = cos(Y)^2
y0 = 0; % Y(x0) = 0
[x,y] = euler_for(x0,y0,x_end,h,fcn);
hold on
plot(x,y),grid
xlabel('x'),ylabel('y')
syms y(x)
Ysol(x) = dsolve(diff(y,x)==cos(y)^2,y(0)==0);
fplot(Ysol,[x0 x_end])
hold off
function [x,y] = euler_for(x0,y0,x_end,h,fcn)
n = fix((x_end-x0)/h)+1;
x = linspace(x0,x_end,n);
y = zeros(n,1);
y(1) = y0;
for i = 2:n
y(i) = y(i-1)+h*fcn(x(i-1),y(i-1));
end
end

Sign in to comment.

Community Treasure Hunt

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

Start Hunting!