Sparse Matrix
Sparse Matrix
「稀疏矩陣」。矩陣元素幾乎都是零。非零元素相當稀少。
矩陣的資料結構=圖的資料結構。稠密矩陣適合Adjacency Matrix,稀疏矩陣適合Adjacency Lists。
實務上採用Adjacency Lists與Edge List的融合版本,化作三條陣列,稱作Compressed Sparse Row。MATLAB即採用此方式。
Sparse Eigenvalue
Sparse Eigenvalue
經典的矩陣運算,諸如反矩陣、QR分解、特徵分解,時間複雜度O(N³)以上,速度太慢不實用。因此數學家旁敲側擊,鑽研稀疏矩陣,精簡計算步驟,以投入實用。
專著《Numerical Methods for Large Eigenvalue Problems》。
幾乎上三角分解演算法(Restarted Arnoldi Iteration)
一口氣求得大量特徵值,其中一種策略是:先變成幾乎上三角矩陣、再變成上三角矩陣。
幾乎上三角矩陣分解,三種演算法Arnoldi Iteration、Householder Reflection、Jacobi Rotation,其中只有Arnoldi Iteration適用於稀疏矩陣。
Arnoldi Iteration的時間複雜度,分成兩個部分。
一、矩陣求值:稀疏矩陣的情況下,需時較少。
二、向量扳正:無法保持稀疏,需時相同。
數學家想到了歪招:restart。被扳正的向量,每回合越來越多,計算量越來越大。因此每隔幾回合就砍掉重練,初始向量換成新的再來一遍,直到總共湊滿N個向量。
restart也能解決「初始向量不足N種特徵分量,無法生成整個空間」的情況。
Large Scale Eigenvalue
Large Scale Eigenvalue
採用平行計算。
特徵向量暨特徵值修正演算法(Olsen's Method)
特徵向量暨特徵值,寫成方程組。
{ Ax = λx { ‖x‖² = 1
移項整理,求解改成求根。
{ (A-λI)x = 0 { ½xᵀx - 1 = 0 為了方便微分,特徵向量長度改成2。
牛頓法求根。
xnext = x - F′(x)ᵀ⁺ F(x) [ xnext ] = [ x ] - [ df₁/dx df₂/dx ]ᵀ⁻¹[ f₁(x,λ) ] [ λnext ] [ λ ] [ df₁/dλ df₂/dλ ] [ f₂(x,λ) ] [ xnext ] = [ x ] - [ (A-λI)ᵀ x ]ᵀ⁻¹[ r ] [ λnext ] [ λ ] [ -xᵀ 0 ] [ 0 ] [ xnext ] = [ x ] - [ A-λI -x ]⁻¹[ r ] where r = (A-λI)x = Ax - λx [ λnext ] [ λ ] [ xᵀ 0 ] [ 0 ]
簡易方式,結果相同。
A(x + Δx) = (λ + Δλ)(x + Δx) (A-λI)Δx - xΔλ = r where r = Ax - λx and ΔλΔx ≈ 0 [ A-λI -x ] [ Δx ] = [ r ] where Δx = xnext - x [ xᵀ 0 ] [ Δλ ] [ 0 ] Δλ = λnext - λ
每回合都要矩陣求解,時間複雜度相當高,似乎沒有用處。
特徵向量暨特徵值修正演算法(Rayleigh Quotient Iteration)
Rayleigh Quotient Iteration是牛頓法的變種。
A(x + Δx) = (λ + Δλ)(x + Δx) A(x + Δx) = λ(x + Δx) + Δλ(x + Δx) (A - λI)(x + Δx) = Δλ(x + Δx) (A - λI)(x + Δx) = Δλx where ΔλΔx ≈ 0 (x + Δx) = (A - λI)⁻¹Δλx xnext = Δλ(A - λI)⁻¹x let λ = xᵀAx / xᵀx is Rayleigh quotient
「向量長度縮放倍率Δλ」改成了「向量長度縮放為1」。
「λnext = λ + Δλ」改成了「λnext是虛擬特徵值」。
每回合都要矩陣求解,時間複雜度相當高,似乎沒有用處。
特徵向量暨特徵值修正演算法(Jacobi–Davidson Method)
追加垂直限制。
A(x + Δx) = (λ + Δλ)(x + Δx) where x ⟂ Δx
移項整理。
(A-λI)Δx - Δλx = r where x ⟂ Δx, r = Ax - λx
利用投影矩陣,消掉-Δλx。
let P = I - xxᵀ such that Px = 0 and Pr = r
P((A-λI)Δx - Δλx) = Pr P(A-λI)Δx - PΔλx = Pr P(A-λI)Δx = r
補上右邊,成為相似變換,保留特徵值。
Pᵀ(A-λI)PΔx = r (P = Pᵀ and PΔx = Δx)
問題變成一次方程組求解。求得特徵向量誤差Δx。
BΔx = r where B = Pᵀ(A-λI)P and r = (A-λI)x = Ax - λx
每回合都要矩陣求解,時間複雜度相當高,似乎沒有用處。