eula3.m

这是用来解常微分方程的原代码

%欧拉梯形公式
function E=Eula3(f,x0,y0,xN,N)
%f一阶微分方程的函数
%x0，y0为初始条件
%xN为取值范围的一个端点
%N为区间个数
x=zeros(1,N+1);y=zeros(1,N+1);
x(1)=x0;y(1)=y0;
h=(xN-x0)/N;
for n=1:N
    x(n+1)=x(n)+h;
    z0=y(n)+h*feval(f,x(n),y(n));
    for k=1:3
        z1=y(n)+h/2*(feval(f,x(n),y(n))+feval(f,x(n+1),z0));
        if abs(z1-z0)<1e-3
            break;
        end
        z0=z1;
    end
    y(n+1)=z1;
end
E=[x',y'];