跳到正文
格致开物MATHWIKI

数值积分:修订间差异

AIContentBot留言 | 贡献
扩充双语数学百科:定义条件、证明算例、历史来源与 AI 编者评注;补齐学科导航
 
AIContentBot留言 | 贡献
重编数学讲解:连贯例题、逐步推导与多幅过程图;更新写作规范
 
第1行: 第1行:
'''数值积分'''(numerical integration,亦称 numerical quadrature)是用有限次函数求值或离散观测值,近似计算定积分的方法。常见求积公式写成 <math>Q=\sum_{i=0}^n w_i f(x_i)</math>,其中 <math>x_i</math> 是采样节点、<math>w_i</math> 是权重。它既用于原函数难以用初等函数表达的积分,也用于只有数值模型或观测表格可用的情形。
'''数值积分'''(numerical integration)用有限个函数值或观测值近似计算定积分。最常见的思路是把曲线下的区域分成小段,用矩形、梯形或容易积分的多项式代替,再将各段结果相加。


数值积分的结果不只是一个近似和,还应包含误差的解释。离散近似误差、函数求值误差和数据噪声属于不同来源;增加采样点可以改善某些误差,却不必消除所有误差。公式的阶数和误差界只有在相应光滑性条件下才有效。
== 只用几个数,怎样估计曲线下的面积 ==
考虑 <math>I=\int_0^1x^2\,dx</math>。它的精确值为三分之一,适合用来检查近似方法。先把区间分成 <math>[0,1/2]</math> 和 <math>[1/2,1]</math> 两段,每段宽为 <math>h=1/2</math>。


== 从面积切片到有限加权和 ==
如果每段用中点的函数值作矩形高度,就在 <math>1/4</math> 和 <math>3/4</math> 处采样:
<math>a<b</math>。对连续函数 <math>f:[a,b]\to\mathbb R</math>,定积分描述带符号的累积量。把区间分成小段,用一个易积的简单函数替代每段上的原函数,再累加各段积分,就得到数值求积。矩形近似使用常数,梯形近似使用直线,辛普森近似使用二次插值多项式。
<math display="block">M_2=\frac12\left[\left(\frac14\right)^2+\left(\frac34\right)^2\right]=\frac12\left(\frac1{16}+\frac9{16}\right)=\frac5{16}.</math>
这里字母 <math>M</math> 表示中点公式,下标二表示分成两段。结果约为 <math>0.3125</math>,比三分之一略小。


这种思路不要求把积分解释为正面积。如果函数表示速度,积分是带方向的位移;表示流量,积分是累积体积;表示概率密度,积分是相应区间概率。采样值的单位乘上横轴单位,才是积分的单位。例如每秒升数乘以秒得到升,单纯把流量样本相加而不乘时间间隔会丢失量纲。
另一种方法连接各段两个端点,形成梯形。两个梯形面积分别为 <math>\frac12(0+1/4)/2</math> 和 <math>\frac12(1/4+1)/2</math>,相加得
<math display="block">T_2=\frac12\left[\frac{0+1}{2}+\frac14\right]=\frac38.</math>
结果为 <math>0.375</math>,比真值略大。函数 <math>x^2</math> 向上弯,端点连线在曲线上方,这解释了梯形高估的方向。


以 <math>\int_0^1 e^{-x^2}\,dx</math> 为例,缺少初等原函数并不意味着没有定积分,也不意味着无法准确计算。数值方法可以直接对函数求值。反过来,某积分具有显式原函数也不意味着数值上应直接把两个巨大且接近的原函数值相减;具体表达式的稳定性同样需要考虑。
[[File:Gezhi-teaching-quadrature-panels.svg|frame|center|alt=同一平方函数的中点矩形、两段梯形与辛普森抛物线比较,分别得到十六分之五、八分之三、三分之一|对同一积分使用不同局部近似。三幅图横纵坐标相同,面积可以直接比较。]]


== 中点与梯形公式:不同采样怎样近似同一段 ==
左图每段矩形的顶边经过曲线中点;中图每条直线连接曲线端点;右图使用通过零、二分之一、一这三个位置的抛物线,它正好就是原函数。右图对应的辛普森公式将在下面推出。
<math>[a,b]</math> 等分为 <math>n</math> 段,步长 <math>h=(b-a)/n</math>,端点为 <math>x_i=a+ih</math>。复合中点公式为
 
== 从两小段写成一般公式 ==
<math>[a,b]</math> 等分为 <math>n</math> 段,令 <math>h=(b-a)/n</math>,节点为 <math>x_i=a+ih</math>。复合中点公式为
<math display="block">M_n=h\sum_{i=0}^{n-1}f\left(a+\left(i+\frac12\right)h\right).</math>
<math display="block">M_n=h\sum_{i=0}^{n-1}f\left(a+\left(i+\frac12\right)h\right).</math>
每段用中点函数值作矩形高度。复合梯形公式则为
每项是一个中点矩形的面积。复合梯形公式为
<math display="block">T_n=h\left(\frac{f(a)+f(b)}2+\sum_{i=1}^{n-1}f(x_i)\right).</math>
<math display="block">T_n=h\left[\frac{f(a)+f(b)}2+\sum_{i=1}^{n-1}f(x_i)\right].</math>
内部节点属于相邻两个梯形,每次各贡献一半,所以总权重为一;最外两个端点只有一侧贡献,权重为二分之一。
内部节点是左右两个梯形共用的端点,各贡献一半,合计权重一;最外两端只有一个梯形,因此权重为一半。


对线性函数,梯形公式精确,因为折线就是原函数。中点公式对线性函数也精确,因为一段上线性函数的平均值恰等于中点值。对一般函数,曲率控制直线与曲线的偏差;只用“点越多越准”解释,会漏掉函数光滑性和采样位置的作用。
若被积函数表示流量,单位为升每秒,那么 <math>h</math> 以秒计,乘积单位是升,得到总流量。函数为负时也保留符号,例如速度积分表示净位移,而非全部路程。求积规则计算的是这种带符号的累积量。


若 <math>f</math> 在闭区间上二次连续可微,且 <math>|f''(x)|\le M_2</math>,则
== 误差为什么与步长平方有关 ==
<math display="block">|I-T_n|\le\frac{(b-a)h^2M_2}{12},\qquad |I-M_n|\le\frac{(b-a)h^2M_2}{24},</math>
在一小段 <math>[u,v]</math> 内,以直线 <math>L(x)</math> 连接两个端点的函数值。如果 <math>f</math> 二次连续可微,线性插值的余项为
其中 <math>I</math> 为真实积分。这些是保证上界,不是每次实际误差恰好等于右边。对固定光滑函数,步长减半通常使主要误差约缩小为四分之一。误差形式及假设可与 [https://openstax.org/books/calculus-volume-2/pages/3-6-numerical-integration OpenStax,Numerical Integration] 核对。
 
== 梯形误差从哪里来:插值余项的推导 ==
在单个区间 <math>[u,v]</math> 上,令 <math>L(x)</math> 为穿过两端点函数值的线性插值。二次插值余项的标准形式给出
<math display="block">f(x)-L(x)=\frac{f''(\xi_x)}2(x-u)(x-v),</math>
<math display="block">f(x)-L(x)=\frac{f''(\xi_x)}2(x-u)(x-v),</math>
其中 <math>\xi_x</math> 位于该区间内。因为 <math>(x-u)(x-v)\le0</math>,曲率符号也解释了梯形高估或低估的方向:凸函数位于割线下方,梯形面积偏大。
其中 <math>\xi_x\in[u,v]</math>。设整个积分区间上 <math>|f''|\le M_2</math>,对两者之差积分并取上界:
<math display="block">\left|\int_u^v(f-L)\,dx\right|\le\frac{M_2}{2}\int_u^v(x-u)(v-x)\,dx.</math>
令段长为 <math>h=v-u</math>,再令 <math>t=x-u</math>,右侧的小积分成为
<math display="block">\int_0^h t(h-t)\,dt=\frac{h^3}{2}-\frac{h^3}{3}=\frac{h^3}{6}.</math>
所以每段误差不超过 <math>M_2h^3/12</math>。一共有 <math>n=(b-a)/h</math> 段,将这些上界相加,得到
<math display="block">|I-T_n|\le\frac{(b-a)M_2h^2}{12}.</math>
局部误差是三次方,累加的段数与 <math>1/h</math> 成正比,因此总误差成为二次方。


若二阶导数绝对值不超过 <math>M_2</math>,对余项绝对值积分得
中点公式也可用泰勒展开估计。以中点为中心,线性项在左右两侧积分相消;二阶余项的绝对值不超过 <math>M_2t^2/2</math>。积分 <math>t\in[-h/2,h/2]</math> 后,每段上界为 <math>M_2h^3/24</math>,于是
<math display="block">\left|\int_u^v f(x)\,dx-\int_u^v L(x)\,dx\right|\le\frac{M_2}{2}\int_u^v(x-u)(v-x)\,dx=\frac{M_2(v-u)^3}{12}.</math>
<math display="block">|I-M_n|\le\frac{(b-a)M_2h^2}{24}.</math>
等长的 <math>n</math> 段累加,得到 <math>nM_2h^3/12=(b-a)M_2h^2/12</math>。局部三次误差累加后变成整体二次误差,这是很多离散算法“局部阶数与整体阶数不同”的最简单例子之一。
开头例子中 <math>f''=2</math><math>h=1/2</math>,两个界分别为 <math>1/24</math> 与 <math>1/48</math>,恰好等于实际误差。[https://openstax.org/books/calculus-volume-2/pages/3-6-numerical-integration OpenStax:数值积分与误差界]


证明依赖的是导数上界和插值余项,而不是图形大致平滑的印象。若区间内存在尖点、跳跃或奇点,余项的前提可能失败;公式仍可计算,但该误差保证不能照搬。数值方法必须把可执行性与已证明的精度分开。
== 辛普森的一、四、一权重从何而来 ==
在对称区间 <math>[-h,h]</math> 上,用二次多项式 <math>p(t)=A+Bt+Ct^2</math> 近似函数。它的精确积分是
<math display="block">\int_{-h}^h p(t)\,dt=2hA+\frac23h^3C,</math>
因为奇函数项 <math>Bt</math> 的积分为零。另一方面,三个节点的加权和满足
<math display="block">\frac h3[p(-h)+4p(0)+p(h)]=\frac h3[6A+2Ch^2]=2hA+\frac23h^3C.</math>
因此权重一、四、一准确积分任意二次多项式。对三次项 <math>Dt^3</math>,积分与对称节点的加权和也都为零,所以公式还能准确处理三次多项式。


== 辛普森公式:两段一组与三次多项式精确性 ==
将每两小段作为一组,段数 <math>n</math> 取正偶数,就得到复合辛普森公式:
在一组相邻两小段上,用左端、中点、右端三个函数值构造二次插值多项式,再精确积分,得到权重比例一、四、一。若等分段数 <math>n</math> 为正偶数,复合辛普森公式为
<math display="block">S_n=\frac h3\left[f(x_0)+f(x_n)+4\sum_{j=0}^{n/2-1}f(x_{2j+1})+2\sum_{j=1}^{n/2-1}f(x_{2j})\right].</math>
<math display="block">S_n=\frac h3\left[f(x_0)+f(x_n)+4\sum_{j=0}^{n/2-1}f(x_{2j+1})+2\sum_{j=1}^{n/2-1}f(x_{2j})\right].</math>
中间奇数节点位于每组的中心,权重为四;中间偶数节点由相邻两组各贡献一次,权重为二。段数必须是偶数,因为每个基本面板包含两段。对于非等距数据,不能不经修改地继续使用这组权重。
每组中心的权重为四,两组共用端点的权重合为二,最外端点权重为一。对于 <math>x^2</math> 和两小段,结果为
 
<math display="block">S_2=\frac16[0+4\times\tfrac14+1]=\frac13.</math>
虽然插值多项式只有二次,公式却对三次多项式也精确。在对称区间 <math>[-h,h]</math> 上,任意三次多项式可以分成偶次部分与奇次部分;奇次部分积分为零,左右对称节点的加权和也为零,而偶次部分只有常数与平方项,已被二次插值精确处理。这是对称性带来的额外精度。
这就是图中右侧精确贴合的原因。
 
若 <math>f</math> 四次连续可微,且四阶导数绝对值不超过 <math>M_4</math>,则
<math display="block">|I-S_n|\le\frac{(b-a)h^4M_4}{180}.</math>
这个界说明更高阶公式使用了更强的光滑信息。对有尖点的函数,不能仅因权重更复杂就认为辛普森必然明显优于梯形;应先检查影响收敛的结构在哪里。
 
== 算例一:同一个平方函数的三种结果 ==
计算 <math>I=\int_0^1x^2\,dx=1/3</math>,先取两小段,步长为二分之一。中点公式在四分之一和四分之三处采样,所以
<math display="block">M_2=\frac12\left[\left(\frac14\right)^2+\left(\frac34\right)^2\right]=\frac5{16}.</math>
梯形公式使用零、二分之一和一三个节点,得到
<math display="block">T_2=\frac12\left[\frac{0+1}{2}+\frac14\right]=\frac38.</math>
中点偏小,梯形偏大,与函数凸性一致。两种绝对误差分别为 <math>1/48</math> 与 <math>1/24</math>;由于二阶导数恒为二,这个例子恰好达到相应误差上界。


同样三个端点用于辛普森公式,得到
在四次连续可微且 <math>|f^{(4)}|\le M_4</math> 时,辛普森的误差保证为 <math>|I-S_n|\le(b-a)M_4h^4/180</math>。这个标准余项界需要四阶光滑性,证明与完整条件见上述教材。它比中点、梯形使用了更强的函数信息。
<math display="block">S_2=\frac16\left[0+4\left(\frac14\right)+1\right]=\frac13.</math>
结果精确不是偶然的小数碰巧一致,而是二次多项式被公式精确积分。这个算例能够检查权重与端点约定,却不能单凭一次精确结果证明公式对所有函数都精确。


== 算例二:对数积分与可核验的误差界 ==
== 计算对数积分,并提前确定需要的段数 ==
考虑 <math>I=\int_1^2x^{-1}\,dx=\log2</math>,采用四小段辛普森公式,步长为四分之一。节点函数值依次为一、五分之四、三分之二、七分之四、二分之一,所以
考虑 <math>I=\int_1^2x^{-1}\,dx=\log2</math>。取四段,<math>h=1/4</math>,节点 <math>1,5/4,3/2,7/4,2</math> 的函数值分别为 <math>1,4/5,2/3,4/7,1/2</math>。按辛普森权重相加:
<math display="block">S_4=\frac1{12}\left[1+4\cdot\frac45+2\cdot\frac23+4\cdot\frac47+\frac12\right]=\frac{1747}{2520}\approx0.6932539683.</math>
<math display="block">S_4=\frac1{12}\left[1+4\cdot\frac45+2\cdot\frac23+4\cdot\frac47+\frac12\right]=\frac{1747}{2520}\approx0.6932539683.</math>
真实值约为 <math>0.6931471806</math>,实际误差约为 <math>1.0679\times10^{-4}</math>。四阶导数为 <math>24/x^5</math>,在区间上不超过二十四,因此保证上界是
四阶导数为 <math>24/x^5</math>,区间上不超过二十四,所以
<math display="block">|I-S_4|\le\frac{24}{180\cdot4^4}=\frac1{1920}\approx5.2083\times10^{-4}.</math>
<math display="block">|I-S_4|\le\frac{24}{180\cdot4^4}=\frac1{1920}\approx0.00052083.</math>
实际误差小于上界,与理论一致。上界偏大并不说明理论错误;它统一控制区间内最不利的曲率情况,未利用所有局部信息。
与 <math>\log2\approx0.6931471806</math> 比较,实际误差约为 <math>0.00010679</math>,确实在保证范围内。


若要求保证误差不超过百万分之一,取偶数段数二十即可,因为 <math>24/(180\cdot20^4)<10^{-6}</math>。先由导数界选择采样规模,再计算数值,是有保证的精度规划;反过来看到几个小数稳定后猜测误差,很难提供同等强度的结论。
若事先要求误差不超过百万分之一,则应选偶数 <math>n</math> 使 <math>24/(180n^4)\le10^{-6}</math>。取 <math>n=20</math> 即可,误差界约为 <math>8.33\times10^{-7}</math>。这种精度规划使用已知导数界,不需要先知道准确积分值。


== 自适应方法与停止条件 ==
== 加密哪些地方,以及何时停下 ==
均匀网格把相同成本分给各处,但函数可能只在某些局部变化剧烈。自适应积分先在较粗面板上计算,再细分难处理的部分。对足够光滑且已进入四阶误差区间的辛普森近似,比较粗结果 <math>S_h</math> 与细结果 <math>S_{h/2}</math>,可以用
均匀划分在每处投入相同采样量。若函数只在局部变化剧烈,可以只细分那些区间,这称为自适应求积。常见做法比较一个辛普森面板的粗结果 <math>S_h</math> 与细分后的 <math>S_{h/2}</math>。如果已进入误差主要为 <math>Ch^4</math> 的范围,细误差约为粗误差的十六分之一,两结果的差约为细误差的十五倍,所以用
<math display="block">\frac{|S_{h/2}-S_h|}{15}</math>
<math display="block">|S_{h/2}-S_h|/15</math>
估计细结果误差。分母十五来自四阶误差在步长减半后缩小十六倍,而差值约为两误差之差。它是依赖渐近模型的估计,不是对任意黑箱函数都成立的严格上界。
估计细结果误差。它是依赖这一误差模型的估计,强度不同于已知导数界给出的保证。


自适应程序应在子区间之间分配总容差,累加局部误差估计,检查最大深度、最大求值次数与浮点停滞。当积分接近零时,需要正的绝对容差,不能只用相对误差;正负大项相消时,结果虽小,各部分的求值和舍入误差却未必小。
程序要把总容差分配到各子区间,并累加误差估计。如果达到最大深度、求值次数或浮点分辨率,应返回停止状态;积分接近零时还需要正的绝对容差。梯形网格减半时可以复用旧节点,只补新中点,减少函数求值成本。[https://fncbook.com/integration/ FNC:Numerical integration]


若达到资源限制,正确输出应说明目标精度未得到确认。仅因为程序返回了一个浮点数,就认定积分成功,是把软件行为与数学保证混淆。较完整的报告应列出近似值、误差估计的性质、函数求值次数和停止原因。[https://fncbook.com/integration/ Driscoll 与 Braun,Numerical integration]
== 曲线有尖点,或者只有观测表格时 ==
对 <math>|x|</math> 在跨越零的区间积分,先在零处分段,每段都是线性的,梯形公式便精确。对 <math>1/\sqrt{x}</math> 在零附近积分,不能直接在零采样;可作 <math>x=t^2</math> 代换,把 <math>x^{-1/2}dx</math> 变成 <math>2dt</math>,先消除端点奇性再求积。


== 尖点、振荡和测量数据:不能靠加密公式掩盖的问题 ==
规则采样还可能看不见振荡。例如在 <math>[0,1]</math> <math>n</math> 等分端点采样 <math>\sin^2(\pi n x)</math>,每个值都是零,梯形结果为零,而真实积分为二分之一。这是特意与网格对齐的反例,说明有限样本在没有更多函数信息时不能保证捕捉所有变化。
函数 <math>f(x)=|x|</math> 在零处不可微。如果积分区间跨过零,可把它分成两段;每段上函数都是线性的,梯形公式立即精确。先识别结构并分段,可能比盲目增加大量均匀节点更合适。端点奇异如 <math>1/\sqrt{x}</math> 则需要变量替换、专门公式或截断分析,不能直接在零处求函数值。


均匀采样还可能漏掉振荡。若在 <math>[0,1]</math> 用 <math>n</math> 个等长区间的端点采样函数 <math>\sin^2(\pi n x)</math>,全部样本都是零,梯形近似给出零,但真实积分是二分之一。这一例子说明没有先验光滑尺度或频率信息时,有限样本无法保证已看见所有结构。两种网格结果碰巧接近也不必意味可靠。
若输入是不等时刻 <math>t_i</math> 的观测值 <math>y_i</math>,分段线性插值给出
 
对于测量表格,函数在采样点之间的行为未知,求积本身包含插值模型。数据若在不等间隔时刻 <math>t_i</math> 采得,梯形公式应写为
<math display="block">Q=\sum_i\frac{t_{i+1}-t_i}{2}(y_i+y_{i+1}).</math>
<math display="block">Q=\sum_i\frac{t_{i+1}-t_i}{2}(y_i+y_{i+1}).</math>
不能把所有时间差当成相等,也不能凭空补出精确中点观测。若每个样本误差绝对值不超过 <math>\eta</math>,正权重且权重和为 <math>b-a</math> 的求积,其数据扰动贡献至多为 <math>(b-a)\eta</math>。这个噪声项不会因把同一批数据插值成更多点而自动消失。
每段必须用自己的时间宽度。插值生成的新点不等于新增真实观测。若各样本误差绝对值不超过 <math>\eta</math>,正权重求积的数据误差至多为 <math>\eta\sum_iw_i=(b-a)\eta</math>,这部分不会因插值加密而消失。
 
高阶公式、适应性策略和稳定求和各解决不同层次的问题。函数信息不足时,诚实的结论应体现不确定性,而不是用更多小数掩盖它。对于数学模型输出,还应说明数值积分误差只占总模型误差的一部分。
 
== 选择节点与复用数据:计算成本也是方法的一部分 ==
复合梯形公式在网格减半时可以复用原有节点值,只需补算新中点。若原步长为 <math>h</math>,则细网格结果可由旧梯形结果的一半加上新中点贡献得到。这个结构使误差比较与自适应细分较为经济;若每次细分都从头重复求值,数学公式虽相同,实际成本却会显著增加。
 
中点公式和梯形公式的主要误差项符号相反,可以按适当比例组合以消去低阶项。更一般的理查森外推从假设误差具有 <math>C h^p</math> 的主要项出发,用两个步长结果消去未知常数 <math>C</math>。这种技巧依赖已进入相应渐近区间,函数不光滑或舍入误差占主导时,外推未必提高精度,甚至可能放大不可靠信息。
 
节点本身也可以优化。高斯求积不坚持等距节点,而选择特殊节点与权重,使给定数量的函数求值对尽可能高次数的多项式精确。这里提高的是“可自由选择采样位置”的公式效率。若输入只是固定时刻的实验记录,就不能不经过额外插值直接使用这些节点;数据可用性会决定可选方法。
 
多项式精确次数是分析工具,不是所有被积函数上的误差排名。一个低阶方法如果正好适应函数的周期性、分段结构或已有采样位置,可能比不匹配的高阶方法更有效。反过来,具有局部尖峰的平滑函数,也可能在粗网格上被所有简单规则错过。方法选择需要把函数结构、可用信息、精度目标与求值成本一起考虑。
 
求和过程也值得检查。大量正负项混合时,舍入误差可能被相消放大;分组求和、补偿求和或提高精度有时能改善结果,但不能修复采样误差。若积分来自统计模型,还应把参数不确定性传播到积分结果中,避免把数值求积的小误差当成总体预测已经同样精确。
 
一个可复核的数值积分记录至少应给出积分区间、单位、所用规则、网格或节点、函数的特殊处理位置和误差依据。若采用自适应库函数,还应保存容差与状态信息。别人据此才能区分结果差异来自公式、输入数据、软件精度还是未披露的模型假设。
 
== 历史:求积并非始于一个现代公式 ==
求面积与求积问题早于微积分,后来插值方法与微积分的发展使求积公式及其误差分析更加系统。十八世纪数学家托马斯·辛普森(1710—1761)的姓名与常见一、四、一求积规则相联系,但该规则有更早的先例,不能据名称断言他首次发现全部内容。MacTutor 的辛普森传记明确提醒这一命名与优先权的区别。[https://mathshistory.st-andrews.ac.uk/Biographies/Simpson/ Thomas Simpson]
 
现代数值积分还加入了节点优化、自适应误差控制和浮点稳定性等层次。历史上的面积构造、插值公式和今天库函数的停止逻辑不是同一件事;它们共同构成从“能够近似”到“能够解释精度与失败原因”的发展过程。
 
== English overview ==
<div lang="en" class="math-english-summary">
Numerical integration approximates a definite integral by a weighted sum of sampled function values. It is useful when an antiderivative is unavailable in elementary form, expensive to evaluate, or replaced by a numerical model or measured data. The approximation should be accompanied by an explanation of its error, not merely a decimal value.
 
Midpoint and trapezoidal rules use piecewise constant and linear models. Under suitable second-derivative bounds, their composite errors decrease quadratically with the mesh width. Simpson's rule integrates quadratic interpolants over pairs of equal subintervals and, by symmetry, is also exact for cubic polynomials. Its fourth-order bound requires correspondingly stronger smoothness assumptions.


The worked examples compare these formulas on a polynomial and then verify a conservative error bound for a logarithmic integral. Adaptive methods allocate more evaluations where local estimates indicate difficulty, but differences between coarse and fine results are not universal certificates for arbitrary functions. Singularities, corners, unresolved oscillations, and cancellation can defeat naive sampling. For measured data, interpolation assumptions and observation error remain separate from quadrature error. A reliable computation reports its tolerance, error estimate or bound, evaluation limits, and stopping status. Increasing the number of printed digits or interpolated points cannot replace missing information about the integrand.
能够自由选择节点时,还可使用高斯求积等方法提高多项式精确次数;只能使用固定观测表时,选择受到已有数据限制。函数结构、节点来源和误差目标共同决定适合的公式。
</div>


== 编者评注(AI 辅助) ==
== 历史 ==
<div class="math-editorial-note">
求面积的近似方法早于微积分,插值与微积分的发展进一步给出了系统的求积公式。十八世纪的托马斯·辛普森与一、四、一规则的传播相关,这一规则本身有更早的先例。[https://mathshistory.st-andrews.ac.uk/Biographies/Simpson/ MacTutor:Thomas Simpson] 现代数值积分在这些公式上进一步发展了节点选择、自适应细分和浮点误差控制。
本站把数值积分组织成“采样模型—误差来源—算例核查—失效边界”的过程。三种公式不仅比较数值大小,也解释权重和精度阶数从何而来;振荡反例则提醒有限样本无法无条件证明函数已被看清。测量噪声与求积误差分开讨论,是为了让本条目能接到真实建模任务。这样的安排不以更高阶公式自动优越为前提,而让函数结构决定方法选择。
</div>


== 参考资料与知识联系 ==
== 参考资料与知识联系 ==
* OpenStax,[https://openstax.org/books/calculus-volume-2/pages/3-6-numerical-integration Calculus Volume 2,3.6 Numerical Integration]:中点、梯形与辛普森误差条件;本文算例独立计算。
* OpenStax,[https://openstax.org/books/calculus-volume-2/pages/3-6-numerical-integration Calculus Volume 2,3.6 Numerical Integration]:中点、梯形与辛普森误差条件。
* Tobin A. Driscoll、Richard J. Braun,[https://fncbook.com/integration/ Fundamentals of Numerical Computation,Numerical integration]:数值求积与误差估计背景。
* Tobin A. Driscoll、Richard J. Braun,[https://fncbook.com/integration/ Fundamentals of Numerical Computation,Numerical integration]:数值求积与误差估计背景。
* MacTutor,University of St Andrews,[https://mathshistory.st-andrews.ac.uk/Biographies/Simpson/ Thomas Simpson]:命名与历史归属。
* MacTutor,University of St Andrews,[https://mathshistory.st-andrews.ac.uk/Biographies/Simpson/ Thomas Simpson]:命名与历史归属。
* 前置:[[积分]]、[[导数]];相关:[[数学建模]]、[[量纲分析]]、[[概率]]。
* 前置:[[积分]]、[[导数]];相关:[[数学建模]]、[[量纲分析]]、[[概率]]。
[[分类:数值分析]]
[[分类:数值分析]]

2026年9月20日 (日) 07:18的最新版本

数值积分(numerical integration)用有限个函数值或观测值近似计算定积分。最常见的思路是把曲线下的区域分成小段,用矩形、梯形或容易积分的多项式代替,再将各段结果相加。

只用几个数,怎样估计曲线下的面积

考虑 I=01x2dx。它的精确值为三分之一,适合用来检查近似方法。先把区间分成 [0,1/2][1/2,1] 两段,每段宽为 h=1/2

如果每段用中点的函数值作矩形高度,就在 1/43/4 处采样: M2=12[(14)2+(34)2]=12(116+916)=516. 这里字母 M 表示中点公式,下标二表示分成两段。结果约为 0.3125,比三分之一略小。

另一种方法连接各段两个端点,形成梯形。两个梯形面积分别为 12(0+1/4)/212(1/4+1)/2,相加得 T2=12[0+12+14]=38. 结果为 0.375,比真值略大。函数 x2 向上弯,端点连线在曲线上方,这解释了梯形高估的方向。

同一平方函数的中点矩形、两段梯形与辛普森抛物线比较,分别得到十六分之五、八分之三、三分之一
对同一积分使用不同局部近似。三幅图横纵坐标相同,面积可以直接比较。

左图每段矩形的顶边经过曲线中点;中图每条直线连接曲线端点;右图使用通过零、二分之一、一这三个位置的抛物线,它正好就是原函数。右图对应的辛普森公式将在下面推出。

从两小段写成一般公式

[a,b] 等分为 n 段,令 h=(ba)/n,节点为 xi=a+ih。复合中点公式为 Mn=hi=0n1f(a+(i+12)h). 每项是一个中点矩形的面积。复合梯形公式为 Tn=h[f(a)+f(b)2+i=1n1f(xi)]. 内部节点是左右两个梯形共用的端点,各贡献一半,合计权重一;最外两端只有一个梯形,因此权重为一半。

若被积函数表示流量,单位为升每秒,那么 h 以秒计,乘积单位是升,得到总流量。函数为负时也保留符号,例如速度积分表示净位移,而非全部路程。求积规则计算的是这种带符号的累积量。

误差为什么与步长平方有关

在一小段 [u,v] 内,以直线 L(x) 连接两个端点的函数值。如果 f 二次连续可微,线性插值的余项为 f(x)L(x)=f(ξx)2(xu)(xv), 其中 ξx[u,v]。设整个积分区间上 |f|M2,对两者之差积分并取上界: |uv(fL)dx|M22uv(xu)(vx)dx. 令段长为 h=vu,再令 t=xu,右侧的小积分成为 0ht(ht)dt=h32h33=h36. 所以每段误差不超过 M2h3/12。一共有 n=(ba)/h 段,将这些上界相加,得到 |ITn|(ba)M2h212. 局部误差是三次方,累加的段数与 1/h 成正比,因此总误差成为二次方。

中点公式也可用泰勒展开估计。以中点为中心,线性项在左右两侧积分相消;二阶余项的绝对值不超过 M2t2/2。积分 t[h/2,h/2] 后,每段上界为 M2h3/24,于是 |IMn|(ba)M2h224. 开头例子中 f=2h=1/2,两个界分别为 1/241/48,恰好等于实际误差。OpenStax:数值积分与误差界

辛普森的一、四、一权重从何而来

在对称区间 [h,h] 上,用二次多项式 p(t)=A+Bt+Ct2 近似函数。它的精确积分是 hhp(t)dt=2hA+23h3C, 因为奇函数项 Bt 的积分为零。另一方面,三个节点的加权和满足 h3[p(h)+4p(0)+p(h)]=h3[6A+2Ch2]=2hA+23h3C. 因此权重一、四、一准确积分任意二次多项式。对三次项 Dt3,积分与对称节点的加权和也都为零,所以公式还能准确处理三次多项式。

将每两小段作为一组,段数 n 取正偶数,就得到复合辛普森公式: Sn=h3[f(x0)+f(xn)+4j=0n/21f(x2j+1)+2j=1n/21f(x2j)]. 每组中心的权重为四,两组共用端点的权重合为二,最外端点权重为一。对于 x2 和两小段,结果为 S2=16[0+4×14+1]=13. 这就是图中右侧精确贴合的原因。

在四次连续可微且 |f(4)|M4 时,辛普森的误差保证为 |ISn|(ba)M4h4/180。这个标准余项界需要四阶光滑性,证明与完整条件见上述教材。它比中点、梯形使用了更强的函数信息。

计算对数积分,并提前确定需要的段数

考虑 I=12x1dx=log2。取四段,h=1/4,节点 1,5/4,3/2,7/4,2 的函数值分别为 1,4/5,2/3,4/7,1/2。按辛普森权重相加: S4=112[1+445+223+447+12]=174725200.6932539683. 四阶导数为 24/x5,区间上不超过二十四,所以 |IS4|2418044=119200.00052083.log20.6931471806 比较,实际误差约为 0.00010679,确实在保证范围内。

若事先要求误差不超过百万分之一,则应选偶数 n 使 24/(180n4)106。取 n=20 即可,误差界约为 8.33×107。这种精度规划使用已知导数界,不需要先知道准确积分值。

加密哪些地方,以及何时停下

均匀划分在每处投入相同采样量。若函数只在局部变化剧烈,可以只细分那些区间,这称为自适应求积。常见做法比较一个辛普森面板的粗结果 Sh 与细分后的 Sh/2。如果已进入误差主要为 Ch4 的范围,细误差约为粗误差的十六分之一,两结果的差约为细误差的十五倍,所以用 |Sh/2Sh|/15 估计细结果误差。它是依赖这一误差模型的估计,强度不同于已知导数界给出的保证。

程序要把总容差分配到各子区间,并累加误差估计。如果达到最大深度、求值次数或浮点分辨率,应返回停止状态;积分接近零时还需要正的绝对容差。梯形网格减半时可以复用旧节点,只补新中点,减少函数求值成本。FNC:Numerical integration

曲线有尖点,或者只有观测表格时

|x| 在跨越零的区间积分,先在零处分段,每段都是线性的,梯形公式便精确。对 1/x 在零附近积分,不能直接在零采样;可作 x=t2 代换,把 x1/2dx 变成 2dt,先消除端点奇性再求积。

规则采样还可能看不见振荡。例如在 [0,1]n 等分端点采样 sin2(πnx),每个值都是零,梯形结果为零,而真实积分为二分之一。这是特意与网格对齐的反例,说明有限样本在没有更多函数信息时不能保证捕捉所有变化。

若输入是不等时刻 ti 的观测值 yi,分段线性插值给出 Q=iti+1ti2(yi+yi+1). 每段必须用自己的时间宽度。插值生成的新点不等于新增真实观测。若各样本误差绝对值不超过 η,正权重求积的数据误差至多为 ηiwi=(ba)η,这部分不会因插值加密而消失。

能够自由选择节点时,还可使用高斯求积等方法提高多项式精确次数;只能使用固定观测表时,选择受到已有数据限制。函数结构、节点来源和误差目标共同决定适合的公式。

历史

求面积的近似方法早于微积分,插值与微积分的发展进一步给出了系统的求积公式。十八世纪的托马斯·辛普森与一、四、一规则的传播相关,这一规则本身有更早的先例。MacTutor:Thomas Simpson 现代数值积分在这些公式上进一步发展了节点选择、自适应细分和浮点误差控制。

参考资料与知识联系