Transcript 第三章

第3章
函数逼近与快速傅里叶变换
1
3.1
函数逼近的基本概念
3.2
正交多项式
3.3
最佳平方逼近
3.4
曲线拟合的最小二乘法
3.5
有理逼近
3.6
三角多项式与快速傅里叶变换
3.1
3.1.1
问题
函数逼近的基本概念
函数逼近与函数空间
1、数值计算中经常要计算函数值,如计算机中计算
基本初等函数及其他特殊函数;
2、当函数只在有限点集上给定函数值,要在包含该
点集的区间上用公式给出函数的简单表达式.
这些都涉及到在区间 [a, b] 上用简单函数逼近已知复杂
函数的问题,这就是函数逼近问题.
2
插值法就是函数逼近问题的一种.
本章讨论的函数逼近,是指“对函数类 A 中给定的函数
f (x), 记作 f ( x)  A , 要在另一类简单的便于计算的函数类
B 中求函数 p( x)  B ,使 p (x) 与 f (x) 的误差在某种度量
意义下最小”.
函数类 A通常是区间 [a, b] 上的连续函数,记作 C[a, b] ,
称为连续函数空间.
3
函数类 B通常为 n 次多项式,有理函数或分段低次多项
式等.
数学上常把在各种集合中引入某些不同的确定关系称为
赋予集合以某种空间结构,并将这样的集合称为空间.
例如将所有实 n维向量组成的集合,按向量加法及向量
与数的乘法构成实数域上的线性空间, 记作 R n, 称为 n 维
向量空间.
4
对次数不超过 n( n为正整数)的实系数多项式全体,
按通常多项式与多项式加法及数与多项式乘法也构成数域
R上一个线性空间,用 H n表示,称为多项式空间.
所有定义在 [a, b] 上的连续函数集合,按函数加法和
数与函数乘法构成数域 R 上的线性空间,记作 C[a, b].
类似地, 记 C p [a, b]为具有 p 阶连续导数的函数空间.
5
定义1
设集合 S 是数域 P上的线性空间,元素
x1 ,, xn  S , 如果存在不全为零的数 1 ,,  n  P ,
使得
1 x1     n xn  0,
(1.1)
则称 x1 ,, xn 线性相关.
否则,若等式(1.1)只对 1   2     n  0成立,
则称 x1 , , xn线性无关.
6
若线性空间 S 是由 n 个线性无关元素 x1 ,, xn 生成的,
x  S 都有
x  1 x1     n xn
则 x1 ,, xn 称为空间 S 的一组基,记为
S  span{x1 ,, xn }
并称空间 S 为 n 维空间,系数 1 ,,  n 称为 x 在基
x1 ,, xn下的坐标, 记作 (1 ,,  n ).
如果 S 中有无限个线性无关元素 x1 , , xn ,, 则称 S
为无限维线性空间.
7
考察次数不超过 n次的多项式集合 H n, 其元素
p( x)  H n 表示为
p( x)  a0  a1 x    an x n ,
(1.2)
它由 n  1个系数 (a0 , a1 ,, an ) 唯一确定.
1, x, , x n是线性无关的,它是 H n 的一组基,故
H n  span{1, x,, x n },
且 (a0 , a1 ,, an ) 是 p (x) 的坐标向量,H n 是 n  1维的.
8
对连续函数 f ( x)  C[a, b],它不能用有限个线性无关的
函数表示,故 C[a, b] 是无限维的,但它的任一元素 f (x)
均可用有限维的 p( x)  H n 逼近,使误差
max f ( x)  p( x)  
a  x b
( 为任给的小正数),这就是著名的魏尔斯特拉斯定理.
9
定理1
设 f ( x)  C[a, b] , 则对任何   0 ,总存在一
个代数多项式 p (x) , 使
f ( x)  p ( x )


在 [a, b] 上一致成立.
伯恩斯坦1912年给出的证明是一种构造性证明.
他根据函数整体逼近的特性构造出伯恩斯坦多项式
n
k
Bn ( f , x)   f ( ) Pk ( x)
n
k 0
(1.3)
10
其中
n k
Pk ( x)    x (1  x) n  k ,
k 
 n  n(n  1)  (n  k  1)
  
k!
k 
为二项式展开系数,并证明了
lim Bn ( f , x)  f ( x) 在 [0,1] 上一致成立;
n 
若 f (x) 在 [0,1]上 m阶导数连续,则
lim Bn( m ) ( f , x)  f
n 
( m)
( x).
这个结果不但证明了定理1,而且由(1.3)给出了 f (x)
的一个逼近多项式,但它收敛太慢,实际中很少使用.
11
更一般地,可用一组在 C[a, b] 上线性无关的函数集合
i ( x) in0 来逼近 f ( x)  C[a, b],此时元素
 ( x)    span{0 ( x), 1 ( x), , n ( x)}  C[a, b]
可表示为
 ( x)  a00 ( x)  a11 ( x)    ann ( x).
(1.4)
函数逼近问题就是对任何 f ( x)  C[a, b],在子空间Φ
*
*
找一个元素  ( x)  , 使 f ( x)   ( x)在某种意义下最小.
12
3.1.2
范数与赋范线性空间
为了对线性空间中元素大小进行衡量,需要引进范数
定义,它是 R n空间中向量长度概念的直接推广.
13
定义2
设 S为线性空间,x  S ,若存在唯一实数‖·‖,
满足条件:
(1)
x  0, 当且仅当 x  0 时, x  0; (正定性)
(2)
x   x ,
(3)
x y  x  y ,
  R;
(齐次性)
x, y  S . (三角不等式)
则称‖·‖为线性空间 S上的范数,S 与‖·‖一起称为赋范
线性空间,记为 X .
14
例如,在 R n上的向量 x  ( x1 , , xn )T  R n , 三种常
用范数为
x
x

1
 max xi ,
1i  n

n

i 1
2
称为 1-范数,
xi ,
n
x
称为  范数或最大范数,
1
2
 ( xi2 ) ,
称为 2-范数.
i 1
15
实际上任何向量的实值函数,只要满足上述三个条件,
就可以定义成一种向量范数.
在 R 中,满足‖·‖2 =1 ,即 x1  x2
2
2
2
 1的向量
为单位圆,
满足‖·‖∞ =1 ,即 max{ x1 , x2 }  1 的向量为单位正
方形,
而满足‖·‖1 =1 的向量 x1  x2  1 则为对角线长
度为1的菱形.
16
所以说,范数是对向量长度的度量,度量方式不同,
结果也不一样,但不同范数之间是存在等价关系的.
17
类似地,对连续函数空间 C[a, b] ,若 f ( x)  C[a, b] ,
f
f
f

1
a  x b


b
称为 1-范数,
f ( x) dx,
a
b
2
称为  范数,
 max f ( x) ,
1
2
 (  f 2 ( x) dx) ,
称为 2-范数.
a
可以验证这样定义的范数均满足定义2中的三个条件.
18
3.1.3
内积与内积空间
在线性代数中, R n中两个向量 x  ( x1 , , xn )T 及
y  ( y1 , , yn )T 的内积定义为
( x, y)  x1 y1    xn yn .
(1.5)
若将它推广到一般的线性空间 X,则有下面的定义.
19
定义3 X 是数域K(R或C)上的线性空间,对
u, v  X, 有K中一个数与之对应,记为 (u, v) ,它满足
以下条件:
(1)
(u, v)  (v, u ),
u, v  X;
(2)
(u, v)   (u, v),
(3)
(u  v, w)  (u, w)  (v, w),
(4)
(u, u)  0,
  K, u, v  X;
u, v, w  X;
当且仅当u  0 时,
(u, u)  0.
则称 (u, v)为X上 u 与 v 的内积.
20
定义了内积的线性空间称为内积空间.
定义中(1)的右端 (u , v) 称为 (u, v) 的共轭,
当K为实数域R时 (u, v)  (v, u ).
如果 (u, v)  0 ,则称 u与 v 正交,这是向量相互垂
直概念的推广.
21
定理2 设X为一个内积空间, 对 u , v  X , 有
(u, v)  (u, u )(v, v).
2
(1.6)
称为柯西-施瓦茨(Cauchy-Schwarz)不等式.
证明
当 v  0 时(1.6)式显然成立.
现设 v  0 ,则 (v, v)  0, 且对任何数  有
0  (u  v, u  v)  (u, u )  2 (u, v)  2 (v, v).
取   (u, v) /( v, v) ,代入上式右端,得
22
(u, u )  2
(u , v)
2
(v, v)

(u, v)
2
(v, v)
0
即得 v  0时
(u, v)  (u, u )(v, v).
2
23
定理3
设X为一个内积空间,u1 , u2 , , un  X,
矩阵
 (u1 , u1 )
(u , u )
G 1 2



(u1 , un )
(u2 , u1 )

(u2 , u2 )


(u 2 , un )

(u n , u1 ) 
(u n , u2 ) 



(u n , u n )
(1.7)
称为格拉姆(Gram)
则 G 非奇异的充分必要条件是 u1 , u2 ,, un 线性无关.
24
G非奇异等价于 det G  0,其充要条件是齐次
证明
方程组
n
n
j 1
j 1
(  j u j , uk )   (u j , uk ) j  0, k  1,2,, n(1.8)
只有零解; 而
n

j 1
j
u j  1u1   2u2     nun  0
n
n
j 1
j 0
(1.9)
 (  j u j ,   j u j )  0
n
 (  j u j , uk )  0,
j 0
k  1,2,, n
25
从以上等价关系知,det G  0 等价于从(1.8)推出
1   2     n  0,
n
n
而后者等价于从(1.9)推出
   n 
0,
(  j u j , uk )   (u j , uk ) j 
1 0, k2 
(1.8)
1,2,, n
j 1
j 1
即 u1 , u2 ,, un 线性无关.
n
在内积空间X
 u  u

j 1
j
j
1 1
  2u 2     n u n  0
u  X, 记
u 
(u, u )
(1.10)
26
利用
(u  v)  u
2
2
2 u v  v
2
 (u, u )  2(u, v)  (v, v)
 (u  v, u  v)  u  v
2
两端开方即得三角不等式
uv  u  v
(1.11)
27
例1
设
R n 与 C n的内积.
T
x, y  R n , x  ( x1 , , xn )T , y  ( y1 ,, yn ) ,
n
( x, y )   xi yi
(1.12)
i 1
向量2-范数为
x
1
2
2
n
 ( x, x)  ( xi2 )
1
2
i 1
28
若给定实数 i  0(i  1,2,, n), 称{i }为权系数,
R n 上的加权内积为
n
( x, y )   i xi yi
(1.13)
i 1
相应的范数为
n
x
2
1
2
 ( i xi2 ) .
i 1
当 i  1(i  1,2,, n) 时,
(1.13)就是前面定义的内积.
29
如果 x, y  C n,带权内积定义为
n
( x, y )   i xi yi
(1.14)
i 1
这里 {i } 仍为正实数序列,yi 为 yi 的共轭.
在 C[a, b]上也可以类似定义带权内积.
30
定义4
设 [a, b] 是有限或无限区间,在 [a, b] 上的非负
函数  (x)满足条件:
(1) a x k  ( x)dx 存在且为有限值 (k  0,1, ),
b
(2) 对 [a, b] 上的非负连续函数 g (x) ,如果

b
a
g ( x)  ( x)dx  0,
g ( x)  0.
则称  (x) 为[a, b]上的一个权函数.
31
例2 C[a, b] 上的内积.
设 f ( x), g ( x)  C[a, b],  (x) 是 [a, b] 上给定的权函数,
则可定义内积
b
( f ( x), g ( x))    ( x) f ( x) g ( x)dx.
(1.15)
a
由此内积导出的范数为
f ( x)
1
2
2
1
2
b

 ( f ( x), f ( x))    ( x) f 2 ( x)dx  (1.16)
 a

称(1.15)和(1.16)为带权  (x) 的内积和范数.
32
常用的是  ( x)  1 的情形,即
b
( f ( x), g ( x))   f ( x) g ( x)dx.
a
f ( x)
2
b

  f 2 ( x)dx 
 a

1
2
33
若 0 , 1 ,,  n是 C[a, b]中的线性无关函数族,记
  span{0 , 1 ,,  n }, 它的格拉姆矩阵为
G  G(0 , 1 , ,  n )
( 0 ,  0 )
 ( ,  )
 1 0



( n ,  0 )
( 0 , 1 )

(1 , 1 )


( n , 1 )

( 0 ,  n ) 
(1 ,  n ) 




( n ,  n ) 
(1.17)
根据定理3可知 0 , 1 ,,  n 线性无关的充要条件是
det G(0 , 1 , ,  n )  0.
34
3.1.4
最佳逼近
函数逼近主要讨论给定 f  C[a, b],求它的最佳逼近
多项式的问题.
若 P* ( x)  H n使误差
f ( x)  P* ( x)  min f ( x)  P( x)
PH n
则称 P* ( x) 是 f (x)在 [ a, b]上的最佳逼近多项式.
若 P( x)    span{0 , 1 ,, n }则称相应的 P* ( x)
为最佳逼近函数.
通常将范数  取为 

或  2.
35
若取  ,即
f ( x)  P * ( x)

 min f ( x)  P( x)
PH n

(1.18)
 min max f ( x)  P( x) ,
PH n a  x b
则称 P* ( x) 是 f (x)在 [ a, b]上的最优一致逼近多项式.
求 P* ( x)就是求 [ a, b]上使最大误差 max f ( x)  P( x)
a  x b
最小的多项式.
36
若取  ,即
2
f ( x)  P * ( x)
2
2
 min f ( x)  P( x)
PH n
 min
PH n

b
a
2
2
(1.19)
[ f ( x)  P( x)]2 dx,
则称 P* ( x) 是 f (x)在 [ a, b]上的最佳平方逼近多项式.
若 f (x)是 [ a, b]上的一个列表函数,在
a  x0  x1    xm  b
上给出 f ( xi )(i  0,1,, m) ,要求 P*   使
f  P * 2  min f  P
P
2
 min
P
m
2
[
f
(
x
)

P
(
x
)]

i
i
(1.20)
i 0
则称 P* ( x) 为 f (x)的最小二乘拟合.
37
3.2
正交多项式
3.2.1
定义5
正交函数族与正交多项式
若 f ( x), g ( x)  C[a, b],  ( x )为 [ a, b]
上的权函数且满足
( f ( x), g ( x)) 
则称 f

b
a
 ( x) f ( x) g ( x)dx  0.
(2.1)
( x )与 g ( x )在 [ a, b]上带权  ( x )正交.
38
若函数族 0 ( x), 1 ( x), , n ( x),  满足关系
( j ,  k ) 

b
a
j  k.
 0,
 ( x) j ( x) k ( x)dx 
j  k.
 Ak  0,
(2.2)
则称{ k ( x)}是 [a, b] 上带权  (x)的正交函数族.
若 Ak  1,则称之为标准正交函数族.
39
三角函数族
1, cos x, sin x, cos 2 x, sin 2 x,
就是在区间 [ ,  ] 上的正交函数族.
定义6
设  n ( x)是 [a, b] 上首项系数 an  0 的 n次多
项式, (x) 为 [ a, b] 上权函数,如果多项式序列 { n ( x)}
0
满足关系式(2.2),则称多项式序列 { n ( x)}0 为在 [a, b] 上
带权  (x) 正交,称  n ( x) 为 [a, b] 上带权  (x) 的 n次正
交多项式.
40
只要给定区间 [a, b] 及权函数  (x) ,均可由一族线性
无关的幂函数 {1, x,, x n ,}, 利用逐个正交化手续构造
出正交多项式序列 { n ( x)}0 :
0 ( x)  1,
n 1
 n ( x)  x n  
j 0
( x n ,  j ( x))
 j ( x)
( j ( x),  j ( x))
(n  1,2, ).
(2.3)
41
正交多项式  n ( x)的最高次项系数为1. 而若 { n ( x)}0
是正交多项式,则 0 ( x), 1 ( x), , n ( x) 在 [ a, b]上是
线性无关的.
事实上,若
c00 ( x)  c11 ( x)    cn n ( x)  0,
用  ( x) j ( x)( j  0,1, , n)乘上式并积分得
b
b
c0   ( x ) 0 ( x ) j ( x )dx  c1   ( x )1 ( x ) j ( x )dx  
a
a
b
 c j   ( x ) j ( x ) j ( x )dx  
a
b
 cn   ( x ) n ( x ) j ( x )dx  0,
a
42
利用正交性有
b
c j   ( x) j ( x) j ( x)dx  0.
a
b
由于 a  ( x) j ( x) j ( x)dx  0,故 c j  0( j  0,1, , n)
即0 ( x), 1 ( x), , n ( x) 线性无关.
关于正交多项式,有
(1) 任何 P( x)  H n 均可表示为 0 ( x), 1 ( x), , n ( x)
的线性组合. 即
n
P( x)   c j j ( x).
j 0
43
(2)  n ( x)与任一次数小于 n的多项式 P( x)  H n 1
正交. 即
( n ( x), P) 

b
a
 ( x) n ( x) P( x)dx  0.
除此之外,还有
44
定理4
设 { n ( x)}0 是 [ a, b]上带权  ( x )的正交多项式,
对 n  0成立关系
n1 ( x)  ( x   n )n ( x)   nn1 ( x)
(n  0,1, ).
(2.4)
其中
0 ( x)  1,
1 ( x)  0,
 n  ( x n ( x),  n ( x)) /( n ( x), n ( x)),
 n  ( n ( x),  n ( x)) /( n 1 ( x),  n 1 ( x)),
(n  1,2, ).
b
2
这里 ( x n ( x), n ( x))  a xn ( x)  ( x)dx.
45
设 { n ( x)}0 是 [ a, b]上带权  ( x )的正交多项式,
则  n ( x)( n  1)在区间 (a, b)内有n 个不同的零点.
定理5
证明
假定  n (x) 在(a, b) 内的零点都是偶数重的,则
 n (x)在 [ a, b]符号保持不变,这与
( n ( x), 0 ( x)) 

b
a
 ( x) n ( x)0 ( x)dx  0.
矛盾.故  n (x) 在 (a, b)内的零点不可能全是偶数重的,现设
xi (i  1,2,, l )为 n (x)在 (a, b)内的奇数重零点,不妨设
a  x1  x2    xl  b,
则  n (x) 在 xi (i  1,2,, l ) 处变号.
46
令
q( x)  ( x  x1 )( x  x2 ) ( x  xl ),
于是 n ( x)q( x)在 [ a, b]上不变号,则得
( n , q) 

b
a
 ( x) n ( x)q( x)dx  0.
若 l  n,由 { n ( x)}0 的正交性可知
( n , q) 

b
a
 ( x) n ( x)q( x)dx  0,
这与 ( n , q)  0矛盾,故 l  n.而  n (x)只有 n个零点,
故 l  n ,即 n个零点都是单重的.
47
3.2.2
勒让德多项式
当区间为[1, 1] ,权函数  ( x)  1时,由 {1, x,, x n ,}
正交化得到的多项式就称为勒让德(Legendre)多项式,
并用 P0 ( x), P1 ( x), , Pn ( x),  表示.
罗德利克(Rodrigul)给出了简单的表达式
P0 ( x)  1,
1 dn
2
n
Pn ( x)  n
{(
x

1
)
}
n
2 n! dx
(n  1,2, ),
(2.5)
48
由于 ( x 2  1) n是 2n 次多项式,所以对其求 n 阶导数后得
Pn ( x) 
1
n
n 1
(
2
n
)(
2
n

1
)

(
n

1
)
x

a
x
   a0 ,
n 1
n
2 n!
于是得首项 x 的系数 an 
n
( 2n)!
.
n
2
2 ( n!)
最高项系数为1的勒让德多项式为
~
Pn ( x) 
n! d n
2
n
[(
x

1
)
].
n
(2n)! dx
(2.6)
49
勒让德多项式重要性质:
性质1
正交性
0,


1 Pn ( x) Pm ( x)dx   2 ,
 2n  1
1
证明
m  n;
m  n.
(2.7)
令  ( x)  ( x 2  1) n,则
 ( k ) (1)  0
(k  0,1, , n  1).
设 Q(x)是在区间 [1, 1] 上 n 阶连续可微的函数,由分部
积分知
50
1

1
Pn ( x)Q( x)dx 
1
1
(n)
Q
(
x
)

( x)dx
n


1
2 n!
1
1
  n  Q( x) ( n 1) ( x)dx
2 n! 1

( 1) n
 n
2 n!
1

1
Q ( n ) ( x ) ( x ) dx
下面分两种情况讨论:
(1) 若 Q(x) 是次数小于 n的多项式,则 Q ( n ) ( x)  0,
故得
1

1
Pn ( x) Pm ( x)dx  0, 当n  m.
51
(2) 若
则
于是
Q( x)  Pn ( x) 
(2n)! n
1
(n)

x  ,

(
x
)
n
2
n
2 (n!)
2 n!
Q ( n ) ( x)  Pn( n ) ( x) 
( 2n)!
,
n
2 n!
n
(

1
)
(2n)! 1
2
2
n
P
(
x
)
dx

(
x

1
)
dx
n
2
2
1


1
2 (n!)
1
(2n)!
 2n
2 (n!) 2
1

1
(1  x 2 ) n dx.
由于
1

0
故

(1  x ) dx   cos 2 n 1 tdt
2
n
2
0
1

1
Pn2 ( x) dx 
2
.
2n  1
52
性质2
奇偶性
Pn ( x)  (1) n Pn ( x).
(2.8)
由于  ( x)  ( x 2  1) n 是偶次多项式,经过偶次求导仍为
偶次多项式,经过奇次求导则为奇次多项式,故 n为偶数时
Pn ( x )为偶函数,n为奇数时 Pn (x)为奇函数,于是(2.8)成
立.
53
性质3
递推关系
考虑 n  1次多项式 xPn ( x), 它可表示为
xPn ( x)  a0 P0 ( x)  a1 P1 ( x)    an 1 Pn 1 ( x).
两边乘 Pk ( x), 并从-1到1积分,得
1

1
1
xPn ( x) Pk ( x)dx  ak  Pk2 ( x)dx.
1
当 k  n  2 时, xPk ( x)次数小于等于 n  1,上式左端积分
为0, 故得 ak  0.
54
当 k  n 时,xPn2 ( x)为奇函数,左端积分仍为0, 故
an  0.
于是
其中
xPn ( x)  an 1 Pn 1 ( x)  an 1 Pn 1 ( x)
an 1
an 1
2n  1 1

xPn ( x ) Pn 1 ( x ) dx


1
2
2n  1
2n
n



,
2
2
4n  1
2n  1
2n  3 1

xPn ( x ) Pn 1 ( x) dx


1
2

2n  3
2( n  1)
n 1


,
2
( 2n  1)( 2n  3)
2n  1
55
从而得到以下的递推公式
(n  1) Pn 1 ( x)  (2n  1) xPn ( x)  nPn 1 ( x)
(2.9)
(n  1,2, ),
由 P0 ( x)  1, P1 ( x)  x, 利用上述递推公式就可推出
P2 ( x)  (3x 2  1) / 2,
P3 ( x)  (5 x 3  3 x) / 2,
P4 ( x)  (35 x 4  30 x 2  3) / 8,
P5 ( x)  (63 x 5  70 x 3  15 x) / 8,
P6 ( x)  (231x 6  315 x 4  105 x 2  5) / 16,

56
图3-1给出了P0 ( x), P1 ( x), P2 ( x), P3 ( x) 的图形.
图3-1
57
性质4
Pn ( x )在区间 [1, 1] 内有 n 个不同的实零点.
58
3.2.3
切比雪夫多项式
当权函数  ( x) 
1
1 x2
,区间为 [1, 1]时,由序
列 {1, x,, x n ,} 正交化得到的正交多项式就是切比雪夫
(Chebyshev)多项式.
它可表示为
Tn ( x)  cos( n arccos x),
x  1.
(2.10)
若令 x  cos , 则 Tn ( x)  cos n , 0     .
59
切比雪夫多项式有很多重要性质:
性质1
递推关系
Tn 1 ( x)  2 xTn ( x)  Tn 1 ( x)
T0 ( x)  1,
(n  1,2, ),
T1 ( x)  x.
(2.11)
这只要在三角恒等式
cos( n  1)  2 cos  cos n  cos( n  1)
(n  1)
中,令 x  cos 即得.
60
由(2.11)可推出
T2 ( x)  2 x 2  1,
3
(
x
)

4
x

3 x,
TT
(
x
)

2
xT
3
n 1
n ( x)  Tn 1 ( x)
T4 ( x)  8 x 4  8 x 2  1,
T5 ( x)  16 x 5  20 x 3  5 x,
T6 ( x)  32 x 6  48 x 4  18 x 2  1,

T0 ( x), T1 ( x), T2 ( x), T3 ( x) 的函数图形见图3-2.
61
图3-2
由递推关系(2.11)还可得到 Tn (x) 的最高次项系数是
2 n 1.
Tn 1 ( x)  2 xTn ( x)  Tn 1 ( x)
62
性质2
切比雪夫多项式 {Tk ( x)} 在区间 [1, 1]上带权
 ( x)  1 / 1  x 2 正交,且
 0,
1 T ( x )T ( x ) dx

n
m
 ,
1
2
1 x
2
 ,
n  m;
n  m  0;
(2.12)
n  m  0.
令 x  cos ,则 dx   sin d ,于是
 0, n  m;
1 T ( x )T ( x ) dx


n
m

  cos n cos md
 , n  m  0;
1
2
0
1 x
2
  , n  m  0.
63
性质3
T2 k ( x) 只含 x 的偶次幂,T2 k 1 ( x)只含 x 的奇次幂.
这个性质由递推关系直接得到.
性质4
Tn (x) 在区间 [1, 1]上有 n 个零点
2k  1
xk  cos
,
2n
性质5
k  1,2,, n.
Tn (x) 的首项 x n 的系数为 2 n 1.
~
若令 T~0 ( x)  1, T~n ( x)  1
,则
T
(
x
),
n

1
,
2
,

Tn ( x)
n 1 n
2
是首项系数为1的切比雪夫多项式.
64
~ 为所有次数小于等于 的首项系数为1的多项
若记 H
n
n
式集合,对 T~n ( x)有以下性质:
定理6
设 T~n ( x)是首项系数为1的切比雪夫多项式,则
~
max Tn ( x)  max P( x) ,
1 x 1
1 x 1
~
P( x)  H n .
且
~
max Tn ( x) 
1 x 1
1
.
n 1
2
65
~ 中,
定理6表明在所有首项系数为1的 n次多项式集合 H
n
~
~
~ 中最大值最小的多项
所以
是
Tn ( x)  min
P
(
x
)
,
T
(
x
)
H
~
n
n


PH
n
式,即
~
max Tn ( x)  min
~ max P ( x ) 
1 x 1
PH n 1 x 1
1
.
n 1
2
(2.13)
利用这一结论,可求 P( x)  H n 在 H n 1中的最佳(一致)
逼近多项式.
66
例3 求 f ( x)  2 x3  x 2  2 x  1 在 [1, 1]上的最佳2次逼
近多项式.
解 由题意,所求最佳逼近多项式 P2* ( x) 应满足
max f ( x)  P2* ( x)  min .
1 x 1
由定理6可知, 当
1
3
3
f ( x)  P ( x)  T3 ( x)  2 x  x
2
2
*
2
*
时,多项式 f ( x)  P2 ( x)
故
67
1
7
2
P ( x)  f ( x)  T3 ( x)  x  x  1
2
2
*
2
就是 f (x) 在 [1, 1]上的最佳2次逼近多项式.
由于切比雪夫多项式是在区间 [1, 1]上定义的,对于
一般区间 [ a, b],要通过变量替换变换到 [1, 1] ,可令
1
x  [(b  a) t  a  b],
2
则可将 x  [a, b] 变换到 t  [1, 1].
(2.14)
68
3.2.4
切比雪夫多项式零点插值
切比雪夫多项式 Tn (x)在区间 [1, 1]上有 n 个零点
xk  cos
2k  1
,
2n
k  1,2,, n.
和 n  1个极值点(包括端点)
k
xk  cos
,
n
k  0,1, , n.
这两组点称为切比雪夫点,它们在插值中有重要作用.
69
从图3-3可以看到切比雪夫点恰好是单位圆周上等距分
布点的横坐标,这些点的横坐标在接近区间 [1, 1]的端点
处是密集的.
利用切比雪夫点做插值,可
使插值区间最大误差最小化.
设插值点 x0 , x1 ,, xn  [1,1]
f  C n 1[1,1], Ln ( x)为相应的 n
次拉格朗日插值多项式,则余项
图3-3
f ( n1 ( )
Rn ( x)  f ( x)  Ln ( x) 
n1 ( x),
(n  1)!
70
于是
max f ( x)  Ln ( x)
1 x 1
M n 1

max ( x  x0 )( x  x1 )  ( x  xn ) ,

(n  1)! 1 x 1
其中
M n 1  f ( n 1) ( x)

 max f ( n 1) ( x)
1 x 1
是由被插函数决定的.
如果插值节点为 Tn1 ( x)的零点
xk  cos
2k  1
,
2(n  1)
k  0,1, , n,
71
则由(2.13)可得
1
~
max n 1 ( x)  max Tn 1 ( x)  n .
1 x 1
1 x 1
2
由此可导出插值误差最小化的理论.
定理7
设插值节点 x0 , x1 ,, xn 为切比雪夫多项式 Tn1 ( x)
的零点,被插函数 f  C n 1[1,1], Ln ( x) 为相应的插值多项
式,则
1
max f ( x)  Ln ( x)  n
f ( n 1) ( x) .

1 x 1
2 (n  1)!
(2.15)
72
对于一般区间 [ a, b]上的插值只要利用变换(2.14)则
可得到相应结果,此时插值节点为
ba
2k  1
ab
xk 
cos

,
2
2(n  1)
2
k  0,1, , n,
73
例4
求 f ( x)  e x在 [0, 1]上的四次拉格朗日插值多项式
x
max
e
 L4 ( x) .
,插值节点用
的零点并估计误差
L4 ( x)
T5 ( x)
0 x 1
解
利用 T5 ( x)的零点和区间变换可知节点
xk 
1
2k  1
(1  cos
 ),
2
10
k  0,1,2,3,4,
即
x0  0.97553,
x1  0.79390,
x3  0.20611,
x4  0.02447,
x2  0.5,
对应的拉格朗日插值多项式为
74
L4 ( x)  1.000 022 74  0.998 862 33x  0.509 022 51x 2
 0.141 841 05 x 3  0.068 494 35 x 4
利用(2.15)可得误差估计
M n 1 (b  a) n 1
max e  L4 ( x) 
,
2 n 1
0 x 1
(n  1)! 2
x
n  4,
而
M n1  f (5) ( x)

 ex

 e1  2.72,
于是有
max e x  L4 ( x) 
0 x 1
e 1
2.72 1
. 9 
 4.4 10 5.
5! 2
6 10240
75
在第2章中曾经指出,高次插值会出现龙格现象,一般
Ln (x)不收敛于 f (x),因此并不适用. 但若用切比雪夫多项
式零点插值却可避免龙格现象,可保证整个区间上收敛.
76
1
,在 [5,5]上利用T11 ( x)的零点
2
1 x
~
作插值点,构造10次拉格朗日插值多项式 L
. 与第2章
10 ( x)
例5
设 f ( x) 
得到的等距节点造出的 L10 ( x)近似 f (x) 作比较.
解
在[1, 1] 上的10次切比雪夫多项式 T11 ( x) 的零点为
t k  cos
21  2k
,
22
k  0,1, ,10.
作变换 xk  5tk ,
k  0,1,,10.它们是 (5,5)内的插值点,
~
由此得到 y  f (x) 在 [5,5]上的拉格朗日插值多项式 L
10 ( x ).
77
~
~
的图形见图3-4,从图上可见
L10 ( x)
f ( x), L10 ( x), L10 ( x)
没有出现龙格现象.
图3-4
78
3.2.5 其他常用的正交多项式
区间 [ a, b]及权函数  (x)不同,则得到的正交多项式也
不同.
除上述两种最重要的正交多项式外,下面是三种较常
用的正交多项式.
1. 第二类切比雪夫多项式
在区间[1, 1] 上带权  ( x) 
1  x 2 的正交多项式
称为第二类切比雪夫多项式.
79
表达式为
U n ( x) 
sin[( n  1) arccos x]
1 x2
.
(2.16)
令 x  cos ,可得
1


U n ( x)U m ( x) 1  x dx   sin( n  1) sin( m  1)d
1
2
0
 0, m  n;

 
, m  n,

2
即 {U n ( x)}是 [1, 1]上带权  ( x) 
1  x 2 的正交多项式族.
80
递推关系
U 0 ( x)  1, U1 ( x)  2 x,
U n1 ( x)  2 xUn ( x)  U n 1 ( x),
(n  1,2, ).
81
2. 拉盖尔多项式
在区间 [0,)上带权  ( x)  e  x 的正交多项式称为拉
盖尔(Laguerre)多项式.
其表达式为
n
d
n x
Ln ( x)  e x
(
x
e ).
n
dx
(2.17)
正交性质


0
 0, m  n;
e Ln ( x) Lm ( x)dx  
2
(n!) , m  n,
x
82
递推关系
L0 ( x)  1,
L1 ( x)  1  x,
Ln 1 ( x)  (1  2n  x) Ln ( x)  n 2 Ln 1 ( x),
(n  1,2, ).
83
3. 埃尔米特多项式
在区间 (,) 上带权  ( x)  e
 x2
的正交多项式称
为埃尔米特多项式.
表达式
H n ( x)  (1) e
n
x2
dn
 x2
(e ),
n
dx
(2.18)
正交关系



e
 x2
0, m  n;

H n ( x) H m ( x)dx   n
2 n!  , m  n,
84
递推关系
H 0 ( x)  1,
H1 ( x)  2 x,
H n1 ( x)  2 xHn ( x)  2nH n 1 ( x),
(n  1,2, ).
85
3.3
最佳平方逼近
86
最佳平方逼近及其计算
3.3.1
对 f ( x)  C[a, b] 及 C[a, b] 中的一个子集
  span{0 ( x), 1 ( x), , n ( x)}
若存在 S * ( x)   ,使
f ( x)  S ( x)
*
2
2
 min f ( x)  S ( x)
S ( x )
 min
S ( x )

b
a
2
(3.1)
2
 ( x)[ f ( x)  S ( x)]2 dx.
则称 S * ( x) 是 f (x) 在子集   C[a, b]中的最佳平方逼近
函数.
87
由(3.1)可知该问题等价于求多元函数
n
I (a0 , a1 , , an )    ( x)[  a j j ( x)  f ( x)]2 dx (3.2)
a
b
j 0
的最小值.
I (a0 , a1 , , an ) 是关于 a0 , a1 , , an 的二次函数,
利用多元函数求极值的必要条件
I
0
ak
即
(k  0,1,, n),
n
b
I
 2  ( x)[  a j j ( x)  f ( x)] k ( x)dx
a
ak
j 0
(k  0,1,, n),
88
于是有
n
 (
j 0
k
( x),  j ( x)) a j  ( f ( x),  k ( x)) (k  0,1, , n). (3.3)
这个关于 a0 , a1 ,, an 的线性方程组,称为法方程.
由于 0 ( x), 1 ( x), , n ( x) 线性无关,故
det G(0 , 1 ,, n )  0
于是方程组(4.3)有唯一解 ak  ak*
(k  0,1, , n),
S * ( x)  a0*0 ( x)    an* n ( x).
89
下面证明 S * ( x)满足(3.1),即对任何 S ( x)   , 有

b
a
b
 ( x)[ f ( x)  S ( x)] dx    ( x)[ f ( x)  S ( x)]2 dx
*
2
a
(3.4)
为此只要考虑
2
2
(3.1)
f ( x)  S ( x)  min2 f ( x)  S ( x) 2
2 S (Sx
D    ( x)[ f ( x) 
)])dx    ( x)[ f ( x)  S * ( x)]2 dx
(x
*
b
a

b
a
b
b
a
 *min 2  ( x)[ f ( x)  S ( x)]2 dx.
 ( x)[ f ( x)  S S(( xx ))]
a
 dx
b
 2  ( x)[ S * ( x)  S ( x)][ f ( x)  S * ( x)]dx.
a
90
由于S * ( x)的系数 a k* 是方程(3.3)的解,故

b
a
 ( x)[ f ( x)  S * ( x)]k ( x)dx  0 (k  0,1,, n),
n
从而上式第二个积分为0,
x),  k ( x)) (k  0,1, , n). (3.3)
 (k ( x),  j ( x)) a j  ( f (于是
j 0
b
D    ( x)[ S ( x)  S * ( x)]2 dx  0,
a
故(3.4)成立.
这就证明了S * ( x)是 f (x) 在  中的最佳平方逼近函数.
b
b
*
2
2

(
x
)[
f
(
x
)

S
(
x
)]
dx


(
x
)[
f
(
x
)

S
(
x
)]
dx


a
a
(3.4)
91
若令  ( x)  f ( x)  S * ( x), 则平方误差为
 ( x)
2
2
 ( f ( x)  S * ( x), f ( x)  S * ( x))
 ( f ( x), f ( x))  ( S * ( x), f ( x))
 f ( x)
2
2
n
  ak* ( k ( x), f ( x)).
(3.5)
k 0
若取  k ( x)  x k ,  ( x)  1, f ( x)  C[0, 1], 则要在 H n
中求 n次最佳平方逼近多项式
S * ( x)  a0*  a1* x    an* x n ,
92
此时
( j ( x),  k ( x)) 
( f ( x),  k ( x)) 
1

0
1

0
x k  j dx 
1
,
k  j 1
f ( x) x k dx  d k .
n
若用 H表示 Gn  G (1, x, , x )
1

 1/ 2
H



1 /( n  1)
1/ 2
1/ 3

1 /( n  2)
1 /( n  1) 
 1 /( n  2) 




 1 /( 2n  1)
即

(3.6)
称为希尔伯特(Hilbert)矩阵.
93
T
T
记 a  (a0 , a1 ,, an ) , d  (d 0 , d1 ,, d n ) , 则
Ha  d
(3.7)
*
a

a
的解 k
k ( k  0,1,  , n) 即为所求.
94
设 f ( x)  1  x 2 , 求 [0, 1]上的一次最佳平方
例6
逼近多项式.
解 利用(3.7),得
1
d0  
0
1
1  x dx  ln( 1 
2
2
2
2) 
 1.147,
2
2 2 1
1
2 3/ 2
2
 0.609,

d1   x 1  x dx  (1  x )
0
3
3
0
1
1
得方程组
1 / 2 a0  1.147 
 1
1 / 2 1 / 3  a   0.609,

 1  

95
解之
a0  0.934,
a1  0.426,
故
S1* ( x)  0.934  0.426 x.
平方误差
 ( x)
2
2
 ( f ( x), f ( x))  ( S1* ( x), f ( x))
1
  (1  x 2 )dx  0.426d1  0.934d 0  0.0026.
0
最大误差
 ( x)

 max
0 x 1
1  x 2  S1* ( x)  0.066.
96
3.3.2
用正交函数族作最佳平方逼近
设 f ( x)  C[a, b],   span{0 ( x), 1 ( x), , n ( x)},
若 0 ( x), 1 ( x), , n ( x) 是满足条件(2.2)的正交函数族,
则
(i ( x),  j ( x))  0, i  j
 (k ( x),  j ( x)) a j  ( f ( x), k ( x)) (k  0,1,, n). (3.3)
n
而 j 0
b
j  k.
(

),)
 j (( xx))
 (0x)dx  0,
j ((xx
( j ,  k )   
)


j
k
a
j  k.
 Ak  0,
故法方程(3.3)的系数矩阵
(2.2)
Gn  G(0 ( x), 1 ( x), , n ( x))
97
为非奇异对角阵,且方程(3.3)的解为
ak*  ( f ( x),  k ( x)) /( k (k ),  k ( x))
(3.8)
(k  0,1,, n).
( k ( x),  j ( x)) a j  ( f ( x),  k ( x))(k  0,1, , n). (3.3)

j 0
f ( x)  C[a, b] 在  中的最佳平方逼近函数为
于是
n
n
S ( x)  
*
k 0
( f ( x),  k ( x))
 k ( x)
2
 k ( x).
(3.9)
2
98
由(3.5)可得均方误差为
 n ( x)
2
 f ( x)  S n* ( x)
2
1
2
2


n 

(
f
(
x
),

(
x
))
2

k
  f ( x) 2   
  .
 k ( x) 2  

k 0 

 2

*
*
 ( x) 2  ( f ( x)  S ( x), f ( x)  S ( x))
(3.10)
由此可得贝塞尔(Bessel)不等式
 ( f ( x), f ( x))  ( S * ( x), f ( x))
n
(3.5)
*
2
n
(
a

(
x
)
)

f
(
x
)
. *
2
(3.11)

k
k
2
2
 f ( x) 2   ak ( k ( x), f ( x)).
k 1
2
k 0
99
若 f ( x)  C[a, b],按正交函数族 { k ( x)} 展开,系数
ak* (k  0,1, ) 按(3.8)计算,得级数

*
a
 k k ( x)
(3.12)
k 0
系数 a k*
*
称这个级数为
的广义傅里叶(Foureir)级数,
f
(x
)
ak  ( f ( x),  k ( x)) /( k (k ),  k ( x))
(3.8)
称为广义傅里叶系数. 它是傅里叶级数的直接推广.
(k  0,1,, n).
讨论特殊情况,设 {0 ( x), 1 ( x), , n ( x)} 是正交多
项式,  span{0 ( x), 1 ( x), , n ( x)}, k ( x)( k  0,1,, n)
可由 1, x,, x n 正交化得到,则有下面的收敛定理.
100
定理8
设 f ( x)  C[a, b], S * ( x) 是由(3.9)给出的
f (x)的最佳平方逼近多项式,其中 {k ( x), k  0,1,, n}
是正交多项式族,则有
lim f ( x)  S n* ( x)
n 
2
 0.
( f ( x),  k ( x))
S ( x)考虑函数

 k ( x).按勒让德多项式
2
(3.9)
f ( (xx) )C[1, 1],

k 0
n
*
k
2
{P0 ( x), P1 ( x), , Pn ( x)} 展开,由(3.8),(3.9)可得
S n* ( x)  a0* P0 ( x)  a1* P1 ( x)    an* Pn ( x),
(3.13)
101
 n ( x)
其中
2
 f ( x)  S n* ( x)
2
1
2
 P ( x))
( f ( x),
2nk (1f (1x),  k ( x))  
2

k
(3.1
a ( x) 
..
  f ( x) 2 
    f ( x) Pk ( x)dx

1 ( x )
( Pk ( x), Pk ( x)) k 02


k
2
  (3.14)

2
*
k
根据均方误差公式(3.10),平方误差为
 k ( x)
2
2
1

1
n
2
*2
f ( x)dx  
ak .
k  0 2k  1
2
(3.15)
由定理8可得
lim f ( x)  S n* ( x)
n 
2
 0.
102
*
*
* *一致收敛于
(3.13)
S n* ( xf) (x
 )a满足光滑性条件,还有
P
(
x
)

a
P
(
x
)



a
如果
f (x)
0 0
1 1
S n n Pn ( x),
的结论.
定理9
设 f ( x) C 2 [1, 1], S n* ( x)由(3.13)给出,
则对任意 x  [1, 1]和   0, 当 n充分大时有
f ( x)  S n* ( x) 

n
.
~ ,
公式(2.6)给出了首项系数为1的勒让德多项式 P
n
它具有以下性质.
~
Pn ( x) 
n! d n
2
n
[(
x

1
)
].
n
(2n)! dx
103
(2.6)
定理10 在所有最高次项系数为1的 n次多项式中,
~
勒让德多项式 P
在 [1, 1] 上与零的平方误差最小.
n ( x)
证明 设 Qn (x)是任意一个最高次项系数为1的 n次
多项式,它可表示为
n 1
~
~
Qn ( x)  Pn ( x)   ak Pk ( x),
k 0
104
于是
Qn ( x)
2
2
1
 (Qn ( x), Qn ( x))   Qn2 ( x)dx
1
n 1
~
~
~
2 ~
 ( Pn ( x), Pn ( x))   ak ( Pk ( x), Pk ( x))
k 0
~
~
 ( Pn ( x), Pn ( x))
~
 Pn ( x)
2
2
当且仅当 a0  a1    an1  0 时等号才成立, 即当
~
Qn ( x)  Pn ( x) 时平方误差最小.
105
例7
求 f ( x)  e x 在 [1, 1]上的三次最佳平方逼近多
项式.
解
~
(
f
(
x
),
P
先计算
k ( x))
(k  0,1,2,3).
1
( f ( x), P0 ( x))   e dx  e   2.3504;
1
e
1
x
1
( f ( x), P1 ( x))   xe x dx  2e 1  0.7358;
1
3
1
( f ( x), P2 ( x))   ( x 2  )e x dx
1 2
2
1
5 3 3
( f ( x), P3 ( x))   ( x  x)e x dx
1 2
2
1
e
7
 0.1431;
e
1
 37  5e  0.02013.
e
106
由(3.14) 得
a0*  ( f ( x), P0 ( x)) / 2  1.1752,
*
a
3x())f ( x),
P1( x
))1/ 2  1.1036,
( f ( x),1P
(
2
k
1
*
k
ak ( x ) 

f ( x) Pk ( x)dx.


1
( Pk ( x),*Pk ( x))
2
a2  5( f ( x), P2 ( x)) / 2  0.3578, (3.14)
a3*  7( f ( x), P3 ( x)) / 2  0.07046.
代入(3.13) 得三次最佳平方逼近多项式
S3* ( x)  0.9963  0.9979 x  0.5367 x 2  0.1761x 3 .
S n* ( x)  a0* P0 ( x)  a1* P1 ( x)    an* Pn ( x),
(3.13)
107
均方误差
 n ( x) 2  e  S ( x)
x
*
3
2
3
2
*2
  e dx  
ak
1
k 0 2k 1
1
2x
 0.0084.
最大误差
 n ( x)
x
*

e

S
3 ( x)


 0.0112.
如果 f ( x)  C[a, b], 求 [ a, b]上的最佳平方逼近多项式,
做变换
ba
ba
x
t
2
2
(1  t  1),
108
于是
ba
ba
F (t )  f (
t
),
2
2
*
在 [1, 1]上可用勒让德多项式做最佳平方逼近多项式 S n (t ),
从而得到区间 [ a, b]上的最佳平方逼近多项式
S n* (
1
(2 x  a  b)).
ba
109
由于勒让德多项式 {Pk ( x)}是在区间 [1, 1]上由
{1, x,, x k ,} 正交化得到的,因此利用函数的勒让德展
开部分和得到最佳平方逼近多项式与由
S * ( x)  a0  a1 x    an x n
直接通过解法方程得到 H n 中的最佳平方逼近多项式是一致
的.
只是当 n 较大时法方程出现病态,计算误差较大,不能
使用,而用勒让德展开不用解线性方程组,不存在病态问题,
因此通常都用这种方法求最佳平方逼近多项式.
110
切比雪夫级数
3.3.3
如果 f ( x) C[1,1] ,按 {Tk ( x)}0 展成广义傅里叶级数,
由(3.12)可得级数

C0*
  Ck*Tk ( x),
2
k 1
(3.16)
其中系数根据(3.8)式,由(2.12)得到
C
这里
*
k

2

1
f ( x)Tk ( x) dx
1
1 x2

,
Tk ( x)  cos(k arccos x),
k  0,1, , (3.17)
x  1.
级数(3.16)称为 f (x) 在 [1,1] 上的切比雪夫级数.
111
若令 x  cos  , 0     ,则(3.16)就是 f (cos  ) 的傅里
叶级数,其中
C 
*
k
2

1

1
f (cos  ) cos k d ,
k  0,1, .
(3.18)
根据傅里叶级数理论,只要 f (x) 在 [1,1]上分段连续,则
f (x)在 [1,1]上的切比雪夫级数(3.16)一致收敛于 f (x).
从而

C0*
f ( x) 
  Ck*Tk ( x).
2
k 1
(3.19)
112
取部分和
n
C0*
C ( x) 
  Ck*Tk ( x),
2
k 1
*
n
(3.20)
其误差为
f ( x)  Cn* ( x)  Cn*1Tn 1 ( x).
在 [1,1] 上 Tn1 ( x) 是均匀分布的,它的最大值 max Tn1 ( x)
1 x 1
最小,因此 Cn* ( x) 可作为 f (x) 在 [1,1] 上的近似最佳一致
逼近多项式.
113
例8
解
求 f ( x)  e x 在 [1, 1]上的切比雪夫部分和 C3* ( x) .
由(3.18)得
C 
*
k
2



0
e cos cos k d ,
k  0,1,2,3.
利用第4章的数值积分法得
C0*  2.53213176, C1*  1.13031821,
C2*  0.27149534, C3*  0.04433685.
由(3.20)及 Tk (x) 的表达式有
C3* ( x)  0.994531  0.997308 x  0.542991x 2  0.177347 x 3 ,
及 e x  C3* ( x)

 0.00607.
114
3.4
曲线拟合的最小二乘法
3.4.1
最小二乘法及其计算
在函数的最佳平方逼近中 f ( x)  C[a, b], 如果 f (x)
只在一组离散点集 {xi , i  0,1, , m} 上给定,这就是科
学实验中经常见到的实验数据 {( xi , yi ), i  0,1,, m}的
曲线拟合.
115
问题为利用 yi  f ( xi ), i  0,1,, m, 求出一个函数
y  S * ( x) 与所给数据{( xi , yi ), i  0,1,, m} 拟合.
记误差
 i  S * ( xi )  yi ,
i  0,1,, m,
  ( 0 , 1 ,,  m )T
116
设 0 ( x), 1 ( x), , n ( x)是 C[a, b] 上线性无关函数族,
在   span{0 ( x), 1 ( x), , n ( x)}中找一函数 S * ( x),
使误差平方和

2
2
m
m
     [ S * ( xi )  yi ]2
i 0
 min
2
i
S ( x )
i 0
m
2
[
S
(
x
)

y
]
 i
i
(4.1)
i 0
这里
S ( x)  a00 ( x)  a11 ( x)    ann ( x)
( n  m)
(4.2)
117
这个问题称为最小二乘逼近,几何上称为曲线拟合的
最小二乘法.
用最小二乘求拟合曲线时,首先要确定 S (x ) 的形式.
确定 S (x ) 的形式问题不仅是数学问题, 还与问题的
实际背景有关.
通常要用问题的运动规律及给定的数据进行数据描图,
确定 S (x )的形式, 然后通过实际计算选出较好的结果.
118
S (x ) 的一般表达式为(4.2)表示的线性形式.
若  k (x)是 k 次多项式,S (x ) 就是 n 次多项式.
为了使问题的提法更有一般性,通常在最小二乘法中
S ( x)  a00 ( x)  a11 ( x)    ann ( x) (n  m)
考虑加权平方和
(4.2)

2
2
m
   ( xi )[ S ( xi )  yi ]2 .
(4.3)
i 0
这里  ( x)  0是 [ a, b]上的权函数,它表示不同点 ( xi , f ( xi ))
处的数据比重不同.
119
用最小二乘法求拟合曲线的问题,就是在形如(4.2)的
S (x )中求一函数 y  S * ( x),使(4.3)取得最小.
这样,最小二乘问题就转化为求多元函数
( n  m)
m
n
(4.2)
   ( xi )[ a j j ( xi )  f ( xi )]2
(4.4)
S (Ix()a
a1
) n)a11 ( x)    ann ( x)
,
,a
0 (x
0, 0
i 0
j 0
m
*
2
,
,
a
的极小点 (a , a 

 n()x问题.
i )[ S ( xi )  yi ] .
*
0
2*
1
2
(4.3)
i 0
由求多元函数极值的必要条件,有
m
n
I
 2  ( xi )[  a j j ( xi )  f ( xi )] k ( xi )  0
ak
i 0
j 0
(k  0,1,, n).
120
若记
m
( j ,  k )    ( xi ) j ( xi ) k ( xi ),
(4.5)
i 0
m
( f ,  k )    ( xi ) f ( xi ) k ( xi )  d k
i 0
(k  0,1, , n).
上式可改写为
n
 (
j 0
k
,  j )a j  d k
(k  0,1,, n). (4.6)
这方程称为法方程,可写成矩阵形式
121
Ga  d .
T
T
a

(
a
,
a
,

,
a
)
,
d

(
d
,
d
,

,
d
)
,
其中
0
1
n
0
1
n
( 0 ,  0 )
 ( ,  )
G 1 0



( n ,  0 )
( 0 , 1 )

(1 , 1 )


( n , 1 )

( 0 ,  n ) 
(1 ,  n ) 
.



( n ,  n )
(4.7)
要使法方程(4.6)有唯一解,就要求矩阵 G非奇异,
而 0 ( x), 1 ( x), , n ( x)在 [ a, b]上线性无关不能推出
n
矩阵 G非奇异,必须加上另外的条件.
 (k ,  j )a j  d k (k  0,1,, n). (4.6)
j 0
122
定义7
设 0 ( x), 1 ( x), , n ( x)  [a, b]的任意线
性组合在点集 {xi , i  0,1,, m}(m  n) 上至多只有 n 个
不同的零点,则称 0 ( x), 1 ( x), , n ( x)在点集
{xi , i  0,1,, m}上满足哈尔(Haar)条件.
显然 1, x,, x n在任意 m(m  n)个点上满足哈尔条件.
如果 0 ( x), 1 ( x), , n ( x)  [a, b] 在 {xi }0m 上满足
哈尔条件,则法方程(4.6) 的系数矩阵(4.7) 非奇异,于是
方程(5.6)存在唯一的解 ak  ak* , k  0,1, , n. 从而得到
n
函数 f (x) 的最小二乘解为
 (k ,  j )a j  d k
j 0
(k  0,1,, n). (4.6)
123
S * ( x)  a0*0 ( x)  a1*1 ( x)    an* n ( x).
这样得到的 S * ( x), 对任何形如(4.2)的 S (x ) ,都有
m
m
i 0
i 0
*
2
2

(
x
)[
S
(
x
)

f
(
x
)]


(
x
)[
S
(
x
)

f
(
x
)]
,
 i
 i
i
i
i
i
故 S * ( x)确是所求最小二乘解.
S ( x)  a00 ( x)  a11 ( x)    ann ( x)
( n  m)
(4.2)
124
给定 f (x) 的离散数据 {( xi , yi ), i  0,1,, m},
一般可取   span{1, x,, x n } ,但这样做当 n  3时,
求解法方程(4.6)将出现系数矩阵 G为病态的问题,
通常对 n  1的简单情形都可通过求法方程(4.6)得到
S * ( x).
有时根据给定数据图形,其拟合函数 y  f (x) 表面上
不是(4.2)的形式,但通过变换仍可化为线性模型.
例如, S ( x)  aebx ,若两边取对数得
S ( x)  a00 ( x)  a11 ( x)    ann ( x)
ln S ( x)  ln a  bx,
( n  m)
125
(4.2)
此时,若令 S ( x)  ln S ( x),
A  ln a,
B  b,
则
S ( x)  A  Bx ,
这样就变成了形如(4.2)的线性模型 .
126
例9
解
已知一组实验数据如下,求它的拟合曲线.
xi
1
2
3
4
5
fi
4
4.5
6
8
8.5
i
2
1
3
1
1
将所给数据在坐标纸上标出,见图3-5.
127
图3-5
从图中看到各点在一条直线附近,故可选择线性函数
作拟合曲线,
令 S1 ( x)  a0  a1 x, 这里 m  4, n  1, 0 ( x)  1,
1 ( x)  x, 故
4
( 0 ,  0 )   i  8,
i 0
4
( 0 , 1 )  (1 ,  0 )   i xi  22,
i 0
4
(1 , 1 )   i xi2  74,
i 0
4
( 0 , f )   i f i  47,
i 0
4
(1 , f )   i xi f i  145.5.
i 0
128
由(4.6)得方程组
8a0  22a1  47

22a0  74a1  145.5.
n
 (
,  j )a j  d k
(k  0,1,, n). (4.6)
j 0
解得 a0  2.77, a1  1.13. 于是所求拟合曲线为
k
S1* ( x)  2.77  1.13x.
129
关于多项式拟合,Matlab中有现成的程序
a  polyfit ( x, y, m)
其中输入参数 x, y为要拟合的数据, m为拟合多项式的次数,
输出参数 a 为拟合多项式的系数.
利用下面的程序,可在Matlab中完成上例的多项式拟合.
130
x=[1 2 3 4 5];
f=[4 4.5 6 8 8.5];
aa=polyfit(x,f,1);
y=polyval(aa,x);
plot(x,f,’r+’,x,y,’k’)
xlabel(‘x’);
ylabel(‘y’);
gtext(‘y=s1(x)’)
131
结果如下:
132
例10 设数据 ( xi , yi )(i  0,1,2,3,4) 由表3-2给出,
表中第4行为 ln yi  yi ,通过描点可以看出数学模型为
y  aebx , 用最小二乘法确定 a及 b .
表3  2
i
0
1
2
3
4
xi
1.00
1.25
1.50
1.75
2.00
yi
5.10
5.79
6.53
7.45
8.46
133
解
用给定数据描图可确定拟合曲线方程为 y  aebx ,
它不是线性形式. 两边取对数得
ln y  ln a  bx,
若令 y  ln y, A  ln a, 则得 y  A  bx,   {1, x}.
为确定 A, b ,先将 ( xi , yi )转化为 ( xi , yi ), 数据表见表3-2.
0 ( x)  1, 1 ( x)  x,  ( x)  1,
得
表3  2
(0 , 0 )  5,
41
i
0
2
( 0 , 1 )   xi  7.5,
xi
1.00
1.25
1.50
i 0
3
4
1.75
2.00
5.10
54 .792
6.53
7.45
(1 , 1 )   xi  11.875,
yi 1.629 1i .0756 1.876 2.008
yi
8.46
2.135
134
4
( 0 , y )   yi  9.404,
i 0
4
(1 , y )   xi yi  14.422.
i 0
故有法方程
5 A  7.50b  9.404


7.50 A  11.875b  14.422.
解得
A 1.122, b  0.505, a  e A  3.071.
于是得最小二乘拟合曲线为
y  3.071e0.505x .
135
利用下面的程序,可在Matlab中完成曲线拟合.
x=[1.00 1.25 1.50 1.75 2.00];
y=[5.10 5.79 6.53 7.45 8.46];
y1=log(y);
aa=polyfit(x,y1,1);
a=aa(1); b=exp(aa(2));
y2=b*exp(a*x);
plot(x,y,’r+’,x,y2,’k’)
xlabel(‘x’);
ylabel(‘y’);
gtext(‘y=a*exp(bx))’;
136
结果如下:
137
3.4.2
用正交多项式做最小二乘拟合
用最小二乘法得到的法方程组(4.6),其系数矩阵
G是病态的.
如果
0 ( x), 1 ( x), , n ( x) 是关于点集
n
( k ,  j )a j  d k
(k  0,1,, n). (4.6)

{xi }j 0
(i  0,1, , m) 带权  ( xi ) (i  0,1,, m)
函数族,即
 0,
( j ,  k )    ( xi ) j ( xi ) k ( xi )  
i 0
 Ak  0,
m
j  k,
j  k,
(4.8)
138
则方程(4.6)的解为
m
ak* n

j 0
  ( x ) f ( x )
( xi )
( f ,k )
 i 0 m
(k  0,1,, n).
 k,,
 k ))a  d
((
k
j
j
k  ( x(i k)k2 (0x,i1),, n). (4.6)
i
i
k
i 0
(4.9)
且平方误差为

2
2
 f
2
2
n
  Ak (ak* ) 2 .
k 0
139
接下来根据给定节点 x0 , x1 , , xm 及权函数  ( x)  0,
构造带权  ( x ) 正交的多项式 {Pn ( x)} .
注意 n  m,用递推公式表示 Pk (x) ,即
 P0 ( x)  1,
(4.10)

 P1 ( x)  ( x  1 ) P0 ( x),
 P ( x)  ( x   ) P ( x)   P ( x)
k 1
k
k k 1
 k 1
(k  1,2,, n  1).
这里 Pk (x) 是首项系数为1的 k 次多项式, 根据 Pk (x) 的
正交性,得
140
m
2

(
x
)
x
P
 i i k ( xi )

( xPk ( x), Pk ( x))
  k 1  i  0m

( Pk ( x ), Pk ( x ))

2

(
x
)
P
(
x
)

i
k
i

i 0


( xPk , Pk )


( Pk , Pk )


m

 ( xi ) Pk2 ( xi )


( Pk , Pk )
i 0




k
m
( Pk 1 , Pk 1 )
2


(
x
)
P
(
x
)
 i k 1 i
i 0
(4.11)
(k  1,2,, n  1).
下面用归纳法证明这样给出的 {Pk ( x)}是正交的.
141
由(4.10)第二式及(4.11)中 1 的表达式,有
( P0 , P1 )  ( P0 , xP0 )  1 ( P0 , P0 )
 ( P0 , xP0 ) 
( xP0 , P0 )
( P0 , P0(4.10)
)  0.
( P0 , P0 )
 P0 ( x)  1,

 P1 ( x)  ( x  1 ) P0 ( x),
( Pl , Ps )  0(l  s) 对 s  0,1,, l  1 及 l  0,1, , k ;
假定
 Pk 1 ( x)  ( x   k 1 ) Pk ( x)   k Pk 1 ( x)
k  n 均成立,要证 ( Pk 1 , Ps )  0 对 s  0,1, , k 均成立.
(k  1,2,, n  1).
由(4.10)有
( Pk 1 , Ps )  (( x   k 1 ) Pk , Ps )   k ( Pk 1 , Ps )
(4.12)
 P0 ( x)  1,
(4.10)

0 
x k ,P1s))P
( xk),
1 ( Pk , Ps )   k ( Pk 1 , Ps ).
 P1 ( x)  ( xP
 P ( x)  ( x   ) P ( x)   P ( x)
142
k 1
k
k k 1
 k 1
(k  1,2,, n  1).
由归纳法假定,当 0  s  k  2 时
( Pl , Ps )  0,
( Pk 1 , Ps )  0.
另外,xPs (x) 是首项系数为1的 s  1 次多项式,它可由
P0 , P1 , , Ps 1 的线性组合表示.
而 s  1  k  1,由归纳法假定又有
( xPk , Ps )  ( Pk , xPs )  0,
于是由(4.12),当 s  k  2 时,( Pk 1 , Ps )  0.
( Pk 1 , Ps )  (( x   k 1 ) Pk , Ps )   k ( Pk 1 , Ps )
(4.12)
 ( xPk , Ps )   k 1 ( Pk , Ps )   k ( Pk 1 , Ps ).
143
再考虑
( Pk 1 , Pk 1 )  ( xPk , Pk 1 )   k 1 ( Pk , Pk 1 )   k ( Pk 1 , Pk 1 ),
(4.13)
由假定有
( Pk , Pk 1 )  0,
k 1
( xPk , Pk 1 )  ( Pk , xPk 1 )  ( Pk , Pk   c j Pj )  ( Pk , Pk ).
j 0
利用(4.11)中  k 表达式及以上结果,得
( Pk 1 , Pk 1 )  ( xPk , Pk 1 )   k ( Pk 1 , Pk 1 )
( Pk , Pk )
k 
).
 ( Pk , Pk )  ( P(kP
,kP

1
k 10
k ,)P
144
最后,由  k 和  k 的表达式(4.11)有
( Pk 1 , Pk )  ( xPk , Pk )   k 1 ( Pk , Pk )   k ( Pk , Pk 1 )
( xPk ,(P
xP
P)
k)
k,
 ( xPk , Pk ) 
( Pkk , Pk )  0.
 k 1 
( Pk , P
(kP)k , Pk )
至此已证明了由(4.10)及(4.11)确定的多项式
( Pk , Pk )
{Pk ( x)}
k 
( Pk 1 , Pk 1 )
{
x
}
组成一个关于点集 i 的正交系.
用正交多项式 {Pk ( x)}的线性组合作最小二乘曲线拟合,
只要根据公式(4.10)及(4.11)逐步求 Pk (x) 的同时,
*
相应计算出系数 a k
145
m
( f , Pk )

a 
( Pk , Pk )
*
k
(x ) f (x )P (x )
i
i 0
i
k
m
  ( x )P
i 0
i
2
k
i
(k  0,1, , n),
( xi )
并逐步把 ak* Pk ( x)累加到 S (x ) 中去,最后就可得到所求的
拟合曲线
y  S ( x)  a0* P0 ( x)  a1* P1 ( x)    an* Pn ( x).
这里 n可事先给定或在计算过程中根据误差确定.
用这种方法编程序不用解方程组,只用递推公式,并
且当逼近次数增加一次时,只要把程序中循环数加1,其余
不用改变.
146
147
2020/5/1
147
2020/5/1
148
148
149
2020/5/1
149
150
151
3.5
有 理 逼 近
3.5.1
有理逼近与连分式
有理函数逼近是指用形如
n
Rnm ( x) 
的函数逼近 f (x ).
Pn ( x)

Qm ( x)
a
k
xk
k
xk
b
k 0
与前面讨论一样,如果 f ( x)  Rnm ( x)
最佳有理一致逼近.
(5.1)
k 0
m

最小就可得到
152
如果 f ( x)  Rnm ( x)
2
最小则可得到最佳有理平方逼近
函数.
本节主要讨论利用函数的泰勒展开获得有理逼近函数
的方法.
对函数 ln( 1  x)用泰勒展开得

ln( 1  x)   (1)
k 1
k 1
xk
k
x  [1,1].
(5.2)
取部分和
n
S n ( x)   (1) k 1
k 1
xk
 ln( 1  x).
k
153
另一方面若对(5.2)式用辗转相除可得到 ln( 1  x) 的
一种连分式展开
x
1 x
1
k

1 x
x
k

1
2

ln( 1  x)   (1)
x2
2[x1,1].
3k
k 1
22  x
4
5 
ln( 1  x ) 
(5.2)
x 1 x 1 x 2 2  x 22  x

.
1  2  3  4  5 
(5.3)
154
(5.3)右端为 ln( 1  x) 的无穷连分式的前5项,最后式子
是它的紧凑形式.
若取(5.3)的前2,4,6,8项,则可分别得到 ln( 1  x)
的以下有理逼近
2x
6 x  3x 2
R22 ( x) 
,
 R11 ( x )  2  x ,
2
6  6x  x


60 x  60 x 2  11x 3

 R33 ( x) 
,
2
3

60  90 x  36 x  3 x


2  260 x 3  25 x 4
420
x

630
x
 R44 ( x) 
. (5.4)

2
3
4
420  840 x  540 x  120 x  6 x

155
若用同样多项的泰勒展开部分和 S 2 n ( x) 逼近 ln( 1  x)
并计算 x  1 处的值 S 2n (1)及 Rnn (1),计算结果见表3-3.
表3  3
n
S 2 n (1)
 s  ln( 2)  S 2 n (1)
Rnn (1)
 R  ln( 2)  Rnn (1)
1
0.50
0.19
0.667
0.026
2
0.58
0.11
0.69231
0.00084
3
0.617
0.076
0.693122
0.000025
4
0.634
0.058
0.69314642
0.00000076
ln 2 的准确值为 0.69314718, 从表3-1可以看出,
R44 (1)  0.69314642,
S8 (1)  0.634,
156
由此看出 R44 (1) 的精度比 S8 (1) 高出近10万倍,
但它们的计算量是相当的,这说明用有理逼近比多项
式逼近好得多.
157
例11
给出有理函数
2 x 4  45 x 3  381x 2  1353x  1511
R43 ( x) 
x 3  21x 2  157 x  409
用辗转相除法将它化为连分式并写成紧凑形式.
解
用辗转相除可逐步得到
158
4 x 2  64 x  284
R43 ( x)  2 x  3 
x 3  21x 2  157 x  409
 2x  3 
4
6( x  9)
x5
x 2  16 x  71
4
 2x  3 
x5
6
x7
 2x  3 
8
x9
4
6
8
.
x5 x7  x9
本例中用连分式计算 R43 ( x)的值只需3次除法,1次乘
法和7次加法.
159
若直接用多项式计算的秦九韶算法则需6次乘法和1次
除法及7次加法.
可见将 Rnm (x) 化成连分式可节省计算乘除法次数.
对一般的有理函数,(5.1)可转化为一个连分式
cl
c2
Rnm ( x)  P1 ( x) 
.
x  d1    x  d l
它的乘除法运算只需 max( m, n) 次.
而直接用有理函数(5.1)计算乘除法次数为 n  m 次.
160
3.5.2
帕德逼近
利用函数 f (x) 的泰勒展开可以得到它的有理逼近.
设 f (x)在 x  0的泰勒展开为
N
1
f ( x)  
f
k  0 k!
(k )
f ( N 1) ( ) N 1
(0) x 
x
.
( N  1)!
k
(5.5)
它的部分和记作
N
1
P( x)  
f
k  0 k!
(k )
(0) x k 
N
c
k
xk .
(5.6)
k 0
161
定义8
设 f ( x)  C N 1 (a, a), N  n  m,
如果有理函数
a0  a1 x    an x n
Pn ( x)
Rnm ( x) 

,
m
1  b1 x    bm x
Qm ( x)
(5.7)
其中 Pn ( x), Qm ( x) 无公因式,且满足条件
(k )
Rnm
(0)  f
(k )
(0)
(k  0,1, , N ),
(5.8)
则称 Rnm (x) 为函数 f (x) 在 x  0 处的 (n, m) 阶帕德逼近,
记作 R(n, m) ,简称 R(n, m) 的帕德逼近.
162
根据定义,若令
h( x)  P( x)Qm ( x)  Pn ( x),
则满足条件(5.8)等价于
h ( k ) (0)  0
即
k  0,1,, N .
(k )
(5.8)
Rnm
(0)  f ( k ) (0)
(k  0,1, , N ),
(k )
(k )
h (0)  ( P( x)Qm ( x)  Pn ( x))
 0,
k  0,1,, N .
x 0
由于 Pn( k ) (0)  k!ak , 应用莱布尼兹求导公式得
( P( x)Qm ( x)  Pn ( x)) ( k )
k
x 0
 k! c j bk  j  k!ak  0
j 0
k  0,1, , N , 163
1
这里 c j  f
j!
( j)
(0) 是由(5.6)得到的,上式两端除 k! ,
并由 b0  1, b j  0(当j  m时), 可得
k 1
ak   c j bk  j  ck
(k  0,1, , n)
(5.9)
j 0
N
及
1
P( x)  
f
k  0 k!
k 1
  c j bk  j  ck
(k )
(0) x k 
N
c
k
xk .
(5.6)
k 0
(k  n  1, , n  m)
(5.10)
j 0
注意当 j  m 时 b j  0, 故(5.10)可写成
164
  cn  m 1bm    cn 1b2  cn b1  cn 1 ,
 c

n  m  2 bm    cn b2  cn 1b1  cn  2 ,

    


 cn bm    cn  m  2b2  cn  m 1b1  cn  m ,
(5.11)
, 若记
其中 j  0 时 c j  0,
  cn  m 1
 c
H   nm 2



  cn

 cn 1

 cn



 cn  m  2
b  (bm , bm1 ,, b1 )T ,
 cn

 cn 1 
, (5.12)



 cn  m 1 
c  (cn1 , cn 2 ,, cn m )T .
165
则方程组(5.11)的矩阵形式为
Hb  c .
定理11
设 f ( x)  C N 1 (a, a),
N  n  m, 则形如
(5.7)的有理函数 Rnm (x) 是 f (x) 的 (n, m)阶帕德逼近的
充分必要条件是多项式 Pn ( x), Qm ( x) 的系数 a0 , a1 ,, an
及 b0 , b1 ,, bm 满足方程组(5.9)及(5.11).
166
根据定理11, 求 f (x) 的帕德逼近时,
首先要由(5.11)解出 Qm (x)的系数 b0 , b1 ,, bm ,
再由(5.9)直接算出 Pn (x) 的系数 a0 , a1 ,, an .
f (x)的各阶帕德逼近可列成一张表,称为帕德表(见表
  cn  m 1bm    cn 1b2  cn b1  cn 1 ,
 c
3-4).

n  m  2 bm    cn b2  cn 1b1  cn  2 ,
表3 - 4
 m    
0
1
2
3

n cn bm    cn  m  2b2  cn  m 1b1  cn  m ,
(5.11)
4

0
(0,0)
(0,1)
(0, 2)
(0,3)
(0, 4)

1
(1,0)
(1,1)
(1, 2)
(1,3)
(1, 4)

( 2,0)
( 2,1)
( 2, 2)
( 2,3)
( 2, 4)

( 4,0)
( 4,1)
( 4, 2)
( 4,3)
( 4, 4)







2
k 1
ak 3  c(3
,2
0,)1, (,3n,)3)
jb
)
j  c(k3,1) ( k(3
,k0
4

j 0
((5.9)
3, 4)

167
例12 求 f ( x)  ln( 1  x)的帕德逼近 R (2,2) 及 R(3,3).
解
由 ln( 1  x) 的泰勒展开
1 2 1 3 1 4
ln( 1  x)  x  x  x  x  
2
3
4
1
1
1
得 c0  0, c1  1, c2   , c3  , c4   , .
2
3
4
当 n  m  2时,由(5.11)得
1
1


b

b

,
2
1

2
3
1
1
1
 b2  b1   .
3
4
2
求得 b1  1, b2  1 , 再由(5.9)得
6
168
a0  0,
a1  1,
a2 
1
,
2
于是得
1 2
x
6 x  3x 2
2

.
R22 ( x) 
2
1
6  6x  x
1 x  x2
6
x
当 n  m  3时,由(5.11)得
1
1
1


b

b

b


,
3
2
1

2
3
4

1
1
1
 1
b

b

b

,

3
2
1
3
4
5
 2
 1 b  1 b  1 b   1 ,
3
2
1

4
5
6
 3
169
解得
b1 
3
,
2
b2 
3
,
5
b3 
1
.
20
代入(5.9)得
a0  0,
a1  1,
a2  1,
11
a3 
.
60
于是得
11 3
xx 
x
60 x  60 x 2  11x 3
60
R33 ( x ) 

.
2
3
3
3 2
1 3
60  90 x  36 x  3 x
1 x  x 
x
2
5
20
2
170
可以看到这里得到的 R22 ( x) 及 R33 ( x) 与 ln( 1  x) 的前面
连分式展开得到的有理逼近(5.4)结果一样.
为了求帕德逼近 Rnm (x) 的误差估计,由(5.9)及(5.11)
求得的 Pn ( x), Qm ( x)系数 a0 , a1 ,, an 及 b0 , b1 ,, bm ,直
接代入则得
f ( x)Qm ( x)  Pn ( x)  x
n  m 1

m
( bk cn  m 1l  k )x l ,
l 0 k 0
将 Qm (x)除上式两端,即得
171
x
n  m 1
f ( x)  Rnm ( x) 

l
r
x
l
l 0
(5.13)
,
Qm ( x)
m
其中 rl   bk cn  m l 1 k .
k 0
当 x  1时可得误差近似表达式
f ( x)  Rnm ( x)  r0 x
n  m 1
m
,
r0   bk cn  m 1 k .
k 0
172
3.6
三角多项式逼近与快速傅里叶变换
当模型数据具有周期性时,用三角函数特别是正弦和余
弦函数作为基函数是合适的,这时前面讨论的用多项式、分
段多项式或有理函数作基函数都是不适合的.
用正弦函数和余弦函数级数表示函数称为傅里叶变换
(简称傅氏变换).
173
3.6.1
最佳平方三角逼近与三角插值
设 f (x) 是以 2π
1
S n ( x)  a0  a1 cos x  b1 sin x    an cos nx  bn sin nx
2
(6.1)
做最佳平方逼近函数. 由于三角函数族
1,cos x,sin x, ,cos kx,sin kx,
在 [0,2 π] 上是正交函数族,于是 f (x) 在 [0,2 π] 上的最小
平方三角逼近多项式 S n (x)
174
1
ak 
π
1
bk 
π

2
0
f ( x) cos kxdx
(k  0,1, , n),
(6.2)

2
0
f ( x) sin kxdx
(k  0,1, , n),
ak , bk 称为傅里叶系数.
函数 f (x) 按傅里叶系数展开得到的级数

1
a0   (ak cos kx  bk sin kx)
2
k 1
(6.3)
就称为傅里叶级数.
175
只要 f (x )在 [0,2 π] 上分段连续,则级数(6.3)一致
收敛到 f (x).
对于最佳平方逼近多项式(6.1)有

1
a0   (ak cos2 kx  bk sin
kx)
2
2 f ( x)k
1 S ( x )
 f ( x)  S
n
2
2
2
n
(6.3)
( x) 2 .
1
由此可以得到相应于(3.11)的贝塞尔不等式
S n ( x)  a0  a1 cos x  b1 sin x    an cos nx  bn sin nx
2
n
(6.1)
1 2
1 2
2
2
2
a0   (ak  bk ) 

[ f ( x)] dx.
2
*
2
(
a

(
x
)
)
 f ( x) 2 .

k
k
2
因为右边不依赖于
k 1
n ,左边单调有界,所以级数
n
2
k 1
0
(3.11)
176

1 2
a0   ( ak2  bk2 )
2
k 1
收敛,并有
lim ak  lim bk  0.
k 
k 
当 f (x) 只在给定的离散点集
2π


x

j
,
j

0
,
1
,

,
N

1
 j

N


上已知时,则可类似得到离散点集正交性与相应的离散傅
里叶系数.
下面只给出奇数个点的情形.
177
令
xj 
2 πj
2m  1
( j  0,1,,2m),
可以证明对任何 0  k , l  m 成立
 2m
0,


 sin l x j sin kx j   2m  1
,
 j 0

 2

0,

 2m  1
 2m
 cos l x j cos kx j  
 j 0
 2
2m  1,

 2m
0  k,
 cos l x j sin kx j  0,
 j 0

l  k , l  k  0;
l  k  0;
l  k;
l  k  0;
l  k  0;
j  m.
178
这表明函数族 {1, cos x,sin x, ,cos mx,sin mx} 在点集
{x j 
2 πj
}
2m  1
上正交.
若令 f j  f ( x j )
( j  0,1,,2m), 则 f (x) 的最小二
乘三角逼近为
n
1
S n ( x)  a0   (ak cos kx  bk sin kx),
2
k 1
n  m,
其中
2m
2
2πjk
ak 
f
cos
 j
2m  1 j  0
2m  1
(k  0,1, , m),
179
2m
2
2 πjk
bk 
f
sin
 j 2m  1
2m  1 j  0
(k  1, , n).
(6.4)
当 n  m时
Sm ( x j )  f j
( j  0,1, ,2m) ,
于是
m
1
S m ( x)  a0   (ak cos kx  bk sin kx)
2
k 1
就是三角插值多项式,系数仍由(6.4)表示.
180
一般情形,假定 f (x) 是以 2π为周期的复函数,给定
在 N个等分点 x j 
2π
 2π
j ( j  0,1, , N  1) 上的值 f j  f 
N
N

j ,

由于
eijx  cos( jx)  i sin( jx) ( j  0,1,, N  1, i   1),
所以函数族 {1, eix ,, ei ( N 1) x }在区间 [0,2 π] 上是正交的.
函数 e ijx 在等距点集 xk  2π k (k  0,1, , N  1) 上的值
N
e ijxk 组成的向量记作
 j  (1, e
ij
2π
N
,, e
ij
2π
( N 1)
N
)T .
181
当 j  0,1,, N  1 时,N 个复向量 0 , 1 , ,  N 1 具有
如下正交性:
N 1
(l , s )   e
il
2π
k
N
e
is
2π
k
N
k 0
N 1
 e
i(l  s )
2π
k
N
k 0
 0, l  s;

 N , l  s.
(6.5)
182
事实上,令 r  e
i (l  s )
0  l  N  1,
2π
N
, 若 0  l , s  N  1, 则有
 ( N  1)  s  0,
于是
 ( N  1)  l  s  N  1,
即
N 1 l  s
N 1
1  


 1;
N
N
N
若 l  s  0, 则 r  1, 从而
r N  ei ( l s ) 2 π  1;
183
于是
(l , s ) 
N 1
r
k 0
k
1 r N

1 r
 0.
若 l  s, 则 r  1, 于是
(s , s ) 
N 1
k
r
  N.
k 0
这就证明了(6.5)成立. 即 0 , 1 , ,  N 1 是正交的.
2π
因此, f (x)在 N个点 {x j 
j ( j  0,1, , N  1)} 上的
N
 0, l  s;
最小二乘傅里叶逼近为
(l ,s )  
(6.5)
N
,
l

s
.

184
S ( x) 
n 1
ikx
c
e
 k ,
n  N,
(6.6)
k 0
其中
1
ck 
N
N 1
f
j 0
j
e
-ikj
2π
N
,
(k  0,1, , n  1). (6.7)
在(6.6)中,若 n  N ,则 S (x )为 f (x) 在点
x j ( j  0,1, , N  1) 上的插值函数, 即 S ( x j )  f ( x j ),
于是由(6.6)得
N 1
f j   ck e
k 0
ik
2π
j
N
,
( j  0,1, , N  1).
(6.8)
185
(6.7)是由 { f j } 求 {ck } 的过程,称为 f (x) 的离散
傅里叶变换. 简称DFT,
而(6.8)是由 {ck }求 { f j } 的过程,称为反变换.
1
ck 
N
N 1
f
j 0
N 1
f j   ck e
j
ik
e
-ikj
2π
j
N
2π
N
,
,
(k  0,1, , n  1). (6.7)
( j  0,1, , N  1).
(6.8)
k 0
186
3.6.2
N 点DFT与FFT算法
不论是按(6.7)式由 { f j } 求 {ck } ,或是按(6.8)
由 {ck } 求 { f j },还是由(6.4)计算傅里叶逼近系数 ak , bk ,
都可归结为计算
N 1
2 kj 2 m
2 πjk
c j ak  xk  N
, f (j cos
j  0,1, , N  1),
2m 1 (6.4)
k  0 2m 1 j 0

(6.9)
2m
其中 {xk }0N 1为已知的输入数据,
2
2{
πcjkj }0N 1为输出数据,而
bk 
N e
i 2N
f

2m 1
j
sin
2m 1
2j0
2
 cos
 isin
, i   1.
N
N
187
(6.9)式称为 N点DFT,表面上看计算 c j 需要 N次复
数乘法和 N次复数加法,称为 N个操作,计算全部 c j 共要 N 2
个操作.
N 1
(6.9)
c  当 Nx 较大且处理数据很多时,就是用高速的电子计算
wkj , ( j  0,1, , N  1),
j

k 0
k
机,很多实际问题仍然无法计算,
直到20世纪60年代中期产生了FFT算法,大大提高了
运算速度,才使DFT得到更广泛的应用.
FFT是快速算法的一个典范,其基本思想就是尽量减
少乘法次数.
188
对于任意正整数 k, j ,成立
 Nj   Nk   Nj  k ,  NjN  k   Nk (周期性),
jk
 Njk  N / 2   Njk (对称性),  jN
  Nk .
由周期性可知所有  Njk ( j, k  0,1,, N  1) 中,做多有 N个不
同的值  N0 , 1N ,,  NN 1.
特别地,有
 N0   NN  1,  NN / 2  1.
当 N  2 p 时,  Njk只有 N / 2个不同的值.
利用这些性质可将(6.9)式对半折成两个和式,再将
对应项相加,有
cj 
N / 2 1
x 
k 0
k
jk
N

N / 2 1
x
k 0
N / 2 k 
j ( N / 2 k )
N

N / 2 1
[ x
k 0
k
 (1) j xN / 2 k ] Njk .189
依下标奇、偶分别考察,则
c2 j 
N / 2 1
 (x
c2 j 1 
若令
k 0
k
 x N / 2 k ) Njk/ 2 ,
N / 2 1
jk
(
x

x
)

 k N / 2 k N / 2 .
k 0
yk  ( xk  xN / 2 k ),
y N / 2 k  ( xk  xN / 2 k ) Nk ,
则可将 N 点DFT归结为两个 N / 2点的DFT:
N / 2 1

jk
c

y


2
j
k
N
/2,


k 0

N / 2 1
jk
c

y


2 j 1
N / 2 k N / 2 .

k 0

j  0,1, , N / 2  1,
如此反复施行二分手续即可得到FFT算法.
190
以
N  23 为例说明FFT的计算方法,此时 k , j  0,1,, N  1  7
在(6.9)中将  N  8 记为  ,于是(6.9)式的和为
cj 
7
x 
k 0
k
jk
,
( j  0,1, ,7).
(6.10)
将 k, j 用二进制表示为
k  k 2 2 2  k1 21  k0 20  (k 2 k1k0 ),
j  j2 2 2  j1 21  j0 20  ( j2 j1 j0 );
其中 kr , jr (r  0,1,2) 只能取0或1.
例如 6  22  22  0  20  (110).
根据 k, j 表示法,有
c j  c( j2 j1 j0 ),
xk  x(k 2 k1k0 ).
191
公式(6.10)可表示为
7
c j   xk  kj , (6.10)
2
1
0
c( j2 j1 j0 )     x(k 2kk10k0 ) ( k 2 k1k0 )( j2 2  j1 2  j0 2 )
1
1
1
(6.11)
k 0  0 k1  0 k 2  0

 j1 ( k1k0 0 ) 
 1  1

j0 ( k 2 k1k 0 )
j ( k 00)
      x(k 2 k1k0 )

  2 0


k0 0 

k1 0  k 2 0

1
若引入记号
 A0 (k2 k1k0 )  x(k2 k1k0 ),

1
 A1 (k1k0 j0 )   A0 (k 2 k1k0 ) w j0 ( k 2 k1k0 ) ,

k2 0

1

j1 ( k1k 0 0 )
A
(
k
j
j
)

A
(
k
k
j
)
w
,

1
1 0 0
 2 0 1 0
k1  0

1
j 2 ( k 0 00)
A ( j j j ) 
A
(
k
j
j
)
w
.

3
2 1 0
2
0 1 0

k0 0

(6.12)
192
则(6.11)变成
c( j2 j1 j0 )  A3 ( j2 j1 j0 ).
若注意  j 2   j N / 2  (1) j , 公式(6.12)还可进一步
p1
0
0
0
简化为
A1 (k1k0 j0 ) 
1
 A (k
k 2 0
0
2
k1k0 ) j0 ( k 2 k1k0 )
 A0 (0k1k0 )
j0 ( 0 k1k 0 )
 A0 (1k1k0 )
j0 2 2
j
0 ( 0 k1k 0 )
 [ A0 (0k1k0 )  (1) j0 A0 (1k1k0 )] j0 ( 0 k1k0 ) ,
A1 (k1k0 0)  A0 (0k1k0 )  A0 (1k1k0 ),
A1 (k1k01)  [ A0 (0k1k0 )  A0 (1k1k0 )] ( 0 k1k0 ) .
193
将这表达式中二进制表示还原为十进制表示:
k  (0k1k0 )  k1 21  k0 20 ,
即 k  0,1,2,3, 得
2
 A1 (2k )  A0 (k )  A0 (k  2 ),


2
k
 A1 ( 2k  1)  [ A0 ( k )  A0 ( k  2 )]w


( k  0,1,2,3).

(6.13)
同样(6.12)中的 A2 也可简化为
A2 (k0 j1 j0 )  [ A1 (0k0 j0 )  (1) j1 A1 (1k0 j0 )] j1 ( 0 k0 0 ) ,
即
A2 (k0 0 j0 )  A1 (0k0 j0 )  A1 (1k0 j0 ),
A2 (k01 j0 )  [ A1 (0k0 j0 )  A0 (1k0 j0 )] ( 0 k0 0 ) .
194
把二进制表示还原为十进制表示,得
 A2 (k 22  j )  A1 (2k  j )  A1 (2k  j  22 ),


2
2
2k
 A2 ( k 2  j  2)  [ A1 ( 2k  j )  A1 ( 2k  j  2 )]w


(6.14)
( k  0,1; j  0,1).

同理(6.12)中 A3 可简化为
A3 ( j2 j1 j0 )  A2 (0 j1 j0 )  (1) j2 A2 (1 j1 j0 ),
即
A3 (0 j1 j0 )  A2 (0 j1 j0 )  A2 (1 j1 j0 ),
A3 (1 j1 j0 )  A2 (0 j1 j0 )  A2 (1 j1 j0 ).
195
表示为十进制,有
2
A
(
j
)

A
(
j
)

A
(
j

2
),
3
2
2



2
2
 A3 ( j  2 )  A2 ( j )  A2 ( j  2 )


( j  0,1,2,3).

(6.15)
根据公式(6.13),(6.14),(6.15),由
A0 (k )  x(k )  xk
(k  0,1, ,7)
逐次计算到 A3 ( j )  c j ( j  0,1, ,7), 见表3-5(略).
从表3-5中看到计算全部8个 c j只用8次乘法运算和24次
加法运算.
196
上面推导的 N  23 的计算公式可类似地推广到 N  2 p
的情形.
根据公式(6.13),(6.14),(6.15),一般情况的FFT
计算公式如下:
q
q 1  j )  A
q 1  j  2 p 1 ),
q 1 (k 2
 Aq (k 2  j )  Aq 1 (k 2

 Aq (k 2q  j  2q 1 )
(6.16)


q 1

  [ Aq 1 (k 2 q 1  j )  Aq 1 ( k 2 q 1  j  2 p 1 )]wk 2 ,
其中 q  1,, p; k  0,1,,2 p q  1; j  0,1,,2q 1  1.
Aq 括号内的数代表它的位置,在计算机中代表存放数的地址.
197
一组 Aq占用 N个复数单元,计算时需给出两组单元,
从 A0 (m)( m  0,1,, N  1) 出发, q 由 1到 p算到
Ap ( j )  c j ( j  0,1, , N  1), 即为所求.
计算过程中只要按地址号存放 Aq , 则最后得到的 Ap ( j )
就是所求离散频谱的次序.
这个计算公式除了具有不倒地址的优点外,计算只有两
重循环.
p q
 1,
外循环 q由 1计算到 p,内循环 k由 0计算到 2
j 由 0计算到 2 q 1  1, 更重要的是整个计算过程省计算量.
198
由公式看到算一个 Aq 共做 2 p q 2q /1  N / 2次复数乘法,
而最后一步计算 Ap 时,由于

k 2 p1
 ( N / 2 ) k  (1) k  (1) 0  1
(注意 q  p 时 2 p q  1  0, 故 k  0 ),因此,总共要算
( p  1) N / 2 次复数乘法,它比直接用(6.9)需 N 2次乘法
快得多,计算量比值是 N : ( p  1) / 2.
当 N  210时比值是 1024 : 4.5  230 : 1, 它比一般FFT的
计算量( pN次乘法)也快一倍.
我们称(6.16)的计算公式为改进的FFT算法 .
199
例13
设 f ( x)  x 4  3x 3  2 x 2  tan x( x  2) .给定
数据 {x j , f ( x j )}7j 0 , x j 
解
j
确定三角插值多项式.
4
先将 [0,2] 变换为 [ ,  ] ,可令 y j   ( x j  1) ,
故输入数据为 { y j , f j )}70 ,
f j  f (1 
yj

).
给定8个点可确定8个参数的4次三角插值多项式
3
1
S 4 ( y )  a0   (ak cos ky  bk sin ky)  a4 cos 4 y (6.17)
2
k 1
这里
2 7

2
a

f
cos
k  0,1, ,4,

k
j
8 kj,

8 j 0


7
2
b 
f j sin 28 kj, k  1,2,3,

k

8 j 0

(6.18)
200
与(6.10)式比较可先计算
7
ck   f j jk ,
j 0
这里 { f j }70代替(6.10)式的 {x j }70,
 e
i 82 
e
i 4
 cos

4
 isin

4
.
对每个 k  0,1,,4有
7
7
 kj
1
1
1
1
i
ik (   4 j )
k
 ik
 ik
4
ck ( 1)  ck e
  f je e
  f je
4
4
4 j 0
4 j 0
1 7
j
j
  f j (cosk (   )  isin k (   ))
4 j 0
4
4
1 7
  f j (cosky j  isin ky j ),
4 j 0
201
所以
(1) k
1
ak  ibk 
ck  ck e ik ,
4
4
即
1
ak  Re( ck e ik ),
4
1
bk  Im( ck e ik ).
4
显然 b0  b4  0 ,用FFT算法求出 ck (k  0,1,2,3,4),也就
得到(6.18)式的系数,从而得到(6.17)式的4次三角
插值多项式
S 4 ( y )  0.761979  0.771841cos y
 0.386374 sin y  0.0173037 cos 2 y
 0.0468750 sin 2 y  0.00686304 cos 3 y
 0.0113738 sin 3 y  0.000578545 cos 4 y.
202
将 y   ( x  1) 代入 S4 ( y) 可以得到 [0,2]上的三角多项式
S 4 ( x).
图3-6给出了 y  f (x) 及 y  S4 ( x) 的图形.
表3-6给出了在点 x j  0.125  j 0.25( j  0,1,,7) 处 f ( x j )
与 S 4 ( x j ) 的值.
203
图3-6
随机数据点上的二元拟合
204
205
206
207
208
209
210
211
212
213