Upcoming Workshops/ Travel Plans
I am going to visit:
- ECrypt’s Tools for Cryptanalysis Workshop, Institute of Mathematics (IMPAN), Krakow (Poland), September 24 - 25, 2007.
- SAGE Days 5, Clay Math Institute, Cambridge (USA), September 30 - October 3, 2007.
- SAGE Days 6, Heilbronn Institute, Bristol (UK), November 10 - 14, 2007.
SD5 will be on computational arithmetic geometry and SD6 will focus cryptography and arithmetic geometry. As usually, the prime resource for SAGE related workshops is the SAGE wiki.
The Sorry State of Sparse Linear Algebra over Finite Fields
By sparse linear algebra I actually only mean computing the (reduced) row echelon form. Surprisingly, there aren’t much implementations out there, not much I am aware of at least.
- Apparently, MAGMA cannot compute the reduced row echelon form. Well, you can compute the nullspace, but: “Since the result will be given in the dense representation, both the nullity of A and the number of rows of A must both be reasonably small.” (MAGMA Documentation)
- LinBox actually does Gaussian elimination for you if you compute the rank using the NoReordering method. However, it kills the rows it doesn’t need anymore to be more memory friendly. Also, in my experiments the Gaussian elimination wasn’t much faster than SAGE’s for random sparse matrices over $GF(127)$. However, it can compute the rank more quickly by using “Symbolic Reordering” (paper). Added bonus: LinBox also does $GF(p^n)$.
- SAGE offers a sparse Gaussian elimination: “We use Gauss elimination, in a slightly intelligent way, in that we clear each column using a row with the minimum number of nonzero entries.”.
The other day I kinda liked the idea to apply William Stein’s integer matrix rational-echelonize-via-solve algorithm to this case. The (adapted) algorithm is as follows (most of it is due to William Stein):
- Compute $r = rank(A)$. This is cheaper than Gaussian elimination because we can use “Symbolic Reordering”.
- Compute the pivot columns of the transpose $A^t$ of $A$. We can convince “Symbolic Reordering” to give us these as well.
- Let $B$ be the submatrix of $A$ consisting of the rows corresponding to the pivot columns found in the previous step. Note that, aside from zero rows at the bottom, $B$ and $A$ have the same reduced row echelon form.
- Compute the pivot columns of $B$. Again, we may do this using “Symbolic Reordering”.
- Let $C$ be the submatrix of $B$ of pivot columns. Let $D$ be the complementary submatrix of $B$ of all all non-pivot columns. Use a solver (such as Wiedemann) to find the matrix $X$ such that $C X=D$ . I.e., solve a bunch of linear systems of the form $ Cx = v$ , where the columns of $ X$ are the solutions.
- Return the matrix $I || X$ where $I$ is the identity matrix of rank $r$.
This algorithm has complexity of two “Symbolic Reordering” applications, and $ncols - r$ applications of a matrix-vector solver. If “Symbolic Reordering” significantly outperforms Gaussian elimination (speed and memory-wise) and if $ncols - r$ is small and the solver fast, this might outperform straight-forward Gaussian elimination. The algorithm in SAGE notation is:
def echelon_form_via_solve(A):
r = A.rank() # Step 1: Compute the rank
if r == self.nrows():
B = A
else:
# Steps 2 and 3: Extract out a submatrix of full rank.
P = A.transpose().pivots()
B = A.matrix_from_rows(P)
# Step 4: Now we instead worry about computing the reduced row echelon form of B.
pivots = B.pivots()
# Step 5: Apply solver
C = B.matrix_from_columns(pivots)
pivots_ = set(pivots)
non_pivots = [i for i in range(B.ncols()) if not i in pivots_]
D = B.matrix_from_columns(non_pivots)
X = C.solve_right(D, algorithm="LinBox:Blackbox")
R = self.parent()()
for i in range(len(pivots)): R[i,pivots[i]] = 1
for i in range(X.nrows()):
for j in range(X.ncols()):
R[i,non_pivots[j]] = X[i,j]
return R
However, as we have to call the solver repeatetly (or find a good matrix-matrix solver) I lost interest in implementing this thing. What is missing, is to hack LinBox to return the pivot columns when performing InPlaceLinearPivoting which is “Symbolic Reordering”. As a side product of my attempts some operations on sparse matrices over GF(p) are way faster now in SAGE (See #655). Also, if you - my dear reader - know about any fast implementation for that problem, please let me know.

