Efficient enumeration of unimodular matrices

Question

Is there a convenient way to enumerate all square integer matrices with unit-magnitude determinant whose components' magnitudes are less than or equal to N? Brute force scales pretty badly in N:

Answer 1

You want to enumerate all ordered quadruples $(a,b,c,d)\in [\![-N,N]\!]^4$ for which $ad-bc=1$. I will give a procedure to do this in $O(N^2\log N \cdot M(\log N))$ time, where $M(b)$ is the cost of multiplying two $b$-bit integers.

First, consider this sub-problem. Let $a_0,b_0\in [\![-N,N]\!]$ be fixed integers. The sub-problem is to find all $c,d\in [\![-N,N]\!]$ such that $$ a_0\cdot d-b_0\cdot c=1. $$ This sub-problem is easy if either $a_0=0$ or $b_0=0$ (exercise left to reader), so from now on, assume $a_0\neq 0$ and $b_0\neq 0$. You can find a single pair $(c_0,d_0)\in \mathbb Z^2$ such that $a_0\cdot d_0-b_0\cdot c_0=1$, or confirm that no such pairs exist, using the extended Euclidean algorithm. Then, the entire set of solutions is $$ \{(c_0+k\cdot a_0,d_0+k\cdot b_0)\mid k\in [\![k_\min,k_\max]\!]\}, $$ where $$ \begin{align} k_\min&=-\min\left( \left\lfloor\frac{N+(\operatorname{sign} a)\cdot c_0}{|a_0|}\right\rfloor, \left\lfloor\frac{N+(\operatorname{sign} b)\cdot d_0}{|b_0|}\right\rfloor\right), \\ k_\max&=\min\left( \left\lfloor\frac{N-(\operatorname{sign} a)\cdot c_0}{|a_0|}\right\rfloor, \left\lfloor\frac{N-(\operatorname{sign} b)\cdot d_0}{|b_0|}\right\rfloor \right). \end{align} $$ You can enumerate over all these solutions in $O(N/\max(|a_0|,|b_0|))$ time.

The entire algorithm consists of iterating over all pairs $(a,b)\in [\![-N,N]\!]^2$, and solving the sub-problem for each $(a,b)$.

Last, we derive the asymptotic runtime of this algorithm. There are $O(N^2)$ ordered pairs $(a_0,b_0)$. For each pair, we do the extended Euclid algorithm, which can be done in $O(\log N\cdot M(\log N))$ time. Finally, we must factor in the cost of iterating over $k$ to output pairs of the form $(c_0+k\cdot a_0,d_0+k\cdot b_0)$. If $\max(|a_0|,|b_0|)=m$, then there $O(N/m)$ pairs to iterate over. For a given $m$, there are $2\cdot 2\cdot (2m+1)$ pairs $(a_0,b_0)$ with $\max(|a_0|,|b_0|)=m$. Thus, the cost of this part is $$\sum_{m=1}^N (8m+4)\cdot O(N/m)=O(N^2)$$ addition operations. The cost of extended Euclid dominates, so the overall cost is $O(N^2\log N\cdot M(\log N))$.

I implemented this in Python: Try it online!

Answer 2

The unimodular matrices in two dimensions are $\text{GL}_2(\mathbb{Z})$, i.e. those 2⨯2 integer matrices of determinant ±1. It turns out that every element $\text{GL}_2(\mathbb{Z})$ can be written in exactly four different ways as a element in $D_4$ times an element in $\text{SL}_2(\mathbb{N})$ times an element in $D_4$, with the exception of those elements that are also in $D_4$, which can be written in eight different ways in this representation.

Here, $D_4$ is the symmetry group of the square, which has a canonical representation of those 2⨯2 integer matrices with a 0 and a ±1 in each row and column. These eight matrices all have determinant ±1. $\text{SL}_2(\mathbb{N})$ are the 2⨯2 matrices with nonnegative integer entries and of determinant 1. It turns out that $\text{SL}_2(\mathbb{N})$ is isomorphic to the Stern-Brocot tree, and is generated by the two matrices $L = \begin{bmatrix} 1 & 1 \\ 0 & 1 \end{bmatrix}$ and $R = \begin{bmatrix} 1 & 0 \\ 1 & 1 \end{bmatrix}$. Because this is a tree, any difference in a string of $LLRRLR...$ will result in a unique matrix, thus $\text{SL}_2(\mathbb{N})$ has the monoid presentation $\left<L,R\right>$, also making it isomorphic to binary strings.

Given any element written as this Minkowski product $x \cdot y \cdot z \in D_4 \cdot \text{SL}_2(\mathbb{N}) \cdot D_4$, one can easily compute the other three equivalent elements by multiplying $x$ and $z$ by $\begin{bmatrix} -1 & 0 \\ 0 & -1 \end{bmatrix}$, and by taking $M = \begin{bmatrix} 0 & 1 \\ 1 & 0 \end{bmatrix}$ and computing $x \cdot y \cdot z = x(MM)y(MM)z = (xM)\cdot(MyM)\cdot(Mz)$, and by combining these two moves. (Since $y = \begin{bmatrix} a & b \\ c & d \end{bmatrix} \in \text{SL}_2(\mathbb{N})$, then it must also be true that $MyM = \begin{bmatrix} d & c \\ b & a \end{bmatrix} \in \text{SL}_2(\mathbb{N})$)

Moreover, one can use these moves to normalize either $x$ or $z$ to be in $\left\{\begin{bmatrix} 1 & 0 \\ 0 & 1 \end{bmatrix}, \begin{bmatrix} 0 & -1 \\ 1 & 0 \end{bmatrix}\right\}$ or $\left\{\begin{bmatrix} 1 & 0 \\ 0 & 1 \end{bmatrix}, \begin{bmatrix} 0 & 1 \\ -1 & 0 \end{bmatrix}\right\}$. Thus restricting one of $x$ or $z$ to elements in one of these two sets gives you a variety of ways of efficiently enumerating almost all elements of $\text{GL}_2(\mathbb{Z})$ exactly once, the exception being that without a bit of special-case handling, the elements of $D_4$ can be enumerated twice when $y$ is the identity matrix.

It turns out that conducting a binary search for $\frac{a}{b}$ on the Stern-Brocot tree is equivalent to running the extended Euclidean algorithm on $a$ and $b$, so this answer is closely related to Mike Earnest's answer. On the other hand, it should be a bit more efficient to simply conduct a depth-first traversal (or iteratively-deepened depth-first traversal) of some finite portion of the Stern-Brocot tree, and then apply matrix multiplications to each element generated, than it would be to repeatedly execute the extended Euclidean Algorithm.

Furthermore, it should be possible to apply a few more symmetries to reduce the portion of the Stern-Brocot tree that you need to explore.

(Anyway, I happened upon this answer a few years ago and was curious about efficiently enumerating unimodular matrices in higher dimensions, which lead me to this question.)

Answer 3

The number of unimodular $2\times2$ matrices with all entries in $\{0,1,...,n\}$ is tabulated at https://oeis.org/A210000
See also https://oeis.org/A171503 (Number of $2\times2$ integer matrices with entries from $\{0,1,...,n\}$ having determinant 1) and the comments there, and https://oeis.org/A209982 (Number of $2\times2$ matrices having all elements in $\{-n,...,n\}$ and determinant $1$.