微分方程数值解
微分方程的数值解法
微分方程的数值解法微分方程是自然科学和现代技术领域中一种最基本的数学描述工具，它可以描述物理世界中的各种现象。
微分方程的解析解往往很难求出，因此数值解法成为解决微分方程问题的主要手段之一。
本文将介绍几种常见的微分方程的数值解法。
一、欧拉法欧拉法是微分方程初值问题的最简单的数值方法之一，它是由欧拉提出的。
考虑一阶常微分方程：$y'=f(t,y),y(t_0)=y_0$其中，$f(t,y)$表示$y$对$t$的导数，则$y(t_{i+1})=y(t_i)+hf(t_i,y_i)$其中，$h$为步长，$t_i=t_0+ih$，$y_i$是$y(t_i)$的近似值。
欧拉法的精度较低，误差随着步长的增加而增大，因此不适用于求解精度要求较高的问题。
二、改进欧拉法改进欧拉法又称为Heun方法，它是由Heun提出的。
改进欧拉法是在欧拉法的基础上进行的改进，它在每个步长内提高求解精度。
改进欧拉法的步骤如下：1. 根据当前$t_i$和$y_i$估算$y_{i+1}$:$y^*=y_i+hf(t_i,y_i),t^*=t_i+h$2. 利用$y^*$和$t^*$估算$f(t^*,y^*)$:$f^*=f(t^*,y^*)$3. 利用$y_i$、$f(t_i,y_i)$和$f^*$估算$y_{i+1}$:$y_{i+1}=y_i+\frac{h}{2}(f(t_i,y_i)+f^*)$改进欧拉法具有比欧拉法更高的精度，但是相较于其他更高精度的数值方法，它的精度仍然较低。
三、龙格-库塔法龙格-库塔法是一种广泛使用的高精度数值方法，它不仅能够求解一阶和二阶常微分方程，还能够求解高阶常微分方程和偏微分方程。
其中，经典的四阶龙格-库塔法是最常用的数值方法之一。
四阶龙格-库塔法的步骤如下：1. 根据当前$t_i$和$y_i$估算$k_1$:$k_1=f(t_i,y_i)$2. 根据$k_1$和$y_i$估算$k_2$:$k_2=f(t_i+\frac{h}{2},y_i+\frac{h}{2}k_1)$3. 根据$k_2$和$y_i$估算$k_3$:$k_3=f(t_i+\frac{h}{2},y_i+\frac{h}{2}k_2)$4. 根据$k_3$和$y_i$估算$k_4$:$k_4=f(t_i+h,y_i+hk_3)$5. 根据$k_1$、$k_2$、$k_3$和$k_4$计算$y_{i+1}$:$y_{i+1}=y_i+\frac{h}{6}(k_1+2k_2+2k_3+k_4)$龙格-库塔法的精度较高，在求解一些对精度要求较高的问题时，龙格-库塔法是一个比较好的选择。
微分方程初值问题的数值解法
积分法:
yk 1 yk h f ( xk , yk ) y ( x0 ) y0
积分项利用矩形公式计算
(1) y( xk 1 ) y( xk )
xk 1
xk
f (t , y(t ))dt
(★)
xk 1
xk
f (t , y(t ))dt h f ( xk , yk ) y( xk 1 ) y( xk ) h f ( xk , yk )
引言
初值问题的数值解法:求初值问题的解在一系列节点的值 y ( xn )的近似值 yn 的方法.本章数值解法的特点:都是采用“步进 式”,即求解过程顺着节点排列的次序一步步向前推进. 常微分方程初值问题: dy f ( x, y ), x [a, b] dx y ( x0 ) y0
替 f (x , y)关于 y 满足Lipschitz条件. 除了要保证(1)有唯一解外,还需保证微分方程本身是稳定的,即 (1)的解连续依赖于初始值和函数 f (x , y). 也就是说, 当初始值 y0 及函数 f (x , y)有微小变化时, 只能引起解的微小变化.
注: 如无特别说明,总假设(1)的解存在唯一且足够光滑. 在 f 连续有界, 则 f (x , y)对变量 y 可微的情形下, 若偏导数 y 可取L为
也称折线法 x
2. 梯形法
若采用梯形公式计算(★)中的积分项,则有 h y ( xk 1 ) y ( xk ) [ f ( xk , y ( xk )) f ( xk 1 , y ( xk 1 ))] 2 h yk 1 yk [ f ( xk , yk ) f ( xk 1 , yk 1 )] 2 称之为梯形公式.这是一个隐式公式,通常用迭代法求解.具体做 法: (0) (0) 先用Euler法求出初值 yk ,1 即 ,将其代入梯形公式 yk 1 yk h f ( xk , yk ) 的右端,使之转化为显式公式,即 h ( l 1) (l ) yk 1 yk [ f ( xk , yk ) f ( xk 1 , yk (☆ ) 1 )] 2
微分方程的数值解法
微分方程的数值解法微分方程是描述自然界中众多现象和规律的重要数学工具。
然而，许多微分方程是很难或者无法直接求解的，因此需要使用数值解法来近似求解。
本文将介绍几种常见的微分方程数值解法。
1. 欧拉方法欧拉方法是最简单的数值解法之一。
它将微分方程转化为差分方程，通过计算离散点上的导数来逼近原方程的解。
欧拉方法的基本思想是利用当前点的导数值来估计下一个点的函数值。
具体步骤如下：首先，将自变量区间等分为一系列的小区间。
然后，根据微分方程的初始条件，在起始点确定初始函数值。
接下来，根据导数的定义，计算每个小区间上函数值的斜率。
最后，根据初始函数值和斜率，递推计算得到每个小区间上的函数值。
2. 龙格-库塔方法龙格-库塔方法是一种常用的高阶精度数值解法。
它通过进行多次逼近和修正来提高近似解的准确性。
相比于欧拉方法，龙格-库塔方法在同样的步长下可以获得更精确的解。
具体步骤如下：首先，确定在每个小区间上的步长。
然后，根据微分方程的初始条件，在起始点确定初始函数值。
接下来，根据当前点的导数值，使用权重系数计算多个中间点的函数值。
最后，根据所有中间点的函数值，计算出当前点的函数值。
3. 改进欧拉方法（改进的欧拉-克罗默法）改进欧拉方法是一种中阶精度数值解法，介于欧拉方法和龙格-库塔方法之间。
它通过使用两公式递推来提高精度，并减少计算量。
改进欧拉方法相对于欧拉方法而言，增加了一个估计项，从而减小了局部截断误差。
具体步骤如下：首先，确定在每个小区间上的步长。
然后，根据微分方程的初始条件，在起始点确定初始函数值。
接下来，利用欧拉方法计算出中间点的函数值。
最后，利用中间点的函数值和斜率，计算出当前点的函数值。
总结：微分方程的数值解法为我们研究和解决实际问题提供了有力的工具。
本文介绍了欧拉方法、龙格-库塔方法和改进欧拉方法这几种常见的数值解法。
选择合适的数值解法取决于微分方程的性质以及对解的精确性要求。
在实际应用中，我们应该根据具体情况选择最合适的数值解法，并注意控制步长以尽可能减小误差。
微分方程数值解使用数值方法求解微分方程
微分方程数值解使用数值方法求解微分方程微分方程是描述自然现象中变化的数学模型，它是数学和科学研究中的重要工具。
然而，许多微分方程并没有精确的解析解，因此需要使用数值方法来近似求解。
本文将介绍一些常用的数值方法来求解微分方程，包括欧拉方法、改进的欧拉方法和龙格-库塔方法。
一、欧拉方法欧拉方法是最简单、最基础的数值方法之一。
它基于微分方程解的定义，通过离散化自变量和因变量来逼近解析解。
假设我们要求解的微分方程为dy/dx = f(x, y)，初始条件为y(x0) = y0。
将自变量x分割成若干个小区间，步长为h，得到x0, x1, x2, ..., xn。
根据微分方程的定义，我们可以得到递推公式 yn+1 = yn + h*f(xn, yn)。
用代码表示即为：```def euler_method(f, x0, y0, h, n):x = [x0]y = [y0]for i in range(n):xn = x[i]yn = y[i]fn = f(xn, yn)xn1 = xn + hyn1 = yn + h*fnx.append(xn1)y.append(yn1)return x, y```二、改进的欧拉方法欧拉方法存在着局部截断误差，即在每个小区间上的误差。
改进的欧拉方法是对欧拉方法的改进，可以减小截断误差。
它的递推公式为yn+1 = yn + h*(f(xn, yn) + f(xn+1, yn+1))/2。
用代码表示即为：```def improved_euler_method(f, x0, y0, h, n):x = [x0]y = [y0]for i in range(n):xn = x[i]yn = y[i]fn = f(xn, yn)xn1 = xn + hyn1 = yn + h*(fn + f(xn1, yn + h*fn))/2x.append(xn1)y.append(yn1)return x, y```三、龙格-库塔方法龙格-库塔方法是一种更加精确的数值方法，它通过计算多个递推式的加权平均值来逼近解析解。
微分方程的数值解法
微分方程的数值解法微分方程（Differential Equation）是描述自然界中变化的现象的重要工具，具有广泛的应用范围。
对于一般的微分方程，往往很难找到解析解，这时候就需要使用数值解法来近似求解微分方程。
本文将介绍几种常见的微分方程数值解法及其原理。
一、欧拉方法（Euler's Method）欧拉方法是最基本也是最容易理解的数值解法之一。
它的基本思想是将微分方程转化为差分方程，通过给定的初始条件，在离散的点上逐步计算出函数的近似值。
对于一阶常微分方程dy/dx = f(x, y)，利用欧拉方法可以得到近似解：y_n+1 = y_n + h * f(x_n, y_n)其中，h是步长，x_n和y_n是已知点的坐标。
欧拉方法的优点在于简单易懂，但是由于是一阶方法，误差较大，对于复杂的微分方程可能不够准确。
二、改进的欧拉方法（Improved Euler's Method）改进的欧拉方法又称为改进的欧拉-柯西方法，是对欧拉方法的一种改进。
它通过在每一步计算中利用两个不同点的斜率来更准确地逼近函数的值。
对于一阶常微分方程dy/dx = f(x, y)，改进的欧拉方法的迭代公式为：y_n+1 = y_n + (h/2) * [f(x_n, y_n) + f(x_n+1, y_n + h * f(x_n, y_n))]相较于欧拉方法，改进的欧拉方法具有更高的精度，在同样的步长下得到的结果更接近真实解。
三、四阶龙格-库塔方法（Fourth-Order Runge-Kutta Method）四阶龙格-库塔方法是一种更高阶的数值解法，通过计算多个点的斜率进行加权平均，得到更为准确的解。
对于一阶常微分方程dy/dx = f(x, y)，四阶龙格-库塔方法的迭代公式为：k1 = h * f(x_n, y_n)k2 = h * f(x_n + h/2, y_n + k1/2)k3 = h * f(x_n + h/2, y_n + k2/2)k4 = h * f(x_n + h, y_n + k3)y_n+1 = y_n + (k1 + 2k2 + 2k3 + k4)/6四阶龙格-库塔方法是数值解法中精度最高的方法之一，它的计算复杂度较高，但是能够提供更为准确的结果。
微分方程的数值解法
微分方程是数学中的一种重要的方程类型，广泛应用于物理、工程、经济等领域。
解微分方程有各种方法，其中数值解法是一种重要而实用的方法。
微分方程的数值解法是通过数值计算来求解微分方程的近似解。
它的基本思想是将微分方程转化为差分方程，并用计算机进行迭代计算，从而求得微分方程的数值解。
数值解法的关键在于如何将微分方程转化为差分方程。
常见的方法有欧拉方法、改进欧拉方法、龙格－库塔方法等。
这些方法都是基于泰勒级数展开的原理进行推导的。
以欧拉方法为例，其基本思路是将微分方程中的导数用差商的方式近似表示，然后通过迭代计算，逐步逼近微分方程的解。
欧拉方法的具体步骤如下：首先确定微分方程的初始条件，即给定t0时刻的函数值y0，然后选取一定的步长ℎ，利用微分方程的导数计算差商y′=dy，进而根据差商dt得到下一个时刻的函数值y n+1=y n+ℎy′。
通过不断迭代计算，即可得到微分方程在一定时间区间内的数值解。
数值解法的另一个重要问题是误差控制。
由于数值计算本身的误差以及近似方法的误差，数值解法所得到的结果通常与真实解存在误差。
为了控制误差，常用的方法有缩小步长ℎ、提高近似方法的阶数等。
此外，还可以通过与解析解进行比较，评估数值解的准确性。
微分方程的数值解法具有以下几点优势。
首先，微分方程的解析解通常较难求得，而数值解法可以给出一个近似解，提供了一种有效的解决方案。
其次，数值解法可以利用计算机的高速运算能力，进行大规模复杂微分方程的求解。
此外，数值解法还可以在实际问题中进行仿真和优化，即通过调整参数来求解微分方程，从而得到最优解。
尽管微分方程的数值解法具有广泛的应用前景，但也存在一些问题和挑战。
首先，数值解法的稳定性和收敛性需要深入研究和分析。
其次，数值解法的计算量通常较大，对计算机运算能力和存储空间的要求较高。
此外，数值解法还需要对问题进行适当的离散化处理，从而可能引入一定的误差。
综上所述，“微分方程的数值解法”是一种重要而实用的方法，可以有效地求解微分方程的近似解。
微分方程数值解法
微分方程数值解法微分方程是数学中的重要概念，它描述了物理系统中变量之间的关系。
解微分方程是许多科学领域中常见的问题，其中又可以分为解析解和数值解两种方法。
本文将重点介绍微分方程的数值解法，并详细讨论其中的常用方法和应用。
一、微分方程的数值解法概述微分方程的解析解往往较为复杂，难以直接求解。
在实际问题中，我们通常利用计算机进行数值计算，以获得方程的数值解。
数值解法的基本思想是将微分方程转化为一组离散的数值问题，通过逼近连续函数来获得数值解。
二、常见的数值解法1. 欧拉法欧拉法是最基础的数值解法之一，其核心思想是将微分方程转化为差分方程，通过逼近连续函数来获得数值解。
欧拉法的基本形式为：yn+1 = yn + h·f(xn, yn)其中，yn表示第n个时间步的数值解，h为时间步长，f为微分方程右端的函数。
欧拉法的精度较低，但计算简单，适用于初步估计或简单系统的求解。
2. 改进的欧拉法（Heun法）改进的欧拉法（Heun法）是对欧拉法的改进，其关键在于求解下一个时间步的近似值时，利用了两个斜率的平均值。
Heun法的基本形式为：yn+1 = yn + (h/2)·(k1 + k2)k1 = f(xn, yn)，k2 = f(xn+h, yn+h·k1)Heun法较欧拉法的精度更高，但计算量较大。
3. 龙格-库塔法（RK方法）龙格-库塔法是一类常用的数值解法，包含了多个不同阶数的方法。
其中，最常用的是经典四阶龙格-库塔法（RK4法），其基本形式为：k1 = f(xn, yn)k2 = f(xn + h/2, yn + (h/2)·k1)k3 = f(xn + h/2, yn + (h/2)·k2)k4 = f(xn + h, yn + h·k3)yn+1 = yn + (h/6)·(k1 + 2k2 + 2k3 + k4)RK4法实现较为复杂，但精度较高，适用于解决大多数常微分方程问题。
求微分方程数值解
求微分方程数值解
微分方程数值解是一种数学方法，用于解决一些复杂的微分方程，特别是那些无法通过解析方法求解的微分方程。
通过数值解法，我们可以得到微分方程的近似解，并且可以在计算机上进行实现，以便更好地理解和分析问题。
我们需要将微分方程转化为差分方程，这样就可以利用数值方法进行求解。
差分方程是一种以离散形式表示微分方程的方法，通过近似替代微分表达式，将连续问题转化为离散问题，从而实现计算机求解。
常见的数值方法包括欧拉方法、龙格-库塔方法等，它们通过不断迭代求解差分方程，逼近微分方程的解。
在应用数值解法求解微分方程时，需要注意选择合适的步长和迭代次数，以确保数值解的准确性和稳定性。
步长过大会导致数值误差增大，步长过小则会增加计算量，影响计算效率。
因此，需要在准确性和效率之间寻找平衡点，选择合适的参数进行计算。
在使用数值解法时，还需要考虑边界条件和初值条件的设定。
这些条件对于微分方程的求解至关重要，不同的条件设定可能会导致不同的数值解，甚至无法得到有效的解。
因此，在进行数值计算之前，需要对问题进行充分的分析和理解，确定合适的条件，以确保数值解的准确性和可靠性。
总的来说，微分方程数值解是一种强大的工具，可以帮助我们解决
复杂的微分方程，探索未知的领域。
通过合理的数值方法和参数选择，我们可以得到准确的数值解，从而更好地理解和应用微分方程的理论。
希望通过不断的探索和实践，我们可以更深入地理解微分方程数值解的原理和方法，为科学研究和工程实践提供更多有益的帮助。
- 1、下载文档前请自行甄别文档内容的完整性，平台不提供额外的编辑、内容补充、找答案等附加服务。
- 2、"仅部分预览"的文档,不可在线预览部分如存在完整性等问题,可反馈申请退款(可完整预览的文档不适用该条件!)。
- 3、如文档侵犯您的权益，请联系客服反馈,我们会尽快为您处理(人工客服工作时间：9:00-18:30)。
微分方程数值解及其应用绪论 自然界中的许多事物的运动和变化规律都可以用微分方程来描述，因此对工程和科学技术中的实际问题的研究中, 常常需要求解微分方程．但往往只有少数较简单和典型的微分方程可求出其解析解，在大多数情况下,只能用近似法求解，数值解法是一类重要的近似方法．本文主要讨论一阶常微分方程的初值问题的数值解法，探讨这些算法在处理来自生活实际问题中的应用，并结合MATLAB 软件，动手编程予以解决． １ 微分方程的初值问题[1]1.1 预备知识在对生活实际问题的研究中，通常需要考虑一阶微分方程的初值问题00(,)()dy f x y dx y x y ⎧=⎪⎨⎪=⎩ （1）这里(),f x y 是矩形区域R ：00,x x a y y b -≤-≤上的连续函数．对初值问题（1）需要考虑以下问题：方程是否一定有解呢?若有解，有多少个解呢?下面给出相关的概念与定理．定义1 Lipschitz 条件[1][2]：矩形区域R ：00,x x a y y b -≤-≤上的连续函数(),f x y 若满足：存在常数0L >，使得不等式()()1212,,f x y f x y L y y -≤-对所有()()12,,,x y x y R ∈都成立，则称(),f x y 在R 上关于y 满足Lipschitz 条件. 定理 1 解的存在唯一性定理[1][3]：设f 在区域()}{,,D x y a x b y R =≤≤∈上连续，关于y 满足Lipschitz 条件，则对任意的[]00,,∈∈x a b y R ，常微分方程初值问题（1）当[],x a b ∈时存在唯一的连续解()y x ．该定理保证若一个函数(),f x y 关于y 满足Lipschitz 条件，它所对应的微分方程的初值问题就有唯一解.在解的存在唯一性得到保证的前提下，自然要考虑方程的求解问题．求解微分方程虽然有多种解析方法，但根据工程和科学实践问题所得到的微分方程往往很复杂，在很多情况下不能或很难给出解析解，有时即使能求出形式解，也往往因形式过于复杂或计算量太大而不实用，因此从实际问题中归结出来的微分方程主要依靠数值解法．定义 2 微分方程数值解：对初值问题（1）寻求数值解就是寻求解()y x 在一系列离散节点上的近似解0121,,,,,,n n y y y y y +L L ，相邻两个节点的间距1n n n h x x +=-称为步长.在一般情况下假定()0,1,i h h i ==L 为常数，这时节点为0,0,1,2,n x x nh n =+=L ．要求微分方程数值解，首先要建立数值算法，即对初值问题（1）中的方程离散化，建立求解数值解法的递推公式．一类是计算1n y +时只用到前一点的值n y ，称为单步法；另一类是用到1n y +前面k 点的值11,,,n n n k y y y --+L 称为k 步法.对初值问题（1）式的单步法可用一般形式表示为11(,,,)n n n n n y y h x y y h ϕ++=+⋅，其中多元函数ϕ与(),f x y 有关，当ϕ含有1n y +时，方法是隐式的；若ϕ中不含1n y +，则为显式方法，所以显式单步法可表示为1(,,)n n n n y y h x y h ϕ+=+⋅.（2）设()y x 是初值问题（1）的准确解，称()()()11(,,)n n n n n T y x y x h x y x h ϕ++=--⋅为显式单步法（2）的局部截断误差. 若存在最大正整数p ,使显式单步法（2）式的局部截断误差满足()()()()11,,p n T y x h y x h x y h O h ϕ++=+--=，则称（2）式有p 阶精度.1.2几种常用的数值解法及其分析、比较1.2.1欧拉法与后退欧拉法1）欧拉法：欧拉曾简单地用差分代替微分，即利用公式将初值问题（1）离散化，则问题（1）可化为1(,),n n n n y y h f x y +=+⋅0n x x n h =+⋅， （3）此方法称为欧拉法.欧拉方法的几何意义在数值计算公式中体现了出来.在xy 平面上，一阶微分方程的解()y y x =称作它的积分曲线.积分曲线上一点(),x y 的切线斜率等于函数(),f x y ，按函数(),f x y 在xy 平面上建立一个方向场，那么，积分曲线上每一点的切线方向均与方向场在该点的方向相一致.基于上述几何解释，从初始点000(,)P x y 出发，先依方向场在该点的方向上推进到1x x =上一点1P ，再从1P 依方向场的方向推进到2x x =上一点2P ，循环前进便作出一条折线012P PP L ，因此欧拉方法又称为折线法.若初值0y 已知，则由（3）式可逐步算出为了分析计算公式的精确度，通常可用泰勒展开将()1n y x +在n x 处展开，则有()()()()()()2''11,,.2++'=+=++∈n n n n n n n n h y x y x h y x y x h y x x ξξ 在()n n y y x =的前提下，()()()(),,.n n n n n f x y f x y x y x '==可得欧拉法（3）的误差为 容易看出，欧拉法（3）式具有一阶精度.2）向后欧拉方法：如果对微分方程（1）从n x 到1n x +积分，得()()()()11,n n x n n x y x y x f t y t dx ++=+⎰，（4） 如果（4）式右端积分用右矩形公式()()11,n n h f x y x ++⋅近似，则得到另一个公式 ()111,n n n n y y hf x y +++=+， （5） 称为后退欧拉法.值得一提的是:后退欧拉法与欧拉公式有着本质的区别，后者是关于1n y +的直接计算公式，它是显式的，而（5）式的右端含有关于1n y +的表达式，它是隐式的.在利用后退欧拉法时，我们通常利用迭代法求解，实质就是逐步显示化.具体迭代过程如下：首先利用欧拉公式(0)1(,)+=+⋅n n n n y y h f x y 给出迭代初值(0)1+n y ，把它代入（5）式的右端，使之转化为显式，直接计算得11(1)(0)1(,)+++=+⋅n n n n y y h f x y .如此反复进行，得 (1)()111(,)++++=+⋅k k n n n n y y h f x y 0,1,k =L ，则得到后退欧拉法的迭代公式(0)1(1)()111(,)(,)+++++⎧=+⋅⎨=+⋅⎩n n n n k k n n n n y y h f x y y y h f x y ， 可以看出，后退欧拉法具有一阶精度，且计算比较麻烦.1.2.2梯形方法为得到比欧拉法精确度高的计算公式，在等式（4）式右端积分中若用梯形求积公式近似，并用n y 代替()n y x ，1n y +代替()1n y x +，则得()()111,,2n n n n n n h y y f x y f x y +++=++⎡⎤⎣⎦， （6） 称其为梯形方法.梯形方法与后退欧拉法一样，都是隐式单步法，可用迭代法求解，其迭代公式为 ()()(0)1(1)()111(,),,2+++++⎧=+⋅⎪⎨⎡⎤=++⎪⎣⎦⎩n n n n k k n n n n n n y y h f x y h y y f x y f x y . （7） 为了分析梯形公式的收敛性，将（6）与（7）式相减，得()()(1)()111111,,2k k n n n n n n h y y f x y f x y +++++++⎡⎤-=-⎣⎦，0,1,2,k =L 因为(),f x y 满足Lipschitz 条件，于是有(1)()11112+++++-≤-k k n n n n hL y y y y ，其中L 为(),f x y 关于y 的Lipschitz 常数.如果选取h 充分小，使得12hL <，则当k →∞时有(1)11+++→k n n y y ，这说明迭代过程（7）式是收敛的[4].容易推导得出梯形法（7）式是二阶方法.经分析，梯形方法虽然提高了精度，但是以增加计算量为代价的.从上述的迭代公式可以看出，每迭代一次都要重新计算(),f x y 的值，而且迭代又要进行若干次，计算相当的复杂.为此，有没有比较简便的计算方法呢？下面给出改进的欧拉方法.1.2.3改进的欧拉方法由前面的讨论可知，梯形法计算相对复杂，现对上面的梯形法进行简化，具体方法是只计算一两次就转入下一步的计算，先用欧拉公式（3）求得一个初步的近似解1n y +，称为预测值，再利用公式（6）把它校正一次，这样建立的预测-校正系统通常称为改进的欧拉公式.具体公式如下()()()11,,,2n n n n n n n n h y y f x y f x y hf x y ++⎡⎤=+++⎣⎦ （8）改进的欧拉法与梯形法一样，是二阶方法.1.2.4 Runge-Kutta 方法由前面讨论可知，从（4）式可以看出，若要使公式阶数提高，就必须使右端积分的数值求积公式精度提高，它必然要增加求积积点，为此将（4）式的右端用求积公式表示为()()()()11,,n n rx i n i n i x i f x y x dx h c f x h y x h λλ+=≈++∑⎰， （9） 一般来说，点数r 越多，精度越高，上式右端相当于增量函数(),,x y h ϕ，为得到便于计算的显式方法，将公式（9）表示为：()1,,,n n n n y y h x y h ϕ+=+（10） 其中()()1111,,,,,2,r n n i i i n n i i n i n ij j j x y h c K K f x y K f x h y h K i r ϕλμ=-=⎧⎪=⎪⎪=⎨⎪⎛⎫⎪=++= ⎪⎪⎝⎭⎩∑∑L （11） 这里,,i i ij c λμ均为常数. i c 为加权因子，i K 为第i 段斜率，共有r 段.我们把（10）和（11）称为r 级显式Runge-Kutta 法，简称为R-K 方法.下面给出其中最经典最常用的一个公式：()()()11234121324322,6,,,22,,22,.n n n n n n n n n n h y y K K K K K f x y h h K f x y K h h K f x y K K f x h y hK +⎧=++++⎪⎪=⎪⎪⎛⎫⎪=++ ⎪⎨⎝⎭⎪⎪⎛⎫=++ ⎪⎪⎝⎭⎪⎪=++⎩（12） Runge-Kutta 方法作为一种重要的单步方法，具有很高的实用价值，它关于初值是稳定的，其解连续地依赖于初值，是一类便于应用的单步法，为了计算1n y +，只用到前面一步的值n y 即可，因此每步的步长可以独立取定.常用的Runge-Kutta 方法精度较高，为了达到预定的精度，与欧拉方法与梯形法相比，步长h 可取得大些，求解区间上的总步数可以少些.但Runge-Kutta 方法也有些缺点，比如四阶Runge-Kutta 方法每算一步需要四次计算(),f x y 的值，计算量较大（对于复杂的(),f x y 而言）． 2 数值方法的应用实例[5-9]例1 对于初值问题()1001y y y '=-⎧⎪⎨=⎪⎩，分别用欧拉法、改进的欧拉法，梯形法求()1y 的近似值.解：易得该方程的解析解()10x y x e -=，()1 4.5400e-005y =，为比较，将按不同数值计算方法所得结果列表如下：表 1 三种不同方法的数值结果图 1 三种不同方法数值解与精确解的误差曲线从表1中可以看出：当0.2h=时，三种方法均不稳定，计算结果严重偏离精确值；h=时，改进后的欧拉和梯形法均稳定，但欧拉法效果很差；当0.01h≤时，三种0.1方法均稳定，但精确度有区别．可以看出，h越小，计算结果越好，要想计算结果充分接近于解析解还须取较小的h值．图1反映的步长0.01h=时，三种数值方法的所得数值解与解析解在[0,1]区间的误差曲线，由图可知，在步长相同的情况下，梯形法的精确度略高于改进的欧拉法；改进的欧拉法和梯形法精确度都明显高于欧拉法．例2用欧拉法、改进的欧拉法和Runge-Kutta法求解初值问题并比较三种方法的结果．解：方程为1n=-的伯努利方程，可求得解析解为现用MATLAB软件编程，用题目要求的方法求解，可得如下图示结果：图2 （a）步长为0.2时R-K法和解析解比较图2 （b）步长为0.2时改进的Euler法和解析解比较图2 （c）步长为0.2时欧拉法和解析解比较上图2(a)，(b)，(c)描述的是步长为0.2时，用欧拉法、改进的欧拉法，Runge-Kutta 法求解方程所得的数值解与解析解之间的对比图．由图可知，Runge-Kutta法所得数值解曲线和解析解曲线吻合的很好，改进的欧拉法和欧拉法随着计算的进行，数值解和解析解之间误差逐步增大，但改进的欧拉法效果要好于欧拉法．图3 (a) 步长为0.1时Euler法和解析解比较图3 (b) 步长为0.1时改进的Euler法和解析解比较图3 (c) 步长为0.1时Runge-Kutta法和解析解比较上图3 (a)，(b)，(c)描述的是步长为0.1时，用欧拉法、改进的欧拉法，Runge-Kutta法求解方程所得的数值解与解析解之间的对比图．由图可知，改进的欧拉法和Runge-Kutta法所得数值解曲线和解析解曲线吻合的很好，而欧拉法随着计算的进行数值解和解析解之间误差逐步增大．相应的程序如下：主程序x=0:0.2:1;jxj=exp(2*x).*(1./exp(4*x) + (2*x)./exp(4*x)).^(1/2);y=Euler(@ff,0,1,0.2,1);gy=geuler(@ff,0,1,0.2,1);Ry=RK(@ff,0,1,0.2,1); figure(1);plot(x,jxj,x,Ry,'*');figure(2);plot(x,jxj,x,gy,'*');figure(3);p lot(x,jxj,x,y,'*')欧拉法程序function y=Euler(f,a,b,h,y0)n=(b-a)/h;x=a:h:b;y=zeros(n+1,1);y(1)=y0;for i=1:ny(i+1)=y(i)+h*feval(f,x(i),y(i));end改进的欧拉法程序function gy=geuler(f,a,b,h,y0)n=(b-a)/h;x=a:h:b;y=zeros(n+1,1);y(1)=y0;for i=1:nyp=y(i)+h*feval(f,x(i),y(i));yc=y(i)+h*feval(f,x(i+1),yp);y(i+1)=(yp+yc)/2;endgy=y;Runge-Kutta法程序function Ry=RK(f,a,b,h,y0)n=(b-a)/h;x=a:h:b;y=zeros(n+1,1);y(1)=y0;for i=1:nk1=feval(f,x(i),y(i));k2=feval(f,x(i)+h/2,y(i)+h*k1/2);k3=feval(f,x(i)+h/2,y(i)+h*k2/2);k4=feval(f,x(i+1),y(i)+h*k3);y(i+1)=y(i)+h*(k1+2*k2+2*k3+k4)/6;endRy=y;3 微分方程数值解法在实际生活中的应用3.1应用实例:耐用消费新产品的销售规律模型一种新产品进入市场以后，常常会经历销售量首先慢慢增加然后逐渐慢慢下降的一个过程，由此给出的时间—销售坐标系下的曲线称为产品的生命曲线,它的形状是钟形的.不过对于较耐用的消费品，情况会有所不一样，它的生命曲线会在开始有个小小的高峰，之后是段比较平坦的曲线，先下降，之后再上升，然后达到高峰，因此它是双峰型的曲线.如何理解这种和传统产品的生命曲线理论相冲突的现象？澳大利亚学者斯蒂芬斯与莫赛观察到的购买耐用消费产品的人大概可分为两类：其一是非常善于接受新的事物，称其为“创新型”消费者，他们会经常从产品广告，制造商给出产品的说明书与商店样品中了解了产品功能与性能之后，再决定其否购买；其二是消费者相对保守，称其为“模仿型”消费者，他们往往会根据大部分已经购买产品的消费者实际使用的经验而提供的信息来决定是否购买产品，下面的例子经过分析，建立相应的数学模型，对这种现象给出了科学解释.3.1.1 模型的建立将顾客获得信息大致可分成两类，其一称之为“搜集型”，源自于产告说明、广告，“创新型”顾客获得这类消息后就可作出其是否购买；第二类称之为“体验型”，即消费者使用之后会获得真实体验，常常以口头的形式散播，“模仿型”类顾客会在获得这种信息之后才能决定是否购买.设K 是潜在用户的总数，1K 与2K 分别为 “创新型”与“模仿型”的人数.再设()N t 是时刻t 已经购买的商品顾客的人数，而()1N t 与()2N t 分别为“创新型”与“模仿型”的顾客的人数,再设()1A t 是时刻t 中已获得“搜集型”信息的人数，由于该部分的信息可直接由外部获得，或者可从已获得该类信息的人群中获得，因此类似巴斯模型，从而建立如下方程：()()()()()()11112112,00,,0dA t K A t A t A dtαααα=-+=> ， 获得“搜集型”信息的“创新型”消费者决定其是否购买的行为，满足如下方程： 对于“模仿型”的顾客，可从已经购买该产品“创新型”或者“模仿型”的顾客中获取信息，于是有在这里，忽略消费者购买产品后需一段短暂的使用后才会散播体验信息的滞后作用． 综上，斯蒂芬斯—莫赛模型是一常微分方程组的初值问题：而()()()12N t N t N t =+为时刻t 购买该商品的总人数.3.1.2 模型的求解对于斯蒂芬斯—莫赛模型中()2N t 的解析解则不能求出，于是可以用四阶Runge Kutta -公式求得，且在它的精度要求达到很高情形下求出()2N t .用MATLAB 软件求解，相应的程序及结果如下function RK=RKFC(fc,a,b,h,y0)n=(b-a)/h;x=a:h:b;m=length(y0);y=zeros(n+1,m);y(1,:)=y0;for i=1:nk1=feval(fc,x(i),y(i,:));k2=feval(fc,x(i)+h/2,y(i,:)+h*k1/2);k3=feval(fc,x(i)+h/2,y(i,:)+h*k2/2);k4=feval(fc,x(i+1),y(i,:)+h*k3); y(i+1,:)=y(i,:)+h*(k1+2*k2+2*k3+k4)/6;endRK=y;function f=FC(x,y)k1=50; k2=70;al=0.0013;be=0.0013;ga=0.0015;f=[(k1-y(1))*(al+be*y(1)),ga*(k2-y(2))*(y(1)+y(2))];x=0:0.3:24;RK=RKFC(@FC,0,24,0.3,[0,0])figure(1)；plot(x,RK(:,1),'+',x,RK(:,2),'*')；legend('N1(t)','N2(t)',2) figure(2)；plot(x,RK(:,1)+RK(:,2),'+')；legend('N(t)',2)图 4 ()()12,N t N t 与时间关系图图 5 ()N t 与[0,25]时间段关系图 由此例可以看出，微分方程数值解在实际生活有着广泛的应用，而数值解法中的Runge-Kutta 方法是一种很重要且应用很广泛的算法.结语微分方程数值解是求解微分方程的一种很重要且应用范围很广的方法，而常用的数值解法如欧拉法、改进的欧拉法、梯形法和Runge-Kutta 方法在一定程度上都有自己的优缺点，理解和掌握各种方法的应用范围，用MATLAB 编制各种方法求解实际问题的通用程序，对用微分方程数值解理论解决现实生活中的实际问题有很重要的意义．参考文献[1] 李庆扬,王能超,易大义.数值分析[M].北京:清华大学出版社,清华大学出版社,2001.[2] 冯康.数值计算方法[M].杭州:浙江大学出版社,2003.[3] 封建湖.计算方法典型题分析解集[M].西安:西北工业大学出版社,2000.[4] 胡建伟,汤怀民.微分方程数值解法(第二版)[M].北京:科学出版社,2007.[5] 王能超.计算方法简明教程[M].北京:高等教育出版社,2004.[6] Nash S G.A history of scientific computing.[M] New York:ACM Press,1990.[7] 于丽妮. ODE 问题解析解及数值解的matlab 实现[J]. 电脑知识与技术.2012,8(14):3303-3305.[8] 霍晓成. 常微分方程数值解法的研究[J]. 临沂师范学院学报. 2011,(6):19-23.[6] 王国立, 陈 瑛.非线性微分方程迭代算法及其在物理学中的应用[J]. 长春师范学院学报( 自然科学版). 2006, 25(2):10-12.。