数值微积分
目录
1 导数与差分算子 3
1.1 导数与差商 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 3
1.2 差分算子 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 3
2 插值型数值导数 5
2.1 一阶数值导数 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 5
2.2 一般的数值导数 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 5
2.3 三点公式 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 6
2.4 基于三次样条的数值导数 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 9
3 Richardson 外推法 9
4 插值型数值积分 10
4.1 梯形公式 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 10
4.2 辛普森公式 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 10
4.3 数值积分的代数精度 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 11
5 Newton-Ctoes 公式 11
5.1 端点为积分节点的 Newton-Ctoes 公式 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 11
5.2 端点不为积分节点的 Newton-Ctoes 公式 . . . . . . . . . . . . . . . . . . . . . . . . . . . . 12
5.3 复化积分公式 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 13
5.3.1 复化 Simpson 公式 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 13
1
5.3.2 复化梯形公式 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 14
5.3.3 复化中点公式 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 14
5.4 自适应求解数值积分 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 15
6 Gauss 求积公式 15
6.1 Gauss 求积公式的基本思想 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 15
6.2 权函数与正交函数组 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 16
6.2.1 权函数 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 16
6.2.2 正交函数组 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 17
6.2.3 正交多项式组 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 17
6.3 Gauss 求积公式的构造 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 18
6.4 Gauss-Legendre 求积公式 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 19
6.4.1 Legendre 多项式 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 19
6.4.2 Gauss-Legendre 求积公式 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 19
1 导数与差分算子
1.1 导数与差商
有三种方法可以计算导数: 前差, 后差, 中心差
𝑓 [𝑥
0
, 𝑥
0
+ ℎ] =
𝑓 (𝑥
0
+ ℎ) − 𝑓 (𝑥
0
)
ℎ
= 𝑓
′
(𝑥
0
) + 𝑂 (ℎ)
𝑓 [𝑥
0
, 𝑥
0
− ℎ] =
𝑓 (𝑥
0
) − 𝑓 (𝑥
0
− ℎ)
ℎ
= 𝑓
′
(𝑥
0
) + 𝑂 (ℎ)
𝑓 [𝑥
0
− ℎ, 𝑥
0
+ ℎ] =
𝑓 (𝑥
0
+ ℎ) − 𝑓 (𝑥
0
− ℎ)
2ℎ
= 𝑓
′
(𝑥
0
) + 𝑂 (ℎ
2
)
1.2 差分算子
平移算子 𝐸 = 𝐸
1
:
(𝐸𝑣)
𝑗
≡ 𝑣
𝑗+1
平移算子 𝐸 是线性算子, 满足以下性质：
1. 递推关系
𝐸
𝑝
𝑣 = 𝐸
𝑝−1
(𝐸𝑣)
2. 平移作用
(𝐸
𝑝
𝑣)
𝑗
= 𝑣
𝑗+ 𝑝
3. 逆运算
(𝐸
−1
𝑣)
𝑗
= 𝑣
𝑗− 1
4. 恒等算子
(𝐸
0
𝑣)
𝑗
= 𝑣
𝑗
且有线性性
(𝛼𝐸
𝑝
+ 𝛽𝐸
𝑞
)𝑣 = 𝛼𝐸
𝑝
𝑣 + 𝛽𝐸
𝑞
𝑣
那么前差算子 𝐷
+
, 后差算子 𝐷
−
, 中心算子 𝐷
0
可以表示为：
𝐷
+
=
𝐸
1
− 𝐸
0
ℎ
, 𝐷
−
=
𝐸
0
− 𝐸
−1
ℎ
, 𝐷
0
=
𝐸
1
− 𝐸
−1
2ℎ
前差和后差的误差为 𝑂(ℎ), 中心差的误差为 𝑂 (ℎ
2
) , 可以如下证明
对于前差有
𝐷
+
𝑓 (𝑥
𝑗
) =
𝑓 (𝑥
𝑗+1
) − 𝑓 (𝑥
𝑗
)
ℎ
在 𝑥
𝑗
展开 𝑓 (𝑥)
𝑓 (𝑥
𝑗+1
) = 𝑓 (𝑥
𝑗
+ ℎ) = 𝑓 (𝑥
𝑗
) + ℎ 𝑓
′
(𝑥
𝑗
) +
ℎ
2
2
𝑓
′′
(𝑥
𝑗
) +
ℎ
3
6
𝑓
′′′
(𝑥
𝑗
) + 𝑂 (ℎ
4
)
代入 𝐷
+
𝑓 (𝑥
𝑗
)
𝐷
+
𝑓 (𝑥
𝑗
) =
[ 𝑓 (𝑥
𝑗
) + ℎ 𝑓
′
(𝑥
𝑗
) +
ℎ
2
2
𝑓
′′
(𝑥
𝑗
) +
ℎ
3
6
𝑓
′′′
(𝑥
𝑗
)] − 𝑓 (𝑥
𝑗
)
ℎ
= 𝑓
′
(𝑥
𝑗
) +
ℎ
2
𝑓
′′
(𝑥
𝑗
) +
ℎ
2
6
𝑓
′′′
(𝑥
𝑗
) + 𝑂 (ℎ
3
)
故
𝐷
+
𝑓 (𝑥
𝑗
) − 𝑓
′
(𝑥
𝑗
) = 𝑂(ℎ)
对于后差有
𝐷
−
𝑓 (𝑥
𝑗
) =
𝑓 (𝑥
𝑗
) − 𝑓 (𝑥
𝑗− 1
)
ℎ
在 𝑥
𝑗
展开 𝑓 (𝑥)
𝑓 (𝑥
𝑗− 1
) = 𝑓 (𝑥
𝑗
− ℎ) = 𝑓 (𝑥
𝑗
) − ℎ 𝑓
′
(𝑥
𝑗
) +
ℎ
2
2
𝑓
′′
(𝑥
𝑗
) −
ℎ
3
6
𝑓
′′′
(𝑥
𝑗
) + 𝑂 (ℎ
4
)
代入 𝐷
−
𝑓 (𝑥
𝑗
) 有
𝐷
−
𝑓 (𝑥
𝑗
) =
𝑓 (𝑥
𝑗
) − [ 𝑓 (𝑥
𝑗
) − ℎ 𝑓
′
(𝑥
𝑗
) +
ℎ
2
2
𝑓
′′
(𝑥
𝑗
) −
ℎ
3
6
𝑓
′′′
(𝑥
𝑗
)]
ℎ
= 𝑓
′
(𝑥
𝑗
) +
ℎ
2
𝑓
′′
(𝑥
𝑗
) −
ℎ
2
6
𝑓
′′′
(𝑥
𝑗
) + 𝑂 (ℎ
3
)
故
𝐷
−
𝑓 (𝑥
𝑗
) − 𝑓
′
(𝑥
𝑗
) = 𝑂(ℎ)
对于中心差有
𝐷
0
𝑓 (𝑥
𝑗
) =
𝑓 (𝑥
𝑗+1
) − 𝑓 (𝑥
𝑗− 1
)
2ℎ
前向展开 𝑥
𝑗+1
= 𝑥
𝑗
+ ℎ
𝑓 (𝑥
𝑗+1
) = 𝑓 (𝑥
𝑗
) + ℎ 𝑓
′
(𝑥
𝑗
) +
ℎ
2
2
𝑓
′′
(𝑥
𝑗
) +
ℎ
3
6
𝑓
′′′
(𝑥
𝑗
) + 𝑂 (ℎ
4
)
后向展开 𝑥
𝑗− 1
= 𝑥
𝑗
− ℎ
𝑓 (𝑥
𝑗− 1
) = 𝑓 (𝑥
𝑗
) − ℎ 𝑓
′
(𝑥
𝑗
) +
ℎ
2
2
𝑓
′′
(𝑥
𝑗
) −
ℎ
3
6
𝑓
′′′
(𝑥
𝑗
) + 𝑂 (ℎ
4
)
代入 𝐷
0
𝑓 (𝑥
𝑗
) 有
𝐷
0
𝑓 (𝑥
𝑗
) =
[ 𝑓 (𝑥
𝑗
) + ℎ 𝑓
′
(𝑥
𝑗
) +
ℎ
2
2
𝑓
′′
(𝑥
𝑗
) +
ℎ
3
6
𝑓
′′′
(𝑥
𝑗
)] − [ 𝑓 (𝑥
𝑗
) − ℎ 𝑓
′
(𝑥
𝑗
) +
ℎ
2
2
𝑓
′′
(𝑥
𝑗
) −
ℎ
3
6
𝑓
′′′
(𝑥
𝑗
)]
2ℎ
=
2ℎ 𝑓
′
(𝑥
𝑗
) +
2ℎ
3
6
𝑓
′′′
(𝑥
𝑗
)
2ℎ
= 𝑓
′
(𝑥
𝑗
) +
ℎ
2
6
𝑓
′′′
(𝑥
𝑗
) + 𝑂 (ℎ
4
)
故
𝐷
0
𝑓 (𝑥
𝑗
) − 𝑓
′
(𝑥
𝑗
) = 𝑂(ℎ
2
)
换句话说
𝐷
+
, 𝐷
−
是
𝜕
𝜕𝑥
的一阶近似,𝐷
0
是
𝜕
𝜕𝑥
的二阶近似
高阶导数可以用差分算子的乘积近似
2 插值型数值导数
2.1 一阶数值导数
由上面的差商公式可以得到
𝐷
+
𝑓 (𝑥
𝑗
) =
𝑓 (𝑥
𝑗+1
) − 𝑓 (𝑥
𝑗
)
ℎ
= 𝑓
′
(𝑥
𝑗
) +
ℎ
2
𝑓
′′
(𝑥
𝑗
) + 𝑂 (ℎ
2
)
因而可以得到导数的近似
𝑓
′
(𝑥
0
) =
1
ℎ
(
𝑓 (𝑥
0
+ ℎ) − 𝑓 (𝑥
0
)
)
−
1
2
ℎ 𝑓
′′
(𝜉
1
)
2.2 一般的数值导数
一般的数值导数可以通过拉格朗日插值多项式得到
𝑓
′
(𝑥
𝑗
) =
𝑛
Õ
𝑘=0
𝑓 (𝑥
𝑘
)𝑙
′
𝑘
(𝑥
𝑗
) +
1
(𝑛 + 1)!
𝑓
(𝑛+1)
(𝜉 (𝑥))
𝑛
Ö
𝑘=0,𝑘≠ 𝑗
(𝑥
𝑗
− 𝑥
𝑘
)
考虑 𝑓 (𝑥) 的拉格朗日插值多项式, 有 𝑛 + 1 个插值点 𝑥
0
, 𝑥
1
, · · · , 𝑥
𝑛
𝑃
𝑛
(𝑥) =
𝑛
Õ
𝑘=0
𝑓 (𝑥
𝑘
)𝑙
𝑘
(𝑥)
其中, 𝑙
𝑘
(𝑥) 为拉格朗日基函数
𝑙
𝑘
(𝑥) =
𝑛
Ö
𝑚=0
𝑚≠𝑘
𝑥 − 𝑥
𝑚
𝑥
𝑘
− 𝑥
𝑚
为了得到导数, 对 𝑃
𝑛
(𝑥) 求导
𝑃
′
𝑛
(𝑥) =
𝑛
Õ
𝑘=0
𝑓 (𝑥
𝑘
)𝑙
′
𝑘
(𝑥) ⇒ 𝑃
′
𝑛
(𝑥
𝑗
) =
𝑛
Õ
𝑘=0
𝑓 (𝑥
𝑘
)𝑙
′
𝑘
(𝑥
𝑗
)
这就得到了前半部分. 为了得到误差, 对拉格朗日插值多项式的插值余项求导
𝑅
′
𝑛
(𝑥) =
𝑑
𝑑𝑥
𝑓
(𝑛+1)
(𝜉 (𝑥))
(𝑛 + 1)!
𝜔(𝑥)
=
𝑓
(𝑛+1)
(𝜉 (𝑥))
(𝑛 + 1)!
𝜔
′
(𝑥)
其中
𝜔(𝑥) =
𝑛
Ö
𝑘=0
(𝑥 − 𝑥
𝑘
)
注意此处未考虑 𝑓
(𝑛+1)
(𝜉 (𝑥)) 的导数, 对于 𝑥
𝑗
点有 𝜔(𝑥
𝑗
) = 0, 认为导数的主要贡献来自于 𝜔
′
(𝑥
𝑗
). 对 𝜔(𝑥)
求导有
𝜔
′
(𝑥) =
𝑛
Õ
𝑘=0
©
«
𝑛
Ö
𝑚=0
𝑚≠𝑘
(𝑥 − 𝑥
𝑚
)
ª
®
®
¬
代入 𝑥
𝑗
, 则含有 𝑥 − 𝑥
𝑗
因子的项全部为零, 因而
𝜔
′
(𝑥
𝑗
) =
𝑛
Ö
𝑚=0
𝑚≠ 𝑗
(𝑥
𝑗
− 𝑥
𝑚
)
故得到余项的导数
𝑅
′
𝑛
(𝑥
𝑗
) =
𝑓
(
𝑛
+
1
)
(𝜉 (𝑥
𝑗
))
(𝑛 + 1)!
𝑛
Ö
𝑘=0
𝑘≠ 𝑗
(𝑥
𝑗
− 𝑥
𝑘
)
即完成了证明
2.3 三点公式
最常使用的是三点公式, 即取 𝑥
𝑗
, 𝑥
𝑗
+ ℎ, 𝑥
𝑗
+ 2ℎ 作为插值点
前差公式
𝑓
′
(𝑥
0
) =
1
2ℎ
[
−3 𝑓 (𝑥
0
) + 4 𝑓 (𝑥
0
+ ℎ) − 𝑓 (𝑥
0
+ 2ℎ)
]
+
1
3
ℎ
2
𝑓
′′′
(𝜉
0
)
中心差公式
𝑓
′
(𝑥
0
) =
1
2ℎ
[
− 𝑓 (𝑥
0
− ℎ) + 𝑓 (𝑥
0
+ ℎ)
]
−
1
6
ℎ
2
𝑓
′′′
(𝜉
1
)
后差公式
𝑓
′
(𝑥
0
) =
1
2ℎ
[
𝑓 (𝑥
0
) − 4 𝑓 (𝑥
0
− ℎ) + 3 𝑓 (𝑥
0
− 2ℎ)
]
+
1
3
ℎ
2
𝑓
′′′
(𝜉
2
)
考察前差, 选取 𝑥
0
, 𝑥
0
+ ℎ, 𝑥
0
+ 2ℎ 作为插值点, 考察 𝑥
0
处的导数. 对 𝑓 (𝑥
0
+ ℎ) 和 𝑓 (𝑥
0
+ 2ℎ) 在 𝑥
0
处进行
泰勒展开
𝑓 (𝑥
0
+ ℎ) = 𝑓 (𝑥
0
) + ℎ 𝑓
′
(𝑥
0
) +
ℎ
2
2
𝑓
′′
(𝑥
0
) +
ℎ
3
6
𝑓
′′′
(𝑥
0
) +
ℎ
4
24
𝑓
′′′′
(𝜉
1
)
𝑓 (𝑥
0
+ 2ℎ) = 𝑓 (𝑥
0
) + 2ℎ 𝑓
′
(𝑥
0
) + 2ℎ
2
𝑓
′′
(𝑥
0
) +
4ℎ
3
3
𝑓
′′′
(𝑥
0
) +
16ℎ
4
24
𝑓
′′′′
(𝜉
2
)
其中,𝜉
1
∈ (𝑥
0
, 𝑥
0
+ ℎ), 𝜉
2
∈ (𝑥
0
, 𝑥
0
+ 2ℎ)
希望消去常数项和尽可能多的高阶导数项, 考察线性组合
𝑎 𝑓 (𝑥
0
) + 𝑏 𝑓 (𝑥
0
+ ℎ) + 𝑐 𝑓 (𝑥
0
+ 2ℎ)
代入展开式, 有
𝑎 𝑓 (𝑥
0
)+𝑏
𝑓 (𝑥
0
) + ℎ 𝑓
′
(𝑥
0
) +
ℎ
2
2
𝑓
′′
(𝑥
0
) +
ℎ
3
6
𝑓
′′′
(𝑥
0
) + · · ·
+𝑐
𝑓 (𝑥
0
) + 2ℎ 𝑓
′
(𝑥
0
) + 2ℎ
2
𝑓
′′
(𝑥
0
) +
4ℎ
3
3
𝑓
′′′
(𝑥
0
) + · · ·
整理得到
(𝑎 + 𝑏 + 𝑐) 𝑓 (𝑥
0
) + (𝑏 · ℎ + 𝑐 · 2ℎ) 𝑓
′
(𝑥
0
) +
𝑏 ·
ℎ
2
2
+ 𝑐 · 2ℎ
2
𝑓
′′
(𝑥
0
) +
𝑏 ·
ℎ
3
6
+ 𝑐 ·
4ℎ
3
3
𝑓
′′′
(𝑥
0
) + · · ·
不可能同时使得三个系数为零, 因为 𝑎, 𝑏, 𝑐 = 0 显然是这样的方程组的解, 而且这样的解是唯一的. 这会
使得 𝑓
′
(𝑥
0
) 项消失. 因而希望消去常数项和二阶导数项. 对于第三个方程, 则添加 𝑓
′
(𝑥) 的归一化条件
𝑎 + 𝑏 + 𝑐 = 0 消去常数项
𝑏ℎ + 2𝑐ℎ = 1 𝑓
′
(𝑥
0
) 系数归一
𝑏
2
ℎ
2
+ 2𝑐ℎ
2
= 0 消去 𝑓
′′
(𝑥
0
)
解得
𝑏 =
2
ℎ
, 𝑐 = −
1
2ℎ
, 𝑎 = −
3
2ℎ
代入即可得到
𝑓
′
(𝑥
0
) ≈
1
2ℎ
[
−3 𝑓 (𝑥
0
) + 4 𝑓 (𝑥
0
+ ℎ) − 𝑓 (𝑥
0
+ 2ℎ)
]
+
1
3
ℎ
2
𝑓
′′′
(𝜉
0
)
对于中心差, 选取三个插值点点 𝑥
0
− ℎ, 𝑥
0
, 𝑥
0
+ ℎ 希望计算 𝑥
0
处的导数. 分别对 𝑓 (𝑥
0
+ ℎ) 和 𝑓 (𝑥
0
− ℎ) 在
𝑥
0
处展开
𝑓 (𝑥
0
+ ℎ) = 𝑓 (𝑥
0
) + ℎ 𝑓
′
(𝑥
0
) +
ℎ
2
2
𝑓
′′
(𝑥
0
) +
ℎ
3
6
𝑓
′′′
(𝑥
0
) +
ℎ
4
24
𝑓
′′′′
(𝜉
3
)
𝑓 (𝑥
0
− ℎ) = 𝑓 (𝑥
0
) − ℎ 𝑓
′
(𝑥
0
) +
ℎ
2
2
𝑓
′′
(𝑥
0
) −
ℎ
3
6
𝑓
′′′
(𝑥
0
) +
ℎ
4
24
𝑓
′′′′
(𝜉
4
)
其中,𝜉
3
∈ (𝑥
0
, 𝑥
0
+ ℎ), 𝜉
4
∈ (𝑥
0
− ℎ, 𝑥
0
)
注意到
𝑓 (𝑥
0
+ ℎ) − 𝑓 (𝑥
0
− ℎ) = 2ℎ 𝑓
′
(𝑥
0
) +
2ℎ
3
6
𝑓
′′′
(𝑥
0
) + 𝑂 (ℎ
4
)
即：
𝑓 (𝑥
0
+ ℎ) − 𝑓 (𝑥
0
− ℎ) = 2ℎ 𝑓
′
(𝑥
0
) +
ℎ
3
3
𝑓
′′′
(𝜉
1
)
得到
𝑓
′
(𝑥
0
) =
𝑓 (𝑥
0
+ ℎ) − 𝑓 (𝑥
0
− ℎ)
2ℎ
−
ℎ
2
6
𝑓
′′′
(𝜉
1
)
其中, 误差项来源于 ℎ
3
𝑓
′′′
(𝑥
0
) 项, 级别为 𝑂(ℎ
2
)
对于后差, 选取三个后向点 𝑥
0
、𝑥
0
− ℎ、𝑥
0
− 2 ℎ 作为插值点, 希望计算 𝑥
0
处的导数. 𝑓 (𝑥
0
− ℎ) 和 𝑓 (𝑥
0
− 2 ℎ)
在 𝑥
0
处展开
𝑓 (𝑥
0
− ℎ) = 𝑓 (𝑥
0
) − ℎ 𝑓
′
(𝑥
0
) +
ℎ
2
2
𝑓
′′
(𝑥
0
) −
ℎ
3
6
𝑓
′′′
(𝑥
0
) +
ℎ
4
24
𝑓
′′′′
(𝜉
1
)
𝑓 (𝑥
0
− 2ℎ) = 𝑓 (𝑥
0
) − 2ℎ 𝑓
′
(𝑥
0
) + 2ℎ
2
𝑓
′′
(𝑥
0
) −
4ℎ
3
6
𝑓
′′′
(𝑥
0
) +
16ℎ
4
24
𝑓
′′′′
(𝜉
2
)
其中,𝜉
1
∈ (𝑥
0
− ℎ, 𝑥
0
), 𝜉
2
∈ (𝑥
0
− 2ℎ, 𝑥
0
)
期望的形式为
𝑎 𝑓 (𝑥
0
) + 𝑏 𝑓 (𝑥
0
− ℎ) + 𝑐 𝑓 (𝑥
0
− 2ℎ) = 𝑓
′
(𝑥
0
) + 高阶项
代入得到
𝑎 𝑓 (𝑥
0
)+𝑏
𝑓 (𝑥
0
) − ℎ 𝑓
′
(𝑥
0
) +
ℎ
2
2
𝑓
′′
(𝑥
0
) −
ℎ
3
6
𝑓
′′′
(𝑥
0
) + · · ·
+𝑐
𝑓 (𝑥
0
) − 2ℎ 𝑓
′
(𝑥
0
) + 2ℎ
2
𝑓
′′
(𝑥
0
) −
4ℎ
3
6
𝑓
′′′
(𝑥
0
) + · · ·
整理有
(𝑎 + 𝑏 + 𝑐) 𝑓 (𝑥
0
) + (−𝑏ℎ − 2𝑐ℎ) 𝑓
′
(𝑥
0
) +
𝑏ℎ
2
2
+ 2𝑐ℎ
2
𝑓
′′
(𝑥
0
) +
−
𝑏ℎ
3
6
−
4𝑐ℎ
3
6
𝑓
′′′
(𝑥
0
) + · · ·
与前差类似, 不可能消去三阶导数项, 考虑消去常数项和二阶导数项, 并添加归一化条件
𝑎 + 𝑏 + 𝑐 = 0 消去常数项
− 𝑏ℎ − 2𝑐ℎ = 1 𝑓
′
(𝑥
0
) 系数归一化
𝑏
2
ℎ
2
+ 2𝑐ℎ
2
= 0 消去 𝑓
′′
(𝑥
0
) 项
解得
𝑎 =
3
2ℎ
, 𝑏 = −
2
ℎ
, 𝑐 =
1
2ℎ
代入整理得到后向差分公式
𝑓
′
(𝑥
0
) ≈
𝑓 (𝑥
0
) − 4 𝑓 (𝑥
0
− ℎ) + 3 𝑓 (𝑥
0
− 2ℎ)
2ℎ
+
1
3
ℎ
2
𝑓
′′′
(𝜉
2
)
其中, 误差项来源于泰勒展开中的 𝑓
′′′
(𝑥
0
) 项, 级别为 𝑂(ℎ
2
)
2.4 基于三次样条的数值导数
可以利用三次样条插值的导数性质来计算导数. 设有插值点
𝑎 = 𝑥
0
< 𝑥
1
< · · · < 𝑥
𝑛
= 𝑏
对于每个区间 [𝑥
𝑖−1
, 𝑥
𝑖
], 有三次多项式
𝑆
𝑖
(𝑥) = 𝑎
𝑖
+ 𝑏
𝑖
(𝑥 − 𝑥
𝑖−1
) + 𝑐
𝑖
(𝑥 − 𝑥
𝑖−1
)
2
+ 𝑑
𝑖
(𝑥 − 𝑥
𝑖−1
)
3
认为
𝑓
′
(𝑥
𝑗
) = 𝑆
′
𝑗
(𝑥
𝑗
)
若对于非插值点处的 𝑥 ∈ [𝑥
𝑖−1
, 𝑥
𝑖
], 认为
𝑓
′
(𝑥) = 𝑆
′
𝑖
(𝑥)
3 Richardson 外推法
假定有一个数值方法, 其采用的步长为 ℎ, 其误差为 𝑂 (ℎ
𝑝
), 即
𝐴(ℎ) = 𝐴 + 𝑐
1
ℎ
𝑝
+ 𝑐
2
ℎ
𝑝+1
+ · · ·
其中 𝐴 为真实值, 而 𝐴(ℎ) 为数值方法计算得到的值. 那么可以取 ℎ →
ℎ
2
, 得到
𝐴
ℎ
2
= 𝐴 + 𝑐
1
1
2
𝑝
ℎ
𝑝
+ 𝑐
2
1
2
𝑝
+
1
ℎ
𝑝+1
+ · · ·
利用线性组合即可消去 𝑐
1
项, 得到
𝐴 =
2
𝑝
𝐴
ℎ
2
− 𝐴(ℎ)
2
𝑝
− 1
+ 𝑂(ℎ
𝑝+1
)
使得误差提高了一阶
4 插值型数值积分
考虑对函数 𝑓 (𝑥) 在区间 [𝑎, 𝑏] 上使用拉格朗日插值多项式
𝑃
𝑛
(𝑥) =
𝑛
Õ
𝑘=0
𝑓 (𝑥
𝑘
)𝑙
𝑘
(𝑥) + 𝐸
𝑛
(𝑥)
其中 𝑙
𝑘
(𝑥) 为拉格朗日基函数
𝑙
𝑘
(𝑥) =
𝑛
Ö
𝑚=0
𝑚≠𝑘
𝑥 − 𝑥
𝑚
𝑥
𝑘
− 𝑥
𝑚
𝐸
𝑛
(𝑥) 为插值余项
𝐸
𝑛
(
𝑥
)
=
𝑓
(𝑛+1)
(𝜉 (𝑥))
(𝑛 + 1)!
𝑛
Ö
𝑘=0
(
𝑥
−
𝑥
𝑘
)
认为 𝑓 (𝑥) 在 [𝑎, 𝑏] 上的积分为 𝑃
𝑛
(𝑥) 的积分, 即
∫
𝑏
𝑎
𝑓 (𝑥)𝑑𝑥 =
𝑛
Õ
𝑘=0
𝑓 (𝑥
𝑘
)
∫
𝑏
𝑎
𝑙
𝑘
(𝑥)𝑑𝑥 +
∫
𝑏
𝑎
𝐸
𝑛
(𝑥)𝑑𝑥
它可以简记为
∫
𝑏
𝑎
𝑓 (𝑥)𝑑𝑥 =
𝑛
Õ
𝑘=0
𝑎
𝑘
𝑓 (𝑥
𝑘
) + 𝑅
𝑛
其中
𝑎
𝑘
=
∫
𝑏
𝑎
𝑙
𝑘
(𝑥)𝑑𝑥, 𝑅
𝑛
=
∫
𝑏
𝑎
𝐸
𝑛
(𝑥)𝑑𝑥
4.1 梯形公式
取 𝑛 = 1, 即取两个插值点 𝑎, 𝑏, 有
∫
𝑏
𝑎
𝑓 (𝑥)𝑑𝑥 =
𝑏 − 𝑎
2
[ 𝑓 (𝑎) + 𝑓 (𝑏)] + 𝑅
1
其中
𝑅
1
=
∫
𝑏
𝑎
𝐸
1
(𝑥)𝑑𝑥 = −
1
12
(𝑏 − 𝑎)
3
𝑓
′′
(𝜉)
4.2 辛普森公式
取 𝑛 = 2, 即取三个插值点 𝑎,
𝑎 + 𝑏
2
, 𝑏, 有
∫
𝑏
𝑎
𝑓 (𝑥)𝑑𝑥 =
𝑏 − 𝑎
6
[ 𝑓 (𝑎) + 4 𝑓
𝑎 + 𝑏
2
+ 𝑓 (𝑏)] + 𝑅
2
其中
𝑅
2
=
∫
𝑏
𝑎
𝐸
2
(𝑥)𝑑𝑥 = −
1
90
𝑏 − 𝑎
2
5
𝑓
(4)
(𝜉)
4.3 数值积分的代数精度
考虑利用多项式函数衡量积分的精度. 对于某个特定的数值积分公式, 令 𝑓 = 𝑥
𝑘
, 则有
∫
𝑏
𝑎
𝑥
𝑘
𝑑𝑥 =
𝑛
Õ
𝑖=0
𝑎
𝑖
𝑥
𝑘
𝑖
+ 𝑅
𝑛
若有 𝑅
𝑛
= 0, 则说明该积分公式对于 𝑘 次多项式能够精确计算, 认为其代数精度至少为 𝑘 阶. 若是在
𝑘 = 𝑀 时 𝑅
𝑛
= 0,𝑘 = 𝑀 + 1 时 𝑅
𝑛
≠ 0, 则其代数精度为 𝑀 阶
5 Newton-Ctoes 公式
5.1 端点为积分节点的 Newton-Ctoes 公式
考虑将积分区间 [𝑎, 𝑏] 等分为 𝑛 份, 即
𝑥
𝑖
= 𝑎 + 𝑖ℎ, 𝑖 = 0, 1, · · · , 𝑛, ℎ =
𝑏 − 𝑎
𝑛
利用拉格朗日插值多项式构造积分
∫
𝑏
𝑎
𝑓 (𝑥)𝑑𝑥 =
𝑛
Õ
𝑖=0
𝑓 (𝑥
𝑖
)
∫
𝑏
𝑎
𝑙
𝑖
(𝑥)𝑑𝑥 + 𝑅
𝑛
设其有一般形式
∫
𝑏
𝑎
𝑓 (𝑥)𝑑𝑥 = (𝑏 − 𝑎)
𝑛
Õ
𝑖=0
𝐶
(𝑛)
𝑖
𝑓 (𝑥
𝑖
) + 𝐸
𝑛
(𝑥)
其中 𝐶
(𝑛)
𝑖
为 Newton-Ctoes 系数, 满足
𝐶
(𝑛)
𝑖
=
1
𝑏 − 𝑎
∫
𝑏
𝑎
𝑙
𝑖
(𝑥)𝑑𝑥
那么首先计算
𝛼
𝑖
=
∫
𝑏
𝑎
𝑙
𝑖
(𝑥)d𝑥
代入 𝑥 = 𝑎 + 𝑡ℎ, 则 𝑥
𝑖
= 𝑎 + 𝑖ℎ 计算得到
𝛼
𝑖
=
∫
𝑛
0
(𝑎 + 𝑡ℎ − 𝑥
0
)(𝑎 + 𝑡ℎ − 𝑥
1
) · · · (𝑎 + 𝑡ℎ − 𝑥
𝑖−1
)(𝑎 + 𝑡ℎ − 𝑥
𝑖+1
) · · · (𝑎 + 𝑡ℎ − 𝑥
𝑛
)
(𝑥
𝑖
− 𝑥
0
)(𝑥
𝑖
− 𝑥
1
) · · · (𝑥
𝑖
− 𝑥
𝑖−1
)(𝑥
𝑖
− 𝑥
𝑖+1
) · · · (𝑥
𝑖
− 𝑥
𝑛
)
ℎd𝑡
简化后得到
𝛼
𝑖
=
ℎ
𝑖!(𝑛 − 𝑖)!
(−1)
𝑛−𝑖
∫
𝑛
0
𝑡(𝑡 − 1) · · · (𝑡 − 𝑖 + 1)(𝑡 − 𝑖 − 1) · · · (𝑡 − 𝑛)d𝑡
因而
系数 𝐶
(𝑛)
𝑖
与积分区间 [𝑎, 𝑏] 的取值无关
𝐶
(𝑛)
𝑖
=
(−1)
𝑛−𝑖
𝑖!(𝑛 − 𝑖)!𝑛
∫
𝑛
0
𝑡(𝑡 − 1) · · · (𝑡 − 𝑖 + 1)(𝑡 − 𝑖 − 1) · · · (𝑡 − 𝑛)d𝑡
𝐸
𝑛
(𝑥) 为积分误差, 需要奇偶讨论. 当 𝑛 为奇数时
𝐸
𝑛
( 𝑓 ) =
𝑓
(𝑛+1)
(𝜂)
(𝑛 + 1)!
∫
𝑏
𝑎
(𝑥 − 𝑥
0
)(𝑥 − 𝑥
1
) · · · (𝑥 − 𝑥
𝑛
)d𝑥
此时积分具有 𝑛 阶精度; 当 𝑛 为偶数时
𝐸
𝑛
( 𝑓 ) =
𝑓
(𝑛+2)
(𝜂)
(𝑛 + 2)!
∫
𝑏
𝑎
𝑥(𝑥 − 𝑥
0
)(𝑥 − 𝑥
1
) · · · (𝑥 − 𝑥
𝑛
)d𝑥
此时积分具有 𝑛 + 1 阶精度
令 ℎ 为区间长度
ℎ =
𝑏 − 𝑎
𝑛
则 𝑛 = 1 时有梯形公式
∫
𝑏
𝑎
𝑓 (𝑥)d𝑥 =
ℎ
2
[
𝑓 (𝑎) + 𝑓 (𝑏)
]
−
ℎ
3
12
𝑓
′′
(𝜂)
𝑛 = 2 时有辛普森公式
∫
𝑏
𝑎
𝑓 (𝑥)d𝑥 =
ℎ
3
[
𝑓 (𝑎) + 4 𝑓 (𝑎 + ℎ) + 𝑓 (𝑏)
]
−
ℎ
5
90
𝑓
(4)
(𝜂)
𝑛 = 3 时有 Simpson 3/8 公式
∫
𝑏
𝑎
𝑓 (𝑥)d𝑥 =
3ℎ
8
[
𝑓 (𝑎) + 3 𝑓 (𝑎 + ℎ) + 3 𝑓 (𝑎 + 2ℎ) + 𝑓 (𝑏)
]
−
3ℎ
5
80
𝑓
(4)
(𝜂)
5.2 端点不为积分节点的 Newton-Ctoes 公式
考虑将积分区间 [𝑎, 𝑏] 等分为 𝑛 + 2 份, 即
𝑥
𝑖
= 𝑎 + 𝑖ℎ, 𝑖 = −1, 0, 1, · · · , 𝑛 + 1, ℎ =
𝑏 − 𝑎
𝑛 + 2
那么
𝑥
−1
= 𝑎, 𝑥
𝑛+1
= 𝑏
同样可以使用 𝑛 阶拉格朗日插值多项式构造积分
∫
𝑎
𝑏
𝑓 (𝑥)𝑑𝑥 =
𝑛
Õ
𝑖=0
𝑓 (𝑥
𝑖
)
∫
𝑏
𝑎
𝑙
𝑖
(𝑥)𝑑𝑥 + 𝑅
𝑛
需要注意的是, 虽然拉格朗日插值多项式形式上不变, 但是其中的插值点发生了变化, 因而需要重新计算
Newton-Ctoes 系数. 积分误差也需要重新计算. 当 𝑛 为偶数时
∫
𝑏
𝑎
𝑓 (𝑥)𝑑𝑥 =
𝑛
Õ
𝑗=0
𝑎
𝑗
𝑓 (𝑥
𝑗
) + ℎ
𝑛+3
𝑓
(𝑛+2)
(𝜉)
(𝑛 + 2)!
∫
𝑛+1
−1
𝑡
2
(𝑡 − 1) · · · (𝑡 − 𝑛)𝑑𝑡
当 𝑛 为奇数时
∫
𝑏
𝑎
𝑓 (𝑥)𝑑𝑥 =
𝑛
Õ
𝑗=0
𝑎
𝑗
𝑓 (𝑥
𝑗
) + ℎ
𝑛+2
𝑓
(𝑛+1)
(𝜉)
(𝑛 + 1)!
∫
𝑛+1
−1
𝑡(𝑡 − 1) · · · (𝑡 − 𝑛)𝑑𝑡
当 𝑛 = 0 时有中点公式
∫
𝑏
𝑎
𝑓 (𝑥)𝑑𝑥 = 2ℎ 𝑓 (𝑥
0
) +
ℎ
3
3
𝑓
′′
(𝜉); 𝜉 ∈ (𝑥
−1
, 𝑥
1
)
𝑛 = 1 有
∫
𝑏
𝑎
𝑓 (𝑥)𝑑𝑥 =
3ℎ
2
[ 𝑓 (𝑥
0
) + 𝑓 (𝑥
1
)] +
3ℎ
3
4
𝑓
′′
(𝜉); 𝜉 ∈ (𝑥
−1
, 𝑥
2
)
𝑛 = 2 有
∫
𝑏
𝑎
𝑓 (𝑥)𝑑𝑥 =
3ℎ
4
[2 𝑓 (𝑥
0
) − 𝑓 (𝑥
1
) + 2 𝑓 (𝑥
2
)] +
14ℎ
5
45
𝑓
(4)
(𝜉); 𝜉 ∈ (𝑥
−1
, 𝑥
3
)
5.3 复化积分公式
希望将大的积分区域划分为小的积分区域, 进而使用低阶的积分公式计算
5.3.1 复化 Simpson 公式
若 𝑓 ∈ 𝐶
4
[𝑎,𝑏]
, 𝑛 为偶数, ℎ =
𝑏−𝑎
𝑛
, 𝑥
𝑗
= 𝑎 + 𝑗 · ℎ, 𝑗 = 0, · · · , 𝑛; 则有复化 Simpson 公式:
∫
𝑏
𝑎
𝑓 (𝑥)𝑑𝑥 =
ℎ
3
𝑓 (𝑎) + 𝑓 (𝑏) + 2
𝑛
2
−1
Õ
𝑗=1
𝑓 (𝑥
2 𝑗
) + 4
𝑛
2
Õ
𝑗=1
𝑓 (𝑥
2 𝑗 −1
)
−
𝑏 − 𝑎
180
ℎ
4
𝑓
(4)
(𝜇),
其中 𝜇 ∈ (𝑎, 𝑏)
将区间 [𝑎, 𝑏] 划分为 𝑛 个小区间, 每个小区间的长度为 ℎ =
𝑏−𝑎
𝑛
. 对应的节点为 𝑥
𝑗
= 𝑎 + 𝑗 ℎ, 其中
𝑗 = 0, 1, 2, . . . , 𝑛
由于 𝑛 为偶数, 可以将整个区间分为
𝑛
2
个子区间, 每个子区间包含两个小区间例如, 第 𝑘 个子区间为
[𝑥
2𝑘
, 𝑥
2𝑘+2
], 其中 𝑘 = 0, 1, 2, . . . ,
𝑛
2
− 1
在每个子区间 [𝑥
2𝑘
, 𝑥
2𝑘+2
] 上应用 Simpson 公式：
∫
𝑥
2𝑘 +2
𝑥
2𝑘
𝑓 (𝑥)𝑑𝑥 =
ℎ
3
[
𝑓 (𝑥
2𝑘
) + 4 𝑓 (𝑥
2𝑘+1
) + 𝑓 (𝑥
2𝑘+2
)
]
−
ℎ
5
90
𝑓
(4)
(𝜉
𝑘
),
其中 𝜉
𝑘
∈ (𝑥
2𝑘
, 𝑥
2𝑘+2
)
将所有
𝑛
2
个子区间的 Simpson 公式相加：
∫
𝑏
𝑎
𝑓 (𝑥)𝑑𝑥 =
𝑛
2
−1
Õ
𝑘=0
∫
𝑥
2𝑘 +2
𝑥
2𝑘
𝑓 (𝑥)𝑑𝑥 =
ℎ
3
𝑛
2
−1
Õ
𝑘=0
(
𝑓 (𝑥
2𝑘
) + 4 𝑓 (𝑥
2𝑘+1
) + 𝑓 (𝑥
2𝑘+2
)
)
−
𝑛
2
−1
Õ
𝑘=0
ℎ
5
90
𝑓
(4)
(𝜉
𝑘
).
整理求和项：
𝑛
2
−1
Õ
𝑘=0
(
𝑓 (𝑥
2𝑘
) + 4 𝑓 (𝑥
2𝑘+1
) + 𝑓 (𝑥
2𝑘+2
)
)
= 𝑓 (𝑎) + 𝑓 (𝑏) + 2
𝑛
2
−1
Õ
𝑗=1
𝑓 (𝑥
2 𝑗
) + 4
𝑛
2
Õ
𝑗=1
𝑓 (𝑥
2 𝑗 −1
).
因此,
∫
𝑏
𝑎
𝑓 (𝑥)𝑑𝑥 =
ℎ
3
𝑓 (𝑎) + 𝑓 (𝑏) + 2
𝑛
2
−1
Õ
𝑗=1
𝑓 (𝑥
2 𝑗
) + 4
𝑛
2
Õ
𝑗=1
𝑓 (𝑥
2 𝑗 −1
)
−
ℎ
5
90
𝑛
2
−1
Õ
𝑘=0
𝑓
(4)
(𝜉
𝑘
).
由于 𝑓 ∈ 𝐶
4
[𝑎,𝑏]
, 且 𝜉
𝑘
∈ (𝑥
2𝑘
, 𝑥
2𝑘+2
) ⊂ (𝑎, 𝑏), 根据平均值定理, 存在 𝜇 ∈ (𝑎, 𝑏) 使得
𝑛
2
−1
Õ
𝑘=0
𝑓
(4)
(𝜉
𝑘
) =
𝑛
2
𝑓
(4)
(𝜇).
因此, 误差项可表示为：
𝑛
2
−1
Õ
𝑘=0
ℎ
5
90
𝑓
(4)
(𝜉
𝑘
) =
ℎ
5
90
·
𝑛
2
𝑓
(4)
(𝜇) =
ℎ
4
(𝑏 − 𝑎)
180
𝑓
(4)
(𝜇).
将误差项代入总公式, 得到：
∫
𝑏
𝑎
𝑓 (𝑥)𝑑𝑥 =
ℎ
3
𝑓 (𝑎) + 𝑓 (𝑏) + 2
𝑛
2
−1
Õ
𝑗=1
𝑓 (𝑥
2 𝑗
) + 4
𝑛
2
Õ
𝑗=1
𝑓 (𝑥
2 𝑗 −1
)
−
𝑏 − 𝑎
180
ℎ
4
𝑓
(4)
(𝜇),
其中 𝜇 ∈ (𝑎, 𝑏)
5.3.2 复化梯形公式
若 𝑓 ∈ 𝐶
2
[𝑎,𝑏]
, ℎ =
𝑏−𝑎
𝑛
, 𝑥
𝑗
= 𝑎 + 𝑗 · ℎ, 𝑗 = 0, · · · , 𝑛; 则有复化梯形公式:
∫
𝑏
𝑎
𝑓 (𝑥) 𝑑𝑥 =
ℎ
2
"
𝑓 (𝑎) + 𝑓 (𝑏) + 2
𝑛−1
Õ
𝑗=1
𝑓 (𝑥
𝑗
)
#
−
𝑏 − 𝑎
12
ℎ
2
𝑓
′′
(𝜇),
其中 𝜇 ∈ (𝑎, 𝑏)
5.3.3 复化中点公式
若 𝑓 ∈ 𝐶
2
[𝑎,𝑏]
, 𝑛 为偶数, ℎ =
𝑏−𝑎
𝑛+2
, 𝑥
𝑗
= 𝑎 + ( 𝑗 + 1) · ℎ, 𝑗 = −1, 0, · · · , 𝑛 + 1；则有复化中点公式:
∫
𝑏
𝑎
𝑓 (𝑥) 𝑑𝑥 = 2ℎ
𝑛
2
Õ
𝑗=1
𝑓 (𝑥
2 𝑗
) +
𝑏 − 𝑎
6
ℎ
2
𝑓
′′
(𝜇),
其中 𝜇 ∈ (𝑎, 𝑏)
5.4 自适应求解数值积分
若给定误差范围 𝜖, 则可以使用下面的方法计算数值积分
1. 任意选定 𝑛, 利用复化梯形公式计算积分值 𝐼
𝑛
2. 加密网格, 计算 𝐼
2𝑛
3. 计算
𝐼
2𝑛
−
𝐼
𝑛
.
考虑到
𝐼
2𝑛
的误差
𝐸
2𝑛
=
𝑏 − 𝑎
12
ℎ
2
𝑓
′′
(𝜇) =
(𝑏 − 𝑎)
3
48𝑛
2
𝑓
′′
(𝜇)
𝐸
𝑛
=
(𝑏 − 𝑎)
3
12𝑛
2
𝑓
′′
(𝜇)
则
𝐸
2𝑛
=
1
4
𝐸
𝑛
又有
𝐼
𝑛
+ 𝐸
𝑛
= 𝐼
2𝑛
+ 𝐸
2𝑛
⇒ 𝐸
2𝑛
− 𝐸
𝑛
= 𝐼
2𝑛
− 𝐼
𝑛
得到
𝐸
2𝑛
=
1
3
|
𝐼
2𝑛
− 𝐼
𝑛
|
因而只需要检验
1
3
|
𝐼
2𝑛
− 𝐼
𝑛
|
< 𝜖
如果满足则停止计算, 否则重复步骤 2
6 Gauss 求积公式
6.1 Gauss 求积公式的基本思想
Gaiss 求积公式希望通过选取适当的积分节点, 使得用这些点构造的积分公式代数精度最高. 即
𝐼 ( 𝑓 ) =
∫
𝑏
𝑎
𝑓 (𝑥)𝜔(𝑥)d𝑥 ≈
𝑛
Õ
𝑗=1
𝑐
𝑗
𝑓 (𝑥
𝑗
)
由” 插值型数值积分” 中的讨论,𝑐
𝑗
由 𝑥
𝑗
唯一确定
𝑐
𝑗
=
∫
𝑏
𝑎
𝑙
𝑗
(𝑥)𝜔(𝑥)d𝑥
其中 𝑙
𝑗
(𝑥) 为拉格朗日基函数
那么只需要选取适当的 𝑥
𝑗
即可. 事实上, 不论如何选取这 𝑛 个积分节点, 公式的代数精度都不会超过
2𝑛 − 1 阶, 可以如下简单地证明
只需要证明该公式对于一个 2𝑛 次积分非零多项式计算结果为零即可, 可以如下构造
𝑝
(
𝑥
)
=
(
𝑥
−
𝑥
1
)
2
(𝑥 − 𝑥
2
)
2
· · · (𝑥 − 𝑥
𝑛
)
2
其中 𝑥
1
, 𝑥
2
, · · · , 𝑥
𝑛
为选取的 𝑛 个积分节点, 那么很显然
𝑝(𝑥
1
) = 𝑝(𝑥
2
) = · · · = 𝑝(𝑥
𝑛
) = 0
公式计算的结果自然为零; 而它又是一个 2𝑛 次多项式, 因而该公式的代数精度至多为 2𝑛 − 1 阶
6.2 权函数与正交函数组
6.2.1 权函数
注意高斯积分公式期望的形式为
∫
𝑏
𝑎
𝑓 (𝑥)𝜔(𝑥)d𝑥
其中 𝜔(𝑥) 为权函数, 利用它可以定义 [𝑎, 𝑏] 区间上的内积
( 𝑓 , 𝑔) =
∫
𝑏
𝑎
𝑓 (𝑥)𝑔(𝑥)𝜔(𝑥)d𝑥
这样的权函数应当满足如下条件
1. 𝜔(𝑥) > 0
2. 𝜔(𝑥) 在 [𝑎, 𝑏] 上可积, 并且积分大于零
3. 对于幂函数, 它的积分应当是有界的
∫
𝑏
𝑎
𝜔(𝑥)𝑥
𝑘
d𝑥
< +∞
若给定了权函数 𝜔(𝑥), 对于对权函数 𝑝 次绝对可积 (带权 Lebesgue 可积函数) 的函数 𝑓 (𝑥), 即满足
∫
𝑏
𝑎
𝜔(𝑥)
|
𝑓 (𝑥)
|
𝑝
d𝑥
< +∞, 𝑝 ≥ 1
所有这样的 𝑓 (𝑥) 构成了一个实数域上的线性空间, 不妨记为
𝐿
𝑝
𝜔
[𝑎, 𝑏]
称为 p-阶带权 Lebesgue 空间. 特别地, 对于 𝜔(𝑥) = 1, 则有
𝐿
𝑝
[𝑎, 𝑏]
称为
经典
Lebesgue
空间
6.2.2 正交函数组
给定权函数 𝜔(𝑥), 考察函数 𝑓 (𝑥) 与自己的内积
( 𝑓 , 𝑓 ) =
∫
𝑏
𝑎
𝜔(𝑥)
|
𝑓 (𝑥)
|
2
d𝑥
它应当是可积的, 有限的. 那么 𝑓 (𝑥) 应当属于 𝐿
2
𝜔
[𝑎, 𝑏]. 这是一个线性空间, 假定找到了空间中 𝑛 + 1 个线
性无关的函数
𝑓
0
(𝑥), 𝑓
1
(𝑥), · · · , 𝑓
𝑛
(𝑥)
它们显然可以作为线性空间的一组基 (虽然并不完备). 可以利用它们构造一组正交的函数
𝑔
0
(𝑥) = 𝑓
0
(𝑥)
𝑔
1
(𝑥) = 𝑓
1
(𝑥) −
( 𝑓
1
, 𝑔
0
)
(𝑔
0
, 𝑔
0
)
𝑔
0
(𝑥)
𝑔
2
(𝑥) = 𝑓
2
(𝑥) −
( 𝑓
2
, 𝑔
0
)
(𝑔
0
, 𝑔
0
)
𝑔
0
(𝑥) −
( 𝑓
2
, 𝑔
1
)
(𝑔
1
, 𝑔
1
)
𝑔
1
(𝑥)
.
.
.
𝑔
𝑛
(𝑥) = 𝑓
𝑛
(𝑥) −
𝑛−1
Õ
𝑗=0
( 𝑓
𝑛
, 𝑔
𝑗
)
(𝑔
𝑗
, 𝑔
𝑗
)
𝑔
𝑗
(𝑥)
6.2.3 正交多项式组
假定上一节的 𝑛 + 1 个正交函数
𝑓
0
(𝑥), 𝑓
1
(𝑥), · · · , 𝑓
𝑛
(𝑥)
是不超过 𝑛 次的多项式, 将它们重新记为
𝑝
0
(𝑥), 𝑝
1
(𝑥), · · · , 𝑝
𝑛
(𝑥)
那么它们应当满足如下的性质
1. 它们可以构成 𝑃
𝑛
空间的一组基, 即任意不超过 𝑛 次的多项式都可以用它们线性表示
2. 𝑝
𝑛
(𝑥) 与任意不超过 𝑛 − 1 次的多项式正交
∫
𝑏
𝑎
𝑝
𝑛
(𝑥)𝑃( 𝑥)𝜔(𝑥)d𝑥 = 0
其中 𝑃(𝑥) ∈ 𝑃
𝑛−1
3. 𝑝
𝑛
(𝑥) 在 [𝑎, 𝑏] 上有 𝑛 个不同的单重零点
4. 𝑝
𝑛
(𝑥) 与 𝑃
𝑛−1
(𝑥) 在 [𝑎, 𝑏] 上的零点呈交叉分布
𝑎 < 𝑥
(𝑛)
1
< 𝑥
(𝑛−1)
1
< 𝑥
(𝑛)
2
< 𝑥
(𝑛−1)
2
< · · · < 𝑥
(𝑛)
𝑛
< 𝑥
(𝑛−1)
𝑛
< 𝑏
6.3 Gauss 求积公式的构造
设需要构造的 Gauss 求积公式的形式为
𝐼 ( 𝑓 ) =
∫
𝑏
𝑎
𝑓 (𝑥)𝜔(𝑥)d𝑥 ≈
𝑛
Õ
𝑗=1
𝑐
𝑗
𝑓 (𝑥
𝑗
)
只需要取权函数 𝜔(𝑥) 在 [𝑎, 𝑏] 上的 𝑛 次正交多项式 𝑝
𝑛
(𝑥) 的 𝑛 个零点为积分节点
𝑥
1
, 𝑥
2
, · · · , 𝑥
𝑛
构造的 Gauss 求积公式的代数精度可以达到 2𝑛 − 1 阶
为了证明它, 利用这 𝑛 个点构造 𝑛 − 1 次 Newton 插值多项式
𝑁
𝑛−1
(𝑥) = 𝑓 (𝑥) + 𝑅(𝑥)
Newton 插值多项式的插值误差为
𝑓 (𝑥) − 𝑁
𝑛−1
(𝑥) = 𝑓 [𝑥, 𝑥
1
, 𝑥
2
, · · · , 𝑥
𝑛
, 𝑥] (𝑥 − 𝑥
1
)(𝑥 − 𝑥
2
) · · · (𝑥 − 𝑥
𝑛
)
只需要证明 𝑓 (𝑥) 为不超过 2𝑛 − 1 次多项式时积分误差恒为零, 即
𝐼 ( 𝑓 ) − 𝐼
𝑛
( 𝑓 ) = 0
𝐼
𝑛
即插值多项式的积分, 代入有
𝐼 ( 𝑓 ) − 𝐼
𝑛
( 𝑓 ) =
∫
𝑏
𝑎
[
𝑓 (𝑥) − 𝑁
𝑛−1
(𝑥)
]
𝜔(𝑥)d𝑥
=
∫
𝑏
𝑎
𝑓 [𝑥, 𝑥
1
, 𝑥
2
, · · · , 𝑥
𝑛
, 𝑥] (𝑥 − 𝑥
1
)(𝑥 − 𝑥
2
) · · · (𝑥 − 𝑥
𝑛
)𝜔(𝑥)d𝑥
注意到 𝑥
1
, 𝑥
2
, · · · , 𝑥
𝑛
为 𝑝
𝑛
(𝑥) 的 𝑛 个零点, 因而
(𝑥 − 𝑥
1
)(𝑥 − 𝑥
2
) · · · (𝑥 − 𝑥
𝑛
) = 𝑐 𝑝
𝑛
(𝑥), 𝑐 = 𝑐𝑜𝑛𝑠𝑡
而又注意到由于 𝑓 (𝑥) 是一个不超过 2𝑛 − 1 次的多项式, 差商
𝑓 [𝑥, 𝑥
1
, 𝑥
2
, · · · , 𝑥
𝑛
, 𝑥]
就是一个不超过 𝑛 − 1 次的多项式 (做一次差商次数减一), 根据前面对 𝑝
𝑛
(𝑥) 性质的讨论, 它与 𝑝
𝑛
(𝑥) 正
交, 即
∫
𝑏
𝑎
𝑓 [𝑥, 𝑥
1
, 𝑥
2
, · · · , 𝑥
𝑛
, 𝑥] (𝑥 − 𝑥
1
)(𝑥 − 𝑥
2
) · · · (𝑥 − 𝑥
𝑛
)𝜔(𝑥)d𝑥 = 0
因而积分误差为零
,
即证
对于如下形式的积分
∫
𝑏
𝑎
𝑓 (𝑥)𝜔(𝑥)d𝑥
可以这样构造具有 𝑛 个积分节点的 Gauss 求积公式
∫
𝑏
𝑎
𝑓 (𝑥)𝜔(𝑥)d𝑥 ≈
𝑛
Õ
𝑗=1
𝑐
𝑗
𝑓 (𝑥
𝑗
)
1. 选取 𝑛 个线性无关的多项式函数, 不妨取 1, 𝑥, 𝑥
2
, · · · , 𝑥
𝑛−1
2. 根据权函数 𝜔(𝑥) 将其正交化, 得到 𝑛 个正交多项式 𝑝
0
(𝑥), 𝑝
1
(𝑥), · · · , 𝑝
𝑛−1
(𝑥)
3. 求解 𝑝
𝑛
(𝑥) 的 𝑛 个零点 𝑥
1
, 𝑥
2
, · · · , 𝑥
𝑛
4. 根据插值型积分确定系数 𝑐
𝑗
𝑐
𝑗
=
∫
𝑏
𝑎
𝑙
𝑗
(𝑥)𝜔(𝑥)d𝑥
6.4 Gauss-Legendre 求积公式
6.4.1 Legendre 多项式
Legendre 多项式是 Legendre 方程的解, 它有微分表达式
𝐿
𝑛
(𝑥) =
1
2
𝑛
𝑛!
d
𝑛
d𝑥
𝑛
(𝑥
2
− 1)
𝑛
Legendre 多项式在 (−1, 1) 上有 𝑛 个不同的零点
𝑥
(𝑛)
1
< 𝑥
(𝑛)
2
< · · · < 𝑥
(𝑛)
𝑛
由于 Legendre 多项式的奇偶性, 它们关于原点对称.Legendre 多项式两两正交
∫
1
−1
𝐿
𝑖
(𝑥)𝐿
𝑗
(𝑥)d𝑥 = 0, 𝑖 ≠ 𝑗
因而有
∫
1
−1
𝐿
𝑛
(𝑥)𝑃( 𝑥)d𝑥 = 0, 𝑃(𝑥) ∈ 𝑃
𝑛−1
𝑃
𝑛−1
是任意不超过 𝑛 − 1 次的多项式
6.4.2 Gauss-Legendre 求积公式
考虑特殊情况,权函数 𝜔(𝑥) = 1, 积分区域为 [−1, 1] , 则可以取勒让德多项式为正交多项式组
∫
1
−1
𝑓 (𝑥)d𝑥 ≈
𝑛
Õ
𝑗=1
𝑐
𝑗
𝐿
𝑗
(𝑥
𝑗
)
那么系数 𝑐
𝑗
为
𝑐
𝑗
=
∫
1
−1
𝑛
Ö
𝑘=1
𝑘≠ 𝑗
𝑥 − 𝑥
𝑘
𝑥
𝑗
− 𝑥
𝑘
d𝑥
对于一般的积分
∫
𝑏
𝑎
𝑓 (𝑥)d𝑥
可以通过变量代换将其转化为 [−1, 1] 上的积分
𝑥 =
𝑏 − 𝑎
2
𝑡 +
𝑎 + 𝑏
2
⇒ 𝑡 =
2𝑥 − 𝑎 − 𝑏
𝑏 − 𝑎
得
∫
𝑏
𝑎
𝑓 (𝑥)d𝑥 =
𝑏 − 𝑎
2
∫
1
−1
𝑓
𝑏 − 𝑎
2
𝑡 +
𝑎 + 𝑏
2
d𝑡
