数值积分  在[a, b]上取 a  x0 < x1 <…< xn  b，做 f 的 n 次插值多项式 ，即得到

Presentation on theme: "数值积分  在[a, b]上取 a  x0 < x1 <…< xn  b，做 f 的 n 次插值多项式 ，即得到"— Presentation transcript:

1 数值积分  在[a, b]上取 a  x0 < x1 <…< xn  b，做 f 的 n 次插值多项式 ，即得到
近似计算 插值型积分公式 8.1 插值型求积公式 思路 利用插值多项式 则积分易算。  在[a, b]上取 a  x0 < x1 <…< xn  b，做 f 的 n 次插值多项式 ，即得到 误差 Ak 由 决定， 与 无关。 节点 f (x)

2 8.2 复化求积公式 如果积分区间比较大，直接地使用上述求积公式， 精度难以保证。 高次插值有Runge 现象，故采用分段低次插值
8.2 复化求积公式 如果积分区间比较大，直接地使用上述求积公式， 精度难以保证。 高次插值有Runge 现象，故采用分段低次插值  分段低次合成的 Newton-Cotes 复合求积公式。 （1）等分求积区间，比如取步长 ，分[a, b]为n等分， 分点为 k = 0, 1, 2,…, n （2）在区间 [xk, xk+1]上使用以上求积公式求得Ik （3）取和值 ，作为整个区间上的积分近似值。

3  复化梯形公式： 在每个 上用梯形公式： = Tn /*积分中值定理*/

4  复化 Simpson 公式： 4 4 4 4 4 = Sn 注：为方便编程，可采用另一记法：令 n’ = 2n 为偶数， 这时 ，有

5 例 8．1：利用数据表 计算积分 这个问题有明显的答案 取n = 8用复化梯形公式 xk 1/8 1/4 3/8 1/2 5/8 3/4
1/8 1/4 3/8 1/2 5/8 3/4 7/8 1 f (xk) 4 2 计算积分 这个问题有明显的答案 取n = 8用复化梯形公式

6 取n=4,用辛卜生公式

7 8.3 变步长梯形方法 8.4 求积公式的误差 当 时，不考虑舍入误差，求积公式是精确成立的。 舍入误差： 取f (x)  1，则 若f (xk)的舍入误差小于 ，则

8 1．梯形公式的截断误差 2．辛卜生公式的截断误差

9 8.5 龙贝格求积公式 龙贝格积分法是在计算梯形和序列的基础上应用了 线性外推的加速方法，由此构成的一种具有超线性 收敛的自动积分法 方法思路 : 1.按照区间逐次分半的方法，计算梯形和序列

10

11 由此生成序列 T0, T1, …, Tn ,… 当 时，就可以结束计算。

12 设Tn为梯形和，I为积分真值，由复化梯形公式
2.加速 设Tn为梯形和，I为积分真值，由复化梯形公式 T f(x) Sn Tn I Tn+1 T o h2 h2

13 由解析几何 令h = 0， 则此直线在T 轴上的截距为 由 ，得：

14 用类似方法可推得： 柯特斯序列 龙贝格序列 由此法，可得如下三角形数表 梯形 辛卜生 柯特斯 龙贝格 T0 T3 T2 T1 S0 S2
C0 C1 D0

15 计算方法的实现： 首先构造T数表：

16 计算步骤： 1．取 ，计算 2．对k = 1, 2,… 计算下列各步 3．对n = 0, 1, 2,…, k = n – 1, n – 2, … 4．收敛控制 若 或 则输出积分值 ，否则转3。

17 8.6 高斯型求积公式 问题：在节点个数一定的情况下，是否可以在[a, b]上自由 选择节点的位置，使求积公式的精度提得更高 ？ 代数精确度：称： （8.9） 为一般求积公式。这里Ak为不依赖f (x)的常数若（8.9） 对任意不高于m次的多项式精确成立，而对于xm+1不能 精确成立，就说（8.9）式具有m次代数精确度。

18 例 8.2：求形如 的两点求积公式。 （1）用梯形公式（即以x0 = -1，x1 = 1为节点的插值型 求积公式）立即可得 一次代数精确度。

19 （2）若对求积公式中的四个待定系数A0, A1, x0, x1适当选取，
使求积公式对f (x) = 1，x，x2，x3都准确成立 y f(x) B A a o b x

20

21 数值积分 /* Numerical Integration */ /*interpolatory quadrature*/
近似计算 插值型积分公式 /*interpolatory quadrature*/ §1 Newton-Cotes 公式 思路 利用插值多项式 则积分易算。  在[a, b]上取 a  x0 < x1 <…< xn  b，做 f 的 n 次插值多项式 ，即得到 误差 Ak 由 决定， 与 无关。 节点 f (x)

22 定义 　　　若某个求积公式所对应的误差R[ f ]满足：R[ Pk ]=0 对任意 k  n 阶的多项式成立，且 R[ Pn+1 ]  0 对某个 n+1 阶多项式成立，则称此求积公式的代数精度为 n 。 例：对于[a, b]上1次插值，有 考察其代数精度。 f(x) a b f(b) 解：逐次检查公式是否精确成立 梯形公式 /* trapezoidal rule*/ f(a) 代入 P0 = 1： = 代入 P1 = x ： = 代数精度 = 1 代入 P2 = x2 ： 

23  注：形如 的求积公式至少有 n 次代数精度  该公式为插值型（即： ）  当节点等距分布时： 令
注：Cotes 系数仅取决于 n 和 i，可查表得到。与 f (x) 及区间[a, b]均无关。 Cotes系数 

24 I could only get arbitrarily
Excuses for not doing homework I could only get arbitrarily close to my textbook. I couldn't actually reach it. Trapezoidal Rule 代数精度 = 1 /* 令 x = a+th, h = ba, 用中值定理 */ n = 2: Simpson’s Rule n 为偶数阶的Newton-Cotes 公式至少有 n+1 次代数精度。 代数精度 = 3 n = 3: Simpson’s 3/8-Rule, 代数精度 = 3, n = 4: Cotes Rule, 代数精度 = 5,

25 §2 复合求积 /* Composite Quadrature */
高次插值有Runge 现象，故采用分段低次插值  分段低次合成的 Newton-Cotes 复合求积公式。  复合梯形公式： 在每个 上用梯形公式： Oh come on, you don’t seriously consider h=(ba)/2 acceptable, do you? Don’t you forget the oscillatory nature of high- degree polynomials! Haven’t we had enough formulae? What’s up now? Why can’t you simply refine the partition if you have to be so picky? Uh-oh = Tn /*中值定理*/

26  复化 Simpson 公式： 4 4 4 4 4 = Sn 注：为方便编程，可采用另一记法：令 n’ = 2n 为偶数， 这时 ，有

27 ~ 定义  收敛速度与误差估计： 若一个积分公式的误差满足 且C  0，则称该公式是 p 阶收敛的。 例：计算 其中 解：
运算量基本相同 其中 解： = 其中 =

28 Q: 给定精度 ，如何取 n ? 例如：要求 ，如何判断 n = ？ ? 上例中若要求 ，则 即：取 n = 409 通常采取将区间不断对分的方法，即取 n = 2k 可用来判断迭代 是否停止。 上例中2k  409  k = 9 时，T512 = S4 = 注意到区间再次对分时

29 §3 龙贝格积分 /* Romberg Integration */
例：计算 已知对于 = 106 须将区间对分 9 次，得到 T512 = 考察 由 来计算 I 效果是否好些？ Romberg 序列 = = S4 一般有：  T1 = ) ( T <  ?  Romberg 算法：  T2 = ) 1 ( T  S1 = ) ( 1 T <  ?  T4 = ) 2 ( T  S2 = ) 1 ( T  C1 = ) ( 2 T <  ?  T8 = ) 3 ( T  S4 = ) 2 ( 1 T  C2 = ) 1 ( 2 T  R1 = ) ( 3 T … … … … … …

30  理查德森外推法 /* Richardson’s extrapolation */
i 与 h 无关 利用低阶公式产生高精度的结果。 设对于某一 h  0，有公式 T0(h) 近似计算某一未知值 I。由Taylor展开得到： T0(h)  I = 1 h + 2 h2 + 3 h3 + … 现将 h 对分，得： ( ) ... 3 2 1 + = - h I T a Q：如何将公式精度由 O(h) 提高到 O(h2) ？ ... 4 3 2 1 ) ( - = h I T a 即：

31  §4 高斯型积分 /* Gaussian Quadrature */ 构造具有2n+1次代数精度的求积公式
将节点 x0 … xn 以及系数 A0 … An 都作为待定系数。令 f (x) = 1, x, x2, …, x2n+1 代入可求解，得到的公式具有2n+1 次代数精度。这样的节点称为Gauss 点，公式称为Gauss 型求积公式。 例：求 的 2 点 Gauss 公式。 解：设 ，应有 3 次代数精度。  +  1 ) ( x f A dx 代入 f (x) = 1, x, x2, x3 不是线性方程组，不易求解。

32   定理 x0 … xn 为 Gauss 点  与任意次数不大于n 的多项式 P(x) （带权）正交。
证明： “” 求 Gauss 点  求w(x) 对任意次数不大于n 的多项式 Pm(x)， Pm(x) w(x)的次数不大于2n+1，则代入公式应精确成立：  = 0 “” 要证明 x0 … xn 为 Gauss 点，即要证公式对任意次数不大于2n+1 的多项式 Pm(x) 精确成立，即证明： 设 

33     正交多项式族{ 0, 1, …, n, … }有性质：任意次数不大于n 的多项式 P(x) 必与n+1 正交。
若取 w(x) 为其中的n+1，则n+1的根就是 Gauss 点。 再解上例：  +  1 ) ( x f A dx Step 1：构造正交多项式2 设 c bx x a + = 2 1 ) ( , j 5 3 - = a ) ( 1 = +  dx a x ) , ( 1 = j   21 5 9 10 = - c b  = + -  1 2 ) )( 5 3 ( , dx c bx x j 即：

34 注：构造正交多项式也可以利用 L-S 拟合中介绍过的递推式进行。
Step 2：求2 = 0 的 2 个根，即为 Gauss 点 x0 ，x1 解线性方程组，简单。 Step 3：代入 f (x) = 1, x 以求解 A0 ，A1 结果与前一方法相同：  利用此公式计算 的值 注：构造正交多项式也可以利用 L-S 拟合中介绍过的递推式进行。

35 注意到积分端点 1 可能是积分的奇点，用普通Newton-Cotes公式在端点会出问题。而Gauss公式可能避免此问题的发生。
 特殊正交多项式族： ① Legendre 多项式族： 1 ) (  x r 定义在[1, 1]上， 满足： 注意到积分端点 1 可能是积分的奇点，用普通Newton-Cotes公式在端点会出问题。而Gauss公式可能避免此问题的发生。 由 有递推 以 Pn+1 的根为节点的求积公式称为Gauss-Legendre 公式。 2 1 ) ( x - = r 定义在[1, 1]上， ② Chebyshev 多项式族： Tn+1 的根为 k = 0, …, n 以此为节点构造公式 称为 Gauss-Chebyshev 公式。

36 Q：什么样的插值多项式在 x0 … xn 上有 2n+1 阶？ 插值多项式的余项 A：Hermite 多项式！ 满足
 Gauss 公式的余项： /* 设P为f 的过x0 … xn的插值多项式 */ /*只要P 的阶数不大于2n+1，则下一步等式成立*/ Q：什么样的插值多项式在 x0 … xn 上有 2n+1 阶？ 插值多项式的余项 A：Hermite 多项式！ 满足