常微分方程牡模型及其数值解.ppt

常微分方程数值解;0、导　　言 ;另一个是方程的求解和结果分析。对一些常系数的或特殊函数形式的微分方程，往往能得到解析解，这对实际问题的分析和应用都是有利的,但是大多数变系数的、非线性函数形式的微分方程都是求不出解析解的,此时就需要应用求解微分方程的另一个重要方法──数值解法。 本章简要介绍有关微分方程模型的概念，微分方程的数值解法和图解法，主要介绍若干建模实例，通过它们展示微分方程模型的建模步骤及解决实际问题的全过程。;欧拉方法和龙格—库塔方法; 所谓数值解法,就是设法将常微分方程离散化,建立差分方程,给出解在一些离散点上的近似值.;对右边的积分应用左矩形公式，则有;对右边的积分应用左矩形公式，则有; 利用Euler方法求初值问题;分别取步长h=0.2 ,0.1 ,0.05,计算结果如下;function [outx,outy]=MyEuler(fun,x0,xt,y0,PointNum) %MyEuler 用前向差分的欧拉方法解微分方程 %fun 表示f（x，y） %x0,xt表示自变量的初值和终值 %y0表示函数在x0处的值,其可以为向量形式 %PointNum表示自变量在[x0,xt]上取的点数 if nargin5 | PointNum=0 %如果函数仅输入4个参数值，则PointNum默认值为100 PointNum=100; end if nargin4 %y0默认值为0 y0=0; end;h=(xt-x0)/PointNum;%计算步长h x=x0+[0:PointNum]*h;%自变量数组 y(1,:) = y0(:);%将输入存为行向量，输入为列向量形式 for k = 1:PointNum f=feval(fun,x(k),y(k,:));%计算f(x,y)在每个迭代点的值 f=f(:); y(k + 1,:) =y(k,:) +h*f; %对于所取的点x迭代计算y值 end outy=y; outx=x; %plot(x,y)%画出方程解的函数图 ;function f=myfun011(x,y) f=1/(1+x*x)-2*y*y;;h;Evaluation only. Created with Aspose.Slides for .NET 3.5 Client Profile 5.2.0.0. Copyright 2004-2011 Aspose Pty Ltd.;Euler中点公式则不然, 计算yn+1时需用到前两步的值yn , yn-1 ,称其为两步方法,两步以上的方法统称为多步法.; 从数值积分的角度来看,梯形公式; 由迭代法收敛的角度看，当; y?=y-2x/y , 0?x?1;计算结果如下:;function [Xout,Yout]=MyEulerPro(fun,x0,xt,y0,PointNumber) %MyEulerPro 用改进的欧拉法解微分方程 if nargin5 | PointNumber=0 %如果函数仅输入4个参数值，则PointNumer默认值为100 PointNumer=100; end if nargin4 %y0默认值为0 y0=0; end;h=(xt-x0)/PointNumber;%计算所取的两离散点之间的距离 x=x0+[0:PointNumber]*h;%表示出离散的自变量x y(1,:)=y0(:); for i=1:PointNumber %迭代计算过程 f1=h*feval(fun,x(i),y(i,:)); f1=f1(:); f2=h*feval(fun,x(i+1),y(i,:)+f1); f2=f2(:); y(i+1,:)=y(i,:)+1/2*(f1+f2); end Xout=x; Yout=y;;function f=myfun012(x,y) f=y-2*x/y;;n;Evaluation only. Created with Aspose.Slides for .NET 3.5 Client Profile 5.2.0.0. Copyright 2004-2011 Aspose Pty Ltd.; 在节点xn+1的误差y(xn+1)-yn+1 ,不仅与yn+1这一步计算有关,而且与前n步计算值yn,yn-1,…,y1都有关. ;另外,在yn=y(xn)的条件下,考虑到y?(x)=?(x,y(x)),则有; yn+1=yn+h?(xn,yn);所以, 改进的Euler方法是二阶方法.;称之为p阶Taylor展开方法.;所以,此差分公式是p阶方法.