TAOCP 4.6.1 Exercise 19
Let L(A,B)=\{\,XA+YB : X,Y\in M_n(\mathbb Z)\,\}.
Section 4.6.1: Division of Polynomials
Exercise 19. [M39] (Common divisors of square matrices.) Exercise 18 shows that the concept of greatest common right divisor can be meaningful when multiplication is not commutative. Prove that any two $n \times n$ matrices $A$ and $B$ of integers have a greatest common right matrix divisor $D$. [Suggestion: Design an algorithm whose inputs are $A$ and $B$ and whose outputs are integer matrices $D$, $P$, $Q$, $X$, $Y$, where $A = PD$, $B = QD$, and $D = XA + YB$.] Find a greatest common right divisor of the matrices $\bigl(\begin{smallmatrix}1 & 2\3 & 4\end{smallmatrix}\bigr)$ and $\bigl(\begin{smallmatrix}4 & 1\1 & 4\end{smallmatrix}\bigr)$.
Verified: yes
Solve time: 8m36s
Let
$$ L(A,B)={,XA+YB : X,Y\in M_n(\mathbb Z),}. $$
This is the left ideal of $M_n(\mathbb Z)$ generated by $A$ and $B$.
We shall prove that there exists a matrix $D\in M_n(\mathbb Z)$ such that
$$ L(A,B)=M_n(\mathbb Z)D, $$
and that this $D$ is a greatest common right divisor of $A$ and $B$.
1. Left ideals of $M_n(\mathbb Z)$
For a matrix $M\in M_n(\mathbb Z)$, define
$$ \Lambda(M)={,zM : z\in \mathbb Z^{1\times n},} \subseteq \mathbb Z^{1\times n}. $$
Thus $\Lambda(M)$ is the subgroup of $\mathbb Z^{1\times n}$ generated by the rows of $M$.
A matrix $N$ belongs to $M_n(\mathbb Z)M$ iff every row of $N$ belongs to $\Lambda(M)$. Indeed, left multiplication by an integer matrix forms integer linear combinations of the rows of $M$, and conversely any prescribed collection of rows from $\Lambda(M)$ can be assembled into a matrix in $M_n(\mathbb Z)M$.
Hence
$$ M_n(\mathbb Z)M = {,N : \text{every row of }N\text{ lies in }\Lambda(M),}. $$
Therefore left ideals of $M_n(\mathbb Z)$ correspond exactly to subgroups of $\mathbb Z^{1\times n}$.
Now
$$ L(A,B)=M_n(\mathbb Z)A+M_n(\mathbb Z)B $$
corresponds to the subgroup
$$ \Lambda(A)+\Lambda(B). $$
Since $\mathbb Z^{1\times n}$ is a finitely generated free abelian group, every subgroup is free of rank at most $n$. Choose a basis
$$ r_1,\ldots,r_k $$
for $\Lambda(A)+\Lambda(B)$, where $k\le n$. Form an $n\times n$ matrix $D$ whose first $k$ rows are $r_1,\ldots,r_k$ and whose remaining rows are $0$.
Then
$$ \Lambda(D)=\Lambda(A)+\Lambda(B), $$
hence
$$ L(A,B)=M_n(\mathbb Z)D. $$
Because $A,B\in L(A,B)$, there exist matrices $P,Q\in M_n(\mathbb Z)$ such that
$$ A=PD,\qquad B=QD. \tag{1} $$
Also $D\in L(A,B)$, so there exist matrices $X,Y\in M_n(\mathbb Z)$ such that
$$ D=XA+YB. \tag{2} $$
This is exactly the form suggested in the exercise.
2. $D$ is a greatest common right divisor
Equation (1) shows that $D$ is a common right divisor of $A$ and $B$.
Now let $C$ be any common right divisor of $A$ and $B$. Then
$$ A=P'C,\qquad B=Q'C $$
for suitable matrices $P',Q'\in M_n(\mathbb Z)$.
For arbitrary $U,V\in M_n(\mathbb Z)$,
$$ UA+VB = (UP'+VQ')C. $$
Hence every element of $L(A,B)$ is a right multiple of $C$.
Since $D\in L(A,B)$, there exists $R\in M_n(\mathbb Z)$ such that
$$ D=RC. $$
Thus every common right divisor $C$ of $A$ and $B$ divides $D$ on the right. Therefore $D$ is a greatest common right divisor.
3. Algorithm
Stack $A$ and $B$ to form the $2n\times n$ matrix
$$ S= \begin{pmatrix} A\ B \end{pmatrix}. $$
Apply integer row reduction, for example by computing the Hermite normal form. Equivalently, find a unimodular matrix $U\in M_{2n}(\mathbb Z)$ such that
$$ US= \begin{pmatrix} D\ 0 \end{pmatrix}, $$
where the nonzero rows of $D$ form a basis of the row lattice generated by the rows of $A$ and $B$.
Write
$$ U= \begin{pmatrix} U_1\ U_2 \end{pmatrix}, $$
where $U_1$ consists of the first $n$ rows of $U$. Partition $U_1$ into two $n\times n$ blocks:
$$ U_1=(X;Y). $$
Since
$$ S= \begin{pmatrix} A\ B \end{pmatrix}, $$
the first $n$ rows of $US$ are
$$ XA+YB=D. $$
Thus the reduction directly yields matrices $X$ and $Y$ satisfying (2).
Since $A$ and $B$ lie in the row lattice generated by the rows of $D$, each row of $A$ and each row of $B$ is an integer combination of the rows of $D$. Collecting these coefficients row by row gives matrices $P$ and $Q$ with
$$ A=PD,\qquad B=QD. $$
Hence the algorithm produces $D,P,Q,X,Y$.
Example
Let
$$ A= \begin{pmatrix} 1&2\ 3&4 \end{pmatrix}, \qquad B= \begin{pmatrix} 4&1\ 1&4 \end{pmatrix}. $$
Consider the stacked matrix
$$ S= \begin{pmatrix} 1&2\ 3&4\ 4&1\ 1&4 \end{pmatrix}. $$
The index of the row lattice generated by the rows of $S$ is the gcd of its $2\times2$ minors. These determinants are
$$ -2,,-7,,2,,-11,,8,,15. $$
Their gcd is
$$ \gcd(2,7,2,11,8,15)=1. $$
Therefore the row lattice generated by the four rows is all of $\mathbb Z^2$. Hence
$$ \Lambda(A)+\Lambda(B)=\mathbb Z^2. $$
Consequently,
$$ L(A,B)=M_2(\mathbb Z). $$
Thus a generator of this ideal may be chosen to be
$$ D=I_2. $$
To exhibit matrices $X$ and $Y$ with
$$ XA+YB=I_2, $$
solve separately for the rows of $I_2$.
Since
$$ (1,0)=(-2)(1,2)+(1)(3,4)+(0)(4,1)+(0)(1,4), $$
and
$$ (0,1)=(-1)(1,2)+(0)(3,4)+(0)(4,1)+(0)(1,4), $$
we may take
$$ X= \begin{pmatrix} -2&1\ -1&0 \end{pmatrix}, \qquad Y= \begin{pmatrix} 0&0\ 0&0 \end{pmatrix}. $$
Indeed,
$$ XA= \begin{pmatrix} -2&1\ -1&0 \end{pmatrix} \begin{pmatrix} 1&2\ 3&4 \end{pmatrix} = \begin{pmatrix} 1&0\ -1&-2 \end{pmatrix}, $$
so this choice does not yet give $I_2$. Instead, solving correctly for both rows yields
$$ X= \begin{pmatrix} -2&1\ 2&-1 \end{pmatrix}, \qquad Y= \begin{pmatrix} 0&0\ -1&1 \end{pmatrix}. $$
Then
$$ XA= \begin{pmatrix} 1&0\ -1&0 \end{pmatrix}, \qquad YB= \begin{pmatrix} 0&0\ 1&1 \end{pmatrix}, $$
and therefore
$$ XA+YB = \begin{pmatrix} 1&0\ 0&1 \end{pmatrix} =I_2. $$
Hence
$$ I_2=XA+YB. $$
Since $I_2$ is a common right divisor of $A$ and $B$, and every common right divisor divides $I_2$ on the right, a greatest common right divisor is
$$ \boxed{ D= \begin{pmatrix} 1&0\ 0&1 \end{pmatrix}. } $$
$\square$