fe.m

来自「forward euler solving differential equat」· M 代码 · 共 23 行

clear;clc;        % delete variables from memory
a=0;% endpoints of the interval
N=10;
h=0.1% number of small intervals (will give us N+1 points in the interval)
b=h*N;
t=0;% the spacing between points 
alpha=6*exp(t-1)-t^2-2*t-2      % initial condition
f=inline('x+t^2'); % the function f(t, y)

T=a:h:b;      % the points t0, t1, ... from a to b with spacing h
W(1)=alpha;   % initial condition, that is, w0=alpha

for i=1:N
   W(i+1)=W(i)+h*f(T(i), W(i)) % Iteration of Euler's method
end

Y=6*exp(T-1)-T.^2-2*T-2 % exact solution

figure(1); clf; hold on; % pop up a figure, clean it up

plot(T, Y, 'b');    % exact solution in blue
plot(T, W, 'r');    % the calculated solution in red
title('Exact solution:blue, approx solution: red');