matlab：欧拉方法求解微分方程

参考书籍：常用数值算法及其matlab实现，作者：夏省祥
%第10章 常微分方程初值问题的数值解法
%欧拉方法，例10.1

%书籍：常用数值算法及其matlab实现
%第10章 常微分方程初值问题的数值解法
%欧拉方法
function E = eulerdif(fun, a, b, y0, h)
%fun:微分方程的右表达式
%a, b 为区间
%y0 为初值
M = floor(b-a)/h ;      %离散点的个数M+1
T =zeros(1, M+1); Y =zeros(1, M+1);         %行向量
T = a:h:b;
Y(1) = y0;
for i = 1:M
    Y(i+1) = Y(i) +h *feval(fun, T(i), Y(i));
end
E = [T' Y'];            %E是M+1行，2列

% function E = eulerdif(fun, a, b, y0, h)
% %fun:微分方程的右表达式
% %a, b 为区间
% %y0 为初值
% M = floor(b-a)/h ;      %离散点的个数M+1
% T =zeros( M+1,1); Y =zeros( M+1,1);         %列向量
% T = a:h:b;
% Y(1) = y0;
% for i = 1:M
%     Y(i+1) = Y(i) +h *feval(fun, T(i), Y(i));
% end
% E = [T  Y];            %E是M+1行，2列

%书籍：常用数值算法及其matlab实现
%第10章 常微分方程初值问题的数值解法
%欧拉方法，例10.1
%function E = eulerdif(fun, a, b, y0, h)
clear all;clc;close  all;
%fun =@(x,y)-3*x*y;
fun = inline('-3*x*y');
a = 0; b =2;
y0 = 1;
h = 0.1
E = eulerdif(fun, a, b, y0, h);
plot(E(:,1),E(:,2) ,'r*-');    hold on;

exa10_1 = dsolve('Dy = -3*x*y', 'y(0) = 1', 'x');       %求出解析解
ezplot(exa10_1, [0 2]);       %画出解析解的图像
legend('数值解','解析解' )

运行结果如下：