书目:《数值方法(matlab版)》,作者周璐等译
%欧拉方法求解微分方程function [t,y] = my_euler(f, t0, tf , y0, h)%f-函数; t0,tf:区间; y0,初值;h 步长M = floor((tf - t0)/h); % 离散点的个数t = zeros(M +1, 1);y = zeros(M +1, 1);t = [t0 : h :tf]';y(1) = y0;for i =1:My(i+1) = y(i) + h*feval(f, t(i), y(i));end
function s = myfun(t,y)s = 1/2*(t - y);end
clear all;clc%调用欧拉方法求解微分方程-%函数原型 function [t,y] = my_euler(f, t0, tf , y0, h)t0 = 0;tf = 3;y0 = 1;h1 = 1; h2 = 0.5; h3 = 0.25; h4 = 1/8;[t1,y1] = my_euler(@myfun, t0, tf , y0, h1)[t2,y2] = my_euler(@myfun, t0, tf , y0, h2)[t3,y3] = my_euler(@myfun, t0, tf , y0, h3)[t4,y4] = my_euler(@myfun, t0, tf , y0, h4)plot(t1,y1,'b>',t2,y2,'y*' ,t3,y3,'ro',t4,y4,'g--');hold on;xlabel('t');ylabel('y');grid on;legend('h=1', 'h=1/2','h=1/4','h=1/8');plot(t4, 3*exp(-t4/2) - 2 +t4);%第二种写法% myfun1 = @(t,y)1/2*(t - y)% [t,y1] = my_euler(myfun1, t0, tf , y0, h1)% [t,y2] = my_euler(myfun1, t0, tf , y0, h2)% [t,y3] = my_euler(myfun1, t0, tf , y0, h3)% [t,y4] = my_euler(myfun1, t0, tf , y0, h4)
结果如下:
结论:可以看出,随着步长的减小,通过欧拉方法得到的微分方程的解,更加接近精确解。