欧拉法、梯形法和龙格-库塔解微分方程.docx

欧拉法、梯形法和龙格-库塔一、解方程： y = 8x(2-y)y(0)=1二、算出方程的解析解为：y= 2 - e-4x2三、实验原理：1.欧拉法原理：将区间a,b分成n段，那么方程在第xi点有y(xi) = f(xi,y(xi)，再用向前差商近似代替导数则为：，在这里，h是步长，即相邻两个结点间的距离。因此可以根据xi点和yi点的数值计算出yi+1来：，i=0,1,2,n这就是欧拉格式，若初值yi + 1是已知的，则可依据上式逐步算出数值解y1,y2.yn2.梯形法原理：将向前欧拉公式中的导数f(xi,yi)改为微元两端f(xi,yi)和f(xi+1,yi+1)的平均，即梯形公式。3. 龙格-库塔方法的基本思想：在区间xn,xn+1内多取几个点，将他们的斜率加权平均，作为导数的近似。令初值问题表述如下。则，对于该问题的RK4由如下方程给出：其中这样，下一个值(yn+1)由现在的值(yn)加上时间间隔(h)和一个估算的斜率的乘积决定。该斜率是以下斜率的加权平均： k1是时间段开始时的斜率； k2是时间段中点的斜率，通过欧拉法采用斜率k1来决定y在点tn + h/2的值； k3也是中点的斜率，但是这次采用斜率k2决定y值； k4是时间段终点的斜率，其y值用k3决定。当四个斜率取平均时，中点的斜率有更大的权值：RK4法是四阶方法，也就是说每步的误差是h5阶，而总积累误差为h4阶。四、欧拉法、梯形法和龙格-库塔的实现代码：h=0.1;x=0:h:1;y1=zeros(size(x);y1(1)=1;y2=zeros(size(x);y2(1)=1;y3=zeros(size(x);y3(1)=1; for i1=2:length(x) y1(i1) = y1(i1-1) + h*8*x(i1-1)*(2-y1(i1-1);%欧拉法 m1= 8*x(i1-1)*(2-y2(i1-1);%梯形法 m2 =8*x(i1)*(2-y2(i1-1)+h*m1); y2(i1)=y2(i1-1)+h*(m1+m2)/2;%梯形法公式 k1=8*x(i1-1)*(2-y2(i1-1); %龙格-库塔 k2=8*(x(i1-1)+h/2)*(2-(y2(i1-1)+h*k1/2); k3=8*(x(i1-1)+h/2)*(2-(y2(i1-1)+h*k2/2); k4=8*(x(i1-1)+h)*(2-(y2(i1-1)+h*k3);y3(i1)=y3(i1-1)+(k1+2*k2+2*k3+k4)*h/6; %龙格-库塔公式end y4=2-exp(-4*(x.2);%解析解plot(x,y1,x,y2,x,y3,x,y4)%解析解与数值解图像legend(y1,y2,y3,y4)plot(x,y4-y1,x,y4-y2,x,y4-y3)% 解析解与数值解误差图像legend(y4-y1,y4-y2,y4-y3)五、图像：1.解析解y4与各个数值解y1,y2,y3的图像:(1)当步长h=0.1时：（2）当步长h=0.05时：（3）当步长h=0.01时： 结论：三个图中，方程的解析解都是y4。Y1表示欧拉法, y2表示梯形法,y3表示龙格-库塔法。从上面三个图像可以看出，y3与y4的拟合度比其他两个都好，即龙格-库塔法得到的数值解与真实解的拟合度比欧拉法和梯形法都高。同时，从上图可以看出，步长h=0.01比h=0.05和h=0.1所得到的的解更接近于真实解，这也说明了，步长越小，区间内的点越多，估计解也就越接近与真实解。2.解析解y4与各个数值解y1,y2,y3的误差图像：（1） 当步长h=0.1时 （2）当步长h=0.05(3)当步长h=0.01时结论：同样的，上面三个误差图可以看出：y4-y3表示龙格-库塔法与真实解的误差，y4-y2表示梯形法与真实解的误差，y4-y1表示欧拉法与真实解的误差。从三条曲线