<figure><img src="image_1.png" /></figure>

HAL
open science

# Pivotal condensation and chemical balancing

Hans-Christian Herbig

► To cite this version:

Hans-Christian Herbig. Pivotal condensation and chemical balancing. 2023. hal-04091017

HAL Id: hal-04091017

https://hal.science/hal-04091017

Preprint submitted on 7 May 2023

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_2.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
# PIVOTAL CONDENSATION AND CHEMICAL BALANCING

HANS-CHRISTIAN HERBIG

ABSTRACT. I present a universal method, called pivotal condensation, for calculating stoichiometric factors of chemical reactions. It can be done by hand, even for rather complicated reactions. The main trick, which I call kernel pivotal condensation (ker pc), to calculate the kernel of a matrix might be of independent interest. The discussion is elaborated for matrices with entries in a principal ideal domain $R$. The ker pc calculates a basis with coefficients in $R$ for the kernel of a matrix, seen as a $Q$-vector space, where $Q$ is the quotient field of $R$. If $W$ is a free saturated $R$-submodule of $R^n$ I address the question how to modify an $R$-basis of the $Q$-vector subspace $Q \otimes_R W$ over the quotient field $Q$ to obtain a basis of the $R$-module $W$. I also indicate how one can solve inhomogeneous linear systems, invert matrices and determine the four subspaces using pivotal condensation. I formulate the balancing by inspection method that is widely used to reduce the size of a linear system arising in chemical balancing in mathematical language.

## CONTENTS

<table><tr><td>1. Introduction</td><td>1</td></tr><tr><td>2. Pivotal condensation</td><td>2</td></tr><tr><td>3. Determinant</td><td>3</td></tr><tr><td>4. Kernel</td><td>3</td></tr><tr><td>5. Consistency checks</td><td>7</td></tr><tr><td>6. Chemical balancing</td><td>7</td></tr><tr><td>7. Quivered ker pc</td><td>12</td></tr><tr><td>8. Miscellaneous investigations</td><td>14</td></tr><tr><td>8.1. Inhomogeneous systems</td><td>14</td></tr><tr><td>8.2. Matrix inversion</td><td>15</td></tr><tr><td>8.3. The four subspaces</td><td>16</td></tr><tr><td>9. Proof of ker pc and 4pc</td><td>18</td></tr><tr><td>10. The Kernel as a free R-module</td><td>20</td></tr><tr><td>10.1. A practical guide to Smith normal form</td><td>21</td></tr><tr><td>10.2. Saturated submodules of R<sup>n</sup></td><td>26</td></tr><tr><td>10.3. Turning a Q-basis into an R-basis</td><td>26</td></tr><tr><td>References</td><td>28</td></tr></table>

## 1. INTRODUCTION

While teaching the basic linear algebra course to chemical engineering students I noticed that for chemical balancing (see, e.g., [Fin09, MKAM15] and references therein) the naive methods, in particular balancing by inspection (see [Fin09, Subsection 13.2]), acquired in the chemistry classes are for all practical purposes way more efficient than applying Gaussian elimination to the corresponding homogeneous linear system. For complicated reactions the naive methods can become cumbersome and require skill. The question arose if there is a universal method that is relatively quick and safe for a human computer. I am presenting such a method that works over the integers. It is inspired by the pivotal condensation of Felice Chiò [Chi53] for calculating determinants of square matrices.

<footnote n="1991">1991 Mathematics Subject Classification. Primary 15A06; Secondary 11D04, 80A32.</footnote>

<footnote>Key words and phrases. Chemical balancing, homogeneous linear systems, pivotal condensation, systems of linear Diophantine equations, matrix inversion, four subspaces.</footnote>
In chemical balancing it sometimes makes sense to balance chemical reactions that depend on a parameter, say $x$. This is for example the case if one has families of compounds. In this situation one has to do linear algebra over the gcd domain $R = \mathbb{Z}[x]$, which is not even a principal ideal domain (PID). The suggested algorithms for the kernel of a matrix are elaborated for the case when $R$ is a gcd domain (e.g. a PID). Readers unfamiliar with abstract algebra can assume that the base ring is $R = \mathbb{Z}$ with quotient field $Q = \mathbb{Q}$. I emphasize that the pivotal condensation provides a basis for the kernel of a matrix with entries in $R$ that has coefficients in $R$ seen as a vector space over the quotient field $Q$. But if the dimension of the kernel is $> 1$ it may happen that the basis is not a basis of the kernel seen as a free $R$-module. This subtlety is typically swept under the rug in the literature on chemical balancing and for the convenience of the reader I spell out in the more mathematical Section 10 how this can be handled in practical terms when $R$ is a PID. The only other work that I am aware of that views the chemical balancing problem as a system of linear Diophantine equations is [PV06]. Their empirical analysis relies on computer implementations of certain algorithms. However, I am trying here to convince the reader that all practically relevant examples can be elaborated by a patient human on a sheet of paper.

The plan of the paper is as follows. In Section 2 I review the ideas of pivotal condensation. In Subsection 3 I recall the pivotal condensation method to calculate the determinant $|A|$ of a square matrix $A \in R^{n \times n}$. In Subsection 4 I propose a variation of pivotal condensation that allows to re-purpose parts the findings of 3 to determine a basis for the kernel of the transposed matrix $A^T$ in case $A$ turns out to be singular. For non-square matrices $A^T$ I simplify the scheme notationally and calculate merely $\ker(A)$. In Subsection 5 I propose consistency checks that can be utilized to localize computational errors. With this preparation we are ready to turn in Section 6 to chemical balancing and work out some examples. In Section 7 I formulate the balancing by inspection method in mathematical language. In Section 8 I show how to solve inhomogeneous systems, invert matrices and determine the four subspaces using pivotal condensation. In the Section 9 I give a formal justification of the pivotal condensation scheme of Subsection 4. In Section 10 I review how $R$-module bases can be calculated by hand. This includes a discussion of the saturation of an $R$-module, a practical guideline to Smith normal form and a quick method to post-process the result of the pivotal condensation the obtain an $R$-module basis for the kernel.

The practitioners among the readership may focus on Sections 5, 6 and, possibly, 8 and ignore the rest of the paper. To do kernel pivotal condensation (ker pc) and, for example, matrix inversion with ker pc, by hand requires only the knowledge of arithmetic acquired in 5th grade and some patience. One does not need to use the calculus of fractions. To have a completely satisfactory answer the results of ker pc have sometimes to be post-processed by a Smith normal form calculation (see Subsection 10.3). This requires division with remainder. I give a detailed account of the example calculations for those readers with a moderate mathematical background who just want to learn the method.

*Acknowledgements.* I would like to thank the undergraduate crowd at the Escola de Química da UFRJ for helping me to think about the matter, mostly by not paying attention to what I said. The paper is an attempt towards a more just and equitable chemistry. Let us bring the wonders of stoichiometry even to the remotest places of the Amazon rain forest without polluting them with consumer electronics. Those who have seen the movie *The Great Race* know that 'Push the button, Max!' is an elitist and dangerous attitude. I would like to thank Benjamin Briggs for pointing out to me the notion of saturatedness. The whole project was triggered by a chat [Gro] that I had with Ben Grossman on stackexchange.com. Thanks to Daniel Herden and Chris Seaton for saving me from embarrassments and Ihsen Yengui for a clarification.

## 2. PIVOTAL CONDENSATION

Chiò's pivotal condensation (pc) is based on the following determinantal identity (for some historical account of the subject see [Abe14]).
**Theorem 2.1.** Let $n \ge 2$ be an integer, $R$ an integral domain with quotient field $Q$ and $A = [a_{ij}]_{i,j=1,...,n} \in R^{n \times n}$ with $a_{11} \ne 0$. Then

$$ (2.1) \qquad |A| = \frac{1}{a_{11}^{n-2}} \left| \begin{array}{c|c} \begin{array}{cc} a_{11} & a_{12} \\ a_{21} & a_{22} \end{array} & \begin{array}{cc} a_{11} & a_{13} \\ a_{21} & a_{23} \end{array} \\ \hline \begin{array}{cc} a_{11} & a_{12} \\ a_{31} & a_{32} \end{array} & \begin{array}{cc} a_{11} & a_{13} \\ a_{31} & a_{33} \end{array} \\ \hline \vdots & \vdots \\ \begin{array}{cc} a_{11} & a_{12} \\ a_{n2} & a_{n2} \end{array} & \begin{array}{cc} a_{11} & a_{13} \\ a_{n2} & a_{n3} \end{array} \end{array} \right| \begin{array}{c} a_{11} \\ a_{21} \\ \vdots \\ a_{n1} \\ a_{n2} \end{array} $$

Formally this expression makes only sense in $Q$, but of course $|A| \in R$. A straight forward proof of Equation (2.1) in terms of row operations can be found in [Eve12, Paragraph 3.6.1]. I am reproducing here an argument from [Abe14] generalizing Equation (2.1). The generalization is based on a block matrix decomposition of Sylvester:

$$ A = \begin{bmatrix} A_{11} & A_{12} \\ A_{21} & A_{22} \end{bmatrix} = \begin{bmatrix} A_{11} & \mathbf{0} \\ A_{21} & \mathbf{1} \end{bmatrix} \begin{bmatrix} \mathbf{1} & A_{11}^{-1}A_{12} \\ \mathbf{0} & A_{22} - A_{21}A_{11}^{-1}A_{12} \end{bmatrix}. $$

Here $n = k+l$, $A_{11} \in R^{k \times k}$ is assumed to be invertible, $A_{12} \in R^{k \times l}$, $A_{21} \in R^{l \times k}$ and $A_{22} \in R^{l \times l}$. Taking determinants we deduce $|A| = |A_{11}| |A_{22} - A_{21}A_{11}^{-1}A_{12}|$. For $k=1$ this means $|A| = A_{11}^{-n+2} |A_{11}A_{22} - A_{21}A_{12}|$ since $A_{11}$ is central, proving (2.1).

### 3. DETERMINANT

The pc for evaluating determinants, referred to as *Chiò's pc*, consists of repeated applications of the following loop in order to reduce the size of the determinant.

**Step 1 (preparation).** If $a_{11} = 0$ search for an $i$ with $a_{i1} \ne 0$. Switch the $i$th row and the 1st row and keep track of the sign. If $a_{i1} = 0$ for all $i = 1, ..., n$ then stop; in fact, $|A| = 0$.

**Step 2 (pc).** Apply formula (2.1) to write $|A|$ as a scalar multiple of a $(n-1) \times (n-1)$-determinant $|A'|$. **Step 3 (cleaning up).** Pull out of the determinant $|A'|$ common factors of rows or columns (step 3 is not compulsory).

Continue the process recursively with the $(n-1) \times (n-1)$-determinant etc., i.e., apply repeatedly Steps 1,2,3, until obtaining a $2 \times 2$-determinant, which in turn is to be evaluated. For clarity I put in Step 1 a box around the pivot $a_{i1}$.

**Example 3.1.** To give an example, let us calculate the area $A$ of the triangle with the sides 2, 3, 4 using the Cayley-Menger determinant (see, e.g., [Fie11, Appendix A.1])

$$ \begin{vmatrix} 0 & 1 & 1 & 1 \\ 1 & 0 & 2^2 & 3^3 \\ 1 & 2^2 & 0 & 4^2 \\ 1 & 3^2 & 4^2 & 0 \end{vmatrix} \stackrel{\text{Step 1}}{=} - \begin{vmatrix} \begin{vmatrix} \mathbf{1} & 9 & 16 & 0 \\ 1 & 0 & 4 & 9 \\ 1 & 4 & 0 & 16 \\ 0 & 1 & 1 & 1 \end{vmatrix} & \begin{vmatrix} -9 & -12 & 9 \\ -5 & -16 & 16 \\ 1 & 1 & 1 \end{vmatrix} & \begin{vmatrix} -3 & -4 & 3 \\ -5 & -16 & 16 \\ 1 & 1 & 1 \end{vmatrix} \\ \hline \begin{vmatrix} 3 & -33 & 28 & -6 \\ -(-3)^{3-2} & 1 & 1 & -11 \end{vmatrix} & \begin{vmatrix} 3 & 28 & -11 & 1 \\ -3 & 1 & 1 & -2 \end{vmatrix} & \begin{vmatrix} 3(-45) & -135 & -2^{2} & 2^{2}A^{2} \end{vmatrix} \end{vmatrix} \stackrel{\text{Step 2}}{=} - \frac{3}{(-3)^{3-2}} \begin{vmatrix} 28 & -33 & 28 & -11 \\ 1 & -6 & 1 & -2 \end{vmatrix} $$

Here the ring $R$ can be taken to be the integers $\mathbb{Z}$. We conclude that $A = \sqrt{135}/4$.

The scheme provides a compact and safe means to evaluate determinants by hand. Typically the symmetry of the matrix can be preserved in the process, cutting the workload almost in half. It combines excellently with Cramer's rule as most of the minors that appear in the denominator show up in the numerator as well.

### 4. KERNEL

When one is confronted with a quadratic matrix $A \in R^{n \times n}$ it is oftentimes unclear if it is regular or singular. In the latter case one is possibly expected to work out the kernel $ker_Q(A)$ as well. Throughout the paper we use the notation

$$ \ker_R(A) = \{ \vec{v} \in \mathbb{R}^n \mid A\vec{v} = \vec{0} \}, \quad \mathrm{im}_R(A) = \{ A\vec{v} \in \mathbb{R}^n \mid \vec{v} \in \mathbb{R}^m \} = AR^m, $$

$$ \ker_Q(A) = \{ \vec{v} \in \mathbb{Q}^n \mid A\vec{v} = \vec{0} \} \simeq Q \otimes_R \ker_R(A), \quad \mathrm{im}_Q(A) = \{ A\vec{v} \in \mathbb{Q}^n \mid \vec{v} \in \mathbb{Q}^m \} = AQ^m. $$
I present here a version of the pc, referred to as the *detker pc*, that allows to re-purpose findings of Chiò's pc for $|A^T|$ to determine the kernel $ker_Q(A)$ in case the matrix A turns out to be singular.

At the start of the $(l+1)$st condensation, $l = 0, 1, \dots, \text{rank}(A) - 2$, one has the pattern

$$ \Delta_l |X_l| Y_l, $$

where $\Delta_l \in Q$, $X_l \in R^{(n-l) \times (n-l)}$ and $Y_l \in R^{(n-l) \times n}$. Here $\Delta_l$ is the product of all the prefactors occuring in the previous steps of the pc, i.e., the signs arising from switching rows, the fractions from Equation (2.1), and the common factors pulled out from rows of the matrix. The matrix $X_l$ is the one arising at Chiò's $l$th condensation. The matrix $Y_l$ is capturing the relevant information of the elementary row operations that are involved. At the start of the recursion we put $\Delta_0 = 1$, $X_0 := A^T$ and $Y_0 := \mathbf{1}_n$. At the $(l+1)$st condensation the procedure consists of the following.

**Step 1 (preparation).** If $X_l = \mathbf{0}$ then stop; $\text{rank}(A^T) = l$ and $\text{ker}_Q(A^T) = Y_l^T R^{n-l}$. Let $k_{l+1} + 1$ be the smallest column index of a nonzero column in $X_l$. This means that the first $k_{l+1}$ columns of $X$ are zero. Pick a row index $i$ with $x_{ik_{l+1}} \neq 0$ and replace $\Delta_l \mapsto \Delta'_l = (-1)^{|\sigma|} \Delta_l$, $X_l \mapsto X'_l$, $Y_l \mapsto Y'_l$ where $X'_l$ and $Y'_l$ are obtained from $X_l$ and $Y_l$, respectively, by applying a row permutation $\sigma$ to move the $i$th row to the first row. The resulting pattern is $\Delta'_l |X'_l| Y'_l$. For clarity I put a box around the pivot $x_{ik_{l+1}}$.

**Step 2 (pc).** The matrix $Z_l'' := [X_l''|Y_l'']$ is obtained from $Z_l' := [X_l'|Y_l']$ as follows. The entry of $Z_l'' \in R^{(n-l-k_{l-1}) \times (2n-l-1)}$ at the $r$th row and $s$th column is

$$ (4.1) \qquad \begin{vmatrix} z_1 k_{l+1} + 1 & z_1 k_{l+1} + s \\ z_{r+1} k_{l+1} + 1 & z_{r+1} k_{l+1} + s \end{vmatrix}, $$

where $Z_l' = [z_{ij}]_{\substack{i=1,\dots,n-l \\ j=1,\dots,2n-l}}$. The resulting pattern is $\Delta_l'' | X_l'' | Y_l''$ with $\Delta_l'' := \Delta_l' z_1^{m+l+1} z_{k_{l+1}+1}^{-m+l+1}$.

**Step 3 (cleaning up).**<fnref n="1" /> Let $L_i$ be the $i$th row of $Z_l'' := [X_l''|Y_l'']$ and $c_i$ be the gcd of the entries of $L_i$ if $L_i$ is nonzero and put $c_i = 1$ otherwise. Note that for a general integral domain $R$ the existence of the gcd is not guaranteed. If $R$ is a PID it can be calculated using Euclid's algorithm. Let $Z_{l+1} := [X_{l+1}|Y_{l+1}]$ be the matrix whose rows are $L_i/c_i$ and put $\Delta_{l+1} := \Delta_l'' \prod_i c_i$. The resulting pattern is $\Delta_{l+1} | X_{l+1} | Y_{l+1}$, completing the $(l+1)$th condensation.

If $X_n \in R$ is nonzero the determinant is $|A| = \Delta_n |X_n|$ and $\text{ker}_Q(A^T)$ is trivial. In practical terms it seems reasonable to work out the $\Delta_l | X_l|$'s first and leave some blank space to fill in the $Y_l$'s later in case the matrix turns out to be singular.

**Example 4.1.** Let $R = \mathbb{Z}$ and consider the following pc

$$ \left| \begin{array}{c|cccc|cccc} \boxed{1} & 2 & 3 & 4 & 1 & 0 & 0 & 0 \\ 5 & 6 & 7 & 8 & 0 & 1 & 0 & 0 \\ 9 & 10 & 11 & 12 & 0 & 0 & 1 & 0 \\ \hline 13 & 14 & 15 & 16 & 0 & 0 & 0 & 1 \end{array} \right\} \mapsto \left| \begin{array}{c|cccc} \boxed{-4} & -8 & -12 & -5 & 1 & 0 & 0 \\ -8 & -16 & -24 & -9 & 0 & 1 & 0 \\ -12 & -24 & -36 & -13 & 0 & 0 & 1 \end{array} \right\} \mapsto \left| \begin{array}{c|cccc} 0 & 0 & -4 & 8 & -4 & 0 & 0 \\ 0 & 0 & -8 & 12 & 0 & -4 & 0 \end{array} \right\} \mapsto \left| \begin{array}{c|cccc} 0 & 0 & -1 & 2 & -1 & 0 & 0 \\ 0 & 0 & -2 & 3 & 0 & -1 & 0 \end{array} \right\} $$

It shows that $\text{ker}_Q \left( \begin{pmatrix} 1 & 5 & 9 & 13 \\ 2 & 6 & 10 & 14 \\ 3 & 7 & 11 & 15 \\ 4 & 8 & 12 & 16 \end{pmatrix} \right) = \mathbb{Q} \begin{bmatrix} -1 \\ 2 \\ -1 \\ 0 \end{bmatrix} \oplus \mathbb{Q} \begin{bmatrix} -2 \\ 3 \\ 0 \\ -1 \end{bmatrix}$, which by accident coincides with the kernel of the transposed matrix. As the gcd $\Delta_2$ of the $2 \times 2$-minors of

$$ \begin{bmatrix} -1 & -2 \\ 2 & 3 \\ -1 & 0 \\ 0 & -1 \end{bmatrix} $$

is a unit in $\mathbb{Z}$ the above generators also form a basis of the $\mathbb{Z}$-module $\text{ker}_Z(A)$ (see Theorem 10.9).

<footnote n="1">¹This is optional up to the last iteration.</footnote>
**Example 4.2.** Let $R$ be the ring of convergent power series $\mathbb{C}\{z\}$ in the variable $z$. The following pc

$$
\begin{align*}
& \left| \begin{array}{ccc|ccc}
\boxed{1} & \cos(z) & \cos(2z) & 1 & 0 & 0 \\
\cos(z) & 1 & \cos(z) & 0 & 1 & 0 \\
\cos(2z) & \cos(z) & 1 & 0 & 0 & 1
\end{array} \right| \\
\mapsto & \left| \begin{array}{ccc|ccc}
\sin^2(z) & \cos(z)(1-\cos(2z)) & -\cos(z) & 1 & 0 & 0 \\
\cos(z)(1-\cos(2z)) & \sin^2(2z) & -\cos(2z) & 0 & 1 & 1
\end{array} \right| \\
= & \left| \begin{array}{ccc|ccc}
\boxed{\sin^2(z)} & 2\cos(z)\sin^2(z) & -\cos(z) & 1 & 0 & \boxed{1} \\
2\cos(z)\sin^2(z) & 4\sin^2(z)\cos^2(z) & -\cos(2z) & 0 & 1 & \boxed{1}
\end{array} \right| \\
= & \left| \begin{array}{ccc|ccc}
\sin^{-2}(z) & 4\sin^4(z)\cos^2(z) - 4\sin^4(z)\cos^2(z) & -\sin^2(z)\cos(2z) + 2\cos^2(z)\sin^2(z) & -2\cos(z)\sin^2(z) & \sin^2(z) & \sin^2(z) \\
\boxed{0} & \boxed{1} + \sin^2(z) & -2\cos(z) & \boxed{1} & \boxed{1} & \boxed{1}
\end{array} \right|
\end{align*}
$$

shows that $\ker_{\mathbb{C}\{z\}} \left( \begin{bmatrix} 1 & \cos(z) & \cos(2z) \\ \cos(z) & 1 & \cos(z) \\ \cos(2z) & \cos(z) & 1 \end{bmatrix} \right) = \mathbb{C}\{z\} \begin{bmatrix} 1 + \sin^2(z) \\ -2\cos(z) \\ 1 \end{bmatrix} \subset \mathbb{C}\{z\}^3$. Here we use the fact that since $\mathbb{C}\{z\}$ is a PID our vector generates $\ker_{\mathbb{C}\{z\}}$ if and only if it generates $\ker_Q$, where $Q$ is the field of convergent Laurent series (see Remark 10.10).

If the matrix $A^T \in R^{n \times m}$ is not quadratic I suggest to apply the same procedure to calculate $\ker_Q(A)$ while not annotating $\Delta_l$ and the first $|$. The scheme is referred to as the *ker pc*. The input of the $(l+1)$th condensation, $l = 0, 1, \dots, \text{rank}(A) - 2$, one is the pattern

$$X_l|Y_l,$$

where $X_l \in R^{(m-l) \times (n-l)}$ and $Y_l \in R^{(n-l) \times n}$. At the start of the recursion we put, $X_0 := A^T$ and $Y_0 := \mathbf{1}_n$. The $(l+1)$th condensation consists of the following procedure.

**Step 1 (preparation).** If $X_l = 0$ then stop; $\text{rank}(A^T) = l$ and $\ker_Q(A^T) = Y_l^T R^{n-l}$. Let $k_{l+1} + 1$ be the smallest column index of a nonzero column in $X_l$. This means that the first $k_{l+1}$ columns of $X$ are zero. Pick a row index $i$ with $x_{ik_{l+1}} \neq 0$ and replace $X_l \mapsto X'_l$, $Y_l \mapsto Y'_l$, where $X'_l$ and $Y'_l$ are obtained from $X_l$ and $Y_l$, respectively, by applying a row permutation $\sigma$ to move the $i$th row to the first row. For clarity I put a box around the pivot $x_{ik_{l+1}}$.

**Step 2 (pc).** The matrix $Z_l'' := [X_l''|Y_l'']$ is obtained from $Z_l' := [X_l'|Y_l']$ as follows. The entry of $Z_l'' \in R^{(m-l-k_{l-1}) \times (n+m-l-1)}$ at the $r$th row and $s$th column is

$$ (4.2) \qquad \left| \begin{array}{cc} z_1 k_{l+1} + 1 & z_1 k_{l+1} + s \\ z_{r+1} k_{l+1} + 1 & z_{r+1} k_{l+1} + s \end{array} \right|, $$

where $Z_l' = [z_{ij}]_{\substack{i=1,\dots,n-l \\ j=1,\dots,2n-l}}$. The resulting pattern is $X_l''|Y_l''$.

**Step 3 (cleaning up).**<fnref n="2" /> Let $L_i$ be the $i$th row of $Z_l'' := [X_l''|Y_l'']$ and $c_i$ be the gcd of the entries of $L_i$ if $L_i$ is nonzero and put $c_i = 1$ otherwise. Let $Z_{l+1} := [X_{l+1}|Y_{l+1}]$ be the matrix whose rows are $L_i/c_i$. The resulting pattern is $X_{l+1}|Y_{l+1}$, completing the $l$th condensation.

If $X_m$ is nonzero (and non-empty) then $\ker_Q(A^T)$ is trivial. Let us look into two examples; more examples are elaborated in Section 6.

<footnote n="2">²This is optional up to the last iteration.</footnote>
**Example 4.3.** Let $R$ be the polynomial ring $\mathbb{Z}[x, y, z]$ with quotient field $\mathbb{Q}(x, y, z)$ the field of rational functions in $x, y, z$ and coefficients in $\mathbb{Q}$. From the pc

$$
\begin{array}{c@{\hspace{1em}}c@{\hspace{1em}}c@{\hspace{1em}}c@{\hspace{1em}}c@{\hspace{1em}}c@{\hspace{1em}}c}
\begin{array}{c} \boxed{1} \\ x \\ x^2 \\ x^3 \end{array} & \begin{array}{c} 1 \\ y \\ y^2 \\ y^3 \end{array} & \begin{array}{c} 1 \\ z \\ z^2 \\ z^3 \end{array} & \begin{array}{cccc} 1 & 0 & 0 & 0 \\ 0 & 1 & 0 & 0 \\ 0 & 0 & 1 & 0 \\ 0 & 0 & 0 & 1 \end{array} & \begin{array}{c} \boxed{y-x} \\ y^2-x^2 \\ y^3-x^3 \end{array} & \begin{array}{c} z-x \\ z^2-x^2 \\ z^3-x^3 \end{array} & \begin{array}{cccc} -x & 1 & 0 & 0 \\ -x^2 & 0 & 1 & 0 \\ -x^3 & 0 & 0 & 1 \end{array} \\[1.5ex]
\mapsto & & & & & & \\[1.5ex]
& & & & & & \\[1.5ex]
\mapsto & (y-x)(z-x)(z-y) & \begin{array}{c} (y-x)xy \\ (y^2-x^2)xy \end{array} & \begin{array}{ccc} -(y^2-x^2) & y-x & 0 \\ -(y^3-x^3) & 0 & y-x \end{array} & \\[1.5ex]
= & & & & & & \\[1.5ex]
& & & & & & \\[1.5ex]
= & \boxed{(z-x)(z-y)} & \begin{array}{c} xy \\ (y+x)xy \end{array} & \begin{array}{ccc} -y-x & 1 & 0 \\ -y^2-xy-x^2 & 0 & 1 \end{array} & \\[1.5ex]
& & & & & & \\[1.5ex]
= & -(z-x)(z-y)xyz & (z-x)(z-y)(xy+yz+xz) & -(z-x)(z-y)(x+y+z) & (z-x)(z-y) \\[1.5ex]
= & -xyz & xy+yz+xz & -x-y-z & 1
\end{array}
$$

we see that $\ker_{\mathbb{Q}(x,y,z)} \left( \begin{pmatrix} 1 & x & x^2 & x^3 \\ 1 & y & y^2 & y^3 \\ 1 & z & z^2 & z^3 \end{pmatrix} \right) = \mathbb{Q}(x,y,z) \begin{bmatrix} -xyz \\ xy+yz+xz \\ -x-y-z \\ 1 \end{bmatrix} \subset \mathbb{Q}(x,y,z)^4$. This is a special case of the identities $\sum_{i=0}^n (-x_j)^{n-i} e_i(x_1, \dots, x_n) = 0$ for elementary symmetric functions $e_i(x_1, \dots, x_n)$, which can in turn easily be checked from the generating function $\prod_{i=1}^n (y-x_i) = \sum_{i=0}^n (-1)^i y^{n-i} e_i(x_1, \dots, x_n)$. I would like to point out some subtleties here. The ring $R = \mathbb{Z}[x, y, z]$ is not a PID and it is not guaranteed that the method produces a generating set for the whole kernel as an $R$-module. It turns out that for this it is sufficient to have a block of a permutation matrix in the matrix whose columns are our generators (more on this in Section 10). In the example at hand this is just the 1 in the last entry of the generator.

**Example 4.4.** To see that something can go wrong when $R$ is not a PID let us look into an example from physics. Let us consider $R = \mathbb{R}[q_1, q_2, q_3, p_1, p_2, p_3]$ and define

$$ \vec{J} = \vec{q} \times \vec{p} = \begin{bmatrix} q_2 p_3 - q_3 p_2 \\ q_3 p_1 - q_1 p_3 \\ q_1 p_2 - q_2 p_1 \end{bmatrix}. $$

In physics the vectors $\vec{q} = [q_1 \ q_2 \ q_3]^T$, $\vec{p} = [p_1 \ p_2 \ p_3]^T$ and $\vec{J}$ are called position, momentum and angular momentum, respectively. Let us apply $\ker \text{pc}$ to the $1 \times 3$-matrix $\vec{J}^T$:

$$ \boxed{\begin{array}{c|ccc} J_1 & 1 & \\ J_2 & 1 & \mapsto \begin{array}[b]{c} -J_2 & J_1 \\ -J_3 & J_1 \end{array} \\ J_3 & 1 & \end{array}} $$

The vectors $\vec{v}_1 := \begin{bmatrix} -J_2 \\ J_1 \\ 0 \end{bmatrix}$, $\vec{v}_2 := \begin{bmatrix} -J_3 \\ 0 \\ J_1 \end{bmatrix}$ generate a rank two free submodule of $\ker_R(\vec{J}^T)$. But evidently $\vec{q}$ and $\vec{p}$ generate another rank two free submodule of $\ker_R(\vec{J}^T)$, and for degree reasons $\vec{q}, \vec{p}$ cannot be in $R\vec{v}_1 + R\vec{v}_2$. People typically employ Gröbner bases for determining such kernels (see e.g. [GP08]). We will not discuss this topic here.

**Remark 4.5.** In this paper, the focus is on the case when $R = \mathbb{Z}$. Let me make some observations that arise when comparing the $\ker \text{pc}$ to Gaussian elimination of $A \in \mathbb{Z}^{n \times m}$:

(1) The $\ker \text{pc}$ of $[A^T|\mathbf{1}_m]$, being equivalent to the first half of Gaussian elimination, requires only half as many steps as the Gaussian elimination of $A$ for determining $\ker_Q(A)$. On the other hand, $[A^T|\mathbf{1}_m]$ is considerably bigger than $A$.

(2) In $\ker \text{pc}$ one can effortlessly avoid fractions<fnref n="3" />.

(3) During $\text{pc}$ the size of the matrices shrinks step-wise.

<footnote n="3">³According to L. Kronecker rational numbers are *Menschenwerk* and require as such a lot of attention, i.e., brain capacity.</footnote>
In fact, the problem is visible in the gcd $\Delta_2$ of the $2 \times 2$-minors of $X := \begin{bmatrix} -2 & -4 \\ -1 & -3 \\ 2 & 0 \\ 0 & 2 \end{bmatrix}$, which is 2. The

image $\mathrm{im}_Z(X) = XZ^2$ of $X$ cannot be saturated (see Theorem 10.9), and the columns of this matrix cannot form an $\mathbb{Z}$-basis for the kernel of $A$. If $R$ is a PID, $A \in R^{n \times m}$ is of rank $p$, $AX = 0$ and $X$ has rank $n - p$, then $\mathrm{ker}_R(A) = \mathrm{im}_R(X)$ if and only if the gcd $\Delta_{n-p}$ of the $(n-p) \times (n-p)$-minors is a unit. Sometimes $\mathrm{ker pc}$ produces generators whose $\Delta_{n-p}$ is not a unit. For details on how one can repair this see Section 10.

**Example 6.3.** Out of a head of a chemist come mostly incidence matrices that have just one chemical equation. For example, the Fischer-Tropsch reaction for the synthesis of paraffin involves the following incidence matrix

$$ (6.1) \qquad \begin{array}{c|cccc} & \mathrm{CO} & \mathrm{H}_2 & \mathrm{C}_n\mathrm{H}_{2n+2} & \mathrm{H}_2\mathrm{O} \\ \mathrm{C} & 1 & & n & \\ \mathrm{H} & & 2 & 2n+2 & 2 \\ \mathrm{O} & 1 & & 1 & \\ \end{array} $$

We work symbolically, i.e., over the ring $\mathbb{Z}[n]$. The $\mathrm{ker pc}$ is

$$ \begin{array}{c|cccccc} \boxed{1} & 1 & 1 & & & & \\ n & 2n+2 & & & & & \\ 2 & 1 & & & & & \\ & & & & & & \\ & & & & & & \\ & & & & & & \\ & & & & & & \\ \end{array} \quad \begin{array}{c|cccccc} \boxed{2} & & & & & & \\ 2n+2 & -n & & & & & \\ 1 & 2 & -n & & & & \\ & & & & & & \\ & & & & & & \\ & & & & & & \\ & & & & & & \\ \end{array} \quad \begin{array}{c|cccccc} 1 & & & & & & \\ 1 & 1 & -n & & & & \\ 2 & 1 & 1 & -n & & & \\ & & & & & & \\ & & & & & & \\ & & & & & & \\ & & & & & & \\ \end{array} \quad \begin{array}{c|cccccc} -2n & -2n-2 & 2 & & & & \\ -2 & 2 & & & & & \\ & & & & & & \\ & & & & & & \\ & & & & & & \\ & & & & & & \\ & & & & & & \\ \end{array} \\ = \boxed{\frac{n}{1}} \begin{array}{c|cccccc} n & n+1 & -1 & & & & \\ -1 & 1 & & & & & \\ & & & & & & \\ & & & & & & \\ & & & & & & \\ & & & & & & \\ & & & & & & \\ \end{array} \quad \begin{array}{c|cccccc} -n & -2n-1 & 1 & n & & & \\ 1 & 2 & & & & & \\ & & & & & & \\ & & & & & & \\ & & & & & & \\ & & & & & & \\ & & & & & & \\ \end{array} $$

and hence $n\mathrm{CO} + (2n+1)\mathrm{H}_2 \to \mathrm{C}_n\mathrm{H}_{2n+2} + n\mathrm{H}_2\mathrm{O}$. As $\begin{bmatrix} 0 & -n & -2n-1 & 1 & n \end{bmatrix}^T$ has a unit entry it actually generates the kernel as a $\mathbb{Z}[n]$-module (see Section 10)$^6$. We will come back to this example in Section 7.

**Example 6.4.** In ionic reactions I treat the electric charge as an atom, as in the following example:

$$ \begin{array}{c|ccccc} & \mathrm{H}_2\mathrm{O} & \mathrm{MnO}_4^- & \mathrm{SO}_3^{2-} & \mathrm{OH}^- & \mathrm{MnO}_2 \\ \mathrm{H} & 2 & & & 1 & \\ \mathrm{Mn} & & 1 & & & \\ \mathrm{S} & & & 1 & & \\ \mathrm{O} & 1 & 4 & 3 & 1 & 2 \\ \mathrm{charge} & -1 & -2 & -1 & -2 & -2 \\ \end{array} $$

If during $\mathrm{ker pc}$ in a row under the pivot there is a zero the corresponding row (with first zero deleted) is simply multiplied with the pivot. This in turn is then divided out. I take the liberty to combine the two steps and simply copy over the row, omitting the first entry. This occurs in the first step of the following $\mathrm{ker pc}$:

$$ \begin{array}{c|cccccc} \boxed{2} & 1 & 1 & & & & \\ 1 & 4 & -1 & 1 & & & \\ & 1 & 3 & -2 & 1 & & \\ 1 & 1 & -1 & 1 & & & \\ & 1 & 2 & & 1 & & \\ & 1 & 4 & -2 & & & \\ & & & & & & \\ & & & & & & \\ \end{array} \quad \begin{array}{c|cccccc} \boxed{1} & 4 & -1 & & & & \\ 1 & 3 & -2 & & & & \\ & 1 & 3 & -2 & 1 & & \\ & & & & & & \\ & & & & & & \\ & & & & & & \\ & & & & & & \\ \end{array} \quad \begin{array}{c|cccccc} 1 & & & & & & \\ 1 & 1 & -n & & & & \\ 2 & 1 & 1 & -n & & & \\ & & & & & & \\ & & & & & & \\ & & & & & & \\ & & & & & & \\ \end{array} \quad \begin{array}{c|cccccc} -2n & -2n-2 & 2 & & & & \\ -2 & 2 & & & & & \\ & & & & & & \\ & & & & & & \\ & & & & & & \\ & & & & & & \\ & & & & & & \\ \end{array} \\ = \boxed{\frac{n}{1}} \begin{array}{c|cccccc} n & n+1 & -1 & & & & \\ -1 & 1 & & & & & \\ & & & & & & \\ & & & & & & \\ & & & & & & \\ & & & & & & \\ & & & & & & \\ \end{array} \quad \begin{array}{c|cccccc} -n & -2n-1 & 1 & n & & & \\ 1 & 2 & & & & & \\ & & & & & & \\ & & & & & & \\ & & & & & & \\ & & & & & & \\ & & & & & & \\ \end{array} $$

$$ \begin{array}{c|cccccc} \boxed{1} & 3 & -2 & 1 & & & \\ 1 & -2 & -1 & 2 & & & \\ & -2 & 1 & -1 & 1 & & \\ 1 & 4 & -2 & 1 & & & \\ & & & & & & \\ & & & & & & \\ \end{array} \quad \begin{array}{c|cccccc} \boxed{1} & -2 & -1 & -1 & & & \\ -2 & 1 & & & & & \\ & 1 & -1 & -1 & 1 & & \\ & & & & & & \\ & & & & & & \\ & & & & & & \\ \end{array} \quad \begin{array}{c|cccccc} 1 & & & & & & \\ 1 & 1 & -n & & & & \\ & & & & & & \\ & & & & & & \\ & & & & & & \\ & & & & & & \\ & & & & & & \\ \end{array} \quad \begin{array}{c|cccccc} -2 & -2 & -3 & -1 & 1 & & \\ -1 & 1 & & & & & \\ & & & & & & \\ & & & & & & \\ & & & & & & \\ & & & & & & \\ & & & & & & \\ \end{array} \\ = \boxed{\frac{n}{2}} \begin{array}{c|cccccc} -2 & -1 & 4 & 1 & & & \\ 1 & -1 & -2 & 1 & & & \\ & & & & & & \\ & & & & & & \\ & & & & & & \\ & & & & & & \\ & & & & & & \\ \end{array} $$

$$ \begin{array}{c|cccccc} \boxed{-3} & -2 & -1 & 4 & 1 & & \\ 2 & 1 & -1 & -2 & 1 & & \\ & & & & & & \\ & & & & & & \\ & & & & & & \\ & & & & & & \\ \end{array} \quad \begin{array}{c|cccccc} 1 & 2 & 3 & -2 & -2 & -3 & \\ -1 & -2 & 1 & -1 & -1 & -3 & \\ & & & & & & \\ & & & & & & \\ & & & & & & \\ & & & & & & \\ & & & & & & \\ \end{array} $$

<footnote>$^6$The ring $\mathbb{Z}[n]$ is not a PID and hence Theorem 10.9 cannot be applied.</footnote>
$$
\begin{array}{r@{\;}c@{\;}l@{\quad}l@{\quad}l@{\quad}l}
 & \left[ \begin{array}{@{}c@{\;}c@{}}
 & -2 \\
 \hline
 3 & 4 \\
 -46 & 565
\end{array} \right] & -1 & & & & -1 & -1 & -1 \\
\hline
\mapsto & & & & & & & & \\
 & & 1 & -1 & & & & & \\
 & & & & 7 & 132 & 96 & -96 & -2 \\
 & & & & & & & & \\
\hline
\mapsto & -5 & & & -2 & 2 & & -3 & 1 & 3 \\
 & -1176 & & -168 & & & -14 & -264 & -146 & 146 & -46 & 4 \\
\hline
\mapsto & \left[ \begin{array}{@{}c@{}} -5 \\ -588 \end{array} \right] & & & -2 & 2 & & -3 & 1 & 3 \\
 & & -84 & & & -7 & 132 & -73 & 73 & -23 & 2 \\
\hline
\mapsto & & 420 & -1176 & 1176 & 35 & 660 & -1399 & 223 & 1879 & -10
\end{array}
$$

The balanced reaction is<fnref n="7" />

$$
\begin{align*}
& 10\text{Cr}_7\text{N}_{66}\text{H}_{96}\text{C}_{42}\text{O}_{24} + 1176\text{KMnO}_4 + 1399\text{H}_2\text{SO}_4 \\
& \rightarrow 420\text{CO}_2 + 1176\text{MnSO}_4 + 35\text{K}_2\text{Cr}_2\text{O}_7 + 660\text{KNO}_3 + 223\text{K}_2\text{SO}_4 + 1879\text{H}_2\text{O}.
\end{align*}
$$

As 4 of 8 atoms in this example occur only in 2 compounds respectively it would also make sense to apply
the quivering method of the next section.

7. QUIVERED KER PC

When balancing reactions chemists often exploit substitutions arising from all the atoms that only occur
in just two compounds. In this way it is often possible to avoid the formalism as the resulting linear systems
are quite simple. I am trying here to elaborate the chemist’s method as an algorithmic tool to preprocess
the linear system of the reaction. The preprocessing is referred to as *quivering* as there are a quivers (i.e.
directed graphs) involved. The difficulty is to find all of the obvious reasons to put variables to zero.

We are aiming to solve the system

$$
(7.1) \qquad \sum_{i \in I} a_{ij} v_j = 0
$$

for all $j \in J$ with $A = [a_{ij}] \in R^{I \times J}$. Here $I$ and $J$ are finite sets and $I, J$ are understood to carry a total order. For the quivering only the order on $J$ is relevant. The plan is to construct minimal subsets $\hat{I} \subseteq I$, $\hat{J} \subseteq J$ and the quivered linear system $\sum_{i \in \hat{I}} \hat{a}_{ij}v_j = 0$ for all $j \in \hat{J}$ such that all solutions of (7.1) arise from solutions of the quivered system by substitutions.

To this end I am introducing a process referred to as *pruning* of the subsets $I' \subseteq I$ and $J' \subseteq J$. For each integer $n \ge 0$ we put

$$
I'_{n} := \{i \in I' | n = |\{j \in J' | a_{ij} \neq 0\}|\}.
$$

For each $i \in I'_1$ there is a unique $j \in J'$ such that $a_{ij} \neq 0$ this defines a map $\phi' : I' \to J'$. We put
$I'' := I' \setminus (I'_0 \cup I'_1)$ and $J'' := J' \setminus \phi'(I')$. Now $I''_2 := \{i \in I'' | 2 = |\{j \in J'' | a_{ij} \neq 0\}|\}$ is a quiver. In fact, we
put $s(i)$ for $i \in I''_2$ to be the smaller of the two $j$'s such that $a_{ij} \neq 0$, the other one is called $t(i)$. This provides
source and target the maps $s, t: I''_2 \to J''$. A path $\gamma = (i_1, \dots, i_{|\gamma|})$ of length $|\gamma|$ in $I''_2$ is a sequence of edges
$i_1, \dots, i_{|\gamma|} \in I''_2$ such that for each $p = 1, \dots, |\gamma| - 1$ we have $t(i_p) = s(i_{p+1})$. Two paths $\gamma = (i_1, \dots, i_{|\gamma|})$
and $\gamma' = (i'_1, \dots, i'_{|\gamma'|})$ in $I''_2$ are called equivalent, $\gamma \sim \gamma'$, if

$$
(1) s(i_1) = s(i'_1) \text{ and } t(i_{|\gamma|}) = s(i'_{|\gamma'|}),
$$

(2) and

$$
(7.2) \qquad (-1)^{|\gamma|} \prod_{p=1}^{|\gamma|} \frac{a_{ip}s(i_p)}{a_{ip'}t(i_{p'})} = (-1)^{|\gamma'|} \prod_{q=1}^{|\gamma'|} \frac{a_{i'_q}s(i'_q)}{a_{i'_q'}t(i'_{q'})}.
$$

If (1) holds and (2) is violated they are called *inconsistent*, $\gamma \perp \gamma'$. The set of problematic edges is

$$
(7.3) \qquad I_2'' := \{i \mid i \text{ an edge of } \gamma \text{ or } \gamma', \gamma \perp \gamma'\}.
$$

In the end by Eqn. (7.1) all the $v_j$ for $j \in s(I_2'' \cup t(I_2'' \cup \phi'(I'))$ have to be put zero. The pruned sets are $\text{pr}(I') := I_2'' \setminus I_2''$ and $\text{pr}(J') := J'' \setminus (s(I_2'') \cup t(I_2''))$. The quiver $s,t:(\text{pr}(I'))_2 \to \text{pr}(J)$ has no oriented cycles, since an edge increases the order by construction. All paths in $s,t:(\text{pr}(I'))_2 \to \text{pr}(J')$ with the

<footnote n="7">⁷I now deserve a grade A from Professor Stout (Central College, Pella, IA) for the entire semester!</footnote>
same initial and terminal vertex are equivalent. Unfortunately, the pruning should be iterated in general, to get an optimal result. But starting with $I$, $J$ there must be a smallest $k$ such that $\mathrm{pr}^{k+1}(I) = \mathrm{pr}^k(I)$ and $\mathrm{pr}^{k+1}(J) = \mathrm{pr}^k(J)$; it is referred to as the *pruning depth $k$*.

Let $\gamma$ be a maximal path (possibly trivial) in the quiver $s, t : (\mathrm{pr}^k(I))_2 \to \mathrm{pr}^k(J)$. Let us denote by $P$ a set of representatives of equivalence classes of such paths. The quivered system is obtained as follows. We put $\hat{I} := \mathrm{pr}^k(I)$, $\hat{J} := \{\mathrm{in}(\gamma) | \gamma \in P\}$ and define for $i \in \hat{I}, j \in \hat{J}$:

$$ (7.4) \qquad b_{ij} := a_{ij} + \sum_{\gamma=(i_1, i_2, \dots, i_{|\gamma|}) \in P : |\gamma| \ge 1, \mathrm{in}(\gamma)=j} (-1)^{|\gamma|} a_{it(i_{|\gamma|})} \prod_{p=1}^{|\gamma|} \frac{a_{ip^s(i_p)}}{a_{ip^t(i_p)}} \in Q. $$

We get a matrix $B = [b_{ij}] \in Q^{\hat{I} \times \hat{J}}$. We multiply the rows in $B$ appropriately with scalars in $R$ to obtain a matrix $\hat{A} \in R^{\hat{I} \times \hat{J}}$.

To reconstruct a solution of (7.1) we proceed as follows. By ker pc we find $\mathrm{ker}_Q \hat{A} = \hat{W}R^m$ for some matrix $\hat{W} = [w_{j\mu}] \in R^{\hat{J} \times m}$ with $m := |\hat{J}| - \mathrm{rank}(\hat{A}) = |J| - \mathrm{rank}(A)$. We put for $\mu = \{1, \dots, m\}$

$$ (7.5) \qquad w_{j\mu} := \begin{cases} 0 & \text{if } j \in J \setminus \mathrm{pr}^k(J) \\ w_{\mathrm{in}(\gamma)\mu}(-1)^{|\gamma|} \prod_{p=1}^{|\gamma|} \frac{a_{ip^s(i_p)}}{a_{ip^t(i_p)}} & \text{if there is a } \gamma \in P(j) \end{cases} $$

to obtain $W = [w_{j\mu}] \in Q^{\hat{J} \times m}$. We multiply appropriately columns of $W$ with elements of $R$ to arrive at $V \in R^J$ with coprime columns (this is always possible if $R$ is a gcd domain). The columns of $V$ form a basis of the $R$-module $\mathrm{ker}_R A$.

**Example 7.1.** We look once again into Example 6.3. The pruning depth of the Fischer-Tropsch incidence matrix (6.1) is 0. Reactions that make sense chemically should not need pruning. We have $I_0 = I_1 = \emptyset$ and the quiver $s, t : I_2 \to J$

$$ \mathrm{H}_2\mathrm{O} \xleftarrow{\mathrm{O}} \mathrm{CO} \xrightarrow{\mathrm{C}} \mathrm{C}_n\mathrm{H}_{2n+2}. $$

The quivered linear system is given by $\hat{J} = \{\mathrm{CO}, \mathrm{H}_2\}$ and $\hat{I} = \{\mathrm{H}\}$ with $\hat{A} = [-2n-1, n]$ with quivered solution $[-n, -2n-1]^T$. Using Eqn. (7.5) we reproduce the result of Example 6.3.

**Example 7.2.** To avoid abuse of terrestrial chemistry (it seems important, see [Moo97, MKAM15, Ris09]) I take a look at an incidence matrix proposed by a chemistry student from Betelgeuse:

<table><thead><tr><th></th><th>ACD</th><th>ABDE</th><th>B<sub>2</sub>C<sub>3</sub>DE</th><th>DEFGH</th><th>E<sub>2</sub>H</th><th>E<sub>6</sub>F</th><th>FG<sub>3</sub></th><th>GH</th></tr></thead><tbody><tr><th>A</th><td>1</td><td>1</td><td></td><td></td><td></td><td></td><td></td><td></td></tr><tr><th>B</th><td></td><td>1</td><td>2</td><td></td><td></td><td></td><td></td><td></td></tr><tr><th>C</th><td>1</td><td></td><td>3</td><td></td><td></td><td></td><td></td><td></td></tr><tr><th>D</th><td>1</td><td>1</td><td>1</td><td>1</td><td></td><td></td><td></td><td></td></tr><tr><th>E</th><td></td><td>1</td><td>1</td><td>1</td><td>2</td><td>6</td><td></td><td></td></tr><tr><th>F</th><td></td><td></td><td></td><td>1</td><td></td><td>1</td><td>1</td><td></td></tr><tr><th>G</th><td></td><td></td><td></td><td>1</td><td></td><td></td><td>3</td><td>1</td></tr><tr><th>H</th><td></td><td></td><td></td><td>1</td><td>1</td><td></td><td></td><td>1</td></tr></tbody></table>

Note that $I_0 = I_1 = \emptyset$ and the two paths

<figure><img src="image_3.png" /></figure>

are inconsistent and $v_{ACD} = v_{B_2C_3DE} = v_{ABDE} = 0$. Hence

$$ \mathrm{pr}(I) = \{D, E, F, G\} \text{ and } \mathrm{pr}(J) = \{DEGH, E_2G, E_6F, FG_3, GH\}. $$

While $\mathrm{pr}(I)_0 = \emptyset$ we have $\mathrm{pr}(I)_1 = \{D\}$ and $v_{DEFGH} = 0$. On the other hand, the two paths
<figure><img src="image_4.png" /></figure>

are equivalent. The pruning depth is 2. The diagram represents the quiver $\mathrm{pr}^2(I)_2 \to \mathrm{pr}^2(J)$. We have $\hat{I} = \emptyset$, $\hat{J} = \{E_2H\}$ and the quivered linear system is empty. The solution of the quivered system is without loss of generality $w_{E_2H} = 1$. The coprime integer solution is $v_{E_2H} = -3$, $v_{GH} = 3$, $v_{E_6F} = 1$ and $v_{FG_3} = -1$. To avoid mistakes it is convenient to annotate on paper the $v$'s in the diagram in front of the vertices, the sign has to change along an edge. The balanced reaction is $3E_2H + FG_3 \to E_6F + 3GH$.<fnref n="8" />

Instead of quivering the practitioner of ker pc can do the following: order $I$ such that 1.) for each $p \ge 1$ and $I_p < I_{p+1}$ and 2.) $I_p < I_0$ and order $J$ in some lexicographic order compatible with that of $J$. The order of $I$ ensures that the ker pc start trivially at $I_1$ and the complexity of the condensations grows when moving from $I_p$ to $I_{p+1}$. The order of $J$ ensures that there are no unnecessary row swaps needed. The cost savings of this strategy and quivering seem similar.

## 8. MISCELLANEOUS INVESTIGATIONS

I would like to point out that other fundamental questions in linear algebra can be solved with pc as well. In fact, with pc one can re-purpose the calculations of less specific investigations for the more specific ones. One can, e.g., calculate the rank of a quadratic matrix with pc. If it turns out to be full then one fills in the necessary details to calculate the determinant or the inverse. Or one figures out the rank of an arbitrary matrix with pc and later decides to determine the four subspaces, or to solve a linear system, just by filling out some extra blocks. I call such a type of analysis *modular*: you can always build on what you already have achieved.

**8.1. Inhomogeneous systems.** Here I am indicating how the ker pc can be adapted to solve systems of inhomogeneous linear equations

$$ (8.1) \qquad A\vec{v} = \vec{w} $$

for $A \in R^{n \times m}$ and $\vec{w} \in R^n$. The solution $\vec{v}$ is in $Q^m$ (actually it can be assumed to be in $(S^{-1}R)^m$ for a certain localization $S^{-1}R$ of $R$.) In principle, we only have to add one new row to the pattern of ker pc. One advantage is that we only have to divide at the very last step. The basic idea is that

$$ A\vec{v} = \vec{w} \iff [A|\vec{w}] \begin{bmatrix} \vec{v} \\ -1 \end{bmatrix} = \vec{0}, $$

so that the solution of an inhomogeneous system is equivalent to an appropriately normalized vector in $\ker_Q[A|\vec{w}]$.

In more detail, the pattern at the start of the $(l+1)$st condensation is

$$ \frac{X_l}{(\vec{w}_l)^T} \mid \frac{Y_l}{(\vec{p}_l)^T} \mid \frac{\vec{0}}{\lambda_l}. $$

The initial condition is $X_0 = A^T$, $Y_0 = \mathbf{1}$, $\vec{w}_0 = \vec{w}$, $\vec{p}_0 = \vec{0}$, $\lambda_0 = -1$. The ker pc algorithm of Section 4 is applied to this pattern. The pivots however are only be chosen in the upper left block. The recursion stops at $l = p$ when $X_p$ is zero or empty, which means rank $A = p$. If in this case $(\vec{w}_p)^T$ is nonzero, then the system has no solution. If $(\vec{w}_p)^T$ is zero or empty the general solution $\vec{v} \in Q^m$ of Equation (8.1) is an element of

$$ \lambda_p^{-1} \vec{p}_p + Y_p^T Q^{n-p}. $$

Note that by construction $\lambda_p \ne 0$ and $Y_p^T$ is injective. Similarly, a solution $\vec{v} \in R_{\lambda_p}^m$ of Equation (8.1) is in $\lambda_p^{-1} \vec{p}_p + Y_p^T R_{\lambda_p}^{n-p}$. Here $R_{\lambda_p}$ denotes the localization of $R$ at the multiplicative subset $S = \{\lambda_p^k | k \ge 1\}$.

<footnote n="8">⁸The Betelgeusian civilization has been saved: all the dangerous $FG_3$ could be removed from the atmosphere and replaced by the toxic $E_6F$.</footnote>
During the process there were no row swaps and the symmetry of the matrix could be maintained, which makes the procedure more effective. I suggest to call this procedure *inv pc*.

**8.3. The four subspaces.** Given a matrix $A \in R^{n \times m}$ (or $\in Q^{n \times m}$) we can use a version of ker pc to actually determine bases for $ker_Q(A)$, $im_Q(A)$, $ker_Q(A^T)$ and $im_Q(A^T)$ in just one calculation.

The procedure is similar to that inv pc of the previous subsection and it goes as follows. Write $A = [\vec{v}_1 \cdots \vec{v}_m]$ and $A^T = [\vec{u}_1 \cdots \vec{u}_m]$ and put $p := \text{rank}(A)$. Before step $l+1$ of the condensation the pattern is

$$ \frac{X_l}{W_l^T} \mid \frac{Y_l}{1}. $$

What is in the lower right block is irrelevant for the four subspaces and we leave it blank. The starting values are $X_0 = A^T$, $Y_0 = \mathbf{1}_m$, $W_0^T = \mathbf{1}_n$. The ker pc algorithm of Section 4 is applied to this pattern, while the pivots can only be chosen in the upper left block. (What appears in the lower right block is irrelevant for the recursion; the starting value is by default 0.) If for $l=0$ there is a pivot in the first column add $\vec{u}_1$ to the list of generators of $im_Q(A^T)$ and denote by $\sigma_1$ row permutation used to move the pivot to the first row. Otherwise, there is a maximal zero block of size $m \times k_1$ in the first columns $X_0$. Add the first $k_1$ columns in $W_0^T$ (that is the columns below the zero block) to the list of generators of $ker_Q(A^T)$. There is a pivot in the $(k_1+1)$th column of $X_0$. Apply now **Step 2 (pc)** and **Step 3 (cleaning up)** and put $l=1$. If there is a pivot in the first column add $\vec{u}_2$ to the list of generators of $im_Q(A^T)$ and denote by $\sigma_2$ row permutation used to move the pivot to the first row. Otherwise, there is a maximal zero block of size $(m-1) \times k_2$ in the first columns $X_1$. Add the first $k_2$ columns in $W_1^T$ (that is the columns below the zero block) to the list of generators of $ker_Q(A^T)$. There is a pivot in the $(k_2+1)$th column of $X_1$ in the first row. Apply now **Step 2 (pc)** and **Step 3 (cleaning up)** and put $l=2$. Continue in that manner until $l=p$. Now $X_p$ is zero or empty. Add the columns of $W_p^T$ to the list of generators of $ker_Q(A^T)$. The transposes of the rows of $Y_p$ form a basis for $ker_Q(A)$. Putting $\sigma := \sigma_1 \cdots \sigma_p$, a basis for $im_Q(A)$ is formed by $\vec{v}_{\sigma^{-1}(1)}, \vec{v}_{\sigma^{-1}(2)}, \dots, \vec{v}_{\sigma^{-1}(p)}$. In this way we obtain simultaneously bases for all the four Q-subspaces $ker_Q(A)$, $im_Q(A)$, $ker_Q(A^T)$ and $im_Q(A^T)$.

I suggest for this procedure the name *4pc*. It is actually a pretty quick and clean way to deduce such bases. Let us look at some examples.

**Example 8.3.** Let me put sort of an academic example that hopefully allows to get an impression what is going on with the row permutations:

$$ A = \begin{bmatrix} 0 & 0 & 0 & 1 \\ 0 & 0 & 1 & 1 \\ 0 & 0 & 1 & 1 \\ 1 & 0 & 1 & 1 \end{bmatrix} = [\vec{v}_1 \ \vec{v}_2 \ \vec{v}_3 \ \vec{v}_4], \quad A^T = \begin{bmatrix} 0 & 0 & 0 & 1 \\ 0 & 0 & 0 & 0 \\ 0 & 1 & 1 & 1 \\ 1 & 1 & 1 & 1 \end{bmatrix} = [\vec{u}_1 \ \vec{u}_2 \ \vec{u}_3 \ \vec{u}_4]. $$

In the first step of the 4pc we apply the permutation $\sigma = (1, 3, 2, 4)$ to the rows of $A^T$. Non-mathematicians might prefer a graphical notation for $\sigma$ instead of the cycle notation used above. It makes it easy to keep track of compositions and inversions of permutations.

<figure><img src="image_5.png" /></figure>
Written in cycle form this is $\sigma = (1,3,2,4) = (2,4)(1,4)(1,3)$ and $\sigma^{-1} = (4,2,3,1)$. The 4pc calculation goes as follows

$$
\begin{array}{c}
\left[
  \begin{array}{cccc|c}
    1 &  &  & 1 \\
    & 1 & 1 & 1 \\
    1 & 1 & 1 & 1 \\
    \hline
    1 &  &  & 1
  \end{array}
\right] \\
\begin{array}{c}
\sigma \mapsto \\
\left[
  \begin{array}{cccc|c}
    1 &  & 1 & 1 & 1 \\
    1 & 1 & 1 & 1 & 1 \\
    1 & 1 & 1 & 1 & 1 \\
    \hline
    1 &  &  & 1 & 1 \\
    1 &  &  & 1 & 1 \\
    1 &  &  & 1 & 1
  \end{array}
\right] \\
\left[
  \begin{array}{cccc|c}
    0 &  &  & 0 \\
    0 &  &  & 0 \\
    0 &  &  & 1 \\
    \hline
    1 &  &  & 1
  \end{array}
\right] \to \mathrm{im}_{\mathbb{Q}}(A^T)
\end{array} \\
\begin{array}{c}
\sigma \mapsto \\
\left[
  \begin{array}{cccc|c}
    1 &  & 1 & 1 & 1 \\
    1 & 1 & 1 & 1 & 1 \\
    1 & 1 & 1 & 1 & 1 \\
    \hline
    1 &  &  & 1 & 1 \\
    1 &  &  & 1 & 1 \\
    1 &  &  & 1 & 1
  \end{array}
\right] \\
\left[
  \begin{array}{cccc|c}
    0 &  &  & 0 \\
    0 &  &  & 1 \\
    1 &  &  & 1 \\
    \hline
    1 &  &  & 1
  \end{array}
\right] \to \mathrm{ker}_{\mathbb{Q}}(A^T)
\end{array} \\
\begin{array}{c}
\sigma \mapsto \\
\left[
  \begin{array}{cccc|c}
    0 &  &  & 0 \\
    -1 &  &  & 1 \\
    1 &  &  & 0 \\
    \hline
    1 &  &  & 1
  \end{array}
\right] \to \mathrm{im}_{\mathbb{Q}}(A^T)
\end{array} \\
\begin{array}{c}
\sigma \mapsto \\
\left[
  \begin{array}{cccc|c}
    0 &  &  & 0 \\
    -1 &  &  & 1 \\
    1 &  &  & 0 \\
    \hline
    1 &  &  & 1
  \end{array}
\right] \to \mathrm{im}_{\mathbb{Q}}(A^T)
\end{array}
\end{array}
$$

and hence

$$
\ker_{\mathbb{Q}}(A) = \mathbb{Q} \begin{bmatrix} 0 \\ 1 \\ 0 \\ 0 \end{bmatrix}, \ker_{\mathbb{Q}}(A^T) = \mathbb{Q} \begin{bmatrix} 0 \\ -1 \\ 1 \\ 0 \end{bmatrix}, \quad \mathrm{im}_{\mathbb{Q}}(A) = \mathrm{im}_{\mathbb{Q}}(\begin{bmatrix} \vec{v}_{\sigma^{-1}(1)} & \vec{v}_{\sigma^{-1}(2)} & \vec{v}_{\sigma^{-1}(3)} \end{bmatrix}) = \begin{bmatrix} 1 & 0 & 0 \\ 1 & 1 & 0 \\ 1 & 1 & 0 \\ 1 & 1 & 1 \end{bmatrix} \mathbb{Q}^3, \\
\mathrm{im}_{\mathbb{Q}}(A^T) = \mathrm{im}_{\mathbb{Q}}(\begin{bmatrix} \vec{u}_1 & \vec{u}_2 & \vec{u}_4 \end{bmatrix}) = \begin{bmatrix} 0 & 0 & 1 \\ 0 & 0 & 0 \\ 0 & 1 & 1 \\ 1 & 1 & 1 \end{bmatrix} \mathbb{Q}^3.
$$

**Example 8.4.** Let us consider<fnref n="10" />

$$
A = \begin{bmatrix} 0 & 1 & 1 & 2 \\ 1 & 1 & 2 & 4 \\ 1 & 2 & 4 & 7 \\ 2 & 4 & 7 & 13 \\ 4 & 7 & 13 & 24 \end{bmatrix} = [\vec{v}_1 \ \vec{v}_2 \ \vec{v}_3 \ \vec{v}_4], \quad A^T = \begin{bmatrix} 0 & 1 & 1 & 2 & 4 \\ 1 & 1 & 2 & 4 & 7 \\ 1 & 2 & 4 & 7 & 13 \\ 2 & 4 & 7 & 13 & 24 \end{bmatrix} = [\vec{u}_1 \ \vec{u}_2 \ \vec{u}_3 \ \vec{u}_4 \ \vec{u}_5].
$$

At the first step we apply the transposition $\sigma = (1,2)$ (i.e., a switch of rows):

<footnote n="10">¹⁰The matrix A is a piece of the Hankel matrix of the tribonacci numbers (entry A000073 in [OEI]).</footnote>
<figure><img src="image_6.png" /></figure>

Here I am already trying to write a bit less, when recording the steps of the 4pc:

<figure><img src="image_7.png" /></figure>

We conclude that

$$ \ker_Q(A) = \mathbb{Q} \begin{bmatrix} -1 \\ -1 \\ -1 \\ 1 \end{bmatrix}, \ker_Q(A^T) = \mathbb{Q}^2, $$

$$ \mathrm{im}_{\mathbb{Q}}(A) = \mathrm{im}_{\mathbb{Q}}(\begin{bmatrix} \vec{v}_{\sigma^{-1}(1)} & \vec{v}_{\sigma^{-1}(2)} & \vec{v}_{\sigma^{-1}(3)} \end{bmatrix}) = \begin{bmatrix} 1 & 0 & 1 \\ 1 & 1 & 2 \\ 2 & 1 & 4 \\ 4 & 2 & 7 \\ 7 & 4 & 13 \end{bmatrix} \mathbb{Q}^3, $$

$$ \mathrm{im}_{\mathbb{Q}}(A^T) = \mathrm{im}_{\mathbb{Q}}(\begin{bmatrix} \vec{u}_1 & \vec{u}_2 & \vec{u}_4 \end{bmatrix}) = \begin{bmatrix} 0 & 1 & 1 \\ 1 & 1 & 2 \\ 1 & 2 & 4 \\ 2 & 4 & 7 \end{bmatrix} \mathbb{Q}^3. $$

Combining the ideas of Subsections 8.1, 8.2 and 8.3 one can use pc to calculate the general solution to the problem of finding the one-sided inverse of a full rank matrix. The reader is invited to work out the details.

## 9. PROOF OF KER PC AND 4PC

Let $R$ be a gcd domain (e.g., a PID). The following argument is based on an idea that I learned from stackexchange user Ben Grossman [Gro].

It is enough to prove 4pc as ker pc is part of it. Recall that one says a matrix $A = [a_{ij}]$ is in echelon form if 1) for pivots $a_{ij}$, $a_{i'j'}$ with $i < i'$ it follows that $j < j'$ and 2) there are no zero rows above a pivot. (A pivot is a non-zero entry $a_{ij}$ such that $a_{ij'} = 0$ for all $j' < j$.) I will define an iteration that brings $A^T$ into echelon form. This will be done by the following operations: permuting rows of $A^T$, adding $R$-multiples of rows to other rows and dividing rows by their gcd's.
By the permutation matrix of the permutation $\sigma \in \text{Aut}(\{1, \dots, m\})$ we mean $P_{\sigma} := \sum_{i=1}^{m} E_{\sigma(i)i} \in \mathbb{R}^{m \times m}$.
By a pc-op matrix we mean a matrix of the form

$$M_{i,j}(a,b) = (\mathbf{1} - bE_{ji} + (a-1)E_{jj}) \in \mathbb{R}^{m \times m}$$

for some $a, b \in \mathbb{R}$, $i, j \in \{1, \dots, m\}$. Here $E_{ij}$ is the $m \times m$ matrix whose only nonzero entry is a 1 in line $i$ and column $j$. For $i, j, k$ pairwise distinct and $a, b, a', b' \in \mathbb{R}$ we easily verify $M_{i,j}(a, b)M_{i,k}(a', b') = M_{i,k}(a', b')M_{i,j}(a, b)$.

For $U \in \mathbb{R}^{m \times n}$ and $1 \le j \le n$, $1 \le \lambda < r \le m$ we define for the integer $k \ge m$

$$\Xi_{\lambda,r,j}^{(k)}(U) := M_{\lambda,r}(u_{\lambda j}, u_{rj}) = \mathbf{1} - u_{rj}E_{r\lambda} + (u_{\lambda j} - 1)E_{rr} \in \mathbb{R}^{k \times k}.$$

Let $B \in \mathbb{R}^{n \times m}$, choose a permutation $\sigma \in \text{Aut}(\{1, \dots, m\})$ and put

$$L_{\lambda, B^T, \sigma} := \prod_{r=\lambda+1}^{m} \Xi_{\lambda, r, j}^{(m)} (P_{\sigma} B^T) \in \mathbb{R}^{m \times m},$$

$$R_{\lambda, B^T, \sigma} := \prod_{r=\lambda+1}^{m} \Xi_{\lambda, r, j}^{(n)} (B^T P_{\sigma}) \in \mathbb{R}^{n \times n}.$$

For the main argument we will use the following.

**Lemma 9.1.** With the notations above we have

$$(1) \quad L_{\lambda, B^T, \sigma} P_{\sigma} B^T - P_{\sigma} B^T (R_{\lambda, B, \sigma^{-1}})^T = \begin{bmatrix} * & * \\ * & \mathbf{0} \end{bmatrix}, \quad \text{where } \mathbf{0} \in \mathbb{R}^{(m-\lambda) \times (n-\lambda)},$$

$$(2) \quad (L_{\lambda, B^T, \sigma} P_{\sigma} B^T)_{rj} = 0 \text{ for } r > \lambda \text{ and}$$

$$(3) \quad (P_{\sigma} B^T R_{\lambda, B^T, \sigma})_{jr} = 0 \text{ for } r > \lambda.$$

*Proof.* We evaluate

$$L_{\lambda, B^T, \sigma} = \prod_{r=\lambda+1}^{m} \Xi_{\lambda, r, j}^{(m)} (P_{\sigma} B^T) = \prod_{r=\lambda+1}^{m} \Xi_{\lambda, r, j}^{(m)} \left( \sum_{p,q} b_{\sigma(q)p} E_{pq} \right) \\ = \prod_{r=\lambda+1}^{m} M_{\lambda, r}(u_{\lambda j}, u_{rj}) = \prod_{r=\lambda+1}^{m} (\mathbf{1} + (b_{\sigma(j)\lambda} - 1)E_{rr} - b_{\sigma(j)r}E_{r\lambda})$$

Here $U = P_{\sigma} B^T$ and $u_{\lambda j} = b_{\sigma(j)\lambda}$, $u_{rj} = b_{\sigma(j)r}$. Expanding the product we obtain

$$L_{\lambda, B^T, \sigma} = \mathbf{1} + \sum_{r=\lambda+1}^{m} ((b_{\sigma(j)\lambda} - 1)E_{rr} - b_{\sigma(j)r}E_{r\lambda}) = \sum_{r=1}^{\lambda} E_{rr} + \sum_{r=\lambda+1}^{m} (b_{\sigma(j)\lambda}E_{rr} - b_{\sigma(j)r}E_{r\lambda}).$$

Hence

$$(9.1) \quad L_{\lambda, B^T, \sigma} P_{\sigma} B^T = \sum_{p,q} \left( \sum_{r=1}^{\lambda} E_{rr} + \sum_{r=\lambda+1}^{m} (b_{\sigma(j)\lambda} E_{rr} - b_{\sigma(j)r} E_{r\lambda}) \right) b_{\sigma(q)p} E_{pq} \\ = \sum_{q} \sum_{r=1}^{\lambda} b_{\sigma(q)r} E_{rq} + \sum_{q} \sum_{r=\lambda+1}^{m} (b_{\sigma(j)\lambda} b_{\sigma(q)r} - b_{\sigma(j)r} b_{\sigma(q)\lambda}) E_{rq}.$$

If $q = j$ the second sum is zero proving 2). Next, we look at

$$R_{\lambda, B^T, \sigma} := \prod_{r=\lambda+1}^{m} \Xi_{\lambda, r, j}^{(n)} (B^T P_{\sigma}) = \prod_{r=\lambda+1}^{m} \Xi_{\lambda, r, j}^{(n)} \left( \sum_{p,q} b_{\sigma(q)p} E_{pq} \right) \\ = \prod_{r=\lambda+1}^{m} M_{\lambda, r}(u_{\lambda j}, u_{rj}) = \prod_{r=\lambda+1}^{m} (\mathbf{1} + (b_{\sigma(j)\lambda} - 1)E_{rr} - b_{\sigma(j)r}E_{r\lambda}).$$

Here $U = B^T P_{\sigma}$ and $u_{\lambda j} = b_{\sigma(j)\lambda}$, $u_{rj} = b_{r\sigma^{-1}(j)}$. Expanding the product we obtain

$$R_{\lambda, B^T, \sigma} = \mathbf{1} + \sum_{r=\lambda+1}^{m} ((b_{\sigma(j)\lambda} - 1)E_{rr} - b_{r\sigma^{-1}(j)}E_{r\lambda}) = \sum_{r=1}^{\lambda} E_{rr} + \sum_{r=\lambda+1}^{m} (b_{\sigma(j)\lambda}E_{rr} - b_{\sigma(j)r}E_{r\lambda}).$$
and hence obtain

$$ (9.2) \qquad \begin{aligned} P_{\sigma} B^T (R_{\lambda, B^T, \sigma})^T &= \sum_{p,q} b_{\sigma(q)p} E_{pq} \left( \sum_{r=1}^{\lambda} E_{rr} + \sum_{r=\lambda+1}^{m} (b_{\sigma(j)\lambda} E_{rr} - b_{\sigma(j)r} E_{\lambda r}) \right) \\ &= \sum_{p} b_{\sigma(r)p} E_{pr} + \sum_{p} \sum_{r=\lambda+1}^{m} (b_{\sigma(r)p} b_{\sigma(j)\lambda} - b_{\sigma(\lambda)p} b_{\sigma(j)r}) E_{pr}. \end{aligned} $$

If $p = j$ the second sum is zero proving 3). For $p, q > \lambda$ in both formulas, Eqns. (9.1) and (9.2), the second second terms coincide, which proves 1). □

Note that the $2 \times 2$-minors in Eqns. (9.1) and (9.2) are exactly those of pivotal condensation (see Step 1 and 2). So the scheme of pc can be realized by left multiplication with a product of pc-op matrices, or, alternatively, by right multiplication with a product of transposes of pc-op matrices. The catch is that the other terms in Eqns. (9.1) and (9.2) do not coincide. This does not cause any serious problems, however.

Now for $0 \le l \le p-1$ we set $L_{l+1} := L_{l+1, A_l^T, \sigma_{l+1}}$, $N_{l+1} := R_{l+1, A_l^T, \sigma_{l+1}}$ where $A_l^T$ is defined by the recursion $A_{l+1}^T = L_l P_{\sigma_l} A_l^T$, $A_0^T = A^T$. The permutations $\sigma_{l+1}$ are the ones from Subsection 8.3. We write $M_{l+1} = L_{l+1} P_{\sigma_{l+1}}$. By Lemma 9.1 2) the matrix

$$ E = M_p \cdots M_2 M_1 A^T \in R^{m \times n} $$

is in echelon form and of rank $p$. So the last $n-p$ rows of $E$ are zero. The lower $(n-p) \times (n-p - \sum_{i=1}^p k_i)$-block of $E$ is just $X_p$, which is zero or empty. The lower $(m-p) \times n$-block of $M_p \cdots M_2 M_1$ multiplied from the left by an appropriate diagonal matrix equals $Y_p$. But since the last $m-p$ columns of $E^T$ are zero we must have $AY_p^T = 0$. Note that $Y_p^T$ is by construction of full rank, so that $ker_Q(A) = Y_p^T Q^{m-p}$.

On the other hand, the matrix

$$ F = [f_{ij}] := P_{\sigma_p} \cdots P_{\sigma_2} P_{\sigma_1} A^T N_1 N_2 \cdots N_p \in R^{m \times n} $$

has by Lemma 9.1 3) the property (P): for transposed pivots $f_{ij}$, $f_{kl}$ with $l > j$ it follows $k > i$. Here we say that an entry in $F$ is a transposed pivot if it is a pivot in $F^T$. Let $I_F = \{i_1, \dots, i_{n-p}\}$ be the set of all column indices of zero columns of $F$. We deduce from property (P) that $ker(F) = R\{\vec{e}_i | i \in I_F\}$. By Lemma 9.1 1) and 3) the set $I_F$ consists of exactly those column indices where one encounters zero columns in some $X_l$. Here it is understood that the columns of $X_l$ are counted starting with column index $\sum_{i=1}^{l-1} (k_i + 1)$. With $N = [\vec{v}_1 \dots \vec{v}_n] := N_1 N_2 \cdots N_p$ put $N' := [\vec{v}_{i_1} \dots \vec{v}_{i_{n-p}}]$. It follows that $A^T N' = 0$, where $N'$ is of full rank by construction. Hence $ker_Q(A^T) = N'Q^{n-p}$.

The first $p$ columns of $E^T = A M_1^T \cdots M_p^T$ are linearly independent. Multiplying from the right by transposed pc-op or diagonal matrices does not mess up linear independence of columns. Now from

$$ A = E^T (M_p^T)^{-1} \cdots (M_1^T)^{-1} $$

we see that the column vectors of $A$ with column indices $\sigma^{-1}(1), \dots, \sigma^{-1}(p)$ form a maximally independent set of vectors in $\mathrm{im}_Q(A)$. Here it is important to keep in mind that a permutation $\sigma$ acts on columns by right multiplying with $P_{\sigma^{-1}1}$. The statement about the image of $A^T$ is, in principle, the usual one from Gaussian elimination. In fact, with $A^T = [\vec{u}_1 \dots \vec{u}_n]$ we have $\mathrm{im}_Q(A^T) = R\{\vec{u}_i | \text{ there is a transposed pivot in the } i\text{th column of } F\}$. This finishes the proof of 4pc.

In general, it is not guaranteed that $ker_R(A^T) = Y_p^T R^n$. A necessary and sufficient criterion for the case when $R$ is a PID is provided by Theorem 10.9. If $R$ is a domain but not a gcd domain it is not guaranteed that the normalizations (i.e., the cleaning up of Step 3) of $ker$ pc are theoretically feasible, even though practically they might be. If in this situation $Y_p^T$ contains a $k \times k$-submatrix $DP$, where $D$ is a diagonal matrix of units and $P$ a permutation matrix then $ker_R(A^T) = Y_p^T R^n$ follows. The reason are explained in the next section.

## 10. THE KERNEL AS A FREE $R$-MODULE

For those who would like understand how to calculate a basis for the kernel $ker_Z(A)$ of and integer matrix $A \in \mathbb{Z}^{n \times m}$ I collect here some material from commutative algebra. I give a practical account of Smith normal form and discuss the saturation of a submodule of a finitely generated free module. Without some experience with Smith normal form the reader will not be able to appreciate the main point: Subsection 10.3, where
I indicate how to turn a Q-basis obtained by ker pc into a Z-basis. A clear and accessible reference for the background material is [Jac85, Section 3.7], [Jac89, Sections 7.2,7.3]. Another good reference on Smith normal forms is [New72].

**10.1. A practical guide to Smith normal form.** In Gaussian elimination over a field one uses the following operations

1.) switching rows,

2.) adding a multiple of a row to another row,

3.) and multiplying a row with a nonzero element.

Working over a PID $R$, it is a natural idea to use operation 3.) only for units. Otherwise we are muddying the waters<fnref n="11" />. The idea of Smith normal form is to compensate for the lost freedom by admitting also column versions of the operations 1.) and 2.). The operations 1.) and 2.) for rows can be understood by left multiplications with matrices in $GL_n(R)$ while the operations 1.) and 2.) for columns correspond to right multiplications with matrices in $GL_m(R)$. Obviously, the matrices used for the operations of type 1.) are transpositions $1 - E_{ii} - E_{jj} + E_{ij} + E_{ji}$ while the matrices used for the operations of type 2.) are of the form $1 + rE_{ij}$ with $r \in R$, the latter actually being unimodular.

**Theorem 10.1** (H.J.S. Smith (1861) [Smi94]). Let $n, m \ge 1$ be integers, $R$ be PID and $A \in R^{n \times m}$ a matrix of rank $p$. Then there exist matrices $U \in GL_n(R)$ and $V \in GL_m(R)$ and a diagonal matrix $D = \sum_{i=1}^p d_i E_{ii} \in R^{n \times m}$ such that

a.) $d_i \neq 0$ for $i = 1, \dots, p$,

b.) $d_i | d_{i+1}$ for $i = 1, \dots, p-1$,

c.) $D = UAV$.

The $d_1, \dots, d_p$ are unique up to multiplication by a unit and are called invariant factors.

For simplicity we will only discuss in detail the case when $R$ be a euclidean domain with euclidean function $\delta : R\setminus\{0\} \to \{1, 2, 3, \dots\}$. We need also the following obvious consequence of the Smith normal form.

**Corollary 10.2.** With the notations of Theorem 10.1 we have the following.

(1) A basis for the free $R$-module $\text{im}_R(A)$ is provided by the first $p$ columns of $U^{-1}D$.

(2) If $m > p$, the last $m-p$ columns of $V$ form a basis for the free $R$-module $\text{ker}_R(A)$.

(3) The last $n-p$ rows of $U$ form a $(n-p) \times n$ matrix $B$ of rank $n-p$ such that $BA = 0$.

For a detailed proof of Theorem 10.1 we refer to Jacobson [Jac85, Section 3.7]. For the reader to understand the basic ideas I work out some examples for the case $R = \mathbb{Z}$. Here the euclidean function is $\delta = |\cdot|$. I use again the notational trick of adding block matrices to keep track of the row and column operations. Other people have found this before me (see, e.g., the method of Gauß-Jordan for matrix inversion), but I do not know a reference where it is used for Smith normal form calculations. The trick has also the advantage that one can do check sums for rows and columns to spot errors.

**Example 10.3.** We would like to find a Smith normal form of the matrix

$$A = \begin{bmatrix} 3 & 4 & 5 & 6 \\ 7 & 8 & 9 & 10 \\ 11 & 12 & 13 & 14 \\ 15 & 16 & 17 & 18 \end{bmatrix}.$$

The first step consists of switching a non-zero entry with minimal $\delta$ to position (1, 1) using operations of type 1.). Here it is not necessary since that entry 3 is already there. Next, we want to annihilate everything left of 3 using column operations of type 2.). As $3 \nmid 4, 5$ we cannot do it in one step but use division by 3

<footnote n="11">¹¹This is of course true for ker pc.</footnote>
with remainder:

<figure><img src="image_8.png" /></figure>

The unit matrix in the bottom keeps track of the column operations and, as a matter of principle, we only work in the first block. Now we switch the first two columns so that an entry with minimal $\delta$ appears at position (1, 1). This enable us to finish the clean up of the first row using column operations of type 2.).

<figure><img src="image_9.png" /></figure>

Next we annihilate everything below the 1 at position (1, 1) using row operations of type 2.). The unit matrix on the right keeps track of the row operations.

<figure><img src="image_10.png" /></figure>

Since $1 | 4$ we can proceed with the block
$$
\begin{array}{ccccc}
4 & 0 & -4 & & \\
8 & 0 & -8 & (\text{1st row and column become spectators}). & \\
12 & 0 & -12 & & \\
\end{array}
$$

4 with minimal $\delta$ appear already at the right spot and we can clean what is on the right of it using a column operation of type 2.).

<figure><img src="image_11.png" /></figure>
I am following here essentially [Ash07, Subsection 4.5]. We switch the first and the third row to move a 2 to position (1, 1):

<figure><img src="image_12.png" /></figure>

Next we clean up everything to the right of that 2, and afterwards, what is below of it:

<figure><img src="image_13.png" /></figure>

Since $2 | 4$ we are fine and can move to the block $\begin{array}{c@{\quad}c} 4 & 0 \\ 0 & 10 \end{array}$, starting by cleaning up what is to the right of the 4:

<figure><img src="image_14.png" /></figure>

Now the problem is that $4 \nmid 10$ (this is the sublety). We add the third row to the second and do column operation of type 2.) dividing by 4 with remainder:

<figure><img src="image_15.png" /></figure>

Next we have to move the 2 to position (2, 2) and clean up what is on the right of it:

<figure><img src="image_16.png" /></figure>
**Theorem 10.6.** With the assumption of Theorem 10.1 let $\Delta_i$ be the gcd of the $i \times i$-minors of $A$. Then $d_1 = \Delta_1$ and $d_i = \frac{\Delta_i}{\Delta_{i-1}}$ for $i = 2, \dots, p$.

**10.2. Saturated submodules of $R^n$.** I am pretty sure the following material is written up somewhere. I just could not allocate it.

Let $S$ be a multiplicative subset of the domain $R$ and $S^{-1}R$ the localization of $R$ at $S$. An $R$-module $W$ of $R^n$ is called $S$-saturated if for any $r \in S$ we have that $rw \in W$ implies $w \in W$. Let us denote by $f$ the canonical $R$-linear map $R^n \to (S^{-1}R)^n$. It is easy to verify that an $R$-module $W$ is $S$-saturated if and only if there is an $S^{-1}R$-submodule $M$ of $(S^{-1}R)^n$ such that $W = f^{-1}(M)$. More generally, the $S$-saturation of $W$ is $f^{-1}(f(W))$.

We now specialize to the case $S = R\setminus\{0\}$. Then $S^{-1}R = Q$ is the field of fractions of $R$. An $(R\setminus\{0\})$-saturated submodule $W$ of $R^n$ we simply call saturated. The key observation is that, since $R$ is a domain, the kernel of a matrix $A \in R^{n \times m}$ is saturated. This is because $A(r\vec{w}) = rA(\vec{w}) = 0$ implies $A(\vec{w}) = 0$ for any $r \in R\setminus\{0\}$. Since an inclusion of two equal dimensional $Q$-vector spaces is an equality we have the following.

**Lemma 10.7.** Let $R$ be a domain and $W$ an $R$-submodule of the saturated $R$-submodule $M \subseteq R^n$ such that $\text{rank} W = \text{rank} M$. Then $W = M$ if and only $W$ is saturated.

**Proposition 10.8.** Let $R$ be a domain and $W$ be a $R$-submodule of $R^n$ and $l \ge 0$ and integer. If $W$ saturated then $\wedge_R^l W$ saturated for each $l \ge 0$.

*Proof.* We have $W = f^{-1}(V)$ for a $Q$-vector subspace $V \subseteq Q^n$, i.e., for every $\vec{w} \in W$ exists $v \in V: f(\vec{w}) = \vec{v}$. Every $\omega \in \wedge^l W$ can be written as $\omega = \vec{w}_1 \wedge \dots \wedge \vec{w}_l$ for some $\vec{w}_1, \dots, \vec{w}_l \in W$ and there exist $\vec{v}_1, \dots, \vec{v}_l \in V$ such that $f(\vec{w}_i) = \vec{v}_i$. Hence $(\wedge_R^l f)(\vec{w}_1 \wedge \dots \wedge \vec{w}_l) = \vec{v}_1 \wedge \dots \wedge \vec{v}_l$ and $\wedge_R^l W = (\wedge_R^l f)^{-1}(\wedge_Q^l V)$. $\square$

If $R$ is a PID we can be more specific.

**Theorem 10.9.** Let $R$ be a PID and $X \in R^{n \times p}$ a matrix of rank $p$ with invariant factors $d_1, \dots, d_p$. Then the following statements are equivalent:

i.) $\text{im}_R(X)$ is saturated.

ii.) $d_i \in R^\times$ for all $i = 1, \dots, p$.

iii.) $\Delta_p \in R^\times$.

*Proof.* The implications i.) $\iff$ ii.) $\implies$ iii.) are obvious from Theorems 10.1 and 10.6. To establish the implication iii.) $\implies$ ii.) we show that if $d_p \notin R^\times$ it follows that $\Delta_p \notin R^\times$. If $A \in R^{r \times s}$ is a matrix we denote by $A_{i_1 \dots i_p}^{j_1 \dots j_p}$ the $p \times p$-minor corresponding the row indices $1 \le i_1 < \dots < i_p \le r$ and column indices $1 \le j_1 < \dots < j_p \le s$. From the Smith normal form $X = U^{-1}DV^{-1}$ and the Cauchy-Binet formula we deduce

$$X_{i_1 \dots i_p}^{1 \dots p} = \sum_{1 \le i_1 < \dots < i_p \le n} (U^{-1})_{i_1 \dots i_p}^{1 \dots p} D_{1 \dots p}^{1 \dots p} (V^{-1})_{1 \dots p}^{1 \dots p} = \left(\prod_{i=1}^{p} d_i\right) \sum_{1 \le i_1 < \dots < i_p \le n} (U^{-1})_{i_1 \dots i_p}^{1 \dots p} (V^{-1})_{1 \dots p}^{1 \dots p}.$$

But this means that $d_p | X_{i_1 \dots i_p}^{1 \dots p}$ and the gcd $\Delta_p$ of the $X_{i_1 \dots i_p}^{1 \dots p}$ cannot be a unit. $\square$

**Remark 10.10.** Let $R$ be a gcd domain (e.g., a PID). Let $V$ be a one-dimensional $Q$-vector space. Clearing the denominators we write $V = Q\vec{w}$ for some $\vec{w} \in R^n$. Put $W := R\vec{w}$ and note that $V = Q \otimes_R W$. Now $W$ is saturated if and only if the gcd $\Delta_1$ of the entries of $\vec{w}$ is a unit. Such vectors are also called primitive in the literature [BAY21]. The generators that are calculated by $\text{ker pc}$ are by construction primitive. That is, if $\text{dim}(\text{ker}_Q(A)) = 1$ then $\text{ker pc}$ provides a generator for $\text{ker}_R(A)$.

**10.3. Turning a $Q$-basis into an $R$-basis.** In this section we assume $R$ is a PID. Once we have already calculated $Q$-basis $\vec{v}_1, \dots, \vec{v}_p \in R^n$ for the $Q$-vector space $\text{ker}_Q(A)$ for $A \in R^{n \times m}$ using $\text{ker pc}$ it seems desirable to turn it into an $R$-basis. In particular, if $n, m$ are large but $\text{rank}(\text{ker}(A))$ is small it is not very attractive to go through the whole Smith normal form calculation for the matrix $A$.

Here one can proceed as follows. Take $Y_p^T = [\vec{v}_1 \ \dots \ \vec{v}_p] \in R^{n \times p}$ and calculate its Smith normal form $D' = U'Y_p^T V'$. Form the matrix $B$ by taking the last $n-p$ rows of $U'$ (see Corollary 10.2(3)) and calculate the Smith normal form $D'' = U''B V''$ of $B$. The last $p$ columns of $V''$ form a basis for the $R$-module $\text{ker}_R(A)$. More generally, this recipe can be used to calculate the $R$-saturation of a submodule of $R^n$.
**Example 10.11.** Let us take a final look at Example 6.2. The kernel of

$$A = \begin{bmatrix} 1 & 0 & 1 & 2 \\ 0 & 2 & 1 & 3 \end{bmatrix}$$

was shown by ker pc to be $\begin{bmatrix} -2 & -4 \\ -1 & -3 \\ 2 & 0 \\ 0 & 2 \end{bmatrix} \mathbb{Q}^2$ as a $\mathbb{Q}$-vector space. We now use a Smith normal form to calculate $B$:

<figure><img src="image_17.png" /></figure>

so that $B = \begin{bmatrix} 3 & -4 & 1 & 0 \\ -1 & 2 & 0 & 1 \end{bmatrix}$. Then we calculate a Smith normal form of $B$

<figure><img src="image_18.png" /></figure>

and deduce once again $\ker_Z(\begin{bmatrix} 1 & 0 & 1 & 2 \\ 0 & 2 & 1 & 3 \end{bmatrix}) = \mathbb{Z} \oplus \mathbb{Z} \begin{bmatrix} 1 \\ 0 \\ -3 \\ 1 \end{bmatrix}$.

In the example it looks like that the last method<fnref n="12" /> to calculate the basis for $\ker_Z(A)$ is more labour intensive than calculating a Smith normal form of $A$. But let us say, to be concrete, $A$ is a $100 \times 100$ matrix and $\ker_Z(A)$ is of rank 3. Then ker pc seems faster than Smith normal form as it does not need division with remainder (it works with bigger numbers though). The saturation procedure explained here takes only about 6% of the total effort and if, by accident, $\Delta_3$ of $X$ is in $R^\times$ it is not needed at all.

I do not know a method to obtain a basis for the non-saturated $\mathrm{im}_R(A)$ by post-processing the result of 4pc.

<footnote n="12">¹²How should one call it? I am voting for *smitherate*.</footnote>
REFERENCES

[Abe14] Francine F. Abeles, Chiò's and Dodgson's determinantal identities, Linear Algebra Appl. 454 (2014), 130–137.

[Ash07] Robert B. Ash, Basic abstract algebra. For graduate students and advanced undergraduates, Mineola, NY: Dover Publications, 2007 (English).

[BAY21] Faten Ben Amor and Ihsen Yengui, Saturation of finitely-generated submodules of free modules over Prüfer domains, Armen. J. Math. 13 (2021), 21, Id/No 1.

[Chi53] F. Chiò, Mémoire sur les fonctions connues sous le nom de résultantes ou de déterminants, Turin, 1853.

[Eve12] H. Eves, Elementary matrix theory, Dover Books on Mathematics, Dover Publications, 2012.

[Fie11] Miroslav Fiedler, Matrices and graphs in geometry, Encyclopedia of Mathematics and its Applications, vol. 139, Cambridge University Press, Cambridge, 2011.

[Fin09] J.K. Fink, Physical chemistry in depth, Springer Berlin Heidelberg, 2009.

[GP08] Gert-Martin Greuel and Gerhard Pfister, A singular introduction to commutative algebra, extended ed., Springer, Berlin, 2008.

[Gro] https://math stackexchange com questions 4239718 kernel-of-a-matrix-using-the-rref-of-its-transpose.

[Jac85] Nathan Jacobson, Basic algebra I. 2nd ed, New York: W. H. Freeman and Company. XVIII, 499 p., 1985.

[Jac89] —, Basic algebra II., 2nd ed. ed., New York, NY: W. H. Freeman and Company, 1989 (English).

[MKAM15] Anna Maria Michalowska-Kaczmarczyk, Agustin G. Asuero, and Tadeusz Michalowski, "Why Not Stoichiometry" versus "Stoichiometry—Why Not?" part i: General context, Critical Reviews in Analytical Chemistry 45 (2015), no. 2, 166–188, PMID: 25558777.

[Moo97] John W. Moore, Balancing the forest and the trees, J. Chem. Educ. 74 (1997), no. 11, 1253.

[New72] Morris Newman, Integral matrices, Pure and Applied Mathematics, 45. New York-London: Academic Press. XVII, 224 p. § 14.00 (1972)., 1972.

[OEI] OEIS Foundation Inc., The On-Line Encyclopedia of Integer Sequences, Published electronically at http://oeis.org.

[PV06] Dávid Papp and Béla Vizvári, Effective solution of linear Diophantine equation systems with an application in chemistry, J. Math. Chem. 39 (2006), no. 1, 15–31 (English).

[Ric68] Jeremias Benjamin Richter, Anfangsgründe der Stöchiometrie: oder, Messkunst chemischer Elemente, no. Bd. 1-2, G. Olms, 1968, Nachdruck der Ausgabe Breslau 1792-94.

[Ris09] Ice B. Risteski, A new singular matrix method for balancing chemical equations and their stability, Journal of the Chinese Chemical Society 56 (2009), no. 1, 65–79.

[Smi94] Henry John Stephen Smith, The collected mathematical papers. Edited by J. W. L. Glaisher. With a mathematical introduction by the editor, biographical sketches and a portrait. In two volumes., Oxford: at the Clarendon Press. Volume I: xcv, 603 p. (with portrait); Volume II: vii, 719 p. (1894)., 1894.

[Sto95] Roland Stout, Redox challenges: Good times for puzzle fanatics, Journal of Chemical Education 72 (1995), no. 12, 1125.

[Yen15] Ihsen Yengui, Constructive commutative algebra. Projective modules over polynomial rings and dynamical Gröbner bases, Lect. Notes Math., vol. 2138, Cham: Springer, 2015 (English).

Email address: herbighc@gmail.com

DEPARTAMENTO DE MATEMÁTICA APLICADA, UNIVERSIDADE FEDERAL DO RIO DE JANEIRO, AV. ATHOS DA SILVEIRA RAMOS 149, CENTRO DE TECNOLOGIA - BLOCO C, CEP: 21941-909 - RIO DE JANEIRO, BRAZIL