昆仑山:眼中无形心中有穴之穴人合一

夫君子之行,静以修身,俭以养德;非澹泊无以明志,非宁静无以致远。夫学须静也,才须学也;非学无以广才,非志无以成学。怠慢则不能励精,险躁则不能冶性。年与时驰,意与岁去,遂成枯落,多不接世。悲守穷庐,将复何及!

 

9.2简单的数值方法

欧拉法与后退欧拉法 深度讲解与推导证明

各位同学,今天我们用多年数值分析研究与教学的视角,从问题本源、几何意义、多视角推导、误差分析、收敛性证明、核心对比六个维度,彻底讲透欧拉法与后退欧拉法这两个常微分方程数值解法的基石方法。


一、前置基础:我们要解决的核心问题

我们的研究对象是一阶常微分方程初值问题(柯西问题),标准形式为:

\[\begin{cases} y' = f(x,y), \quad x \in [a,b] \\ y(x_0) = y_0 \end{cases} \]

  • 核心目标:求解析解\(y=y(x)\)在离散节点\(x_0 < x_1 < x_2 < \dots < x_n < \dots\)上的近似值\(y_0,y_1,\dots,y_n,\dots\)
  • 步长定义:等步长下\(h = x_{n+1} - x_n\),即相邻两个节点的横坐标间隔
  • 几何本质:微分方程\(y'=f(x,y)\)\(xy\)平面上定义了一个方向场,解\(y=y(x)\)是一条处处与方向场相切的积分曲线,初值条件\(y(x_0)=y_0\)给定了积分曲线的起点\(P_0(x_0,y_0)\),我们的数值方法就是用折线逼近这条积分曲线。

二、欧拉法(显式欧拉法) 详细讲解与推导

欧拉法是最经典的显式单步数值解法,核心思想是用切线近似代替积分曲线,用向前均差近似导数。

2.1 多视角公式推导

我们从三个完全独立的视角推导欧拉法公式,彻底理解其本质。

视角1:几何直观(折线法)推导

从初始点\(P_0(x_0,y_0)\)出发,沿着方向场在\(P_0\)点的切线方向(斜率为\(f(x_0,y_0)\))前进到\(x=x_1\),得到点\(P_1(x_1,y_1)\)
根据直线斜率的定义:

\[\frac{y_1 - y_0}{x_1 - x_0} = f(x_0,y_0) \]

代入步长\(h=x_1-x_0\),整理得:

\[y_1 = y_0 + h f(x_0,y_0) \]

以此类推,从\(P_n(x_n,y_n)\)沿该点切线方向推进到\(P_{n+1}(x_{n+1},y_{n+1})\),得到通用递推公式:

\[\boldsymbol{y_{n+1} = y_n + h f(x_n,y_n), \quad n=0,1,2,\dots} \]

这就是欧拉法公式,也叫显式欧拉法——\(y_{n+1}\)可由前一步的\(x_n,y_n\)直接计算,无需解方程。

视角2:数值微分(向前均差近似)推导

根据导数的定义,当\(h\)足够小时,可用向前均差近似代替\(x_n\)处的导数:

\[y'(x_n) \approx \frac{y(x_{n+1}) - y(x_n)}{h} \]

根据微分方程,\(y'(x_n) = f(x_n, y(x_n))\),代入上式得:

\[\frac{y(x_{n+1}) - y(x_n)}{h} \approx f(x_n, y(x_n)) \]

用近似值\(y_n\)代替准确值\(y(x_n)\)\(y_{n+1}\)代替\(y(x_{n+1})\),将近似号改为等号,即可得到欧拉法公式,与几何推导结果完全一致。

视角3:数值积分(左矩形公式)推导

对微分方程\(y'=f(t,y(t))\)在区间\([x_n, x_{n+1}]\)上两边同时积分:

\[\int_{x_n}^{x_{n+1}} y'(t) dt = \int_{x_n}^{x_{n+1}} f(t, y(t)) dt \]

左边积分结果为\(y(x_{n+1}) - y(x_n)\),因此得到所有单步解法的核心基础式:

\[y(x_{n+1}) = y(x_n) + \int_{x_n}^{x_{n+1}} f(t, y(t)) dt \]

欧拉法采用左矩形数值积分公式近似右端积分:用区间左端点\(x_n\)的函数值\(f(x_n,y(x_n))\)作为整个区间的近似,即:

\[\int_{x_n}^{x_{n+1}} f(t, y(t)) dt \approx h \cdot f(x_n, y(x_n)) \]

代入基础式,替换近似值后,同样得到欧拉法公式。

2.2 欧拉法的局部截断误差与精度阶数

误差分析是数值方法的核心,我们重点分析局部截断误差——即假设前一步计算值完全准确(\(y_n = y(x_n)\)),单步计算带来的误差,它决定了方法的精度阶数。

详细推导过程

将准确解\(y(x_{n+1})=y(x_n+h)\)\(x=x_n\)处做泰勒展开(带拉格朗日余项),展开到二阶项:

\[y(x_n + h) = y(x_n) + y'(x_n) h + \frac{h^2}{2!} y''(\xi_n), \quad \xi_n \in (x_n, x_{n+1}) \]

根据微分方程,\(y'(x_n) = f(x_n, y(x_n))\);在\(y_n = y(x_n)\)的前提下,\(f(x_n,y_n) = y'(x_n)\)
因此欧拉法的计算值为:

\[y_{n+1} = y_n + h f(x_n,y_n) = y(x_n) + h y'(x_n) \]

用准确值减去计算值,得到局部截断误差\(T_{n+1}\)

\[T_{n+1} = y(x_{n+1}) - y_{n+1} = \left[ y(x_n) + h y'(x_n) + \frac{h^2}{2} y''(\xi_n) \right] - \left[ y(x_n) + h y'(x_n) \right] \]

化简后低阶项完全抵消,最终得到:

\[\boldsymbol{T_{n+1} = \frac{h^2}{2} y''(\xi_n) \approx \frac{h^2}{2} y''(x_n)} \]

精度阶数定义

若数值方法的局部截断误差为\(O(h^{p+1})\),则称该方法为\(p\)阶方法
欧拉法的局部截断误差为\(O(h^2)\),因此欧拉法是一阶方法——步长\(h\)缩小10倍,误差仅缩小10倍,这也是教材案例中欧拉法精度不高的核心原因。


三、后退欧拉法(隐式欧拉法) 详细讲解与推导

后退欧拉法是欧拉法的隐式形式,核心是用向后均差近似导数、右矩形公式近似积分,在稳定性上有本质优势。

3.1 多视角公式推导

同样从三个对应视角推导,与欧拉法形成对照。

视角1:数值微分(向后均差近似)推导

与欧拉法的向前均差相反,后退欧拉法用向后均差近似\(x_{n+1}\)处的导数:

\[y'(x_{n+1}) \approx \frac{y(x_{n+1}) - y(x_n)}{h} \]

根据微分方程,\(y'(x_{n+1}) = f(x_{n+1}, y(x_{n+1}))\),代入上式得:

\[\frac{y(x_{n+1}) - y(x_n)}{h} \approx f(x_{n+1}, y(x_{n+1})) \]

用近似值替换准确值,近似号改等号,得到后退欧拉法的核心公式:

\[\boldsymbol{y_{n+1} = y_n + h f(x_{n+1}, y_{n+1}), \quad n=0,1,2,\dots} \]

视角2:数值积分(右矩形公式)推导

回到单步解法的核心基础式:

\[y(x_{n+1}) = y(x_n) + \int_{x_n}^{x_{n+1}} f(t, y(t)) dt \]

后退欧拉法采用右矩形数值积分公式近似右端积分:用区间右端点\(x_{n+1}\)的函数值\(f(x_{n+1},y(x_{n+1}))\)作为整个区间的近似,即:

\[\int_{x_n}^{x_{n+1}} f(t, y(t)) dt \approx h \cdot f(x_{n+1}, y(x_{n+1})) \]

代入基础式后,即可得到后退欧拉法公式。

核心特点:隐式属性

后退欧拉法的公式中,未知量\(y_{n+1}\)同时出现在等式两端,无法直接计算,必须通过迭代法求解,因此也叫隐式欧拉法,这是它与显式欧拉法最本质的区别。

3.2 后退欧拉法的局部截断误差与精度阶数

同样在\(y_n = y(x_n)\)的前提下,推导局部截断误差。

\(y(x_n) = y(x_{n+1} - h)\)\(x=x_{n+1}\)处做泰勒展开,展开到二阶项:

\[y(x_{n+1} - h) = y(x_{n+1}) - y'(x_{n+1}) h + \frac{h^2}{2!} y''(\eta_n), \quad \eta_n \in (x_n, x_{n+1}) \]

整理得:

\[y(x_{n+1}) = y(x_n) + h y'(x_{n+1}) - \frac{h^2}{2} y''(\eta_n) \]

结合微分方程\(y'(x_{n+1})=f(x_{n+1},y(x_{n+1}))\),与后退欧拉法的计算式对比,利用利普希茨条件忽略高阶小量后,最终得到局部截断误差:

\[\boldsymbol{T_{n+1} \approx - \frac{h^2}{2} y''(x_n)} \]

由此可知:后退欧拉法的局部截断误差同样为\(O(h^2)\),也是一阶方法,与欧拉法精度阶数相同,仅误差主导项的符号相反。

3.3 后退欧拉法的迭代求解与收敛性证明

隐式公式无法直接计算,必须通过迭代法将其逐步显式化,我们需要推导迭代格式,并证明其收敛条件。

迭代格式推导

  1. 预测步(初始值):用显式欧拉法给出迭代的初始近似值

\[y_{n+1}^{(0)} = y_n + h f(x_n, y_n) \]

  1. 校正步(迭代更新):将上一次迭代的结果代入隐式公式右端,计算新的近似值

\[y_{n+1}^{(1)} = y_n + h f(x_{n+1}, y_{n+1}^{(0)}) \\ y_{n+1}^{(2)} = y_n + h f(x_{n+1}, y_{n+1}^{(1)}) \\ \vdots \]

最终得到通用迭代格式:

\[\boldsymbol{y_{n+1}^{(k+1)} = y_n + h f(x_{n+1}, y_{n+1}^{(k)}), \quad k=0,1,2,\dots} \]

其中\(k\)为迭代次数,直到相邻两次迭代结果的差值小于预设精度,迭代终止。

迭代收敛性严格证明

收敛性定义:当\(k \to \infty\)时,\(y_{n+1}^{(k)} \to y_{n+1}\)\(y_{n+1}\)为后退欧拉法的准确解),则迭代收敛。

证明过程

  1. 后退欧拉法的准确解满足:

\[y_{n+1} = y_n + h f(x_{n+1}, y_{n+1}) \]

  1. 迭代格式的第\(k+1\)步为:

\[y_{n+1}^{(k+1)} = y_n + h f(x_{n+1}, y_{n+1}^{(k)}) \]

  1. 两式相减并取绝对值:

\[\left| y_{n+1}^{(k+1)} - y_{n+1} \right| = \left| h \cdot \left[ f(x_{n+1}, y_{n+1}^{(k)}) - f(x_{n+1}, y_{n+1}) \right] \right| \]

  1. 利用利普希茨条件\(f(x,y)\)\(y\)满足利普希茨条件,存在常数\(L>0\),使得

\[\left| f(x,y_1) - f(x,y_2) \right| \leq L \left| y_1 - y_2 \right| \]

代入上式得:

\[\left| y_{n+1}^{(k+1)} - y_{n+1} \right| \leq hL \cdot \left| y_{n+1}^{(k)} - y_{n+1} \right| \]

  1. 递推展开得:

\[\left| y_{n+1}^{(k+1)} - y_{n+1} \right| \leq (hL)^{k+1} \cdot \left| y_{n+1}^{(0)} - y_{n+1} \right| \]

  1. 收敛条件:当\(k \to \infty\)时,要使误差趋近于0,必须满足

\[\boldsymbol{hL < 1} \]

即步长\(h < \frac{1}{L}\)\(L\)为利普希茨常数),此时迭代收敛,且\(hL\)越小,收敛速度越快。


四、欧拉法与后退欧拉法 核心知识点总结表

对比维度 欧拉法(显式欧拉法) 后退欧拉法(隐式欧拉法)
核心递推公式 \(y_{n+1} = y_n + h f(x_n, y_n)\) \(y_{n+1} = y_n + h f(x_{n+1}, y_{n+1})\)
公式类型 显式公式:\(y_{n+1}\)可由前一步\(x_n,y_n\)直接计算,无需求解方程 隐式公式:\(y_{n+1}\)出现在等式两端,必须通过迭代法求解
核心推导来源 1. 向前均差近似导数
2. 左矩形数值积分公式
3. 切线法折线逼近积分曲线
1. 向后均差近似导数
2. 右矩形数值积分公式
3. 反向斜率折线逼近积分曲线
局部截断误差主导项 \(\frac{h^2}{2} y''(x_n)\) \(- \frac{h^2}{2} y''(x_n)\),与欧拉法大小相等、符号相反
精度阶数 一阶方法(局部截断误差\(O(h^2)\) 一阶方法(局部截断误差\(O(h^2)\)
计算复杂度 低,每步仅需计算1次\(f(x,y)\),无迭代过程 较高,每步需多次迭代计算\(f(x,y)\),迭代次数由步长和精度要求决定
数值稳定性 条件稳定:仅当步长\(h\)满足限定条件时稳定,步长过大易发散 绝对稳定(A-稳定):对任意\(h>0\),只要\(\text{Re}(\lambda)<0\)(模型方程\(y'=\lambda y\))均稳定,适合刚性方程
迭代收敛条件 无迭代,直接计算 需满足\(hL < 1\)\(L\)\(f\)\(y\)的利普希茨常数)
核心适用场景 步长较小的非刚性方程、快速估算、教学演示 刚性微分方程、对稳定性要求高的场景,可使用较大步长
单步/多步属性 单步方法:仅需前一步\(y_n\)即可计算\(y_{n+1}\) 单步方法:仅需前一步\(y_n\)即可迭代计算\(y_{n+1}\)

五、补充说明与核心结论

  1. 两种方法均为一阶单步方法,是常微分方程数值解法的基础,后续的梯形法、龙格-库塔法均是在此基础上的改进与拓展。
  2. 显式欧拉法的优势是计算简单、无需迭代,劣势是稳定性差、步长受限;隐式欧拉法的优势是稳定性极强,劣势是计算复杂、需要迭代。
  3. 两种方法的误差符号相反,将两者加权平均即可得到二阶精度的梯形法,兼顾精度与稳定性,这也是预测-校正方法的核心思想。

梯形方法 深度讲解与完整推导证明

承接前文欧拉法与后退欧拉法的内容,我们继续以资深数值分析研究员的视角,从方法本源、公式推导、迭代求解、收敛性证明、误差与精度分析、工程应用拓展六个维度,完整讲透梯形方法这一兼顾精度与稳定性的核心常微分方程数值解法。


一、梯形方法的核心思想与前置基础

1.1 方法本源

前文我们已经知道,欧拉法(左矩形积分)与后退欧拉法(右矩形积分)均为一阶方法,精度较低,核心原因是矩形求积公式仅用区间单点的函数值近似整个区间的积分,截断误差较高。

梯形方法的核心改进是:对微分方程积分形式的右端项,采用梯形数值求积公式近似——用区间两端点函数值的平均作为区间的平均高度,本质是用连接两个端点的直线围成的梯形面积近似曲边梯形面积,大幅降低积分近似的截断误差,从而提升数值解法的精度阶数。

1.2 前置核心恒等式

对于一阶常微分方程初值问题

\[\begin{cases} y' = f(x,y), \quad x \in [a,b] \\ y(x_0) = y_0 \end{cases} \]

在区间\([x_n, x_{n+1}]\)上对等式两边同时积分,结合牛顿-莱布尼茨公式,得到所有单步数值解法的核心基础式:

\[y(x_{n+1}) = y(x_n) + \int_{x_n}^{x_{n+1}} f(t, y(t)) dt \tag{9.7} \]

其中步长\(h = x_{n+1} - x_n\),梯形方法正是基于此式,通过高精度的梯形求积公式推导而来。


二、梯形方法公式的完整推导

2.1 梯形数值求积公式

对于定积分\(\int_a^b g(t)dt\),梯形求积公式的形式为:

\[\int_a^b g(t)dt \approx \frac{b-a}{2} \left[ g(a) + g(b) \right] \]

其几何意义是用连接\((a,g(a))\)\((b,g(b))\)的梯形面积近似曲边梯形面积,截断误差为\(-\frac{(b-a)^3}{12}g''(\xi)\)\(\xi \in (a,b)\)),是二阶精度的求积公式,远高于矩形公式的一阶精度。

2.2 梯形法递推公式推导

将梯形求积公式应用到核心基础式(9.7)中,令\(a=x_n\)\(b=x_{n+1}\)\(g(t)=f(t,y(t))\),步长\(h=x_{n+1}-x_n\),可得积分项的近似:

\[\int_{x_n}^{x_{n+1}} f(t, y(t)) dt \approx \frac{h}{2} \left[ f(x_n, y(x_n)) + f(x_{n+1}, y(x_{n+1})) \right] \]

用数值近似值\(y_n\)代替准确解\(y(x_n)\)\(y_{n+1}\)代替准确解\(y(x_{n+1})\),将近似号改为等号,代入基础式(9.7),最终得到梯形方法的核心递推公式

\[\boldsymbol{y_{n+1} = y_n + \frac{h}{2} \left[ f(x_n, y_n) + f(x_{n+1}, y_{n+1}) \right], \quad n=0,1,2,\dots} \tag{9.10} \]

2.3 公式核心属性

梯形方法是隐式单步法

  • 单步法:仅需前一步的计算值\(y_n\)即可完成当前步计算;
  • 隐式属性:未知量\(y_{n+1}\)同时出现在等式两端,无法通过\(y_n\)直接计算,必须通过迭代法求解,这是它与显式欧拉法的核心区别。

三、梯形方法的迭代求解格式

隐式公式的求解核心是「预测-校正」思路:先用显式公式给出迭代初始值,再用梯形公式反复校正,直到满足精度要求。

3.1 完整迭代格式推导

  1. 预测步(迭代初值):用显式欧拉法给出\(y_{n+1}\)的初始近似值,作为迭代的起点

\[y_{n+1}^{(0)} = y_n + h f(x_n, y_n) \]

其中上标\((k)\)表示第\(k\)次迭代的计算结果。

  1. 校正步(迭代更新):将上一次迭代的结果\(y_{n+1}^{(k)}\)代入梯形公式右端,计算得到第\(k+1\)次迭代的近似值

\[y_{n+1}^{(k+1)} = y_n + \frac{h}{2} \left[ f(x_n, y_n) + f(x_{n+1}, y_{n+1}^{(k)}) \right], \quad k=0,1,2,\dots \]

最终得到梯形法的完整迭代格式:

\[\boldsymbol{ \begin{cases} y_{n+1}^{(0)} = y_n + h f(x_n, y_n) \\ y_{n+1}^{(k+1)} = y_n + \frac{h}{2} \left[ f(x_n, y_n) + f(x_{n+1}, y_{n+1}^{(k)}) \right], \quad k=0,1,2,\dots \end{cases} } \tag{9.11} \]

3.2 迭代终止条件

工程计算中,通常预设精度阈值\(\varepsilon\),当相邻两次迭代结果满足

\[\left| y_{n+1}^{(k+1)} - y_{n+1}^{(k)} \right| < \varepsilon \]

时,迭代终止,取\(y_{n+1} = y_{n+1}^{(k+1)}\)作为当前步的最终计算结果。


四、梯形法迭代的收敛性严格证明

我们需要证明:当步长\(h\)满足特定条件时,迭代次数\(k \to \infty\)时,迭代结果\(y_{n+1}^{(k)}\)会收敛到梯形法隐式方程的准确解\(y_{n+1}\)

4.1 证明前提

假设\(f(x,y)\)\(y\)满足利普希茨条件:存在常数\(L>0\)(利普希茨常数),对任意\(x\)、任意\(y_1,y_2\),均有

\[\left| f(x,y_1) - f(x,y_2) \right| \leq L \left| y_1 - y_2 \right| \]

该条件是常微分方程初值问题解存在唯一的核心条件,也是数值方法收敛性的基础。

4.2 完整证明过程

  1. 梯形法的准确解\(y_{n+1}\)满足隐式方程:

\[y_{n+1} = y_n + \frac{h}{2} \left[ f(x_n, y_n) + f(x_{n+1}, y_{n+1}) \right] \tag{1} \]

  1. 迭代格式的第\(k+1\)步为:

\[y_{n+1}^{(k+1)} = y_n + \frac{h}{2} \left[ f(x_n, y_n) + f(x_{n+1}, y_{n+1}^{(k)}) \right] \tag{2} \]

  1. 用(1)式减去(2)式,消去相同项\(y_n\)\(\frac{h}{2}f(x_n,y_n)\),得到:

\[y_{n+1} - y_{n+1}^{(k+1)} = \frac{h}{2} \left[ f(x_{n+1}, y_{n+1}) - f(x_{n+1}, y_{n+1}^{(k)}) \right] \]

  1. 对等式两边取绝对值,结合利普希茨条件,可得:

\[\left| y_{n+1} - y_{n+1}^{(k+1)} \right| = \frac{h}{2} \left| f(x_{n+1}, y_{n+1}) - f(x_{n+1}, y_{n+1}^{(k)}) \right| \leq \frac{hL}{2} \left| y_{n+1} - y_{n+1}^{(k)} \right| \]

  1. 对不等式进行递推展开:

\[\left| y_{n+1} - y_{n+1}^{(k+1)} \right| \leq \frac{hL}{2} \left| y_{n+1} - y_{n+1}^{(k)} \right| \leq \left( \frac{hL}{2} \right)^2 \left| y_{n+1} - y_{n+1}^{(k-1)} \right| \leq \dots \leq \left( \frac{hL}{2} \right)^{k+1} \left| y_{n+1} - y_{n+1}^{(0)} \right| \]

  1. 收敛条件分析:
    初始迭代值\(y_{n+1}^{(0)}\)是欧拉法的预测值,为有限值,因此\(\left| y_{n+1} - y_{n+1}^{(0)} \right|\)是有限常数。当\(k \to \infty\)时,要使误差趋近于0,必须满足:

\[\boldsymbol{\frac{hL}{2} < 1} \]

即步长\(h < \frac{2}{L}\)

当该条件满足时,\(\left( \frac{hL}{2} \right)^{k+1} \to 0\)\(k \to \infty\)),因此\(\left| y_{n+1} - y_{n+1}^{(k+1)} \right| \to 0\),即\(y_{n+1}^{(k)} \to y_{n+1}\),迭代收敛。

4.3 收敛性补充说明

与后退欧拉法的收敛条件\(hL < 1\)相比,梯形法的收敛条件\(\frac{hL}{2} < 1\)更宽松,允许的最大步长是后退欧拉法的2倍,且每迭代一次误差的衰减速度更快,计算效率更高。


五、梯形法的局部截断误差与精度阶数

这是梯形法相比欧拉法最核心的优势,我们通过泰勒展开严格推导其误差量级与精度阶数。

5.1 局部截断误差定义

局部截断误差是指:假设前一步计算值完全准确(\(y_n = y(x_n)\)),单步计算带来的误差,即

\[T_{n+1} = y(x_{n+1}) - y_{n+1} \]

它决定了数值方法的精度阶数:若局部截断误差为\(O(h^{p+1})\),则称该方法为\(p\)阶方法。

5.2 误差与精度推导

  1. 将准确解\(y(x_{n+1})=y(x_n+h)\)\(x=x_n\)处做泰勒展开,展开到三阶项:

\[y(x_n + h) = y(x_n) + y'(x_n)h + \frac{h^2}{2!}y''(x_n) + \frac{h^3}{3!}y'''(\xi_n), \quad \xi_n \in (x_n, x_{n+1}) \tag{3} \]

根据微分方程,\(y'(x) = f(x,y(x))\),因此\(y'(x_n)=f(x_n,y(x_n))\)\(y'(x_{n+1})=f(x_{n+1},y(x_{n+1}))\)

  1. \(y'(x_{n+1})=y'(x_n+h)\)\(x=x_n\)处做泰勒展开,展开到一阶项:

\[y'(x_n + h) = y'(x_n) + h y''(x_n) + \frac{h^2}{2!}y'''(\eta_n), \quad \eta_n \in (x_n, x_{n+1}) \tag{4} \]

  1. \(y_n = y(x_n)\)的前提下,梯形法的计算式可写为:

\[y_{n+1} = y(x_n) + \frac{h}{2} \left[ y'(x_n) + y'(x_{n+1}) \right] \tag{5} \]

  1. 将(4)式代入(5)式右端,展开化简:

\[\begin{align*} \text{右端} &= y(x_n) + \frac{h}{2} \left[ y'(x_n) + y'(x_n) + h y''(x_n) + \frac{h^2}{2}y'''(\eta_n) \right] \\ &= y(x_n) + h y'(x_n) + \frac{h^2}{2}y''(x_n) + \frac{h^3}{4}y'''(\eta_n) \end{align*} \]

  1. 用准确解(3)式减去计算值,得到局部截断误差:

\[\begin{align*} T_{n+1} &= y(x_{n+1}) - y_{n+1} \\ &= \left[ y(x_n) + h y'(x_n) + \frac{h^2}{2}y''(x_n) + \frac{h^3}{6}y'''(\xi_n) \right] - \left[ y(x_n) + h y'(x_n) + \frac{h^2}{2}y''(x_n) + \frac{h^3}{4}y'''(\eta_n) \right] \end{align*} \]

其中\(h^0, h^1, h^2\)的低阶项完全抵消,仅剩\(h^3\)量级的项。由于\(\xi_n, \eta_n\)均在区间\((x_n, x_{n+1})\)内,当\(h\)足够小时,\(y'''(x)\)在区间内连续,因此存在\(\zeta_n \in (x_n, x_{n+1})\),可将误差合并为:

\[\boldsymbol{T_{n+1} = - \frac{h^3}{12} y'''(\zeta_n)} \]

5.3 精度阶数结论

梯形法的局部截断误差为\(O(h^3)\),因此梯形法是二阶方法,相比一阶的欧拉法、后退欧拉法,精度提升一个量级:当步长\(h\)缩小10倍时,梯形法的误差缩小1000倍,而一阶欧拉法仅缩小10倍,精度优势极其显著。

同时,梯形法是A-稳定的(绝对稳定),对任意步长\(h>0\),对于模型方程\(y'=\lambda y\)\(\text{Re}(\lambda)<0\))均保持稳定,稳定性与后退欧拉法相当,远优于显式欧拉法,是兼顾精度与稳定性的优秀方法。


六、工程拓展:改进欧拉法(预测-校正法)

梯形法需要多次迭代,计算量较大,工程应用中通常采用一次预测-校正的简化形式,称为改进欧拉法,无需迭代,计算量与显式欧拉法接近,同时保持二阶精度,是工程中最常用的形式之一。

改进欧拉法格式

  1. 预测步:用显式欧拉法计算预测值\(y_{n+1}^p\)

\[y_{n+1}^p = y_n + h f(x_n, y_n) \]

  1. 校正步:用梯形公式计算校正值,作为最终结果

\[y_{n+1} = y_n + \frac{h}{2} \left[ f(x_n, y_n) + f(x_{n+1}, y_{n+1}^p) \right] \]


七、核心知识点总结对比表

对比维度 欧拉法(显式) 后退欧拉法(隐式) 梯形法(隐式)
核心递推公式 \(y_{n+1}=y_n+hf(x_n,y_n)\) \(y_{n+1}=y_n+hf(x_{n+1},y_{n+1})\) \(y_{n+1}=y_n+\frac{h}{2}[f(x_n,y_n)+f(x_{n+1},y_{n+1})]\)
公式类型 显式,直接计算 隐式,需迭代求解 隐式,需迭代求解
推导来源 向前均差近似导数、左矩形积分 向后均差近似导数、右矩形积分 梯形数值求积公式
局部截断误差主导项 \(\frac{h^2}{2}y''(x_n)\) \(-\frac{h^2}{2}y''(x_n)\) \(-\frac{h^3}{12}y'''(x_n)\)
精度阶数 一阶(\(O(h^2)\) 一阶(\(O(h^2)\) 二阶(\(O(h^3)\)
迭代收敛条件 无迭代 \(hL < 1\) \(\frac{hL}{2} < 1\)(更宽松)
数值稳定性 条件稳定,步长受限 A-稳定,任意步长稳定 A-稳定,任意步长稳定
计算复杂度 低,每步1次\(f(x,y)\)计算 较高,每步多次迭代计算 较高,每步多次迭代计算
核心优势 计算简单、无需迭代 稳定性极强 兼顾二阶精度与优秀稳定性

八、核心结论

  1. 梯形法是二阶隐式单步法,通过梯形求积公式实现了精度的跃升,是欧拉法体系中精度最高的基础方法;
  2. 迭代收敛条件宽松,稳定性优秀,既可以通过多次迭代实现高精度计算,也可以通过一次预测-校正简化为工程常用的改进欧拉法;
  3. 梯形法是后续龙格-库塔法、线性多步法等高精度方法的核心基础,其「平均斜率提升精度」的思想贯穿了整个常微分方程数值解法体系。

改进欧拉公式(预测-校正法)深度讲解与完整推导

承接前文欧拉法、后退欧拉法与梯形法的内容,我们系统讲解改进欧拉法——这一工程中最常用的、兼顾显式计算便捷性二阶精度的常微分方程数值解法,完整覆盖方法由来、公式推导、精度证明、计算流程、案例解析与核心对比。


一、方法的核心背景与设计思想

1.1 现有方法的痛点

  • 显式欧拉法:公式简单、计算便捷,无需迭代,但仅为一阶精度,误差大,工程应用中精度不足;
  • 梯形法:二阶精度,误差小,稳定性优秀,但属于隐式公式,每一步都需要多次迭代计算函数\(f(x,y)\),计算量极大,且迭代终止条件难以预判,工程落地性差。

1.2 改进欧拉法的核心设计思路

用「一次预测 + 一次校正」的两步计算,替代梯形法的多次迭代,在不增加过多计算量的前提下,保留梯形法的二阶精度:

  1. 预测步:用显式欧拉法快速计算一个粗糙的近似值(预测值),为梯形法提供迭代初值;
  2. 校正步:将预测值代入梯形公式,仅做一次校正计算,得到最终的校正值,无需反复迭代。

这种设计既保留了显式计算的便捷性(每步仅需计算2次\(f(x,y)\)),又将精度从一阶提升到二阶,完美平衡了精度与计算效率,是工程计算中最常用的一阶常微分方程数值解法之一。


二、改进欧拉法的两种等价公式形式

2.1 预测-校正标准形式

这是改进欧拉法的核心形式,直接对应「先预测、后校正」的设计逻辑:

\[\boldsymbol{ \begin{cases} \text{预测} \quad \bar{y}_{n+1} = y_n + h f(x_n, y_n) \\ \text{校正} \quad y_{n+1} = y_n + \frac{h}{2} \left[ f(x_n, y_n) + f(x_{n+1}, \bar{y}_{n+1}) \right] \end{cases} \quad n=0,1,2,\dots } \tag{9.12} \]

  • 预测步:用欧拉法计算\(x_{n+1}\)处的预测值\(\bar{y}_{n+1}\),仅需1次\(f(x_n,y_n)\)计算;
  • 校正步:将预测值代入梯形公式右端,替代隐式的\(y_{n+1}\),计算得到最终校正值\(y_{n+1}\),仅需额外1次\(f(x_{n+1},\bar{y}_{n+1})\)计算;
  • 整体为显式计算:所有右端项均为已知量,无需迭代,直接得到\(y_{n+1}\)

2.2 平均化等价形式

我们可以将预测-校正形式展开,改写为更直观的「平均斜率」形式,更易理解其精度提升的本质:

  1. 第一步(欧拉法预测):\(y_p = y_n + h f(x_n, y_n)\),本质是用\(x_n\)处的斜率\(k_1=f(x_n,y_n)\)推进一步;
  2. 第二步(后退欧拉法形式):\(y_c = y_n + h f(x_{n+1}, y_p)\),本质是用预测值得到\(x_{n+1}\)处的斜率\(k_2=f(x_{n+1},y_p)\),再推进一步;
  3. 第三步(平均):将两个斜率推进的结果取算术平均,得到最终值。

最终平均化形式为:

\[\boldsymbol{ \begin{cases} y_p = y_n + h f(x_n, y_n) \\ y_c = y_n + h f(x_{n+1}, y_p) \\ y_{n+1} = \frac{1}{2} \left( y_p + y_c \right) \end{cases} \quad n=0,1,2,\dots } \]

精度提升的本质:欧拉法仅用区间左端点的斜率\(k_1\)近似整个区间的平均斜率,而改进欧拉法用区间左右两个端点的斜率的平均值作为区间的平均斜率,大幅降低了斜率近似的误差,从而将精度从一阶提升到二阶,这也是二阶龙格-库塔法的核心思想。


三、改进欧拉法的精度阶数严格证明

我们通过泰勒展开,严格证明改进欧拉法是二阶方法,局部截断误差为\(O(h^3)\),与梯形法精度一致,远高于显式欧拉法的一阶精度。

3.1 证明前提

假设前一步计算值完全准确,即\(y_n = y(x_n)\),我们需要证明局部截断误差\(T_{n+1} = y(x_{n+1}) - y_{n+1} = O(h^3)\)

3.2 完整证明过程

  1. 首先,将准确解\(y(x_{n+1})=y(x_n+h)\)\(x=x_n\)处做泰勒展开,展开到三阶项:

\[y(x_{n+1}) = y(x_n) + h y'(x_n) + \frac{h^2}{2!} y''(x_n) + \frac{h^3}{3!} y'''(\xi_n), \quad \xi_n \in (x_n, x_{n+1}) \tag{1} \]

根据微分方程\(y'=f(x,y)\),有\(y'(x_n) = f(x_n, y(x_n)) = f(x_n, y_n)\),记为\(k_1\)

  1. \(y''(x)\)求导,由复合函数求导法则:

\[y''(x) = \frac{d}{dx}f(x,y(x)) = \frac{\partial f}{\partial x} + \frac{\partial f}{\partial y} \cdot y'(x) = \frac{\partial f}{\partial x} + f \cdot \frac{\partial f}{\partial y} \]

因此\(y''(x_n) = \left. \left( \frac{\partial f}{\partial x} + f \cdot \frac{\partial f}{\partial y} \right) \right|_{(x_n, y_n)}\)

  1. 展开改进欧拉法的计算式:
    首先,预测值\(\bar{y}_{n+1} = y_n + h f(x_n, y_n) = y(x_n) + h y'(x_n)\)
    然后,校正步的\(f(x_{n+1}, \bar{y}_{n+1})\)\((x_n, y_n)\)处做二元泰勒展开,展开到一阶项:

\[f(x_{n+1}, \bar{y}_{n+1}) = f(x_n + h, y_n + h k_1) = f(x_n, y_n) + h \frac{\partial f}{\partial x}\bigg|_{(x_n,y_n)} + h k_1 \frac{\partial f}{\partial y}\bigg|_{(x_n,y_n)} + O(h^2) \]

代入\(k_1 = f(x_n,y_n)\),化简得:

\[f(x_{n+1}, \bar{y}_{n+1}) = y'(x_n) + h \left. \left( \frac{\partial f}{\partial x} + f \cdot \frac{\partial f}{\partial y} \right) \right|_{(x_n,y_n)} + O(h^2) = y'(x_n) + h y''(x_n) + O(h^2) \tag{2} \]

  1. 将(2)代入改进欧拉法的校正公式:

\[\begin{align*} y_{n+1} &= y_n + \frac{h}{2} \left[ f(x_n,y_n) + f(x_{n+1}, \bar{y}_{n+1}) \right] \\ &= y(x_n) + \frac{h}{2} \left[ y'(x_n) + y'(x_n) + h y''(x_n) + O(h^2) \right] \\ &= y(x_n) + h y'(x_n) + \frac{h^2}{2} y''(x_n) + O(h^3) \tag{3} \end{align*} \]

  1. 计算局部截断误差:
    用准确解(1)减去计算值(3),低阶项完全抵消,得到:

\[T_{n+1} = y(x_{n+1}) - y_{n+1} = \frac{h^3}{6} y'''(\xi_n) - O(h^3) = O(h^3) \]

3.3 精度结论

改进欧拉法的局部截断误差为\(O(h^3)\),因此是二阶数值方法,与梯形法精度阶数一致。相比一阶的显式欧拉法,精度提升一个量级:步长\(h\)缩小10倍,改进欧拉法的误差缩小1000倍,而显式欧拉法仅缩小10倍。


四、改进欧拉法的标准计算步骤

我们将计算流程拆解为可直接落地的6个步骤,方便工程实现与手动计算:

  1. 初始化:确定初值条件\(y(x_0)=y_0\),求解区间\([x_0, b]\),步长\(h\),计算总步数\(N = \frac{b - x_0}{h}\)
  2. 节点赋值:计算离散节点\(x_n = x_0 + n h\)\(n=0,1,2,\dots,N-1\)
  3. 预测步:对当前节点\(x_n\),用欧拉法计算预测值\(\bar{y}_{n+1} = y_n + h f(x_n, y_n)\)
  4. 校正步:代入预测值,计算校正值\(y_{n+1} = y_n + \frac{h}{2} \left[ f(x_n, y_n) + f(x_{n+1}, \bar{y}_{n+1}) \right]\)
  5. 循环推进:令\(n = n+1\),重复步骤3-4,直到计算到\(x_N = b\)
  6. 结果输出:输出所有节点的近似值\(y_n\),并可与解析解对比计算误差。

五、教材案例完整解析(例9.2)

5.1 案例初值问题

求解伯努利方程初值问题:

\[\begin{cases} y' = y - \frac{2x}{y}, \quad 0 < x < 1 \\ y(0) = 1 \end{cases} \]

该问题的解析解为\(y = \sqrt{1+2x}\),取步长\(h=0.1\),用改进欧拉法计算各节点的近似值。

5.2 改进欧拉法公式代入

\(f(x,y) = y - \frac{2x}{y}\)代入平均化形式,得到该问题的专属计算格式:

\[\begin{cases} y_p = y_n + h \left( y_n - \frac{2x_n}{y_n} \right) \\ y_c = y_n + h \left( y_p - \frac{2x_{n+1}}{y_p} \right) \\ y_{n+1} = \frac{1}{2} \left( y_p + y_c \right) \end{cases} \]

其中\(x_n = 0 + 0.1n\),初值\(y_0 = 1\)

5.3 手动计算演示(以\(n=0\)\(x_0=0\)为例)

  1. 预测步:\(y_p = 1 + 0.1 \times \left( 1 - \frac{2 \times 0}{1} \right) = 1.1\)
  2. 校正步:\(x_1=0.1\)\(y_c = 1 + 0.1 \times \left( 1.1 - \frac{2 \times 0.1}{1.1} \right) \approx 1.09182\)
  3. 平均得到最终值:\(y_1 = \frac{1}{2} \times (1.1 + 1.09182) = 1.09591\),与教材表9-2中的\(y_1=1.0959\)完全一致。
  4. 解析解对比:\(y(0.1) = \sqrt{1+2\times0.1} = \sqrt{1.2} \approx 1.0954\),误差仅为\(0.0005\),而显式欧拉法的计算值为\(1.1\),误差为\(0.0046\),精度提升近10倍。

5.4 计算结果与精度分析

教材表9-2的完整计算结果与误差对比如下:

\(x_n\) 改进欧拉法计算值\(y_n\) 解析解\(y(x_n)\) 改进欧拉法绝对误差 显式欧拉法计算值 显式欧拉法绝对误差
0.1 1.0959 1.0954 0.0005 1.1000 0.0046
0.2 1.1841 1.1832 0.0009 1.1918 0.0086
0.3 1.2662 1.2649 0.0013 1.2774 0.0125
0.4 1.3434 1.3416 0.0018 1.3582 0.0166
0.5 1.4164 1.4142 0.0022 1.4351 0.0209
0.6 1.4860 1.4832 0.0028 1.5090 0.0258
0.7 1.5525 1.5492 0.0033 1.5803 0.0311
0.8 1.6165 1.6125 0.0040 1.6498 0.0373
0.9 1.6782 1.6733 0.0049 1.7178 0.0445
1.0 1.7379 1.7321 0.0058 1.7848 0.0527

核心结论
改进欧拉法的计算值与解析解的误差远小于显式欧拉法,大多数结果的前三位有效数字与解析解一致,精度显著提升;随着\(x\)增大,误差缓慢累积,但始终保持在较小量级,而显式欧拉法的误差增长速度远快于改进欧拉法。


六、三种核心方法的对比总结表

对比维度 显式欧拉法 梯形法(多次迭代) 改进欧拉法(预测-校正)
核心公式 \(y_{n+1}=y_n+hf(x_n,y_n)\) \(y_{n+1}=y_n+\frac{h}{2}[f(x_n,y_n)+f(x_{n+1},y_{n+1})]\) 预测+校正两步显式计算
公式类型 显式,直接计算 隐式,需多次迭代 显式,无需迭代,两步计算
精度阶数 一阶(局部截断误差\(O(h^2)\) 二阶(局部截断误差\(O(h^3)\) 二阶(局部截断误差\(O(h^3)\)
每步\(f(x,y)\)计算次数 1次 迭代\(k\)次则\(k\) 固定2次,与迭代无关
计算复杂度 极低 高,迭代次数不确定 低,固定计算量
数值稳定性 条件稳定,步长受限 A-稳定,任意步长稳定 条件稳定,稳定区间大于显式欧拉法
工程落地性 差,精度不足 差,计算量过大 优秀,平衡精度与效率
核心适用场景 教学演示、快速粗算 高精度理论计算 工程实际计算、通用场景

七、核心结论与拓展

  1. 改进欧拉法是显式欧拉法与梯形法的完美结合,通过「一次预测+一次校正」的设计,在仅增加1次函数计算的前提下,将精度从一阶提升到二阶,是工程中最实用的二阶单步方法;
  2. 改进欧拉法本质是二阶龙格-库塔法的一种形式,其「用多个点的斜率平均近似区间平均斜率」的思想,是后续高阶龙格-库塔法(如经典四阶龙格-库塔法)的核心基础;
  3. 该方法计算流程固定,无需迭代,易于编程实现,是求解一阶常微分方程初值问题的入门首选工程方法。

单步法的局部截断误差与阶 深度讲解与完整推导

本节是常微分方程数值解法的核心理论基础。我们之前学习的欧拉法、后退欧拉法、梯形法、改进欧拉法,本质都是求解初值问题的单步法,而局部截断误差与精度阶数,是量化评判这些方法精度高低、收敛性好坏的核心标准。本节我们将从通用定义、核心概念辨析、经典方法推导、关键结论四个维度,彻底讲透这一知识点。


一、单步法的通用表达形式

我们研究的一阶常微分方程初值问题标准形式为:

\[\begin{cases} y' = f(x,y), \quad x \in [a,b] \\ y(x_0) = y_0 \end{cases} \]

所有求解该问题的单步法,都可以用如下通用形式统一表达:

\[\boldsymbol{y_{n+1} = y_n + h \varphi(x_n, y_n, y_{n+1}, h), \quad n=0,1,2,\dots} \tag{9.13} \]

其中核心参数定义如下:

  1. \(h = x_{n+1} - x_n\)为计算步长;
  2. \(\varphi\)称为增量函数:它是与微分方程右端\(f(x,y)\)相关的多元函数,本质是单步法在区间\([x_n, x_{n+1}]\)上的“平均斜率”,决定了方法的核心形式与精度;
  3. 显式与隐式的核心区分
    • 若增量函数\(\varphi\)不含未知量\(y_{n+1}\),即\(\varphi = \varphi(x_n, y_n, h)\),则方法为显式单步法,通用形式简化为:

      \[\boldsymbol{y_{n+1} = y_n + h \varphi(x_n, y_n, h), \quad n=0,1,2,\dots} \tag{9.14} \]

      特点:\(y_{n+1}\)可由前一步的\(x_n,y_n\)直接计算,无需迭代。
    • 若增量函数\(\varphi\)包含未知量\(y_{n+1}\),则方法为隐式单步法,特点:\(y_{n+1}\)无法直接计算,需通过迭代求解。

经典方法的增量函数对应表(衔接前文内容)

方法名称 增量函数\(\varphi\) 方法类型
显式欧拉法 \(\varphi(x,y,h) = f(x,y)\) 显式单步法
后退欧拉法 \(\varphi(x_n,y_n,y_{n+1},h) = f(x_{n+1}, y_{n+1})\) 隐式单步法
梯形法 \(\varphi(x_n,y_n,y_{n+1},h) = \frac{1}{2}\left[ f(x_n,y_n) + f(x_{n+1}, y_{n+1}) \right]\) 隐式单步法
改进欧拉法 \(\varphi(x,y,h) = \frac{1}{2}\left[ f(x,y) + f(x+h, y+hf(x,y)) \right]\) 显式单步法

二、局部截断误差的严格定义与核心内涵

2.1 局部截断误差的定义(定义9.1)

\(y(x)\)是初值问题的准确解析解,则显式单步法(9.14)的局部截断误差定义为:

\[\boldsymbol{T_{n+1} = y(x_{n+1}) - y(x_n) - h \varphi(x_n, y(x_n), h)} \]

核心内涵:“局部”的本质含义

定义中最关键的“局部”二字,有两个严格的前提假设,也是理解该概念的核心:

  1. 前序步骤无误差:假设在\(x_n\)处,数值解完全等于准确解,即\(y_n = y(x_n)\),不存在之前步骤的误差累积;
  2. 单步公式固有误差:仅衡量「用准确解代入递推公式」时,公式本身带来的截断误差,与数值计算的累积误差无关。

我们可以通过推导更直观地理解:
\(y_n = y(x_n)\)时,显式单步法的计算值为:

\[y_{n+1} = y_n + h \varphi(x_n, y_n, h) = y(x_n) + h \varphi(x_n, y(x_n), h) \]

用准确解\(y(x_{n+1})\)减去该计算值,得到单步计算的真实误差:

\[y(x_{n+1}) - y_{n+1} = y(x_{n+1}) - \left[ y(x_n) + h \varphi(x_n, y(x_n), h) \right] = T_{n+1} \]

因此,局部截断误差,就是假设前一步完全准确时,单步计算产生的理论误差,它完全由递推公式的形式决定,是方法本身的固有属性。

2.2 局部截断误差的主项

局部截断误差通常可以展开为步长\(h\)的幂级数形式:

\[T_{n+1} = \psi(x_n, y(x_n)) h^{p+1} + O(h^{p+2}) \]

其中,阶数最低的非零项\(\psi(x_n, y(x_n)) h^{p+1}\),称为局部截断误差的主项
主项决定了局部截断误差的量级和变化趋势,是衡量方法精度的核心指标。


三、单步法的精度阶数定义与核心意义

3.1 精度阶数的定义(定义9.2)

\(y(x)\)是初值问题的准确解,若存在最大的整数\(p\),使得单步法的局部截断误差满足:

\[\boldsymbol{T_{n+1} = y(x_n + h) - y(x_n) - h \varphi(x_n, y(x_n), h) = O(h^{p+1})} \tag{9.15} \]

则称该单步法具有\(p\)阶精度,简称\(p\)阶方法。

补充说明:该定义对隐式单步法完全适用,只需将隐式增量函数代入定义即可。

3.2 精度阶数的核心意义

精度阶数\(p\),直接决定了方法的精度高低和误差收敛速度:

  1. 误差衰减速度\(p\)越大,步长\(h\)缩小时,局部截断误差下降得越快。例如:
    • 1阶方法(欧拉法):\(h\)缩小10倍,误差缩小\(10^2=100\)倍;
    • 2阶方法(梯形法、改进欧拉法):\(h\)缩小10倍,误差缩小\(10^3=1000\)倍;
    • 4阶方法(经典龙格-库塔法):\(h\)缩小10倍,误差缩小\(10^5=100000\)倍。
  2. 特殊性质:若微分方程的准确解\(y(x)\)\(p\)次多项式,则具有\(p\)阶精度的方法可以求出精确解。
    原理:\(p\)次多项式的\(p+1\)阶及以上导数恒为0,因此局部截断误差\(T_{n+1}=O(h^{p+1})=0\),单步计算无误差,自然能得到精确解。

四、经典单步法的局部截断误差与阶数完整推导

我们用上述定义,完整推导前文所有经典方法的局部截断误差、误差主项与精度阶数,既验证定义的应用,也形成完整的知识闭环。

4.1 显式欧拉法

  1. 增量函数:\(\varphi(x,y,h) = f(x,y)\)
  2. 代入局部截断误差定义:

    \[T_{n+1} = y(x_{n+1}) - y(x_n) - h f(x_n, y(x_n)) \]

  3. 泰勒展开:将\(y(x_{n+1})=y(x_n+h)\)\(x_n\)处泰勒展开到二阶项:

    \[y(x_n + h) = y(x_n) + h y'(x_n) + \frac{h^2}{2!} y''(x_n) + O(h^3) \]

  4. 代入化简:根据微分方程,\(y'(x_n) = f(x_n, y(x_n))\),因此:

    \[\begin{align*} T_{n+1} &= \left[ y(x_n) + h y'(x_n) + \frac{h^2}{2} y''(x_n) + O(h^3) \right] - y(x_n) - h y'(x_n) \\ &= \frac{h^2}{2} y''(x_n) + O(h^3) \end{align*} \]

  5. 结论:局部截断误差为\(O(h^2)\),因此显式欧拉法是1阶精度方法,局部截断误差主项为\(\boldsymbol{\frac{h^2}{2} y''(x_n)}\)

4.2 后退欧拉法(隐式)

  1. 增量函数:\(\varphi = f(x_{n+1}, y_{n+1})\),代入局部截断误差定义:

    \[T_{n+1} = y(x_{n+1}) - y(x_n) - h f(x_{n+1}, y(x_{n+1})) \]

  2. 泰勒展开:
    • \(y(x_{n+1}) = y(x_n) + h y'(x_n) + \frac{h^2}{2} y''(x_n) + O(h^3)\)
    • \(f(x_{n+1}, y(x_{n+1})) = y'(x_{n+1}) = y'(x_n + h) = y'(x_n) + h y''(x_n) + O(h^2)\)
  3. 代入化简:

    \[\begin{align*} T_{n+1} &= \left[ y(x_n) + h y'(x_n) + \frac{h^2}{2} y''(x_n) + O(h^3) \right] - y(x_n) - h \left[ y'(x_n) + h y''(x_n) + O(h^2) \right] \\ &= -\frac{h^2}{2} y''(x_n) + O(h^3) \end{align*} \]

  4. 结论:局部截断误差为\(O(h^2)\),因此后退欧拉法是1阶精度方法,局部截断误差主项为\(\boldsymbol{-\frac{h^2}{2} y''(x_n)}\),与显式欧拉法主项大小相等、符号相反。

4.3 梯形法(隐式)

  1. 增量函数:\(\varphi = \frac{1}{2}\left[ f(x_n,y_n) + f(x_{n+1}, y_{n+1}) \right]\),代入局部截断误差定义,结合\(y'(x)=f(x,y(x))\),公式可改写为:

    \[T_{n+1} = y(x_{n+1}) - y(x_n) - \frac{h}{2} \left[ y'(x_n) + y'(x_{n+1}) \right] \]

  2. 泰勒展开:将\(y(x_{n+1})\)\(y'(x_{n+1})\)都在\(x_n\)处泰勒展开到三阶项:

    \[y(x_n + h) = y(x_n) + h y'(x_n) + \frac{h^2}{2!} y''(x_n) + \frac{h^3}{3!} y'''(x_n) + O(h^4) \]

    \[y'(x_n + h) = y'(x_n) + h y''(x_n) + \frac{h^2}{2!} y'''(x_n) + O(h^3) \]

  3. 代入化简:

    \[\begin{align*} T_{n+1} &= \left[ y(x_n) + h y'(x_n) + \frac{h^2}{2} y''(x_n) + \frac{h^3}{6} y'''(x_n) + O(h^4) \right] - y(x_n) \\ &\quad - \frac{h}{2} \left[ y'(x_n) + y'(x_n) + h y''(x_n) + \frac{h^2}{2} y'''(x_n) + O(h^3) \right] \\ &= -\frac{h^3}{12} y'''(x_n) + O(h^4) \end{align*} \]

  4. 结论:局部截断误差为\(O(h^3)\),因此梯形法是2阶精度方法,局部截断误差主项为\(\boldsymbol{-\frac{h^3}{12} y'''(x_n)}\)

4.4 改进欧拉法(显式)

  1. 增量函数:\(\varphi(x,y,h) = \frac{1}{2}\left[ f(x,y) + f(x+h, y+hf(x,y)) \right]\)
  2. 代入定义后,通过二元泰勒展开可推导得到\(T_{n+1} = O(h^3)\)
  3. 结论:改进欧拉法是2阶精度方法,与梯形法精度阶数一致。

五、关键结论与拓展

  1. 数值积分与精度阶数的对应关系
    所有单步法都可以从微分方程的积分形式推导而来:

    \[y(x_{n+1}) = y(x_n) + \int_{x_n}^{x_{n+1}} f(t,y(t)) dt \]

    右端积分的数值求积公式的代数精度,直接决定了单步法的精度阶数:

    • 左/右矩形求积公式(1阶代数精度)→ 欧拉法/后退欧拉法(1阶精度);
    • 梯形求积公式(2阶代数精度)→ 梯形法(2阶精度);
    • 辛普森求积公式(4阶代数精度)→ 对应4阶精度的单步法。
      这一规律揭示了提升单步法精度的核心路径:采用更高代数精度的数值求积公式,近似右端的积分项。
  2. 局部与全局截断误差的关系
    局部截断误差是单步的理论误差,而实际计算中每一步的误差都会累积,最终的总误差称为全局截断误差。可以严格证明:

    • 若单步法是\(p\)阶精度,且满足利普希茨条件,则全局截断误差为\(O(h^p)\),即全局误差的阶数比局部误差低1阶;
    • 局部截断误差是决定全局误差的核心因素,局部精度越高,全局误差收敛越快。
  3. 显式与隐式方法的精度对比
    同阶的显式与隐式方法,精度量级一致,但隐式方法通常具有更优秀的数值稳定性,适合求解刚性微分方程;显式方法计算便捷,适合非刚性方程的快速求解。


六、经典单步法核心参数总结表

方法名称 方法类型 增量函数\(\varphi\) 局部截断误差主项 精度阶数\(p\) 局部截断误差量级
显式欧拉法 显式单步法 \(f(x,y)\) \(\frac{h^2}{2}y''(x_n)\) 1 \(O(h^2)\)
后退欧拉法 隐式单步法 \(f(x_{n+1},y_{n+1})\) \(-\frac{h^2}{2}y''(x_n)\) 1 \(O(h^2)\)
改进欧拉法 显式单步法 \(\frac{1}{2}\left[f(x,y)+f(x+h,y+hf(x,y))\right]\) \(O(h^3)\)主项 2 \(O(h^3)\)
梯形法 隐式单步法 \(\frac{1}{2}\left[f(x_n,y_n)+f(x_{n+1},y_{n+1})\right]\) \(-\frac{h^3}{12}y'''(x_n)\) 2 \(O(h^3)\)

核心总结

本节建立了单步法的通用精度评判体系,通过局部截断误差的定义,将不同形式的单步法统一到同一个量化标准下。理解了局部截断误差与阶数,就掌握了判断数值方法精度高低的核心能力,也为后续学习更高阶的龙格-库塔法、线性多步法奠定了理论基础。

posted on 2026-03-06 16:34  Indian_Mysore  阅读(0)  评论(0)    收藏  举报

导航