数值分析第八章常微分方程数值解法解析

科学计算—理论、方法及其基于MATLAB 的程序实现与分析微分方程（组）数值解法§1 常微分方程初值问题的数值解法微分方程（组）是科学研究和工程应用中最常用的数学模型之一。

如揭示质点运动规律的Newton 第二定律：()()()⎪⎩⎪⎨⎧'='==000022x t x x t x t F dt xd m （１） 和刻画回路电流或电压变化规律的基尔霍夫回路定律等，但是，只有一些简单的和特殊的常微分方程及常微分方程组，可以求得用公式给出的所谓“解析解”或“公式解”，如一阶线性微分方程的初值问题：()()00y y t f ay dtdy=+= （２） 的解为：()()()τττd f e y e t y tt a at ⎰-+=00 （３）但是，绝大多数在实际中遇到的常微分方程和常微分方程组得不到“解析解”，因此，基于如下的事实：１、绝大多数的常微分方程和常微分方程组得不到（有限形式的）解析解；２、实际应用中往往只需要知道常微分方程（组）的解在（人们所关心的）某些点处的函数值（可以是满足一定精度要求的近似值）；如果只需要常微分方程（组）的解在某些点处的函数值，则没有必要非得通过求得公式解，然后再计算出函数值不可，事实上，我们可以采用下面将介绍的常微分方程（组）的初值问题的数值解法，就可以达到这一目的。

一般的一阶常微分方程（组）的初值问题是指如下的一阶常微分方程（组）的定解问题：()()000,y t y t t t y t F dtdyf=≤≤= （７）其中()()()()⎪⎪⎪⎪⎪⎭⎫ ⎝⎛=t y t y t y t y n 21 （８）()()()()⎪⎪⎪⎪⎪⎭⎫ ⎝⎛=y t f y t f y t f y t F n ,,,,21 （９）常微分方程（组）的初值问题通常是对一动态过程（动态系统、动力系统）演化规律的描述，求解常微分方程（组）的初值问题就是要了解和掌握动态过程演化规律。

许多实际问题的数学模型是微分方程或微分方程的定解问题。

如物体运动、电路振荡、化学反映及生物群体的变化等。

常微分方程可分为线性、非线性、高阶方程与方程组等类；线性方程包含于非线性类中，高阶方程可化为一阶方程组。

若方程组中的所有未知量视作一个向量，则方程组可写成向量形式的单个方程。

因此研究一阶微分方程的初值问题⎪⎩⎪⎨⎧=≤≤=0)(),(y a y bx a y x f dxdy， (9-1) 的数值解法具有典型性。

常微分方程的解能用初等函数、特殊函数或它们的级数与积分表达的很少。

用解析方法只能求出线性常系数等特殊类型的方程的解。

对非线性方程来说，解析方法一般是无能为力的，即使某些解具有解析表达式，这个表达式也可能非常复杂而不便计算。

因此研究微分方程的数值解法是非常必要的。

只有保证问题(9-1)的解存在唯一的前提下，研究其数值解法或者说寻求其数值解才有意义。

由常微分方程的理论知，如果(9-1)中的),(y x f 满足条件（1）),(y x f 在区域} ),({+∞<<∞-≤≤=y b x a y x D ，上连续； （2）),(y x f 在上关于满足Lipschitz 条件，即存在常数，使得y y L y x f y x f -≤-),(),(则初值问题(9-1)在区间],[b a 上存在惟一的连续解)(x y y =。

在下面的讨论中，我们总假定方程满足以上两个条件。

所谓数值解法，就是求问题(9-1)的解)(x y y =在若干点b x x x x a N =<<<<= 210处的近似值),,2,1(N n y n =的方法。

),,2,1(N n y n =称为问题(9-1)的数值解，n n x x h -=+1称为由到1+n x 的步长。

今后如无特别说明，我们总假定步长为常量。

建立数值解法，首先要将微分方程离散化，一般采用以下几种方法： (1) 用差商近似导数在问题(9-1)中，若用向前差商hx y x y n n )()(1-+代替)(n x y '，则得)1,,1,0( ))(,()()(1-=≈-+N n x y x f hx y x y n n n n n)(n x y 用其近似值代替，所得结果作为)(1+n x y 的近似值，记为1+n y ，则有 1(,) (0,1,,1)n n n n y y hf x y n N +=+=-这样，问题(9-1)的近似解可通过求解下述问题100(,) (0,1,,1)()n n n n y y hf x y n N y y x +=+=-⎧⎨=⎩(9-2)得到，按式(9-2)由初值经过步迭代，可逐次算出N y y y ,,21。

一、 常微分方程的解析解常微分方程的解析解也就是常微分方程的精确解，也称为常微分方程的符号解；一般可理解为求微分方程的通解或者特解的解析式或表达式；但只有少数的微分方程存在解析解。

在MA TLAB 中，由函数dsolve()求解常微分方程（组）的解析解，其具体格式如下： X=dsolve(‘方程1’，‘方程2’，…‘方程n ’，‘初始条件’，‘自变量’)函数dsolve 用来解符号常微分方程、方程组，如果没有初始条件，则求出通解，如果有初始条件，则求出特解。

例1：求解常微分方程1dy dx x y =+的MA TLAB 程序为：dsolve('Dy=1/(x+y)','x')，注意，系统缺省的自变量为t ，因此这里要把自变量写明。

结果为：-lambertw(-C1*exp(-x-1))-x-1其中：Y=lambertw(X)表示函数关系Y*exp(Y)=X 。

例2：求解常微分方程2'''0yy y -=的MA TLAB 程序为：Y2=dsolve('y*D2y-Dy^2=0’,’x’) 结果为：Y2 =[ exp((x+C2)/C1)][ C2]我们看到有两个解，其中一个是常数。

例3：求常微分方程组253t tdx x y e dt dy x y e dt ⎧++=⎪⎪⎨⎪--=⎪⎩通解的MA TLAB 程序为：[X,Y]=dsolve('Dx+5*x+y=exp(t),Dy-x-3*y=exp(2*t)','t')例4：求常微分方程组020210cos ,224,0t t t dx dy x t x dt dt dx dy y e y dt dt =-=⎧+-==⎪⎪⎨⎪++==⎪⎩通解的MA TLAB 程序为：[X,Y]=dsolve('Dx+2*x-Dy=10*cos(t),Dx+Dy+2*y=4*exp(-2*t)','x(0)=2','y(0)=0')二、 常微分方程的数值解在生产和科研中所处理的微分方程往往很复杂且大多得不出一般解。

介绍常微分方程数值解法常微分方程（ordinary differential equations，ODE）可用于描述许多日常存在的物理系统。

处理ODE问题常常被称为数值求解法，这指的是找到概括ODE或者其他适用于数学模型的解决方案来模括这些ODE。

这种解决方案可能在一系列不同方案中发挥重要作用，以此来提供更好的解释和预测。

常微分方程与几何图形更为相关，它利用二维或者三维空间中曲线的绘制以及分析。

通过引入一些不同的方法，可以对不同的常微分方程中的量进行描述，使得可以通过数值方法的解析来进行研究。

数值解法可能是时间消耗较多的，但有助于验证几何图形中的某些过程，以此帮助揭示数学模型。

四种常见的常微分方程数值解法四种常见的常微分方程数值解法是：前向差分法、向后差分法、中点法和全分方法。

•前向差分法：前向差分法的基本概念是利用ODE的特定解来表达时间步的影响。

这是一种基本的数值法，可以在ODE中确定任意位置的点作为终点。

在这里，任何这样的点都可以表示为ODE右边的时间步。

•向后差分法：它是反过来基于前向差分法。

它要求对ODE中的时间步进行逆向推导，以获得某一特定点的解。

向后差分法要求推导反向解中点，以便可以从每一步中获取该点的解。

•中点法：这是一种非常基本的数值解法，可以用来求解ODE中的某一步的解，但不具有直观的方法解释。

主要的思想是在每一次时间步中通过求出ODE的中点来寻找解。

•全分方法：这是一种更复杂的数值解法，它要求将ODE中的每一步解细分并解决。

与前面提到的三种解法不同，它首先要求将ODE分解成若干离散区间，然后计算每一段区间中的点。

这种解法可以更准确地进行处理，但时间消耗较多，因此比较少被使用。

优化方案在需要解决常微分方程时，为了得到最佳的结果，有必要考虑一些优化措施。

•首先，应考虑将一个复杂的ODE拆分成一些更易解决的问题。

这样做的结果是，预见到解决此ODR的总耗时将会降低。

•其次，为了加快计算速度，可以考虑使用预解算法。

常微分方程的数值解法常微分方程是研究变量的变化率与其当前状态之间的关系的数学分支。

它在物理、工程、经济等领域有着广泛的应用。

解常微分方程的精确解往往十分困难甚至不可得，因此数值解法在实际问题中起到了重要的作用。

本文将介绍常见的常微分方程的数值解法，并比较其优缺点。

1. 欧拉方法欧拉方法是最简单的数值解法之一。

它基于近似替代的思想，将微分方程中的导数用差商近似表示。

具体步骤如下：（1）确定初始条件，即问题的初值。

（2）选择相应的步长h。

（3）根据微分方程的定义使用近似来计算下一个点的值。

欧拉方法的计算简单，但是由于误差累积，精度较低。

2. 改进欧拉方法为了提高欧拉方法的精度，改进欧拉方法应运而生。

改进欧拉方法通过使用两个点的斜率的平均值来计算下一个点的值。

具体步骤如下：（1）确定初始条件，即问题的初值。

（2）选择相应的步长h。

（3）根据微分方程的定义使用近似来计算下一个点的值。

改进欧拉方法相较于欧拉方法而言，精度更高。

3. 龙格-库塔法龙格-库塔法（Runge-Kutta）是常微分方程数值解法中最常用的方法之一。

它通过迭代逼近精确解，并在每一步中计算出多个斜率的加权平均值。

具体步骤如下：（1）确定初始条件，即问题的初值。

（2）选择相应的步长h。

（3）计算各阶导数的导数值。

（4）根据权重系数计算下一个点的值。

与欧拉方法和改进欧拉方法相比，龙格-库塔法的精度更高，但计算量也更大。

4. 亚当斯法亚当斯法（Adams）是一种多步法，它利用之前的解来近似下一个点的值。

具体步骤如下：（1）确定初始条件，即问题的初值。

（2）选择相应的步长h。

（3）通过隐式或显式的方式计算下一个点的值。

亚当斯法可以提高精度，并且比龙格-库塔法更加高效。

5. 多步法和多级法除了亚当斯法，还有其他的多步法和多级法可以用于解常微分方程。

多步法通过利用多个点的值来逼近解，从而提高精度。

而多级法则将步长进行分割，分别计算每个子问题的解，再进行组合得到整体解。

第8章常微分方程的数值解法8.4单步法的收敛性与稳定性8.4.1相容性与收敛性上面所介绍的方法都是用离散化的方法，将微分方程初值问题化为差分方程初值问题求解的.这些转化是否合理？即当h →∞时，差分方程是否能无限逼近微分方程，差分方程的解n y 是否能无限逼近微分方程初值问题的准确解()n y x ，这就是相容性与收敛性问题.用单步法(8.3.14)求解初值问题(8.1.1)，即用差分方程初值问题100(,,)()n n n n y y h x y h y x y ϕ+=+⎧⎨=⎩(8.4.1)的解作为问题(8.1.1)的近似解，如果近似是合理的，则应有()()(,(),)0 (0)y x h y x x y x h h hϕ+--→→(8.4.2)其中()y x 为问题(8.1.1)的精确解.因为0()()lim ()(,)h y x h y x y x f x y h→+-'==故由(8.4.2)得lim (,,)(,)h x y h f x y ϕ→=如果增量函数(,(),)x y x h ϕ关于h 连续，则有(,,0)(,)x y f x y ϕ=(8.4.3)定义8.3如果单步法的增量函数(,,)x y h ϕ满足条件(8.4.3)，则称单步法(8.3.14)与初值问题(8.1.1)相容.通常称(8.4.3)为单步法的相容条件.满足相容条件(8.4.3)是可以用单步法求解初值问题(8.1.1)的必要条件.容易验证欧拉法和改进欧拉法均满足相容性条件.一般地，如果单步法有p 阶精度（1p ≥），则其局部截断误差为[]1()()(,(),)()p y x h y x h x y x h O h ϕ++-+=上式两端同除以h ，得()()(,,)()p y x h y x x y h O h hϕ+--=令0h →，如果(,(),)x y x h ϕ连续，则有()(,,0)0y x x y ϕ'-=所以1p ≥的单步法均与问题(8.1.1)相容.由此即得各阶龙格－库塔法与初值问题(8.1.1)相容.定义8.4一种数值方法称为是收敛的，如果对于任意初值0y 及任意固定的(,]x a b ∈，都有lim () ()n h y y x x a nh →==+其中()y x 为初值问题(8.1.1)的精确解.如果我们取消局部化假定，使用某单步法公式，从0x 出发，一步一步地推算到1n x +处的近似值1n y +.若不计各步的舍入误差，而每一步都有局部截断误差，这些局部截断误差的积累就是整体截断误差.定义8.5称111()n n n e y x y +++=-为某数值方法的整体截断误差.其中()y x 为初值问题(8.1.1)的精确解，1n y +为不计舍入误差时用某数值方法从0x 开始，逐步得到的在1n x +处的近似值（不考虑舍入误差的情况下，局部截断误差的积累）.定理8.1设单步法(8.3.14)具有p 阶精度，其增量函数(,,)x y h ϕ关于y 满足利普希茨条件，问题(8.1.1)的初值是精确的，即00()y x y =，则单步法的整体截断误差为111()()p n n n e y x y O h +++=-=证明由已知，(,,)x y h ϕ关于y 满足利普希茨条件，故存在0L >，使得对任意的12,y y 及[,]x a b ∈，00h h <≤，都有1212(,,)(,,)x y h x y h L y y ϕϕ-≤-记1()(,(),)n n n n y y x h x y x h ϕ+=+，因为单步法具有p 阶精度，故存在0M >，使得1111()p n n n R y x y Mh ++++=-≤从而有111111111()()()(,(),)(,,)()(,(),)(,,)n n n n n n n p n n n n n n p n n n n n n e y x y y x y y y Mh y x h x y x h y h x y h Mh y x y h x y x h x y h ϕϕϕϕ+++++++++=-≤-+-≤++--≤+-+-1(1)p nMh hL e +≤++反复递推得11111101110(1)(1)1(1)(1)(1)(1)1(1)p p n n n p n n p n e Mh hL Mh hL e hL hL Mh hL e hL Mh hL e hL+++-+++++⎡⎤≤++++⎣⎦⎡⎤≤+++++++⎣⎦+-≤++因为00()y x y =，即00e =，又(1)n h b a +≤-，于是ln(1)1()(1)(1)b a b a hL n L b a h h hL hL e e --++-+≤+=≤所以()11()p L b a p n M e h e O h L -+⎡⎤≤-=⎣⎦推论设单步法具有p (1p ≥)阶精度，增量函数(,,)x y h ϕ在区域G ：, , 0a x b y h h ≤≤-∞<<+∞≤≤上连续，且关于y 满足利普希茨条件，则单步法是收敛的.当(,)f x y 在区域:,D a x b y ≤≤-∞<<+∞上连续，且关于y 满足利普希茨条件时，改进欧拉法，各阶龙格－库塔法的增量函数(,,)x y h ϕ在区域G 上连续，且关于y 满足利普希茨条件，因而它们都是收敛的.关于单步法收敛的一般结果是：定理8.2设增量函数(,,)x y h ϕ在区域G 上连续，且关于y 满足利普希茨条件，则单步法收敛的充分必要条件是相容性条件(8.4.3).8.4.2稳定性稳定性与收敛性是两个不同的概念，收敛性是在假定每一步计算都准确的前提下，讨论当步长0h →时，方法的整体截断误差是否趋于零的问题.而稳定性则是讨论舍入误差的积累能否对计算结果有严重影响的问题.定义8.6若一种数值方法在节点值n y 上有一个大小为δ的扰动，于以后各节点()m y m n >上产生的偏差均不超过δ，则称该方法是稳定的.我们以欧拉法为例进行讨论.假设由于舍入误差，实际得到的不是n y 而是n n n y y δ=+，其中n δ是误差.由此再计算一步，得到1(,)n n n n y y hf x y +=+把它与不考虑舍入误差的欧拉公式相减，并记111n n n y y δ+++=-，就有[]1(,)(,)1(,)n n n n n n y n nh f x y f x y hf x δδηδ+⎡⎤=+-=+⎣⎦其中y f f y∂=∂.如果满足条件1(,)1y n hf x η+≤，(8.4.4)则从n y 到1n y +的计算，误差是不增的，可以认为计算是稳定的.如果条件(8.4.4)不满足，则每步误差将增大.当0y f >时，显然条件(8.4.4)不可能满足，我们认为问题本身具有先天的不稳定性.当0y f <时，为了满足稳定性要求(8.4.4)，有时h 要很小.一般的，稳定性与方法有关，也与步长h 的大小有关，当然也与方程中的(,)f x y 有关.为简单起见，通常只考虑数值方法用于求解模型方程的稳定性，模型方程为y y λ'=(8.4.5)其中λ为复数.一般的方程可以通过局部线性化转化为模型方程，例如在(,)x y 的邻域(,)(,)(,)()(,)()x y y f x y f x y f x y x x f x y y y '==+-+-+略去高阶项，再作变量替换就得到u u λ'=的形式.对于模型方程(8.4.5)，若Re 0λ>，类似以上分析，可以认为方程是不稳定的.所以我们只考虑Re 0λ<的情形，这时不同的数值方法可能是数值稳定的或者是数值不稳定的.当一个单步法用于试验方程y y λ'=，从n y 计算一步得到1()n n y E h y λ+=(8.4.6)其中()E h λ依赖于所选的方法.因为通过点(,)n n x y 试验方程的解曲线（它满足,()n n y y y x y λ'==）为[]exp ()n n y y x x λ=-，而一个p 阶单步法的局部截断误差在()n n y x y =时有1111()()p n n n T y x y O h ++++=-=，所以有1exp()()()p n n y h E h y O h λλ+-=(8.4.7)这样可以看出()E h λ是h e λ的一个近似值.由(8.4.6)可以看到，若n y 计算中有误差ε，则计算1n y +时将产生误差()E h λε，所以有下面定义.定义8.7如果(8.4.6)式中，()1E h λ<，则称单步法(8.3.14)是绝对稳定的.在复平面上复变量h λ满足()1E h λ<的区域，称为方法(8.3.14)的绝对稳定区域，它与实轴的交称为绝对稳定区间.在上述定义中，规定严格不等式成立，是为了和线性多步法的绝对稳定性定义一致.事实上，()1E h λ=时也可以认为误差不增长.(1)欧拉法的稳定性欧拉法用于模型方程(8.4.5)，得1(1)n n y h y λ+=+，所以有()1E h h λλ=+.所以绝对稳定条件是11h λ+<，它的绝对稳定区域是h λ复平面上以(1,0)-为中心的单位圆，见图8.3.而λ为实数时，绝对稳定区间是(2,0)-.Im()h λRe()h λ2-1-O 图8.3欧拉法的绝对稳定区域(2)梯形公式的稳定性对模型方程，梯形公式的具体表达式为11()2n n n n h y y y y λλ++=++，即11212n nh y y h λλ++=-，所以梯形公式的绝对稳定区域为12112h h λλ+<-.化简得Re()0h λ<，因此梯形公式的绝对稳定区域为h λ平面的左半平面，见图8.4.特别地，当λ为负实数时，对任意的0h >，梯形公式都是稳定的.Im()h λRe()h λO 图8.4梯形公式的绝对稳定区域(3)龙格－库塔法的稳定性与前面的讨论相仿，将龙格－库塔法用于模型方程(8.4.5)，可得二、三、四阶龙格－库塔法的绝对稳定区域分别为211()12h h λλ++<23111()()126h h h λλλ+++<2341111()()()12624h h h h λλλλ++++<当λ为实数时，二、三、四阶显式龙格－库塔法的绝对稳定区域分别为20h λ-<<、2.510h λ-<<、 2.780h λ-<<.例8.5设有初值问题21010101(0)0xy y x x y ⎧'=-≤≤⎪+⎨⎪=⎩用四阶经典龙格－库塔公式求解时，从绝对稳定性考虑，对步长h 有何限制？解对于所给的微分方程有2100,(010)1f x x y xλ∂==-<≤≤∂+在区间[0,10]上，有201010max ||max51t x x λ<<==+由于四阶经典龙格－库塔公式的绝对稳定区间为 2.7850h λ-<<，则步长h 应满足00.557h <<.。

数值分析中的常微分方程数值求解常微分方程是自然科学中一类最为普遍的数学模型，涉及到热力学、物理、化工等多个领域。

然而，解常微分方程并非易事。

尤其是当我们面对一些复杂、非线性、多维的方程组时，常微分方程数值求解成为了一个十分关键的问题。

因此，数值求解方法成为了常微分方程研究中的重要组成部分。

本文将介绍一些数值解常微分方程的常见方法和应用。

1. 一般线性方法一般线性方法（general linear methods）是经典的常微分方程数值解法之一。

它以一种特殊的形式给出步进公式：$$ y_{n+1}=\sum_{i=0}^{s-1}\alpha_i y_{n-i}+h\sum_{i=0}^{s-1}\beta_i f(t_{n-i},y_{n-i}) $$ 其中，$y_{n}$为第$n$步的项值，$f(t_n,y_n)$为时间$t_n$处函数$y(t)$的导数。

$\alpha_i$和$\beta_i$是常数，可以通过确定如下特征方程来选择：$$ \sum_{i=0}^{s-1}\alpha_i\ lambda^{i}=0,~(\lambda\in C) $$ 与此同时，也可以通过选择$\beta_i$来使方法达到一定的准确性和稳定性。

2. Runge-Kutta方法比一般线性方法更为流行的方法是Runge-Kutta方法。

通常附加一个或多个修正以获得更好的数值稳定性和误差控制。

第1阶Runge-Kutta方法仅使用导数$f(t_n,y_n)$估算下一个项的值：$$y_{n+1}=y_n+hf(t_n,y_n)$$ 许多高阶方法可以使用中间的“插值”来更准确地估计下一个步骤：$$y_{n+1}=y_n+h\sum_{i=1}^kb_ik_i$$$$k_i=f(t_n+c_ih,y_n+h\sum _{j=1}^{i-1}a_{ij}k_j)$$ $k_i$是第$i$台车的估计值，$a_{ij}$和$b_i$在经典Runge-Kutta方法和其他变体中具有不同的取值。

第八章 常微分方程数值解法教学目的 1. 掌握解常微分方程的单步法：Euler 方法、Taylor 方法和Runge-Kutta 方法；2. 掌握解常微分方程的多步法：Adams 步法、Simpson 方法和Milne 方法等；3. 了解单步法的收敛性、相容性与稳定性；多步法的稳定性。

教学重点及难点 重点是解常微分方程的单步法：Euler 方法、Taylor 方法和Runge-Kutta 方法和解常微分方程的多步法：Adams 步法、Simpson 方法和Milne 方法等；难点是理解单步法的收敛性、相容性与稳定性及多步法的稳定性。

教学时数 20学时教学过程§1基本概念1．1常微分方程初值问题的一般提法常微分方程初值问题的一般提法是求函数b x a x y ≤≤),(，满足⎪⎩⎪⎨⎧=<<=)2.1()()1.1(),,(αa yb x a y x f dx dy其中),(y x f 是已知函数，α是已知值。

假设),(y x f 在区域},),{(+∞<≤≤=y b x a y x D 上满足条件： （1）),(y x f 在D 上连续； （2）),(y x f 在D 上关于变量y 满足Lipschitz 条件：2121),(),(y y L y x f y x f -≤-，21,,y y b x a ∀≤≤ （1.3）其中常数L 称为Lipschitz 常数。

我们简称条件（1）、（2）的基本条件。

由常微分方程的基本理论，我们有：定理 1 当),(y x f 在D 上满足基本条件时，一阶常微分方程初值问题（1.1）、（1.2）对任意给定α存在唯一解)(x y 在],[b a 上连续可微。

定义1 方程（1.1）、（1.2）的解)(x y 称为适定的，若存在常数0>ε和0>K ，对任意满足条件εδ≤及εη≤∞)(x 的δ和)(x η，常微分方程初值问题⎪⎩⎪⎨⎧+=<<+=δηa a z b x a x z x f dx dz)(),(),(（1.4）存在唯一解)(x z ，且}.{)()(δη+≤-∞∞K x z x y适定问题的解)(x y 连续依赖于（1.1）右端的),(y x f 和初值α。

数值解常微分方程的方法和技巧在科学和工程领域，我们经常遇到一些复杂的常微分方程（Ordinary Differential Equations, ODEs），这些方程往往很难用解析方法得到精确解。

而数值解常微分方程的方法和技巧提供了一种有效的途径来近似求解这些方程。

本文将介绍一些常用的数值解ODEs的方法和技巧。

一、欧拉方法（Euler Method）欧拉方法是最简单的数值解ODEs的方法，它利用初始条件和微分方程的导数来计算下一个点的近似值。

具体来说，假设我们要求解的ODE为dy/dx = f(x, y)，其中f(x, y)是已知函数，初始条件为x0 = x(0)，y0 = y(0)。

欧拉方法的迭代公式为：y[i+1] = y[i] + h * f(x[i], y[i])其中，h是步长，x[i]表示第i个点的x坐标，y[i]表示对应的y坐标。

二、龙格-库塔方法（Runge-Kutta Method）龙格-库塔方法是一族常用的数值解ODEs方法，其基本思想是通过计算不同阶数的导数来提高求解的精度。

最常用的龙格-库塔方法是四阶龙格-库塔方法，也称为RK4方法。

它的迭代公式如下：k1 = h * f(x[i], y[i])k2 = h * f(x[i] + h/2, y[i] + k1/2)k3 = h * f(x[i] + h/2, y[i] + k2/2)k4 = h * f(x[i] + h, y[i] + k3)y[i+1] = y[i] + 1/6 * (k1 + 2*k2 + 2*k3 + k4)其中，k1、k2、k3、k4是中间变量，h是步长。

三、改进的欧拉方法（Improved Euler Method）改进的欧拉方法是对欧拉方法的改进，它通过使用导数的平均值来提高求解的精度。

其迭代公式为：k1 = h * f(x[i], y[i])k2 = h * f(x[i] + h, y[i] + k1)y[i+1] = y[i] + 1/2 * (k1 + k2)其中，k1、k2是中间变量，h是步长。

第八章 常微分方程的数值解法一．内容要点考虑一阶常微分方程初值问题:⎪⎩⎪⎨⎧==00)(),(y x y y x f dx dy微分方程的数值解：设微分方程的解y (x )的存在区间是[a,b ]，在[a,b ]内取一系列节点a= x 0< x 1<…< x n =b ,其中h k =x k+1-x k ；(一般采用等距节点，h=(b-a)/n 称为步长)。

在每个节点x k 求解函数y(x)的近似值:y k ≈y(x k )，这样y 0 , y 1 ,...,y n 称为微分方程的数值解。

用数值方法，求得f(x k )的近似值y k ，再用插值或拟合方法就求得y(x)的近似函数。

(一)常微分方程处置问题解得存在唯一性定理对于常微分方程初值问题:⎪⎩⎪⎨⎧==00)(),(y x y y x f dx dy如果：(1) 在B y y A x x 00≤-≤≤,的矩形内),(y x f 是一个二元连续函数。

(2) ),(y x f 对于y 满足利普希茨条件，即2121y y L y x f y x f -≤-),(),(则在C x x 0≤≤上方程⎪⎩⎪⎨⎧==00)(),(y x y y x f dxdy的解存在且唯一，这里C=min((A-x 0),x 0+B/L),L 是利普希茨常数。

定义：任何一个一步方法可以写为),,(h y x h y y k k k 1k Φ+=+，其中),,(h y x k k Φ称为算法的增量函数。

收敛性定理：若一步方法满足： (1)是p 解的.(2) 增量函数),,(h y x k k Φ对于y 满足利普希茨条件.(3) 初始值y 0是精确的。

则),()()(p h O x y kh y =-kh =x -x 0,也就是有0x y y lim k x x kh 0h 0=--=→)((一)、主要算法 1．局部截断误差局部截断误差：当y(x k )是精确解时，由y(x k )按照数值方法计算出来的1~+k y 的误差y (x k+1)- 1~+k y 称为局部截断误差。

常微分方程数值解法常微分方程是研究函数的导数与自变量之间的关系的数学分支，广泛应用于物理、工程、生物等领域的建模与分析。

在实际问题中，我们常常遇到无法通过解析方法求得精确解的常微分方程，因此需要利用数值解法进行求解。

本文将介绍几种常用的常微分方程数值解法。

一、欧拉方法(Euler's Method)欧拉方法是最基本的数值解法之一。

它的思想是将微分方程转化为差分方程，通过逐步逼近解的方式求得数值解。

具体步骤如下：1. 将微分方程转化为差分方程：根据微分方程的定义，可以得到差分方程形式。

2. 选择步长：将自变量范围进行离散化，确定步长h。

3. 迭代计算：根据差分方程递推公式，利用前一步的数值解计算后一步的数值解。

二、改进的欧拉方法(Improved Euler's Method)改进的欧拉方法通过使用欧拉方法中的斜率来进行更准确的数值计算。

具体步骤如下：1. 计算欧拉方法的斜率：根据当前节点的数值解计算斜率。

2. 根据斜率计算改进的数值解：将得到的斜率代入欧拉方法的递推公式中，计算改进的数值解。

三、龙格-库塔方法(Runge-Kutta Method)龙格-库塔方法是一类常微分方程数值解法，其中最著名的是四阶龙格-库塔方法。

它通过计算各阶导数的加权平均值来逼近解，在精度和稳定性方面相对较高。

具体步骤如下：1. 计算每一步的斜率：根据当前节点的数值解计算每一步的斜率。

2. 计算权重：根据斜率计算各个权重。

3. 计算下一步的数值解：根据计算得到的权重，将其代入龙格-库塔方法的递推公式中，计算下一步的数值解。

四、多步法(多步差分法)多步法是需要利用多个前面节点的数值解来计算当前节点的数值解的数值方法。

常见的多步法有Adams-Bashforth法和Adams-Moulton法。

具体步骤如下：1. 选择初始值：根据差分方程的初始条件，确定初始值。

2. 迭代计算：根据递推公式，利用前面节点的数值解计算当前节点的数值解。