Nonsymmetric Gaussian elimination
N
N
o
o
n
n
s
s
y
y
m
m
m
m
e
e
t
t
r
r
i
i
c
c
G
G
a
a
u
u
s
s
s
s
i
i
a
a
n
n
e
e
l
l
i
i
m
m
i
i
n
n
a
a
t
t
i
i
o
o
n
n
•
A = LU:
does not always exist, can be unstable
•
PA = LU:
Partial pivoting
•
At each elimination step, pivot on largest-magnitude element in column
•
“
GEPP
”
is the standard algorithm for dense nonsymmetric systems
•
PAQ = LU:
Complete pivoting
•
Pivot on largest-magnitude element in the entire uneliminated matrix
•
Expensive to search for the pivot
•
No freedom to reorder for sparsity
•
Hardly ever used in practice
•
Conflict between permuting for sparsity and for numerics
•
Lots of different approaches to this tradeoff
L
L
e
e
f
f
t
t
-
-
l
l
o
o
o
o
k
k
i
i
n
n
g
g
C
C
o
o
l
l
u
u
m
m
n
n
L
L
U
U
F
F
a
a
c
c
t
t
o
o
r
r
i
i
z
z
a
a
t
t
i
i
o
o
n
n
for
column j = 1 to n
do
solve
scale
:
l
j
=
l
j
/
u
jj
•
Column j of A becomes column j of L and U
S
y
m
b
o
l
i
c
s
p
a
r
s
e
G
a
u
s
s
i
a
n
e
l
i
m
i
n
a
t
i
o
n
:
A
=
L
U
•
Add fill edge a
->
b if there is a path from a to b
through lower-numbered vertices.
•
But this doesn
’
t work with numerical pivoting!
A
G (A)
C
C
o
o
l
l
u
u
m
m
n
n
P
P
r
r
e
e
o
o
r
r
d
d
e
e
r
r
i
i
n
n
g
g
f
f
o
o
r
r
S
S
p
p
a
a
r
r
s
s
i
i
t
t
y
y
•
P
A
Q
T
= LU
:
Q
preorders columns for sparsity,
P
is row pivoting
•
Column permutation of A
Symmetric permutation of A
T
A (or
G
(A))
•
S
y
m
m
e
t
r
i
c
o
r
d
e
r
i
n
g
:
A
p
p
r
o
x
i
m
a
t
e
m
i
n
i
m
u
m
d
e
g
r
e
e
,
e
t
c
.
•
But, forming A
T
A
is expensive (sometimes bigger than L+U).
P
Q
C
C
o
o
l
l
u
u
m
m
n
n
I
I
n
n
t
t
e
e
r
r
s
s
e
e
c
c
t
t
i
i
o
o
n
n
G
G
r
r
a
a
p
p
h
h
•
S
y
m
b
o
l
i
c
v
e
r
s
i
o
n
o
f
t
h
e
n
o
r
m
a
l
e
q
u
a
t
i
o
n
s
A
T
A
x
=
A
T
b
•
G
(
A
)
=
G
(
A
T
A
)
i
f
n
o
c
a
n
c
e
l
l
a
t
i
o
n
(
o
t
h
e
r
w
i
s
e
)
•
Permuting the rows of A does not change
G
(A)
A
G
(A)
A
T
A
F
F
i
i
l
l
l
l
e
e
d
d
C
C
o
o
l
l
u
u
m
m
n
n
I
I
n
n
t
t
e
e
r
r
s
s
e
e
c
c
t
t
i
i
o
o
n
n
G
G
r
r
a
a
p
p
h
h
•
G
(
A
)
=
s
y
m
b
o
l
i
c
C
h
o
l
e
s
k
y
f
a
c
t
o
r
o
f
A
T
A
•
In
PA=LU
,
G(U)
G
(A)
and
G(L)
G
(A)
•
Tighter bound on
L
from symbolic
QR
•
Bounds are best possible if
A
is strong Hall
A
c
h
o
l
(
A
T
A
)
C
C
o
o
l
l
u
u
m
m
n
n
E
E
l
l
i
i
m
m
i
i
n
n
a
a
t
t
i
i
o
o
n
n
T
T
r
r
e
e
e
e
•
E
l
i
m
i
n
a
t
i
o
n
t
r
e
e
o
f
A
T
A
(
i
f
n
o
c
a
n
c
e
l
l
a
t
i
o
n
)
•
Depth-first spanning tree of
G
(A)
•
Represents column dependencies in various factorizations
A
c
h
o
l
(
A
T
A
)
T
(A)
+
C
C
o
o
l
l
u
u
m
m
n
n
D
D
e
e
p
p
e
e
n
n
d
d
e
e
n
n
c
c
i
i
e
e
s
s
i
i
n
n
P
P
A
A
=
=
L
L
U
U
•
If column j modifies
column k, then j
T
[k].
*
definition of “strong Hall” coming up in a few slides…
E
E
f
f
f
f
i
i
c
c
i
i
e
e
n
n
t
t
S
S
t
t
r
r
u
u
c
c
t
t
u
u
r
r
e
e
P
P
r
r
e
e
d
d
i
i
c
c
t
t
i
i
o
o
n
n
G
i
v
e
n
t
h
e
s
t
r
u
c
t
u
r
e
o
f
(
u
n
s
y
m
m
e
t
r
i
c
)
A
,
o
n
e
c
a
n
f
i
n
d
.
.
.
•
column elimination tree
T
(A)
•
row and column counts for
G
(A)
•
supernodes of
G
(A)
•
nonzero structure of
G
(A)
.
.
.
w
i
t
h
o
u
t
f
o
r
m
i
n
g
G
(
A
)
o
r
A
T
A
+
+
+
S
y
m
m
e
t
r
i
c
A
i
m
p
l
i
e
s
G
+
(
A
)
i
s
c
h
o
r
d
a
l
,
w
i
t
h
l
o
t
s
o
f
s
t
r
u
c
t
u
r
e
a
n
d
e
l
e
g
a
n
t
t
h
e
o
r
y
For unsymmetric A, things are not as nice
•
N
o
k
n
o
w
n
w
a
y
t
o
c
o
m
p
u
t
e
G
+
(
A
)
f
a
s
t
e
r
t
h
a
n
G
a
u
s
s
i
a
n
e
l
i
m
i
n
a
t
i
o
n
•
No fast way to recognize perfect elimination graphs
•
No theory of approximately optimal orderings
•
Directed analogs of elimination tree:
Smaller graphs that preserve path structure
D
D
i
i
r
r
e
e
c
c
t
t
e
e
d
d
g
g
r
r
a
a
p
p
h
h
•
A is square, unsymmetric, nonzero diagonal
•
Edges from rows to columns
•
S
y
m
m
e
t
r
i
c
p
e
r
m
u
t
a
t
i
o
n
s
P
A
P
T
r
e
n
u
m
b
e
r
v
e
r
t
i
c
e
s
A
G(A)
S
S
t
t
r
r
o
o
n
n
g
g
l
l
y
y
c
c
o
o
n
n
n
n
e
e
c
c
t
t
e
e
d
d
c
c
o
o
m
m
p
p
o
o
n
n
e
e
n
n
t
t
s
s
•
Symmetric permutation to block triangular form
•
Diagonal blocks are Strong Hall
(irreducible / strongly connected)
•
Find P in linear time by depth-first search
[Tarjan]
•
Row and column partitions are independent of choice of nonzero diagonal
•
Solve Ax=b by block back substitution
PAP
T
G(A)
S
S
o
o
l
l
v
v
i
i
n
n
g
g
A
A
*
*
x
x
=
=
b
b
i
i
n
n
b
b
l
l
o
o
c
c
k
k
t
t
r
r
i
i
a
a
n
n
g
g
u
u
l
l
a
a
r
r
f
f
o
o
r
r
m
m
% Permute A to block form
[p,q,r] = dmperm(A);
A = A(p,q); x = b(p);
% Block backsolve
nblocks = length(r) – 1;
for k = nblocks : –1 : 1
% Indices above the k-th block
I = 1 : r(k) – 1;
% Indices of the k-th block
J = r(k) : r(k+1) – 1;
x(J) = A(J,J) \ x(J);
x(I) = x(I) – A(I,J) * x(J);
end;
% Undo the permutation of x
x(q) = x;
B
i
p
a
r
t
i
t
e
m
a
t
c
h
i
n
g
:
P
e
r
m
u
t
a
t
i
o
n
t
o
n
o
n
z
e
r
o
d
i
a
g
o
n
a
l
•
Represent A as an undirected bipartite graph (one
node for each row and one node for each column)
•
Find
perfect matching
: set of edges that hits each
vertex exactly once
•
Permute rows to place matching on diagonal
1
5
2
3
4
1
5
2
3
4
A
d
d
m
m
p
p
e
e
r
r
m
m
:
:
M
M
a
a
t
t
c
c
h
h
i
i
n
n
g
g
a
a
n
n
d
d
b
b
l
l
o
o
c
c
k
k
t
t
r
r
i
i
a
a
n
n
g
g
u
u
l
l
a
a
r
r
f
f
o
o
r
r
m
m
•
Dulmage-Mendelsohn decomposition:
•
Bipartite matching followed by strongly connected components
•
Square A with nonzero diagonal:
•
[p, p, r] = dmperm(A);
•
connected components of an undirected graph
•
strongly connected components of a directed graph
•
Square, full rank A:
•
[p, q, r] = dmperm(A);
•
A(p,q) has nonzero diagonal and is in block upper triangular form
•
Arbitrary A:
•
[p, q, r, s] = dmperm(A);
•
maximum-size matching in a bipartite graph
•
minimum-size vertex cover in a bipartite graph
•
decomposition into strong Hall blocks
S
S
t
t
r
r
o
o
n
n
g
g
H
H
a
a
l
l
l
l
c
c
o
o
m
m
p
p
s
s
a
a
r
r
e
e
i
i
n
n
d
d
e
e
p
p
e
e
n
n
d
d
e
e
n
n
t
t
o
o
f
f
m
m
a
a
t
t
c
c
h
h
i
i
n
n
g
g
D
D
u
u
l
l
m
m
a
a
g
g
e
e
-
-
M
M
e
e
n
n
d
d
e
e
l
l
s
s
o
o
h
h
n
n
T
T
h
h
e
e
o
o
r
r
y
y
•
A. L. Dulmage & N. S. Mendelsohn.
“
Coverings of bipartite
graphs.
”
Can. J. Math.
10: 517-534, 1958.
•
A. L. Dulmage & N. S. Mendelsohn.
“
The term and stochastic ranks
of a matrix.
”
Can. J. Math.
11: 269-279, 1959.
•
A. L. Dulmage & N. S. Mendelsohn.
“
A structure theory of bipartite
graphs of finite exterior dimension.
”
Trans. Royal Soc. Can.,
ser. 3,
53: 1-13, 1959.
•
D. M. Johnson, A. L. Dulmage, & N. S. Mendelsohn.
“
Connectivity
and reducibility of graphs.
”
Can. J. Math.
14: 529-539, 1962.
•
A. L. Dulmage & N. S. Mendelsohn.
“
Two algorithms for bipartite
graphs.
”
SIAM J.
11: 183-194, 1963.
•
A. Pothen & C.-J. Fan.
“
Computing the block triangular form of a
sparse matrix.
”
ACM Trans. Math. Software
16: 303-324, 1990.
H
H
a
a
l
l
l
l
a
a
n
n
d
d
s
s
t
t
r
r
o
o
n
n
g
g
H
H
a
a
l
l
l
l
p
p
r
r
o
o
p
p
e
e
r
r
t
t
i
i
e
e
s
s
Let
G
be a bipartite graph with
m
“
row
”
vertices and
n
“
column
”
vertices.
•
A
matching
is a set of edges of
G
with no common endpoints.
•
G
has the
Hall property
if for all
k >= 0
, every set of
k
columns is
adjacent to at least
k
rows.
•
Hall
’
s theorem:
G
has a matching of size
n
iff
G
has the Hall property.
•
G
has the
strong Hall property
if for all
k
with
0 < k < n
, every set
of
k
columns is adjacent to at least
k+1
rows.
A
A
l
l
t
t
e
e
r
r
n
n
a
a
t
t
i
i
n
n
g
g
p
p
a
a
t
t
h
h
s
s
•
Let
M
be a matching. An
alternating walk
is a sequence of edges with
every second edge in
M
. (Vertices or edges may appear more than
once in the walk.) An
alternating tour
is an alternating walk whose
endpoints are the same. An
alternating path
is an alternating walk with
no repeated vertices. An
alternating cycle
is an alternating tour with no
repeated vertices except its endpoint.
•
Lemma.
Let M and
N
be two maximum matchings. Their symmetric
difference
(M
N)
–
(M
N)
consists of vertex-disjoint components,
each of which is either
1.
an alternating cycle in both
M
and
N,
or
2.
an alternating path in both
M
and
N
from an
M
-unmatched column to an
N
-unmatched column, or
3.
same as
2
but for rows.
D
D
u
u
l
l
m
m
a
a
g
g
e
e
-
-
M
M
e
e
n
n
d
d
e
e
l
l
s
s
o
o
h
h
n
n
d
d
e
e
c
c
o
o
m
m
p
p
o
o
s
s
i
i
t
t
i
i
o
o
n
n
(
(
c
c
o
o
a
a
r
r
s
s
e
e
)
)
Let
M
be a maximum-size matching. Define:
•
VR
=
{ rows reachable via alt. path from some unmatched row
}
•
VC
=
{ cols reachable via alt. path from some unmatched row
}
•
HR
=
{ rows reachable via alt. path from some unmatched col
}
•
HC
=
{ cols reachable via alt. path from some unmatched col
}
•
S
R
=
R
–
V
R
–
H
R
•
S
C
=
C
–
V
C
–
H
C
D
D
u
u
l
l
m
m
a
a
g
g
e
e
-
-
M
M
e
e
n
n
d
d
e
e
l
l
s
s
o
o
h
h
n
n
d
d
e
e
c
c
o
o
m
m
p
p
o
o
s
s
i
i
t
t
i
i
o
o
n
n
D
D
u
u
l
l
m
m
a
a
g
g
e
e
-
-
M
M
e
e
n
n
d
d
e
e
l
l
s
s
o
o
h
h
n
n
t
t
h
h
e
e
o
o
r
r
y
y
•
Theorem 1.
VR
,
HR
, and
SR
are pairwise disjoint.
VC
,
HC
, and
SC
are pairwise disjoint.
•
Theorem 2.
No matching edge joins
xR
and
yC
if
x
and
y
are different.
•
Theorem 3.
No edge joins
VR
and
SC
, or
VR
and
HC
, or
SR
and
HC
.
•
Theorem 4.
SR
and
SC
are perfectly matched to each other.
•
Theorem 5.
The subgraph induced by
VR
and
VC
has the
strong Hall property. The transpose of the subgraph
induced by
HR
and
HC
has the strong Hall property.
•
Theorem 6.
The vertex sets
VR
,
HR
,
SR
,
VC
,
HC
,
SC
are
independent of the choice of maximum matching
M
.
D
D
u
u
l
l
m
m
a
a
g
g
e
e
-
-
M
M
e
e
n
n
d
d
e
e
l
l
s
s
o
o
h
h
n
n
d
d
e
e
c
c
o
o
m
m
p
p
o
o
s
s
i
i
t
t
i
i
o
o
n
n
(
(
f
f
i
i
n
n
e
e
)
)
•
Consider the perfectly matched square block induced by
SR
and
SC
. In the sequel we shall ignore
VR
,
VC
,
HR
,
and
HC
. Thus,
G
is a bipartite graph with
n
row vertices
and
n
column vertices, and
G
has a perfect matching
M
.
•
Call two columns
equivalent
if they lie on an alternating
tour. This is an equivalence relation; let the equivalence
classes be
C
1
,
C
2
, . . .,
C
p
. Let
R
i
be the set of rows
matched to
C
i
.
T
T
h
h
e
e
f
f
i
i
n
n
e
e
D
D
u
u
l
l
m
m
a
a
g
g
e
e
-
-
M
M
e
e
n
n
d
d
e
e
l
l
s
s
o
o
h
h
n
n
d
d
e
e
c
c
o
o
m
m
p
p
o
o
s
s
i
i
t
t
i
i
o
o
n
n
Matrix A
Bipartite graph H(A)
Directed graph
G(A)
D
D
u
u
l
l
m
m
a
a
g
g
e
e
-
-
M
M
e
e
n
n
d
d
e
e
l
l
s
s
o
o
h
h
n
n
t
t
h
h
e
e
o
o
r
r
y
y
•
Theorem 7.
The
R
i
’
s and the
C
j
’
s can be renumbered so no edge
joins
R
i
and
C
j
if
i > j
.
•
Theorem 8.
The subgraph induced by
R
i
and
C
i
has the strong Hall property.
•
Theorem 9.
The partition
R
1
C
1
,
R
2
C
2
, . . .,
R
p
C
p
is independent of the
choice of maximum matching.
•
Theorem 10.
If non-matching edges are directed from rows to columns and
matching edges are shrunk into single vertices, the resulting directed graph
G(A)
has strongly connected components
C
1
,
C
2
, . . .,
C
p
.
•
Theorem 11.
A bipartite graph
G
has the strong Hall property
iff
every pair of edges of
G
is on some alternating tour
iff
G
is connected and every edge of
G
is in some perfect matching.
•
Theorem 12.
Given a square matrix
A
, if we permute rows and columns to get
a nonzero diagonal and then do a symmetric permutation to put the strongly
connected components into topological order (i.e. in block triangular form),
then the grouping of rows and columns into diagonal blocks is independent
of the choice of nonzero diagonal.
S
S
t
t
r
r
o
o
n
n
g
g
l
l
y
y
c
c
o
o
n
n
n
n
e
e
c
c
t
t
e
e
d
d
c
c
o
o
m
m
p
p
o
o
n
n
e
e
n
n
t
t
s
s
a
a
r
r
e
e
i
i
n
n
d
d
e
e
p
p
e
e
n
n
d
d
e
e
n
n
t
t
o
o
f
f
c
c
h
h
o
o
i
i
c
c
e
e
o
o
f
f
p
p
e
e
r
r
f
f
e
e
c
c
t
t
m
m
a
a
t
t
c
c
h
h
i
i
n
n
g
g
M
M
a
a
t
t
r
r
i
i
x
x
t
t
e
e
r
r
m
m
i
i
n
n
o
o
l
l
o
o
g
g
y
y
•
S
q
u
a
r
e
m
a
t
r
i
x
A
i
s
i
r
r
e
d
u
c
i
b
l
e
i
f
t
h
e
r
e
d
o
e
s
n
o
t
e
x
i
s
t
a
n
y
p
e
r
m
u
t
a
t
i
o
n
m
a
t
r
i
x
P
s
u
c
h
t
h
a
t
P
A
P
T
h
a
s
a
n
o
n
t
r
i
v
i
a
l
b
l
o
c
k
t
r
i
a
n
g
u
l
a
r
f
o
r
m
[
A
1
1
A
1
2
;
0
A
2
2
]
.
•
S
q
u
a
r
e
m
a
t
r
i
x
A
i
s
f
u
l
l
y
i
n
d
e
c
o
m
p
o
s
a
b
l
e
i
f
t
h
e
r
e
d
o
n
o
t
e
x
i
s
t
a
n
y
p
e
r
m
u
t
a
t
i
o
n
m
a
t
r
i
c
e
s
P
a
n
d
Q
s
u
c
h
t
h
a
t
P
A
Q
T
h
a
s
a
n
o
n
t
r
i
v
i
a
l
b
l
o
c
k
t
r
i
a
n
g
u
l
a
r
f
o
r
m
[
A
1
1
A
1
2
;
0
A
2
2
]
.
•
Fully indecomposable implies irreducible, not vice versa.
•
Fully indecomposable = square and strong Hall.
•
A square matrix with nonzero diagonal is irreducible
iff
fully indecomposable
iff
strong Hall
iff
strongly connected.
A
A
p
p
p
p
l
l
i
i
c
c
a
a
t
t
i
i
o
o
n
n
s
s
o
o
f
f
D
D
-
-
M
M
d
d
e
e
c
c
o
o
m
m
p
p
o
o
s
s
i
i
t
t
i
i
o
o
n
n
•
Permutation to block triangular form for Ax=b
•
Connected components of undirected graphs
•
Strongly connected components of directed graphs
•
Minimum-size vertex cover for bipartite graphs
•
Extracting vertex separators from edge cuts
for arbitrary graphs
•
For strong Hall matrices, several upper bounds in
nonzero structure prediction are best possible:
•
Column intersection graph factor is R in QR
•
Column intersection graph factor is tight bound on U in PA=LU
•
Row merge graph is tight bound on Lbar and U in PA=LU
Intricacies of nonsymmetric Gaussian elimination, LU factorization, partial pivoting, left-looking column LU factorization, symbolic sparse Gaussian elimination, column preordering for sparsity, and more in numerical linear algebra algorithms.
Download Presentation
Please find below an Image/Link to download the presentation.
The content on the website is provided AS IS for your information and personal use only. It may not be sold, licensed, or shared on other websites without obtaining consent from the author.If you encounter any issues during the download, it is possible that the publisher has removed the file from their server.
You are allowed to download the files provided on this website for personal or commercial use, subject to the condition that they are used lawfully. All files are the property of their respective owners.
The content on the website is provided AS IS for your information and personal use only. It may not be sold, licensed, or shared on other websites without obtaining consent from the author.
E N D
Presentation Transcript
Nonsymmetric Gaussian elimination A = LU: does not always exist, can be unstable PA = LU: Partial pivoting At each elimination step, pivot on largest-magnitude element in column GEPP is the standard algorithm for dense nonsymmetric systems PAQ = LU: Complete pivoting Pivot on largest-magnitude element in the entire uneliminated matrix Expensive to search for the pivot No freedom to reorder for sparsity Hardly ever used in practice Conflict between permuting for sparsity and for numerics Lots of different approaches to this tradeoff
Left-looking Column LU Factorization j for column j = 1 to n do solve L 0 L I ( ) U lj( )= aj for uj, lj L uj A scale: lj= lj/ ujj L Column j of A becomes column j of L and U
Symbolic sparse Gaussian elimination: A = LU 2 1 4 5 7 6 3 + A L+U G (A) Add fill edge a -> b if there is a path from a to b through lower-numbered vertices. But this doesn t work with numerical pivoting!
Column Preordering for Sparsity Q P = x PAQT = LU: Q preorders columns for sparsity, P is row pivoting Column permutation of A Symmetric permutation of ATA (or G (A)) Symmetric ordering: Approximate minimum degree, etc. But, forming ATA is expensive (sometimes bigger than L+U).
Column Intersection Graph 3 1 2 3 4 5 1 2 3 4 5 1 2 3 4 4 5 2 5 1 A ATA G (A) Symbolic version of the normal equations ATAx=AT b G (A) = G(ATA) if no cancellation(otherwise ) Permuting the rows of A does not change G (A)
Filled Column Intersection Graph 3 1 2 3 4 5 1 2 3 4 5 1 2 3 4 4 5 2 5 1 + A chol(ATA) G (A) + G (A) = symbolic Cholesky factor of ATA In PA=LU, G(U) G (A) and G(L) G (A) Tighter bound on L from symbolic QR Bounds are best possible if A is strong Hall + +
Column Elimination Tree 5 1 2 3 4 5 1 2 3 4 5 1 2 3 4 4 2 3 5 1 T (A) A chol(ATA) Elimination tree of ATA (if no cancellation) + Depth-first spanning tree of G (A) Represents column dependencies in various factorizations
Column Dependencies in PA = LU If column j modifies column k, then j T [k]. k j T [k] If A is strong Hall* then, for some pivot sequence P, every column modifies its parent in T (A). * definition of strong Hall coming up in a few slides
Efficient Structure Prediction Given the structure of (unsymmetric) A, one can find . . . column elimination tree T (A) row and column counts for G (A) supernodes of G (A) nonzero structure of G (A) + + + . . . without forming G (A) or ATA
Symmetric A implies G+(A) is chordal, with lots of structure and elegant theory For unsymmetric A, things are not as nice No known way to compute G+(A) faster than Gaussian elimination No fast way to recognize perfect elimination graphs No theory of approximately optimal orderings Directed analogs of elimination tree: Smaller graphs that preserve path structure
Directed graph 2 1 4 5 7 6 3 A G(A) A is square, unsymmetric, nonzero diagonal Edges from rows to columns Symmetric permutations PAPT renumber vertices
Strongly connected components 1 2 4 7 5 3 6 2 1 1 2 4 7 4 5 7 5 3 6 6 3 G(A) PAPT Symmetric permutation to block triangular form Diagonal blocks are Strong Hall (irreducible / strongly connected) Find P in linear time by depth-first search [Tarjan] Row and column partitions are independent of choice of nonzero diagonal Solve Ax=b by block back substitution
Solving A*x = b in block triangular form % Permute A to block form [p,q,r] = dmperm(A); A = A(p,q); x = b(p); 1 2 3 4 5 6 7 1 2 3 4 % Block backsolve nblocks = length(r) 1; for k = nblocks : 1 : 1 = % Indices above the k-th block I = 1 : r(k) 1; 5 6 7 % Indices of the k-th block J = r(k) : r(k+1) 1; A x b x(J) = A(J,J) \ x(J); x(I) = x(I) A(I,J) * x(J); end; % Undo the permutation of x x(q) = x;
Bipartite matching: Permutation to nonzero diagonal 1 2 3 4 5 1 2 3 4 5 1 1 1 2 3 4 4 5 3 1 2 2 3 3 4 4 5 2 5 5 A PA Represent A as an undirected bipartite graph (one node for each row and one node for each column) Find perfect matching: set of edges that hits each vertex exactly once Permute rows to place matching on diagonal
dmperm: Matching and block triangular form Dulmage-Mendelsohn decomposition: Bipartite matching followed by strongly connected components Square A with nonzero diagonal: [p, p, r] = dmperm(A); connected components of an undirected graph strongly connected components of a directed graph Square, full rank A: [p, q, r] = dmperm(A); A(p,q) has nonzero diagonal and is in block upper triangular form Arbitrary A: [p, q, r, s] = dmperm(A); maximum-size matching in a bipartite graph minimum-size vertex cover in a bipartite graph decomposition into strong Hall blocks
Strong Hall comps are independent of matching 1 2 4 7 5 3 6 1 1 22 11 1 2 4 7 2 2 3 3 44 55 4 4 77 5 5 5 3 6 6 6 66 33 7 7 1 1 1 2 4 7 5 3 6 12 41 4 1 7 2 2 2 3 3 74 55 4 4 27 5 5 5 6 3 6 6 36 63 7 7
Dulmage-Mendelsohn Theory A. L. Dulmage & N. S. Mendelsohn. Coverings of bipartite graphs. Can. J. Math. 10: 517-534, 1958. A. L. Dulmage & N. S. Mendelsohn. The term and stochastic ranks of a matrix. Can. J. Math. 11: 269-279, 1959. A. L. Dulmage & N. S. Mendelsohn. A structure theory of bipartite graphs of finite exterior dimension. Trans. Royal Soc. Can., ser. 3, 53: 1-13, 1959. D. M. Johnson, A. L. Dulmage, & N. S. Mendelsohn. Connectivity and reducibility of graphs. Can. J. Math. 14: 529-539, 1962. A. L. Dulmage & N. S. Mendelsohn. Two algorithms for bipartite graphs. SIAM J. 11: 183-194, 1963. A. Pothen & C.-J. Fan. Computing the block triangular form of a sparse matrix. ACM Trans. Math. Software 16: 303-324, 1990.
Hall and strong Hall properties Let G be a bipartite graph with m row vertices and n column vertices. A matching is a set of edges of G with no common endpoints. G has the Hall property if for all k >= 0, every set of k columns is adjacent to at least k rows. Hall s theorem: G has a matching of size n iff G has the Hall property. G has the strong Hall property if for all k with 0 < k < n, every set of k columns is adjacent to at least k+1 rows.
Alternating paths Let M be a matching. An alternating walk is a sequence of edges with every second edge in M. (Vertices or edges may appear more than once in the walk.) An alternating tour is an alternating walk whose endpoints are the same. An alternating path is an alternating walk with no repeated vertices. An alternating cycle is an alternating tour with no repeated vertices except its endpoint. Lemma. Let M and N be two maximum matchings. Their symmetric difference (M N) (M N) consists of vertex-disjoint components, each of which is either 1. an alternating cycle in both M and N, or 2. an alternating path in both M and N from an M-unmatched column to an N-unmatched column, or 3. same as 2 but for rows.
Dulmage-Mendelsohn decomposition (coarse) Let M be a maximum-size matching. Define: VR = { rows reachable via alt. path from some unmatched row } VC = { cols reachable via alt. path from some unmatched row } HR = { rows reachable via alt. path from some unmatched col } HC = { cols reachable via alt. path from some unmatched col } SR = R VR HR SC = C VC HC
Dulmage-Mendelsohn decomposition 1 2 1 2 3 4 5 6 7 8 9 10 11 1 2 3 4 1 3 HR HC 2 4 3 5 4 6 5 6 7 8 9 5 7 SR SC 6 8 7 9 8 10 10 11 9 11 VR VC 10 12 11 12
Dulmage-Mendelsohn theory Theorem 1. VR, HR, and SR are pairwise disjoint. VC, HC, and SC are pairwise disjoint. Theorem 2. No matching edge joins xR and yC if x and y are different. Theorem 3. No edge joins VR and SC, or VR and HC, or SR and HC. Theorem 4. SR and SC are perfectly matched to each other. Theorem 5. The subgraph induced by VR and VC has the strong Hall property. The transpose of the subgraph induced by HR and HC has the strong Hall property. Theorem 6. The vertex sets VR, HR, SR, VC, HC, SC are independent of the choice of maximum matching M.
Dulmage-Mendelsohn decomposition (fine) Consider the perfectly matched square block induced by SR and SC. In the sequel we shall ignore VR, VC, HR, and HC. Thus, G is a bipartite graph with n row vertices and n column vertices, and G has a perfect matching M. Call two columns equivalent if they lie on an alternating tour. This is an equivalence relation; let the equivalence classes be C1, C2, . . ., Cp. Let Ri be the set of rows matched to Ci.
The fine Dulmage-Mendelsohn decomposition 1 2 3 4 5 6 7 Matrix A 1 2 3 4 1 1 2 2 R1 C1 3 3 5 6 7 4 4 5 5 R2 C2 Directed graph G(A) 6 6 R3 C3 2 1 7 7 Bipartite graph H(A) 3 5 4 3 6
Dulmage-Mendelsohn theory Theorem 7. The Ri s and the Cj s can be renumbered so no edge joins Ri and Cj if i > j. Theorem 8. The subgraph induced by Ri and Ci has the strong Hall property. Theorem 9. The partition R1 C1 , R2 C2 , . . ., Rp Cp is independent of the choice of maximum matching. Theorem 10. If non-matching edges are directed from rows to columns and matching edges are shrunk into single vertices, the resulting directed graph G(A) has strongly connected components C1 , C2 , . . ., Cp. Theorem 11. A bipartite graph G has the strong Hall property iff every pair of edges of G is on some alternating tour iff G is connected and every edge of G is in some perfect matching. Theorem 12. Given a square matrix A, if we permute rows and columns to get a nonzero diagonal and then do a symmetric permutation to put the strongly connected components into topological order (i.e. in block triangular form), then the grouping of rows and columns into diagonal blocks is independent of the choice of nonzero diagonal.
Strongly connected components are independent of choice of perfect matching 1 2 4 7 5 3 6 1 1 22 11 1 2 4 7 2 2 3 3 44 55 4 4 77 5 5 5 3 6 6 6 66 33 7 7 1 1 1 2 4 7 5 3 6 12 41 4 1 7 2 2 2 3 3 74 55 4 4 27 5 5 5 6 3 6 6 36 63 7 7
Matrix terminology Square matrix A is irreducible if there does not exist any permutation matrix P such that PAPT has a nontrivial block triangular form [A11 A12 ; 0 A22]. Square matrix A is fully indecomposable if there do not exist any permutation matrices P and Q such that PAQT has a nontrivial block triangular form [A11 A12 ; 0 A22]. Fully indecomposable implies irreducible, not vice versa. Fully indecomposable = square and strong Hall. A square matrix with nonzero diagonal is irreducible iff fully indecomposable iff strong Hall iff strongly connected.
Applications of D-M decomposition Permutation to block triangular form for Ax=b Connected components of undirected graphs Strongly connected components of directed graphs Minimum-size vertex cover for bipartite graphs Extracting vertex separators from edge cuts for arbitrary graphs For strong Hall matrices, several upper bounds in nonzero structure prediction are best possible: Column intersection graph factor is R in QR Column intersection graph factor is tight bound on U in PA=LU Row merge graph is tight bound on Lbar and U in PA=LU
