Module ordering 

Chapter 5.10 — Fill-reducing orderings: Natural, Amd, Colamd

New with task L1.3 (PLAN.md trilha L; D-0080, decision D10: the orderings are our own from v0.1, no C, no GPL/LGPL).

Purpose

The cost of a sparse factorization is decided before the first flop, by the order of elimination. An Ordering looks at a Pattern (Chapter 5.9) and returns a Perm; it never sees values. Two targets (OrderingTarget):

• OrderingTarget::Symmetric — a symmetric order P (S + Sᵀ) Pᵀ for Cholesky or LDLᵀ (the KKT systems of the interior-point method, L1.5): the permutation of the rows and columns;
• OrderingTarget::Columns — a column order S Q for an LU with row pivoting (or a QR): it reduces the fill of the Cholesky factor of QᵀSᵀSQ, which bounds the U of every row-pivoted LU of S Q (George & Ng 1987), without forming SᵀS.

Three orderings, stateless policies:

• Natural — the identity (the reference of every table);
• Amd — approximate minimum degree on the graph of S + Sᵀ (Amestoy, Davis & Duff 1996), with aggressive absorption, mass elimination, supervariables and dense-row detection; for the Columns target it orders the graph of SᵀS, formed explicitly;
• Colamd — column approximate minimum degree (Davis, Gilbert, Larimore & Ng 2004): the same process on the quotient graph of SᵀS started from the rows of S, never forming SᵀS; for the Symmetric target, the SYMAMD construction (every off-diagonal pair {i, j} of S + Sᵀ is a row with two entries, whose MᵀM has the pattern of S + Sᵀ).

Both run Algorithm Q of Chapter 5.11; this chapter is their setup.

Notation

S is m × n; r_i and c_j are its row and column counts. κ(α, N) = max(16, ⌊α √N⌋) is the dense threshold (SuiteSparse’s rule); a negative α disables dense detection. Excluded (dense or empty) indices are ordered last, in increasing index order.

Algorithm A (AMD)

A1. [Graph.] G ← S + Sᵀ without its diagonal (Algorithm S of Chapter 5.9; for Columns, the graph of SᵀS by Algorithm G). An error unless S is square for Symmetric. A2. [Dense rows.] With δ = κ(dense, n) (dense default 10), every node of degree > δ is excluded: ordering it in the middle would cost as much as it fills, and its presence would spoil every degree it touches. A3. [Order.] Algorithm Q with A_i = G’s adjacency, no initial elements, the bound min(d̄_i, δ_i), AMD’s conventions (largest index first among equal degrees, the postorder Q8 as output); aggressive absorption as asked (default on). A4. Append the excluded nodes. Return the Perm (fwd = pivot order).

Algorithm C (COLAMD)

C1. [Dense columns.] With κ_c = κ(dense_col, min(m, n)) (default dense_col = 10), exclude every column with c_j > κ_c. C2. [Dense rows.] With κ_r = κ(dense_row, n) (default 10), drop every row with more than κ_r entries in non-excluded columns: a dense row makes SᵀS dense and says nothing useful about the order. C3. [Null columns.] Exclude every column left without a row (they are ordered last, before the dense columns, as COLAMD does). C4. [Order.] Algorithm Q with no A_i, the surviving rows as initial elements, COLAMD’s score δ_i (no min with the previous score) and COLAMD’s conventions (smallest index first among ties, initial rows before pivot rows in Q3, the pivot order Q8′ as output); aggressive absorption on. C5. Append the null columns, then the dense ones.

For Symmetric (SYMAMD): exclude the nodes of S + Sᵀ of degree > κ(dense_col, n); every edge {i, j} (i < j) between the others is an initial element {i, j}; then C4–C5 with the excluded nodes.

Invariants

(O1) Every ordering returns a Perm of n (the number of columns) satisfying (P1) of Chapter 5.9, and fails only on an input error (a non-square pattern for Symmetric, (P2) violated). (O2) Natural is the identity; the others are deterministic functions of the pattern.

Defaults and their source

Amd { dense: None (= 10), aggressive: true } and Colamd { dense_row: None (= 10), dense_col: None (= 10) } are SuiteSparse’s AMD_DENSE, AMD_AGGRESSIVE, COLAMD_DENSE_ROW and COLAMD_DENSE_COL defaults (AMD 2.x and COLAMD 2.x user guides).

References

Amestoy, Davis & Duff, SIAM J. Matrix Anal. Appl. 17(4), 1996, and Algorithm 837 (AMD), ACM TOMS 30(3), 2004; Davis, Gilbert, Larimore & Ng, “A column approximate minimum degree ordering algorithm” and Algorithm 836 (COLAMD), ACM TOMS 30(3), 2004; George & Ng (1987); the contract of the crate (Ordering, OrderingTarget, Natural, Amd, Colamd).

Structs

Amd

Approximate minimum degree (Amestoy, Davis & Duff 1996), Algorithm A: quotient graph with element absorption, approximate external degrees, mass elimination of indistinguishable nodes and dense-row detection. Our own implementation (D10); SuiteSparse AMD is the oracle in tests.

Colamd

Column approximate minimum degree (Davis, Gilbert, Larimore & Ng 2004), Algorithm C: the AMD idea applied to the pattern of SᵀS without forming it, for an LU with row pivoting. Our own implementation (D10).

Natural

The identity ordering.

Enums

OrderingTarget

What an ordering is asked to reduce fill for.

Traits

Ordering

Fill-reducing orderings. Implementations are stateless policies (invariants (O1)–(O2)).