关于Toeplitz矩阵的计算_21_25

j =0
n −1
第二步，计算 µ k = λ ， （ k = 1,2, " , n − 1 ） ；
−1 k
第 三 步 ， 求 出 循 环 Toeplitz 矩 阵 C 的 逆 矩 阵 C −1 = circ(η 0 ,η 1 , "η n −1 ) ，
η k = ∑ µ j ω − kj ， （ k = 1,2, " , n − 1 ） 。
x 0 = 0 ， x1 = 1 ， y 0 = 0 ， y1 = 1 1 x k +1 = − (αx k + γx k −1 ) ， （ k = 1,2, " , n − 1 ）
β
y k +1 = −
1
γ
(αy k + β y k −1 ) ， （ k = 1,2, " , n − 1 ）
p = c T Jx = γx n −1 + αx n ， q = y T Jr = β y n −1 + αy n
15
第三章 特殊 Toeplitz 矩阵的逆矩阵的研究
特殊的 Toeplitz 矩阵有广泛的应用前景，特别在生物学、物理学、数学、社会 科学中的许多问题都和 Toeplitz 矩阵的理论有着密切的关系。在数理统计、石油、 地震物探及其它应用科学中，尤其是图象、数字信号处理中，常会遇到循环矩阵 这类特殊矩阵。带状 Toeplitz 方程组广泛应用于科学计算和工程计算中，尤其是在 双曲偏微分方程的数值求解中，随着科学技术的不断发展，问题的规模越来越大， 适应于大规模计算的并行算法也被陆续提出来。因此，对于这类特殊的 Toeplitz 矩 阵的研究是有必要的。
⎞ ⎛β ⎛γ ⎟ ⎜ ⎜ ⎟ ⎜α β ⎜ ⎟ ⎜ ，N =⎜ M= γ α β ⎟ ⎜ ⎜ % % % ⎟ ⎜ ⎜ ⎟ ⎜ ⎜ γ α β⎠ ⎝ ⎝
α β ⎞ ⎟ γ α % ⎟ γ % β⎟ ⎟ % α⎟ γ⎟ ⎠
分别是 n − 1 阶非奇异下三角 Toeplitz 矩阵和上三角 Toeplitz 矩阵。令 x1 = 1 ， x = ( x 2 , " x n ) T = − M −1 c ；
因此，公式 T −1 = T1U 1 + T2U 2 是向前稳定的。
2.4.3 算法分析
现在我们来考虑一下算法，对于公式 T −1 = T1U 1 + T2U 2 ，关键就是要解线性方 程组 Tx = f 和 Ty = e1 。 如果用 Zohar 算法求解 Toeplitz 线性方程组 Tx = f 和 Ty = e1 ， 则要 8n 2 − 4n − 4 次乘除法运算， 8n 2 − 14n + 6 次加减法运算；如果用快速傅里叶 （FFT）方法，则运算复杂性为 O(n log n) 。这些方法将在第四章里做详细的阐述。 因此， 如果用快速傅里叶方法， 得到分解式 T −1 = T1U 1 + T2U 2 只需要 O(n log n) 的运 算量。
j =0
n −1
如果第一步和第三步借助于离散富氏变换，则整个算法的计算量为
17
O(n log n) 。
下面简单地讨论一下 r-循环 Toeplitz 矩阵的逆矩阵。 设形如 n 阶矩阵称为 r-循环 Toeplitz 矩阵为 ξ 1 " ξ n −1 ⎤ ⎡ ξ0 ⎢ rξ ξ 0 " ξ n−2 ⎥ ⎥ C (r ) = ⎢ n −1 ⎥ ⎢ # # % ⎥ ⎢ ⎣ r ξ 1 rξ 2 " ξ 0 ⎦ 特别地，当 r = 1 时，即得循环 Toeplitz 矩阵；当 r = −1 时，称之为斜循环 Toeplitz 矩阵（或反循环 Toeplitz 矩阵） ；当 r = 0 时，即得上三角 Toeplitz 矩阵。由 r-循环
利用这两个式子，可得 ⎛ I ⎜ ⎜ xT J ⎝
0 ⎞⎛ c M ⎞⎛ 1 0 T ⎞ ⎛ 0 M ⎞ ⎟=⎜ ⎟⎜ ⎟ ⎟， T T ⎟ ⎜ ⎟⎜ ⎜ ⎟ ⎜ 1⎟ ⎠⎝ c J ⎠⎝ x I ⎠ ⎝ p 0 ⎠ ⎛ 1 y T ⎞⎛ r T 0 ⎞⎛ I Jy ⎞ ⎛ 0 T q ⎞ ⎜ ⎟⎜ ⎟⎜ T ⎟ ⎟ ⎟=⎜ ⎜ 0 I ⎟⎜ N Jr ⎟⎜ ⎜ N 0⎟ 。 0 I ⎝ ⎠ ⎝ ⎠⎝ ⎠ ⎝ ⎠ 从而 Tn 可逆的充要条件是 p ≠ 0 和 q ≠ 0 ，且 ⎛ 0T 0⎞ 1 ⎛ 1⎞ t ⎟ T =⎜ ⎜ x⎟ ⎟x J 1 ⎜ M −1 0 ⎟ + p ⎜ ⎝ ⎠ ⎝ ⎠ ⎛ 0 N −1 ⎞ 1 ⎛ Jy ⎞ T ⎟ Tn−1 = ⎜ ⎜1⎟ ⎟1 y ⎜ 0 0T ⎟ + q ⎜ ⎝ ⎠ ⎝ ⎠ −1 −1 因为矩阵 M 和矩阵 N 分别是下三角和上三角矩阵， 由上面的式子可知， Tn−1
18
y1 = 1 ， y = ( y 2 , " y n ) T = −r T N −1 ；
p = c T Jx = γx n −1 + αx n ， q = y T Jr = β y n −1 + αy n 。
于是，由矩阵 M 与 N 的次对称性得
x T J = −c T M −T J = −c T JM −1 ， Jy = − JN −T r = − N −1 Jr 。
(4).循环 Toeplitz 矩阵 C = ∑ ξ j K J ；
j =0
n −1
(5). CK = KC 。 定理 3.1.1 矩阵 C 为 n 阶循环 Toeplitz 矩阵的充要条件是存在数 λ 0 , λ1 , " λ n −1
使
F −1CF = diag (λ 0 , λ1 , " λ n −1 ) 。
， （ i = − 1 ）是 n 次单位根。易知矩阵 F 是酉矩阵，即 F −1 = F T ，且
F −1 KF = diag (1, ω , ω 2 , " , ω n −1 ) ，
可见， ω j ， （ j = 1,2, " , n − 1 ）是矩阵 K 的 n 个特征值，而 F 的各列是对应的 单位特征正交特征向量；
j
K n = I n （规定 K 0 = I n ） ； (3). 取 n 阶 Fourier 变换矩阵 F ：
16
1 ⎛1 1 ⎜ ω2 ⎜1 ω F = ⎜1 ω 2 ω4 ⎜ # # ⎜# ⎜1 ω n −1 ω 2 ( n −1) ⎝
其中 ω = e 有
− 2πi n
" 1 ⎞ ⎟ n −1 ω " ⎟ 2 ( n −1) ⎟ " ω ⎟ % # ⎟ ( n −1)( n −1) ⎟ " ω ⎠
分析：经过简单计算可知，
λ k = ∑ ξ j ω kj ， （ k = 1,2, " , n − 1 ） 。
j =0
n −1
通过这个定理，可以知道循环 Toeplitz 矩阵的逆矩阵仍然为循环 Toeplitz 矩阵，而 且
−1 −1 1 −1 C −1 = Fdiag (λ 0 , λ1 , " λ− ， n −1 ) F
Toeplitz 矩阵的定义和性质可知， r-循环 Toeplitz 矩阵的逆矩阵仍然是 r-循环 Toeplitz
矩阵。因此，求 r-循环 Toeplitz 矩阵的逆矩阵的算法同循环 Toeplitz 矩阵的算法类 似。
3.２ 三对角 Toeplitz 矩阵的逆矩阵
在构造三次样条函数及用差分法求解微分方程等问题中，均须要求解以三对 角 Toeplitz 矩阵为系数矩阵的线性方程组。 考虑如下的 n 阶三对角 Toeplitz 矩阵： ⎞ ⎛α β ⎟ ⎜ ⎟ ⎜γ α % Tn = ⎜ ， βγ ≠ 0 % % β⎟ ⎟ ⎜ ⎜ γ α⎟ ⎠ ⎝ 首先，将 n 阶三对角 Toeplitz 矩阵写成分块形式 ⎛c M ⎞ ⎛rT 0 ⎞ ⎟，n ≥ 3 Tn = ⎜ ⎜ cT J ⎟ ⎟=⎜ ⎜ ⎟ ⎝ ⎠ ⎝ N Jr ⎠ 其中， c = (α , λ ,0, " 0) T ， r = (α , β ,0, " 0) T ， J = (e1 , e 2 , " , e n ) ，而
−1 n
(
)
(
)
的上三角部分由矩阵
1 ⎛1⎞ t ⎜ ⎟ ( x J 1) 的上三角部分给出，而 Tn−1 的下三角部分由矩阵 ⎜ ⎟ p ⎝ x⎠
1 ⎛1⎞ t ⎜ ⎟ ⎟(x J 1) 的下三角部分给出。 p⎜ ⎝ x⎠
综上所述，可得求三对角三角 Toeplitz 矩阵 Tn 的逆矩阵 Tn−1 = (v ij ) in, j =1 的算法：
3.1 循环 Toeplitz 矩阵的逆矩阵
在数理统计、石油、地震物探及其它应用科学中，尤其是在图象、数字信号 处理中，常会遇到循环矩阵这类特殊矩阵。 设 n 阶方阵 C 称为循环 Toeplitz 矩阵为 ⎡ ξ 0 ξ 1 " ξ n −1 ⎤ ⎢ξ ξ 0 " ξ n−2 ⎥ ⎥ C = ⎢ n −1 ⎥ ⎢ # # % ⎥ ⎢ ⎣ ξ1 ξ 2 " ξ 0 ⎦ 显然， 矩阵 C 由其首行元素唯一确定， 简记为 C = circ(ξ 0 , ξ1 ,"ξ n −1 ) ， 特别地，n 阶
2 2
由 Ty = e1 和 Tx = f ，可得 y 2 ≤ T 和 x2≤ T 2 ~ T −1 − T −1 ~ + nε )(1 + 2 T −1 f 2 ≤ n(2ε −1 2 T
2
−1
−1 2
f
2
，因此，相对误差为
2
)+ε n 。
因为 T 是良态的，所以 T −1 是有界的。而 f
2
2
也是有界的。
Toeplitz 矩阵 ⎛ 0 I n −1 ⎞ K = circ(0,1,0, " 0) = ⎜ ⎜1 0 ⎟ ⎟ ⎝ ⎠