Sparse Linear Algebra

189 Pages·2015·2.46 MB·English

Sparse Linear Algebra: Direct Methods, advanced features
P. Amestoy and A. Buttari
(INPT(ENSEEIHT)-IRIT)
A. Guermouche (Univ. Bordeaux-LaBRI),
J.-Y. L’Excellent and B. U¸car (INRIA-CNRS/LIP-ENS Lyon)
F.-H. Rouet (Lawrence Berkeley National Laboratory)
C. Weisbecker INPT(ENSEEIHT)-IRIT
2013-2014
Sparse Matrix Factorizations
m×m T
• A ∈ ℜ , symmetric positive definite → LL = A
Ax = b
m×m T
• A ∈ ℜ , symmetric → LDL = A
Ax = b
m×m
• A ∈ ℜ , unsymmetric → LU = A
Ax = b
m×n
• A ∈ ℜ , m ≠ n → QR = A
min ‖Ax − b‖ if m > n
x
min ‖x‖ such that Ax = b if n > m
Cholesky on a dense matrix
 
a11

 a21 a22 
 a
31 a32 a33 
a41 a42 a43 a44
left-looking Cholesky right-looking Cholesky
for k = 1, ..., n do for k = 1, ..., n do
√
for i = k, ..., n do (k−1)
for j = 1, ..., k − 1 do lkk = akk
for i = k + 1, ..., n do
(k) (k)
endaifkor = aik − lijlkj lik = a(ik−1)/lkk
for j = k + 1, ..., i do
end fo√r (k) (k)
(k−1) aij = aij − likljk
lkk = akk
end for
for i = k + 1, ..., n do
end for
Left−looking (k−1) Right−looking
lik = aik /lkk end for
end for
end for used for modification
modified
138 JOSEPH W. H. LIU
a
d
fa
o
d
I
o o h
0
o oo oj
FIG. 2.1. An example of matrix structures.
Factorization of sparse matrices: problems
The factorization of a sparse matrix is problematic due to the
3 10
presence of fill-in.
The basic LU step:
(k) (k)
8 a a
(k+1) (k) i,k k,j 4
a = a −
i,j i,j (k)
a
k,k
(k) (k+1)
Even if a is null, a can be a nonzero
i,j i,j
G(A
G(F)
G(Ft)= T(A
FIG. 2.2. Graph structures of the example in Fig. 2.1.
• the factorization is more expensive than O(nz)
• higher amount of memory required than O(nz)
• more complicated algorithms to achieve the factorization
Factorization of sparse matrices: problems
Because of the fill-in, the matrix factorization is commonly preceded
by an analysis phase where:
• the fill-in is reduced through matrix permutations
• the fill-in is predicted though the use of (e.g.) elimination graphs
• computations are structured using elimination trees
The analysis must complete much faster than the actual factorization
Analysis
Fill-reducing ordering methods
Three main classes of methods for minimizing fill-in during
factorization
Local approaches : At each step of the factorization, selection of the
pivot that is likely to minimize fill-in.
• Method is characterized by the way pivots are
selected.
• Markowitz criterion (for a general matrix)
• Minimum degree (for symmetric matrices)
Global approaches : The matrix is permuted so as to confine the
fill-in within certain parts of the permuted matrix
• Cuthill-McKee, Reverse Cuthill-McKee
• Nested dissection
Hybrid approaches : First permute the matrix globally to confine the
fill-in, then reorder small parts using local heuristics.
The elimination process in the graphs
GU(V,E) ← undirected graph of A
for k = 1 : n − 1 do
V ← V − {k} {remove vertex k}
E ← E − {(k, ℓ) : ℓ ∈ adj(k)} ∪ {(x, y) : x ∈ adj(k) and y ∈
adj(k)}
Gk ← (V,E) {for definition}
end for
Gk are the so-called elimination graphs (Parter,’61).
1
1 2 3
2
3
4
G0 : H0 = 4
5
6 5
6
138 JOSEPH W. H. LIU
a
d
fa
o
d
I
o o h
0
o oo oj
FIG. 2.1. An example of matrix structures.
The Elimination Tree
The elimination tree is a graph that expresses all the dependencies in
3 10
the factorization in a compact way.
Elimination tree: basic idea
If i depends on j (i.e., lij ≠ 0, i > j)8 and k depends on i (i.e.,
4
lki ≠ 0, k > i), then k depends on j for the transitive relation and
therefore, the direct dependency expressed by lkj ≠ 0, if it exists, is
redundant
For each column j < n of L, remove all the nonzeros in the column j
except the first one below the diagonal.
G(A
Let Lt denote the remaining structure and cGo(nF)sider the matrix
G(Ft)= T(A
T
Ft = Lt + Lt . The graph G(Ft) is a tree called the elimination tree.
FIG. 2.2. Graph structures of the example in Fig. 2.1.
Uses of the elimination tree
The elimination tree has several uses in the factorization of a sparse
matrix:
• it expresses the order in which variables can be eliminated:
◦ the elimination of a variable only affects (directly or indirectly) its
ancestors and ...
◦ ... only depends on its descendants
Therefore, any topological order of the elimination tree leads to a
correct result and to the same fill-in
• it expresses concurrence: because variables in separate subtrees do
not affect each other, they can be eliminated in parallel