${\mathcal{V}}_{m}^{({\mathfrak{g}})}\sim(\theta_{m}^{({\mathfrak{g}})})^{2}\,V_{m}^{({\mathfrak{g}})}$,
where the difference operator $V_{m}^{({\mathfrak{g}})}$ is independent of
$t$. Similarly
${\mathcal{D}}_{m}^{({\mathfrak{g}})}\sim(\theta_{m}^{({\mathfrak{g}})})^{2}\,D_{m}^{({\mathfrak{g}})}$
where $D_{m}^{({\mathfrak{g}})}$ independent of $t$, for
$m=1,2,...,N_{\mathfrak{g}}$. By inspection, we find the same leading behavior
for $N_{\mathfrak{g}}<m\leq N$, using (A.35). This leads to the definitions:

(A.38)
$D_{m}^{({\mathfrak{g}})}(x):=\lim_{t\to\infty}(\theta_{m}^{({\mathfrak{g}})})^{-2}\,{\mathcal{D}}_{m}^{({\mathfrak{g}})}(x),\
\
V_{m}^{({\mathfrak{g}})}:=\lim_{t\to\infty}(\theta_{m}^{({\mathfrak{g}})})^{-2}\,{\mathcal{V}}_{m}^{({\mathfrak{g}})}$

with $\theta_{m}^{({\mathfrak{g}})}$ as in (A.34-A.35).

######  Remark A.17.

The limit $t\to\infty$ of ${\mathcal{D}}^{({\mathfrak{g}})}_{m}$ is the same
as that of $\tilde{\mathcal{D}}^{({\mathfrak{g}})}_{m}$ of Sect. A.3.3: Using
${\mathcal{D}}^{({\mathfrak{g}})}_{m}\simeq(\theta_{m}^{({\mathfrak{g}})})^{2}\,{D}^{({\mathfrak{g}})}_{m}$
at leading order in $t$, the statement follows immediately from (A.36) by
noting that for $m\leq N_{\mathfrak{g}}$:

$\displaystyle\lim_{t\to\infty}\frac{\theta_{m-1}^{({\mathfrak{g}})}}{\theta_{m}^{({\mathfrak{g}})}}$
$\displaystyle=$ $\displaystyle 0\quad(m\geq
1;{\mathfrak{g}}=C_{N}^{(1)},D_{N+1}^{(2)}),$
$\displaystyle\lim_{t\to\infty}\frac{\theta_{m-2}^{({\mathfrak{g}})}}{\theta_{m}^{({\mathfrak{g}})}}$
$\displaystyle=$ $\displaystyle 0\quad(m\geq
2;{\mathfrak{g}}=B_{N}^{(1)},A_{2N-1}^{(2)},A_{2N}^{(2)}).$

For any ${\mathfrak{g}}$, define
$\Pi_{\lambda}^{({\mathfrak{g}})}=\lim_{t\to\infty}t^{-\rho^{*}\cdot\lambda}P_{\lambda}^{({\mathfrak{g}})}$.
We refer to these as $q$-Whittaker functions, although they coincide with the
usual definition only for the case of ${\mathfrak{g}}=D_{N}^{(1)}$ and the
twisted algebras in Table 1. The spectral theorem for
$D_{m}^{({\mathfrak{g}})}$ is as follows.

######  Theorem A.18.

The functions $\Pi^{({\mathfrak{g}})}_{\lambda}$ are common eigenfunctions of
$D_{m}^{({\mathfrak{g}})}$, with

$D_{m}^{({\mathfrak{g}})}(x)\,\Pi^{({\mathfrak{g}})}_{\lambda}(x)=\Lambda^{\omega_{m}^{*}}\,\Pi^{({\mathfrak{g}})}_{\lambda}(x).$

###### Proof.

The eigenvalues are the limits

$\delta_{\lambda;m}^{({\mathfrak{g}})}=\lim_{t\to\infty}(\theta_{m}^{({\mathfrak{g}})})^{-1}{\hat{e}}_{m}^{(R^{*})}(s),\qquad
s_{i}=q^{\lambda_{i}}t^{N+\xi_{\mathfrak{g}}-i}.$

The limit extracts the leading monomial of the $s$ variables in the symmetric
function, which is the term involving $s_{1},s_{2},...,s_{m}$ with positive
powers. For $m=1,2,...,N_{\mathfrak{g}}$ this gives, using (A.30) (A.34) and
(A.12):
$\delta_{\lambda;m}^{({\mathfrak{g}})}=\lim_{t\to\infty}t^{-m(N+\xi_{g}-\frac{m+1}{2})}{\hat{e}}_{m}(s)=q^{\lambda_{1}+\lambda_{2}+\cdots+\lambda_{m}}$.
By inspection, using (A.31), (A.35) and (A.14),

$\displaystyle
D_{N}^{(1)}:\quad\delta_{\lambda;N}^{(D_{N}^{(1)})}=q^{\frac{1}{2}(\lambda_{1}+\lambda_{2}+\cdots+\lambda_{N-1}+\lambda_{N})},\
\delta_{\lambda;N-1}^{(D_{N}^{(1)})}=q^{\frac{1}{2}(\lambda_{1}+\lambda_{2}+\cdots+\lambda_{N-1}-\lambda_{N})},$
$\displaystyle
C_{N}^{(1)}:\quad\delta_{\lambda;N}^{(C_{N}^{(1)})}=q^{\frac{1}{2}(\lambda_{1}+\lambda_{2}+\cdots+\lambda_{N-1}+\lambda_{N})},$
$\displaystyle
D_{N+1}^{(2)}:\delta_{\lambda;N}^{(D_{N+1}^{(2)})}=q^{\frac{1}{2}(\lambda_{1}+\lambda_{2}+\cdots+\lambda_{N-1}+\lambda_{N})}.$

The result can be uniformly written as
$\Lambda^{\omega_{m}^{*}}=q^{\omega_{m}^{*}\cdot\lambda}$, $\omega_{m}^{*}$
the fundamental weights of $R^{*}$. ∎

######  Remark A.19.

Using the convention $V_{0}^{({\mathfrak{g}})}=1$, Def. A.8 implies

(A.39) $D_{m}^{({\mathfrak{g}})}=\sum_{j=0}^{m}V_{j}^{({\mathfrak{g}})},\qquad
m\leq N_{\mathfrak{g}}.$

##  Appendix B Pieri rules

In this section, we list some of the Pieri operators implementing the first
Pieri rule for Koornwinder and ${\mathfrak{g}}$-Macdonald polynomials. Their
$q$-Whittaker limit gives $q$-difference Toda Hamiltonians for certain root
systems.

###  B.1. First Pieri rule for Koornwinder polynomials

The first Pieri operator for generic $(a,b,c,d)$ is the $q$-difference
operator ${\mathcal{H}}_{1}^{(a,b,c,d)}(s)$ of Theorem (3.3).

######  Theorem B.1.

The first Pieri operator for Koornwinder polynomials is

$\displaystyle{\mathcal{H}}_{1}^{(a,b,c,d)}(s;q,t)$ $\displaystyle=$
$\displaystyle
a^{-1}t^{1-N}\varphi^{(a^{*},b^{*},c^{*},d^{*})}(s)+\sum_{i=1}^{N}\left[\left(\prod_{j=1}^{i-1}\frac{s_{i}-t^{-1}s_{j}}{s_{i}-s_{j}}\frac{qs_{i}-ts_{j}}{qs_{i}-s_{j}}\right)\,T_{i}\right.$
$\displaystyle+\frac{\prod_{u\in\\{a^{*},b^{*},c^{*},d^{*}\\}}(1-u^{-1}s_{i})(q-us_{i})}{(1-s_{i}^{2})(q-s_{i}^{2})^{2}(q^{2}-s_{i}^{2})}$
$\displaystyle\times\left.\left(\prod_{j=i+1}^{N}\frac{t^{-1}s_{i}-s_{j}}{s_{i}-s_{j}}\frac{ts_{i}-qs_{j}}{s_{i}-qs_{j}}\prod_{j\neq
i}\frac{q-ts_{i}s_{j}}{q-s_{i}s_{j}}\frac{1-t^{-1}s_{i}s_{j}}{1-s_{i}s_{j}}\right)T_{i}^{-1}\right],$

where

$\displaystyle\varphi^{(a^{*},b^{*},c^{*},d^{*})}(x)$ $\displaystyle=$
$\displaystyle\sum_{\epsilon=\pm
1}\frac{\big{(}1-\frac{\epsilon}{q^{1/2}}a^{*}\big{)}\big{(}1-\frac{\epsilon}{q^{1/2}}b^{*}\big{)}\big{(}1-\frac{\epsilon}{q^{1/2}}c^{*}\big{)}\big{(}1-\frac{\epsilon}{q^{1/2}}d^{*}\big{)}}{2(1-t)(1-q^{-1}t)}$
$\displaystyle\qquad\qquad\qquad\times\left\\{\prod_{i=1}^{N}\frac{1-\frac{\epsilon\,t}{q^{1/2}}x_{i}}{1-\frac{\epsilon}{q^{1/2}}x_{i}}\,\frac{1-\frac{\epsilon\,t}{q^{1/2}}x_{i}^{-1}}{1-\frac{\epsilon}{q^{1/2}}x_{i}^{-1}}-t^{N}\right\\},$

and $(a^{*},b^{*},c^{*},d^{*})$ are the dual Koornwinder parameters (3.11).

###### Proof.

We use the formula (3.18) with $(a,b,c,d)\to(a^{*},b^{*},c^{*},d^{*})$. For
simplicity, we work with the values $(a,b,c,d)$ and interchange them with
their duals at the end. We need the following.

######  Lemma B.2.

$\displaystyle\Delta^{(a,b,c,d)}(x)^{-1}\Gamma_{i}\Delta^{(a,b,c,d)}(x)$
$\displaystyle\qquad=\frac{(1-\frac{1}{x_{i}^{2}})(1-\frac{1}{qx_{i}^{2}})}{(1-\frac{1}{ax_{i}})(1-\frac{1}{bx_{i}})(1-\frac{1}{cx_{i}})(1-\frac{1}{dx_{i}})}\prod_{j=1}^{i-1}\frac{1-\frac{qx_{i}}{tx_{j}}}{1-\frac{qx_{i}}{x_{j}}}\prod_{j=i+1}^{N}\frac{1-\frac{x_{j}}{x_{i}}}{1-\frac{x_{j}}{tx_{i}}}\prod_{j\neq
i}\frac{1-\frac{1}{x_{i}x_{j}}}{1-\frac{1}{tx_{i}x_{j}}}\,\Gamma_{i},$
$\displaystyle\Delta^{(a,b,c,d)}(x)^{-1}\Gamma_{i}^{-1}\Delta^{(a,b,c,d)}(x)$
$\displaystyle\qquad=\frac{(1-\frac{q}{ax_{i}})(1-\frac{q}{bx_{i}})(1-\frac{q}{cx_{i}})(1-\frac{q}{dx_{i}})}{(1-\frac{q}{x_{i}^{2}})(1-\frac{q^{2}}{x_{i}^{2}})}\prod_{j=1}^{i-1}\frac{1-\frac{x_{i}}{x_{j}}}{1-\frac{x_{i}}{tx_{j}}}\prod_{j=i+1}^{N}\frac{1-\frac{qx_{j}}{tx_{i}}}{1-\frac{qx_{j}}{x_{i}}}\prod_{j\neq
i}\frac{1-\frac{q}{tx_{i}x_{j}}}{1-\frac{q}{x_{i}x_{j}}}\,\Gamma_{i}^{-1}.$

Using Lemma A.11, the conjugation of the ${\mathcal{D}}_{1}^{(a,b,c,d)}$
operator boils down to that of the terms
$\Phi_{i,\epsilon}^{(a,b,c,d)}\Gamma_{i}^{\epsilon}$. We have:

$\displaystyle\Delta^{(a,b,c,d)}(x)^{-1}\Phi_{i,+}^{(a,b,c,d)}\Gamma_{i}\Delta^{(a,b,c,d)}(x)=\frac{abcd}{q}t^{2N-2}\left(\prod_{j=1}^{i-1}\frac{1-\frac{x_{j}}{tx_{i}}}{1-\frac{x_{j}}{x_{i}}}\frac{1-\frac{qx_{i}}{tx_{j}}}{1-\frac{qx_{i}}{x_{j}}}\right)\Gamma_{i},$
$\displaystyle\Delta^{(a,b,c,d)}(x)^{-1}\Phi_{i,-}^{(a,b,c,d)}\Gamma_{i}^{-1}\Delta^{(a,b,c,d)}(x)$
$\displaystyle\qquad=t^{2N-2}\frac{(1-\frac{a}{x_{i}})(1-\frac{q}{ax_{i}})(1-\frac{b}{x_{i}})(1-\frac{q}{bx_{i}})(1-\frac{c}{x_{i}})(1-\frac{q}{cx_{i}})(1-\frac{d}{x_{i}})(1-\frac{q}{dx_{i}})}{(1-\frac{1}{x_{i}^{2}})(1-\frac{q}{x_{i}^{2}})^{2}(1-\frac{q^{2}}{x_{i}^{2}})}$
$\displaystyle\qquad\qquad\quad\times\,\left(\prod_{j=i+1}^{N}\frac{1-\frac{x_{i}}{tx_{j}}}{1-\frac{x_{i}}{x_{j}}}\frac{1-\frac{qx_{j}}{tx_{i}}}{1-\frac{qx_{j}}{x_{i}}}\prod_{j\neq
i}\frac{1-\frac{q}{tx_{i}x_{j}}}{1-\frac{q}{x_{i}x_{j}}}\frac{1-\frac{x_{i}x_{j}}{t}}{1-x_{i}x_{j}}\right)\Gamma_{i}^{-1}.$

We finally use the formula (3.18) for $m=1$ and substitute the above results
with the change $(a,b,c,d)\to(a^{*},b^{*},c^{*},d^{*})$, and the Theorem
follows. ∎

###  B.2. First Pieri rules for ${\mathfrak{g}}$

We may now specialize the result of Theorem B.1 to ${\mathfrak{g}}$, according
to Table 1. We first treat the cases where
${\mathfrak{g}}=D_{N}^{(1)},C_{N}^{(1)},A_{2N-1}^{(2)}$ for which the constant
term $\varphi^{({\mathfrak{g}}^{*})}(x)$ vanishes101010 Here we show the
explicit dependence on ${\Lambda}_{i}$ and $t$ (rather than
$s_{i}={\Lambda}_{i}\,t^{\xi_{\mathfrak{g}}+N-i}$), as $\xi_{\mathfrak{g}}$
varies with ${\mathfrak{g}}$. This makes the limit $t\to\infty$ easier to
follow.:

$\displaystyle{\mathcal{H}}_{1}^{(D_{N}^{(1)})}$ $\displaystyle=$
$\displaystyle\sum_{i=1}^{N}\left\\{\prod_{j<i}\frac{t^{i-j-1}{\Lambda}_{j}-{\Lambda}_{i}}{t^{i-j}{\Lambda}_{j}-{\Lambda}_{i}}\frac{t^{i-j+1}{\Lambda}_{j}-q{\Lambda}_{i}}{t^{i-j}{\Lambda}_{j}-q{\Lambda}_{i}}\,T_{i}\right.$
$\displaystyle\left.+\prod_{j>i}\frac{t^{j-i+1}{\Lambda}_{i}-q{\Lambda}_{j}}{t^{j-i}{\Lambda}_{i}-q{\Lambda}_{j}}\frac{t^{j-i-1}{\Lambda}_{i}-{\Lambda}_{j}}{t^{j-i}{\Lambda}_{i}-{\Lambda}_{j}}\prod_{j\neq
i}\frac{1-t^{2N-i-j-1}{\Lambda}_{i}{\Lambda}_{j}}{1-t^{2N-i-j}{\Lambda}_{i}{\Lambda}_{j}}\frac{t^{2N-i-j+1}{\Lambda}_{i}{\Lambda}_{j}-q}{t^{2N-i-j}{\Lambda}_{i}{\Lambda}_{j}-q}\,T_{i}^{-1}\right\\},$
$\displaystyle{\mathcal{H}}_{1}^{(C_{N}^{(1)})}$ $\displaystyle=$
$\displaystyle\sum_{i=1}^{N}\left\\{\prod_{j<i}\frac{t^{j-i+1}{\Lambda}_{i}-{\Lambda}_{j}}{t^{j-i}{\Lambda}_{i}-{\Lambda}_{j}}\frac{{\Lambda}_{j}-qt^{j-i-1}{\Lambda}_{i}}{{\Lambda}_{j}-qt^{j-i}{\Lambda}_{i}}\,T_{i}\right.$
$\displaystyle+\frac{1-t^{N-i}{\Lambda}_{i}}{1-t^{N+1-i}{\Lambda}_{i}}\frac{t^{N+2-i}{\Lambda}_{i}-q}{t^{N+1-i}{\Lambda}_{i}-q}\prod_{j>i}\frac{t^{j-i+1}{\Lambda}_{i}-q{\Lambda}_{j}}{t^{j-i}{\Lambda}_{i}-q{\Lambda}_{j}}\frac{{\Lambda}_{j}-t^{j-i-1}{\Lambda}_{i}}{{\Lambda}_{j}-t^{j-i}{\Lambda}_{i}}$
$\displaystyle\qquad\times\left.\prod_{j\neq
i}\frac{1-t^{2N+1-i-j}{\Lambda}_{i}{\Lambda}_{j}}{1-t^{2N+2-i-j}{\Lambda}_{i}{\Lambda}_{j}}\frac{t^{2N+3-i-j}{\Lambda}_{i}{\Lambda}_{j}-q}{t^{2N+2-i-j}{\Lambda}_{i}{\Lambda}_{j}-q}\,T_{i}^{-1}\right\\}.$
$\displaystyle{\mathcal{H}}_{1}^{(A_{2N-1}^{(2)})}$ $\displaystyle=$
$\displaystyle\sum_{i=1}^{N}\left\\{\prod_{j<i}\frac{t^{j-i+1}\Lambda_{i}-\Lambda_{j}}{t^{j-i}\Lambda_{i}-\Lambda_{j}}\frac{{\Lambda}_{j}-qt^{j-i-1}{\Lambda}_{i}}{{\Lambda}_{j}-qt^{j-i}{\Lambda}_{i}}\,T_{i}\right.$
$\displaystyle\left.+\frac{t^{2N-2i}\Lambda_{i}^{2}-1}{t^{2N+1-2i}\Lambda_{i}^{2}-1}\frac{t^{2N+2-2i}\Lambda_{i}^{2}-q^{2}}{t^{2N+1-2i}\Lambda_{i}^{2}-q^{2}}\prod_{j>i}\frac{t^{j-i+1}\Lambda_{i}-q\Lambda_{j}}{t^{j-i}\Lambda_{i}-q\Lambda_{j}}\frac{t^{j-i-1}\Lambda_{i}-\Lambda_{j}}{t^{j-i}\Lambda_{i}-\Lambda_{j}}\right.$
$\displaystyle\left.\times\prod_{j\neq
i}\frac{t^{2N-i-j}\Lambda_{i}\Lambda_{j}-1}{t^{2N+1-i-j}\Lambda_{i}\Lambda_{j}-1}\frac{t^{2N+2-i-j}\Lambda_{i}\Lambda_{j}-q}{t^{2N+1-i-j}\Lambda_{i}\Lambda_{j}-q}\,T_{i}^{-1}\right\\}.$

We now treat the cases ${\mathfrak{g}}=B_{N}^{(1)},D_{N+1}^{(2)},A_{2N}^{(2)}$
for which $\varphi^{({\mathfrak{g}}^{*})}$ has a non-trivial contribution. We
find:

$\displaystyle{\mathcal{H}}_{1}^{(B_{N}^{(1)})}=G^{(B_{N}^{(1)})}(\Lambda)+\sum_{i=1}^{N}\left\\{\prod_{j<i}\frac{t^{j-i+1}{\Lambda}_{i}-{\Lambda}_{j}}{t^{j-i}{\Lambda}_{i}-{\Lambda}_{j}}\frac{{\Lambda}_{j}-qt^{j-i-1}{\Lambda}_{i}}{{\Lambda}_{j}-qt^{j-i}{\Lambda}_{i}}\,T_{i}\right.$
$\displaystyle+\frac{1-t^{2N-2i}{\Lambda}_{i}^{2}}{1-t^{2N+1-2i}{\Lambda}_{i}^{2}}\,\frac{q-t^{2N-2i}{\Lambda}_{i}^{2}}{q-t^{2N+1-2i}{\Lambda}_{i}^{2}}\,\frac{t^{2N+2-2i}{\Lambda}_{i}^{2}-q^{2}}{t^{2N+1-2i}{\Lambda}_{i}^{2}-q^{2}}\,\frac{t^{2N+2-2i}{\Lambda}_{i}^{2}-q}{t^{2N+1-2i}{\Lambda}_{i}^{2}-q}$
$\displaystyle\times\left.\prod_{j>i}\frac{t^{j-i+1}{\Lambda}_{i}-q{\Lambda}_{j}}{t^{j-i}{\Lambda}_{i}-q{\Lambda}_{j}}\frac{{\Lambda}_{j}-t^{j-i-1}{\Lambda}_{i}}{{\Lambda}_{j}-t^{j-i}{\Lambda}_{i}}\prod_{j\neq
i}\frac{1-t^{2N-i-j}{\Lambda}_{i}{\Lambda}_{j}}{1-t^{2N+1-i-j}{\Lambda}_{i}{\Lambda}_{j}}\frac{t^{2N+2-i-j}{\Lambda}_{i}{\Lambda}_{j}-q}{t^{2N+1-i-j}{\Lambda}_{i}{\Lambda}_{j}-q}T_{i}^{-1}\right\\},$

where

$G^{(B_{N}^{(1)})}({\Lambda})=-1+\frac{1}{2}\sum_{\epsilon=\pm
1}\prod_{i=1}^{N}\frac{t^{N+\frac{3}{2}-i}\Lambda_{i}-\epsilon
q^{1/2}}{t^{N+\frac{1}{2}-i}\Lambda_{i}-\epsilon
q^{1/2}}\,\frac{t^{N-\frac{1}{2}-i}\Lambda_{i}-\epsilon
q^{-1/2}}{t^{N+\frac{1}{2}-i}\Lambda_{i}-\epsilon q^{-1/2}}.$
$\displaystyle{\mathcal{H}}_{1}^{(D_{N+1}^{(2)})}=G^{(D_{N+1}^{(2)})}({\Lambda})+\sum_{i=1}^{N}\left\\{\prod_{j<i}\frac{t^{j-i+1}{\Lambda}_{i}-{\Lambda}_{j}}{t^{j-i}{\Lambda}_{i}-{\Lambda}_{j}}\frac{qt^{j-i-1}{\Lambda}_{i}-{\Lambda}_{j}}{qt^{j-i}{\Lambda}_{i}-{\Lambda}_{j}}\,T_{i}\right.$
$\displaystyle+\frac{1-t^{N-i}\Lambda_{i}}{1-t^{N+1-i}\Lambda_{i}}\frac{q^{1/2}-t^{N-i}\Lambda_{i}}{q^{1/2}-t^{N+1-i}\Lambda_{i}}\frac{q-t^{N+2-i}\Lambda_{i}}{q-t^{N+1-i}\Lambda_{i}}\frac{q^{1/2}-t^{N+2-i}\Lambda_{i}}{q^{1/2}-t^{N+1-i}\Lambda_{i}}$
$\displaystyle\left.\times\prod_{j\neq
i}\frac{1-t^{2N+1-i-j}{\Lambda}_{i}{\Lambda}_{j}}{1-t^{2N+2-i-j}{\Lambda}_{i}{\Lambda}_{j}}\frac{q-t^{2N+3-i-j}{\Lambda}_{i}{\Lambda}_{j}}{q-t^{2N+2-i-j}{\Lambda}_{i}{\Lambda}_{j}}\prod_{j>i}\frac{t^{j-i-1}{\Lambda}_{i}-{\Lambda}_{j}}{t^{j-i}{\Lambda}_{i}-{\Lambda}_{j}}\frac{t^{j-i+1}{\Lambda}_{i}-q{\Lambda}_{j}}{t^{j-i}{\Lambda}_{i}-q{\Lambda}_{j}}\,T_{i}^{-1}\right\\},$

where

$G^{(D_{N+1}^{(2)})}({\Lambda})=\frac{(1-tq^{-1/2})(1+q^{-1/2})}{1-tq^{-1}}\left\\{\prod_{i=1}^{N}\frac{q^{1/2}-{\Lambda}_{i}t^{N+2-i}}{q^{1/2}-{\Lambda}_{i}t^{N+1-i}}\frac{1-q^{1/2}{\Lambda}_{i}t^{N-i}}{1-q^{1/2}{\Lambda}_{i}t^{N+1-i}}-1\right\\},$

and finally

$\displaystyle{\mathcal{H}}_{1}^{(A_{2N}^{(2)})}=G^{(A_{2N}^{(2)})}({\Lambda})+\sum_{i=1}^{N}\left\\{\prod_{j<i}\frac{t^{j-i+1}{\Lambda}_{i}-{\Lambda}_{j}}{t^{j-i}{\Lambda}_{i}-{\Lambda}_{j}}\frac{qt^{j-i-1}{\Lambda}_{i}-{\Lambda}_{j}}{qt^{j-i}{\Lambda}_{i}-{\Lambda}_{j}}\,T_{i}\right.$
$\displaystyle+\frac{1-t^{N-i}\Lambda_{i}}{1-t^{N+1-i}\Lambda_{i}}\frac{q-t^{2N-2i+1}\Lambda_{i}^{2}}{q-t^{2N+2-2i}\Lambda_{i}^{2}}\frac{q^{2}-t^{2N+3-2i}\Lambda_{i}^{2}}{q^{2}-t^{2N+2-2i}\Lambda_{i}^{2}}\frac{q-t^{N+2-i}\Lambda_{i}}{q-t^{N+1-i}\Lambda_{i}}$
$\displaystyle\left.\times\prod_{j>i}\frac{t^{j-i-1}{\Lambda}_{i}-{\Lambda}_{j}}{t^{j-i}{\Lambda}_{i}-{\Lambda}_{j}}\frac{t^{j-i+1}{\Lambda}_{i}-q{\Lambda}_{j}}{t^{j-i}{\Lambda}_{i}-q{\Lambda}_{j}}\prod_{j\neq
i}\frac{1-t^{2N+1-i-j}{\Lambda}_{i}{\Lambda}_{j}}{1-t^{2N+2-i-j}{\Lambda}_{i}{\Lambda}_{j}}\frac{q^{2}-t^{2N+3-i-j}{\Lambda}_{i}{\Lambda}_{j}}{q^{2}-t^{2N+2-i-j}{\Lambda}_{i}{\Lambda}_{j}}\,T_{i}^{-1}\right\\},$

where

$\displaystyle G^{(A_{2N}^{(2)})}({\Lambda})$ $\displaystyle=$
$\displaystyle\sum_{\epsilon=\pm 1}\frac{(1-\epsilon tq^{-1/2})(1+\epsilon
q^{-1/2})}{2(1-tq^{-1})}\left\\{\prod_{i=1}^{N}\frac{q^{1/2}-\epsilon{\Lambda}_{i}t^{N+2-i}}{q^{1/2}-\epsilon{\Lambda}_{i}t^{N+1-i}}\frac{\epsilon-q^{1/2}{\Lambda}_{i}t^{N-i}}{\epsilon-q^{1/2}{\Lambda}_{i}t^{N+1-i}}-1\right\\}.$

######  Remark B.3.

When we consider the $q$-Whittaker limit $t\to\infty$, the constant pieces
$G^{({\mathfrak{g}})}({\Lambda})$ in the latter three cases
${\mathfrak{g}}=B_{N}^{(1)},D_{N+1}^{(2)},A_{2N}^{(2)}$ provide an extra
constant term in the limiting Hamiltonians $H_{1}^{({\mathfrak{g}})}$. These
read respectively:

$\lim_{t\to\infty}G^{(B_{N}^{(1)})}(\Lambda)=\frac{q^{-1}}{{\Lambda}_{N-1}{\Lambda}_{N}}-\frac{1+q^{-1}}{{\Lambda}_{N}^{2}},\
\lim_{t\to\infty}G^{(D_{N+1}^{(2)})}(\Lambda)=-\frac{1+q^{-1/2}}{{\Lambda}_{N}},\
\lim_{t\to\infty}G^{(A_{2N}^{(2)})}(\Lambda)=-\frac{1}{{\Lambda}_{N}}.$

##  Appendix C Proof of Theorem 4.10

Here, we provide the proof of Theorem 4.10, that the time translation
operators $g^{({\mathfrak{g}})}$ of Equation (4.13) commute with the limiting
first Pieri operator/Toda Hamiltonian $H_{1}^{({\mathfrak{g}})}$ of Equation
(3.46):

(C.1)
$g^{({\mathfrak{g}})}H_{1}^{({\mathfrak{g}})}=H_{1}^{({\mathfrak{g}})}g^{({\mathfrak{g}})}.$

We divide each Hamiltonian into bulk and boundary pieces, study how $g$
commutes with each, and then sum up the contributions to prove the commutation
(C.1).

######  Lemma C.1.

For all ${\mathfrak{g}}$,

(C.2) $\displaystyle g_{T}\,\left(1-{\Lambda}^{-{\alpha}_{a}}\right)$
$\displaystyle=$
$\displaystyle\left(1-q{\Lambda}^{-{\alpha}_{a}}\right)T^{-{\alpha}_{a}}g_{T}\quad(a=0,1,...,N-1),$
(C.3) $\displaystyle
g_{\Lambda}\,\left(1-{\Lambda}^{-{\alpha}_{a}}\right)T_{a}$ $\displaystyle=$
$\displaystyle\left(1-q^{-1}{\Lambda}^{-{\alpha}_{a+1}}\right)T_{a+1}\,g_{\Lambda}\quad(a=1,2,...,N-2),$
(C.4) $\displaystyle
g_{\Lambda}\,\left(1-{\Lambda}^{-{\alpha}_{a}}\right)T_{a}^{-1}$
$\displaystyle=$
$\displaystyle\left(1-q^{-1}{\Lambda}^{-{\alpha}_{a-1}}\right)T_{a}^{-1}\,g_{\Lambda}\quad(a=1,2,...,N-1),$

with the convention ${\alpha}_{0}=0$.

###  C.1. Bulk and boundary terms

Each of the first Pieri operators (3.46) is the sum of a bulk and boundary
piece, where the bulk piece has the same form for all ${\mathfrak{g}}$:

$H_{1}^{[m]}:=\sum_{a=1}^{m}\left(1-\frac{{\Lambda}_{a}}{{\Lambda}_{a-1}}\right)\,T_{a}+\left(1-\frac{{\Lambda}_{a+1}}{{\Lambda}_{a}}\right)\,\frac{1}{T_{a}},$

with the convention that $\frac{1}{{\Lambda}_{0}}=0$, for suitable value of
$m=m_{\mathfrak{g}}$: $m_{\mathfrak{g}}=N-2$ for
${\mathfrak{g}}=D_{N}^{(1)},B_{N}^{(1)},A_{2N}^{(2)}$, and
$m_{\mathfrak{g}}=N-1$ for all other cases.

The operators $g^{({\mathfrak{g}})}$ of (4.13) also have some common
structures. For ${\mathfrak{g}}=B_{N}^{(1)},A_{2N-1}^{(2)},D_{N+1}^{(2)}$,
$g=g_{T}^{\frac{1}{2}}{\alpha}g_{T}^{\frac{1}{2}}g_{\Lambda}\beta$ where
${\alpha},\beta$ are functions of ${\Lambda}_{N}$ only. For
${\mathfrak{g}}=D_{N}^{(1)}$ $g=g_{T}g_{\Lambda}\beta$ with $\beta$ a function
of ${\Lambda}_{N-1}{\Lambda}_{N}$ only. The cases
${\mathfrak{g}}=C_{N}^{(1)},A_{2N}^{(2)}$ are different, and have a $g$
operator of the form $g=g_{T}g_{\Lambda}{\alpha}g_{T}g_{\Lambda}\beta$, where
${\alpha},\beta$ are functions of ${\Lambda}_{N}$ only.

######  Lemma C.2.

For all ${\mathfrak{g}}$ and $m\leq m_{\mathfrak{g}}$, and $t_{1}=2$ for
${\mathfrak{g}}=C_{N}^{(1)},A_{2N}^{(2)}$, $t_{1}=1$ otherwise,

$g^{({\mathfrak{g}})}\,H_{1}^{[m]}\,(g^{({\mathfrak{g}})})^{-1}=(g_{T}\,g_{\Lambda})^{t_{1}}\,H_{1}^{[m]}\,(g_{T}\,g_{\Lambda})^{-t_{1}}.$

###### Proof.

The only subtleties arise in the cases
${\mathfrak{g}}=B_{N}^{(1)},A_{2N}^{(2)}$. When ${\mathfrak{g}}=B_{N}^{(1)}$,
with $g=g_{T}^{\frac{1}{2}}{\alpha}g_{T}^{\frac{1}{2}}g_{\Lambda}\beta$, the
action of $g_{T}^{\frac{1}{2}}g_{\Lambda}\beta$ on $H_{1}^{[m]}$ creates a
term proportional to $T_{m+1}$, as a consequence of

(C.5)
$g_{T}\,g_{\Lambda}\,\left(1-\frac{{\Lambda}_{a}}{{\Lambda}_{a-1}}\right)\,T_{a}=\left(T_{a}-\frac{{\Lambda}_{a+1}}{{\Lambda}_{a}}\,T_{a+1}\right)\,g_{T}g_{\Lambda}$

for $a=m$. This term commutes with ${\alpha}$ only if $m+1\leq N-1$, hence we
set $m_{\mathfrak{g}}=N-2$. Similarly, when ${\mathfrak{g}}=A_{2N}^{(2)}$,
with $g=g_{T}g_{\Lambda}{\alpha}g_{T}g_{\Lambda}\beta$, the action of
$g_{T}g_{\Lambda}\beta$ on $H_{1}^{[m]}$ creates a term proportional to
$T_{m+1}$ which commutes with ${\alpha}$ only if $m+1\leq N-1$, hence we set
$m_{\mathfrak{g}}=N-2$ as well. ∎

######  Lemma C.3.

We have the commutation relations for $m=1,2,...,N-1$:

$g_{T}\,g_{\Lambda}\,H_{1}^{[m]}=\left\\{H_{1}^{[m]}+\frac{{\Lambda}_{m+1}}{{\Lambda}_{m}}\,\left(\frac{1}{T_{m}}-T_{m+1}\right)\right\\}\,g_{T}\,g_{\Lambda}.$

###### Proof.

We simply sum over the relevant values of $a$ the commutation relations (C.5),
as well as the following:

(C.6)
$g_{T}\,g_{\Lambda}\,\left(1-\frac{{\Lambda}_{a+1}}{{\Lambda}_{a}}\right)\,\frac{1}{T_{a}}=\left(\frac{1}{T_{a}}-\frac{{\Lambda}_{a}}{{\Lambda}_{a-1}}\,\frac{1}{T_{a-1}}\right)\,g_{T}g_{\Lambda},$

obtained by combining (C.2-C.4). ∎

The cases ${\mathfrak{g}}=C_{N}^{(1)}$ and $A_{2N}^{(1)}$ are special, as we
will have to use Lemma C.3 twice. This is deferred to Section C.6 below.

For all the other cases, we apply Lemmas C.2 and C.3 to commute $g$ through
the relevant bulk contributions to $H_{1}^{({\mathfrak{g}})}$ corresponding to
$a=1,2,..,m_{\mathfrak{g}}$. The remaining terms of $H_{1}^{({\mathfrak{g}})}$
are the “boundary” terms, corresponding to $m_{\mathfrak{g}}<a\leq N$ and
possibly to constant terms independent of $T_{a}^{\pm 1}$. We address them
individually in the following Sections C.2 through C.5. In all cases, the
Lemma 4.10 follows from summing all contributions from the bulk
$a=1,2,...,m_{\mathfrak{g}}$ and the boundary $a=m_{\mathfrak{g}}+1,...,N$
plus constant terms.

###  C.2. Case ${\mathfrak{g}}=D_{N}^{(1)}$

Recall that
$g=g_{T}g_{\Lambda}\prod_{n=0}^{\infty}\left(1-\frac{q^{n}}{{\Lambda}_{N-1}{\Lambda}_{N}}\right)^{-1}$,
and that $m_{\mathfrak{g}}=N-2$. For $a=N-1$, we have:

$\displaystyle
g\left(1-\frac{{\Lambda}_{N-1}}{{\Lambda}_{N-2}}\right)T_{N-1}=\left(T_{N-1}-\frac{\Lambda_{N}}{\Lambda_{N-1}}\,T_{N}-\frac{1}{\Lambda_{N-1}\Lambda_{N}}T_{N}^{-1}+\frac{1}{\Lambda_{N-1}^{2}}T_{N-1}^{-1}\right)\,g,$
$\displaystyle
g\,\left(1-\frac{\Lambda_{N}}{\Lambda_{N-1}}\right)\left(1-\frac{1}{\Lambda_{N-1}\Lambda_{N}}\right)T_{N-1}^{-1}=\left(T_{N-1}^{-1}-\frac{\Lambda_{N-1}}{\Lambda_{N-2}}\,T_{N-2}^{-1}\right)\,g,$

and for $a=N$:

$\displaystyle
g\,\left(1-\frac{\Lambda_{N}}{\Lambda_{N-1}}\right)T_{N}=\left(T_{N}-\frac{1}{\Lambda_{N-1}\Lambda_{N}}T_{N-1}^{-1}\right)\,g,$
$\displaystyle
g\,\left(1-\frac{1}{\Lambda_{N-1}\Lambda_{N}}\right)\,T_{N}^{-1}=\left(T_{N}^{-1}-\frac{\Lambda_{N}}{\Lambda_{N-1}}\,T_{N-1}^{-1}\right)\,g.$

Eq. (C.1) follows by summing all bulk and boundary contributions.

###  C.3. Case ${\mathfrak{g}}=B_{N}^{(1)}$

Recall that $g=g_{1}\,g_{2}$, $g_{1}=g_{T}^{\frac{1}{2}}{\alpha}$,
$g_{2}=g_{1}g_{\Lambda}$ with
${\alpha}=\prod_{n=0}^{\infty}\left(1-\frac{q^{n}}{{\Lambda}_{N}^{2}}\right)^{-1}$,
and that $m_{\mathfrak{g}}=N-2$. The terms $a=N-1,N$ and the constant term
read

$\displaystyle
g\,\left(1-\frac{\Lambda_{N}}{\Lambda_{N-1}}\right)\,T_{N-1}^{-1}=T_{N-1}^{-1}\,g_{T}^{\frac{1}{2}}\,{\alpha}g_{T}^{\frac{1}{2}}\,\left(1-\frac{\Lambda_{N-1}}{\Lambda_{N-2}}\right){\alpha}\,g_{\Lambda}=\left(T_{N-1}^{-1}-\frac{\Lambda_{N-1}}{\Lambda_{N-2}}T_{N-2}^{-1}\right)\,g,$
$\displaystyle
g\,\left(1-\frac{\Lambda_{N-1}}{\Lambda_{N-2}}\right)\,T_{N-1}=T_{N-1}\,g_{1}g_{T}^{\frac{1}{2}}\,\left(1-\frac{\Lambda_{N}}{\Lambda_{N-1}}\right)\,{\alpha}g_{\Lambda}=g_{1}\,\left(T_{N-1}-q^{-1}\frac{\Lambda_{N}}{\Lambda_{N-1}}T_{N}^{\frac{1}{2}}T_{N-1}^{\frac{1}{2}}\right)\,g_{2},$
$\displaystyle
g\,\left(1-\frac{\Lambda_{N}}{\Lambda_{N-1}}\right)\,T_{N}=g_{1}\,\left(1-\frac{1}{\Lambda_{N}^{2}}T_{N}^{-1}\right)\left(1-q^{-2}\frac{1}{\Lambda_{N}^{2}}T_{N}^{-1}\right)\,T_{N}\,g_{2},$
$\displaystyle
g\,\left(1-\frac{1}{\Lambda_{N}^{2}}\right)\left(1-q^{2}\frac{1}{\Lambda_{N}^{2}}\right)\,T_{N}^{-1}=g_{1}\,\left(T_{N}^{-1}-q^{-1}\frac{\Lambda_{N}}{\Lambda_{N-1}}T_{N-1}^{-\frac{1}{2}}T_{N}^{-\frac{1}{2}}\right)\,g,$
$\displaystyle
g\,\left(\frac{q^{-2}}{\Lambda_{N-1}\Lambda_{N}}-\frac{1+q^{-2}}{\Lambda_{N}^{2}}\right)=g_{1}\left(\frac{q^{-1}}{\Lambda_{N-1}\Lambda_{N}}T_{N-1}^{-\frac{1}{2}}T_{N}^{-\frac{1}{2}}-q^{2}\frac{1+q^{-2}}{\Lambda_{N}^{2}}T_{N}^{-1}\right)\,g_{2}.$

Summing these five terms gives

$\displaystyle
g\left\\{\left(1-\frac{\Lambda_{N}}{\Lambda_{N-1}}\right)\,T_{N-1}^{-1}+\left(1-\frac{\Lambda_{N-1}}{\Lambda_{N-2}}\right)\,T_{N-1}+\left(1-\frac{\Lambda_{N}}{\Lambda_{N-1}}\right)\,T_{N}\right.$
$\displaystyle\qquad\qquad+\left.\left(1-\frac{1}{\Lambda_{N}^{2}}\right)\left(1-q^{2}\frac{1}{\Lambda_{N}^{2}}\right)\,T_{N}^{-1}+\frac{q^{-2}}{\Lambda_{N-1}\Lambda_{N}}-\frac{1+q^{-2}}{\Lambda_{N}^{2}}\right\\}$
$\displaystyle=\left\\{T_{N-1}+\left(1-\frac{\Lambda_{N}}{\Lambda_{N-1}}\right)T_{N}-\frac{{\Lambda}_{N-1}}{{\Lambda}_{N-2}}T_{N-2}^{-1}+\left(1-\frac{\Lambda_{N}}{\Lambda_{N-1}}\right)T_{N-1}^{-1}\right.$
$\displaystyle\qquad\qquad+\left.\left(1-\frac{1}{\Lambda_{N}^{2}}\right)\left(1-q^{2}\frac{1}{\Lambda_{N}^{2}}\right)\,T_{N}^{-1}+\frac{q^{-2}}{\Lambda_{N-1}\Lambda_{N}}-\frac{1+q^{-2}}{\Lambda_{N}^{2}}\right\\}\,g.$

Eq. (C.1) follows by summing all bulk and boundary contributions.

###  C.4. Case ${\mathfrak{g}}=A_{2N-1}^{(2)}$

Recall that $g=g_{T}g_{\Lambda}{\alpha}$, with
${\alpha}=\prod_{n=0}^{\infty}\left(1-\frac{q^{2n}}{{\Lambda}_{N}^{2}}\right)^{-1}$,
and $m_{\mathfrak{g}}=N-1$. For the two boundary terms for $a=N$ we have:

$\displaystyle g\,\left(1-\frac{\Lambda_{N}}{\Lambda_{N-1}}\right)T_{N}$
$\displaystyle=$ $\displaystyle
g_{T}\,\left(1-q^{-2}\frac{1}{\Lambda_{N}^{2}}\right)T_{N}\,g_{\Lambda}\,{\alpha}=\left(T_{N}-\frac{1}{\Lambda_{N}^{2}}T_{N}^{-1}\right)\,g,$
$\displaystyle g\,\left(1-\frac{1}{\Lambda_{N}^{2}}\right)T_{N}^{-1}\,g^{-1}$
$\displaystyle=$ $\displaystyle
g_{T}\,\left(1-q^{-1}\frac{\Lambda_{N}}{\Lambda_{N-1}}\right)T_{N}^{-1}\,g_{\Lambda}\,{\alpha}=\left(T_{N}^{-1}-\frac{\Lambda_{N}}{\Lambda_{N-1}}T_{N-1}^{-1}\right)\,g.$

Eq. (C.1) follows by summing all bulk and boundary contributions.

###  C.5. Case ${\mathfrak{g}}=D_{N+1}^{(2)}$

Recall that $g=g_{T}g_{\Lambda}{\alpha}$, with
${\alpha}=\prod_{n=0}^{\infty}\left(1-\frac{q^{\frac{n}{2}}}{{\Lambda}_{N}}\right)^{-1}$,
and $m_{\mathfrak{g}}=N-1$. For the two boundary terms for $a=N$ and the
constant term we have:

$\displaystyle
g\,\left(1-\frac{\Lambda_{N}}{\Lambda_{N-1}}\right)T_{N}=\left(T_{N}-\frac{1+q^{-\frac{1}{2}}}{\Lambda_{N}}+\frac{q}{\Lambda_{N}^{2}}T_{N}^{-1}\right)\,g,$
$\displaystyle
g\,\left(1-\frac{1}{\Lambda_{N}}\right)\left(1-\frac{q^{\frac{1}{2}}}{\Lambda_{N}}\right)T_{N}^{-1}=g_{T}\,\left(1-q^{-1}\frac{\Lambda_{N}}{\Lambda_{N-1}}\right)T_{N}^{-1}\,g_{\Lambda}\,{\alpha}=\left(T_{N}^{-1}-\frac{\Lambda_{N}}{\Lambda_{N-1}}T_{N-1}^{-1}\right)\,g,$
$\displaystyle
g\,\left(-\frac{1+q^{-\frac{1}{2}}}{\Lambda_{N}}\right)=-\frac{1+q^{\frac{1}{2}}}{\Lambda_{N}}T_{N}^{-1}\,g.$

Eq. (C.1) follows by summing all bulk and boundary contributions.

###  C.6. Cases ${\mathfrak{g}}=C_{N}^{(1)},A_{2N}^{(2)}$

Recall the general structure of the operator $g=g_{1}g_{2}$ with
$g_{1}=g_{T}g_{\Lambda}{\alpha}$ and $g_{2}=g_{T}g_{\Lambda}\beta$, where
${\alpha},\beta$ are functions of ${\Lambda}_{N}$ only. We note that $g_{2}$
is the same for both cases, with
$\beta=\prod_{n=0}^{\infty}\left(1-\frac{q^{n}}{{\Lambda}_{N}}\right)^{-1}$.
We note also that $H_{1}^{[N-1]}$ is the same for both cases, and Lemma C.3
gives the commutation of $g_{2}$ through $H_{1}^{[N-1]}$:

(C.7)
$g_{2}\,H_{1}^{[N-1]}=\left\\{H_{1}^{[N-1]}+\frac{{\Lambda}_{N}}{{\Lambda}_{N-1}}\,\left(\frac{1}{T_{N-1}}-T_{N}\right)\right\\}\,g_{2}.$

We now need to add the boundary and constant terms. For
${\mathfrak{g}}=C_{N}^{(1)}$ we have for $a=N$ (and no constant term):

$\displaystyle
g_{2}\,\left(1-\frac{\Lambda_{N}}{\Lambda_{N-1}}\right)T_{N}=\left(T_{N}-\frac{q^{-\frac{1}{2}}}{\Lambda_{N}}\right)\,g_{2},\
\
g_{2}\,\left(1-\frac{1}{\Lambda_{N}}\right)T_{N}^{-1}=\left(T_{N}^{-1}-\frac{\Lambda_{N}}{\Lambda_{N-1}}T_{N-1}^{-1}\right)\,g_{2}.$

For ${\mathfrak{g}}=A_{2N}^{(2)}$ we have for $a=N$ and the constant term:

$\displaystyle
g_{2}\,\left(1-\frac{{\Lambda}_{N}}{{\Lambda}_{N-1}}\right)T_{N}$
$\displaystyle=$ $\displaystyle
g_{T}\,\left(1-\frac{q^{-1}}{{\Lambda}_{N}}\right)T_{N}\,g_{\Lambda}\,\beta=\left(T_{N}-\frac{q^{-\frac{1}{2}}}{{\Lambda}_{N}}\right)\,g_{2},$
$\displaystyle g_{2}\,\left(1-\frac{1}{{\Lambda}_{N}}\right)T_{N}^{-1}$
$\displaystyle=$
$\displaystyle\left(T_{N}^{-1}-\frac{{\Lambda}_{N}}{{\Lambda}_{N-1}}T_{N-1}^{-1}\right)\,g_{2},$
$\displaystyle g_{2}\,\left(1-\frac{1}{{\Lambda}_{N}}\right)$ $\displaystyle=$
$\displaystyle\left(1-\frac{q^{\frac{1}{2}}}{{\Lambda}_{N}}T_{N}^{-1}\right)\,g_{2}.$

Adding these to the bulk contributions (C.7), we deduce the following
commutations:

(C.8) $\displaystyle{\mathfrak{g}}=C_{N}^{(1)}:$ $\displaystyle
g_{2}H_{1}=\left\\{H_{1}+\frac{1}{{\Lambda}_{N}}(T_{N}^{-1}-q^{-\frac{1}{2}})\right\\}g_{2},$
(C.9) $\displaystyle{\mathfrak{g}}=A_{2N}^{(2)}:$ $\displaystyle
g_{2}H_{1}=\left\\{H_{1}+\frac{1-q^{\frac{1}{2}}}{{\Lambda}_{N}}(T_{N}^{-1}-q^{-\frac{1}{2}})\right\\}g_{2}.$

We must now commute $g_{1}$ through this. For both cases, we split again
$H_{1}$ into a bulk piece $H_{1}^{[N-1]}$ and boundary pieces and constant
terms. For the bulk contribution, we use again Lemma C.3 to write in both
cases

(C.10)
$g_{1}\,H_{1}^{[N-1]}=\left\\{H_{1}^{[N-1]}+\frac{{\Lambda}_{N}}{{\Lambda}_{N-1}}\,\left(\frac{1}{T_{N-1}}-T_{N}\right)\right\\}\,g_{1}.$

Then for ${\mathfrak{g}}=C_{N}^{(1)}$ and $a=N$:

$\displaystyle g_{1}\,\left(1-\frac{\Lambda_{N}}{\Lambda_{N-1}}\right)T_{N}$
$\displaystyle=$ $\displaystyle T_{N}\,g_{1},$ $\displaystyle
g_{1}\,\left(1-\frac{1}{\Lambda_{N}}\right)T_{N}^{-1}$ $\displaystyle=$
$\displaystyle\left(1-\frac{q^{\frac{1}{2}}}{\Lambda_{N}}T_{N}^{-1}\right)\left(T_{N}^{-1}-\frac{\Lambda_{N}}{\Lambda_{N-1}}T_{N-1}^{-1}\right)\,g_{1},$

and for ${\mathfrak{g}}=A_{2N}^{(2)}$, $a=N$ and the constant term:

$\displaystyle
g_{1}\,\left(1-\frac{{\Lambda}_{N}}{{\Lambda}_{N-1}}\right)T_{N}$
$\displaystyle=$ $\displaystyle
g_{T}\,\left(1-\frac{q^{-\frac{1}{2}}}{{\Lambda}_{N}}\right)T_{N}\,g_{\Lambda}\,{\alpha}=\left(T_{N}-\frac{1}{{\Lambda}_{N}}\right)\,g_{1},$
$\displaystyle
g_{1}\,\left(1-\frac{q^{\frac{1}{2}}}{{\Lambda}_{N}}\right)T_{N}^{-1}$
$\displaystyle=$
$\displaystyle\left(T_{N}^{-1}-\frac{{\Lambda}_{N}}{{\Lambda}_{N-1}}T_{N-1}^{-1}\right)\,g_{1},$
$\displaystyle g_{1}\,\left(1-\frac{q^{-\frac{1}{2}}}{{\Lambda}_{N}}\right)$
$\displaystyle=$
$\displaystyle\left(1-\frac{1}{{\Lambda}_{N}}T_{N}^{-1}\right)\,g_{1}.$

Summing these with the bulk contributions (C.8)and (C.9) respectively, we
arrive at

(C.11) $\displaystyle{\mathfrak{g}}=C_{N}^{(1)}:$ $\displaystyle
g_{1}\,\left\\{H_{1}+\frac{1}{{\Lambda}_{N}}(T_{N}^{-1}-q^{-\frac{1}{2}})\right\\}=H_{1}\,g_{1},$
(C.12) $\displaystyle{\mathfrak{g}}=A_{2N}^{(2)}:$ $\displaystyle
g_{1}\,\left\\{H_{1}+\frac{1-q^{\frac{1}{2}}}{{\Lambda}_{N}}(T_{N}^{-1}-q^{-\frac{1}{2}})\right\\}=H_{1}\,g_{1}.$

Finally, Lemma 4.10 follows for ${\mathfrak{g}}=C_{N}^{(1)},A_{2N}^{(2)}$ by
combining the commutation relations (C.8-C.9) and (C.11-C.12). In all cases,
Eq. (C.1) follows by summing all bulk and boundary contributions.

## References

  * [BFN19] Alexander Braverman, Michael Finkelberg, and Hiraku Nakajima.  Coulomb branches of $3d$ $\mathcal{N}=4$ quiver gauge theories and slices in the affine Grassmannian.  Adv. Theor. Math. Phys., 23(1):75–166, 2019.  With two appendices by Braverman, Finkelberg, Joel Kamnitzer, Ryosuke Kodera, Nakajima, Ben Webster and Alex Weekes.
  * [BZ05] Arkady Berenstein and Andrei Zelevinsky.  Quantum cluster algebras.  Adv. Math., 195(2):405–455, 2005.
  * [Cha02] Oleg A. Chalykh.  Macdonald polynomials and algebraic integrability.  Adv. Math., 166(2):193–259, 2002.
  * [Che95a] Ivan Cherednik.  Double affine Hecke algebras and Macdonald’s conjectures.  Ann. of Math. (2), 141(1):191–216, 1995.
  * [Che95b] Ivan Cherednik.  Macdonald’s evaluation conjectures and difference Fourier transform.  Invent. Math., 122(1):119–145, 1995.
  * [Che97] Ivan Cherednik.  Difference Macdonald-Mehta conjecture.  Internat. Math. Res. Notices, (10):449–467, 1997.
  * [Che05] Ivan Cherednik.  Double affine Hecke algebras, volume 319 of London Mathematical Society Lecture Note Series.  Cambridge University Press, Cambridge, 2005.
  * [Che09] Ivan Cherednik.  Whittaker limits of difference spherical functions.  Int. Math. Res. Not. IMRN, (20):3793–3842, 2009.
  * [CM06] Vyjayanthi Chari and Adriano Moura.  The restricted Kirillov-Reshetikhin modules for the current and twisted current algebras.  Comm. Math. Phys., 266(2):431–454, 2006.
  * [CNV] S. Cecotti, A. Neitzke, and C. Vafa.  R-twisting and 4d/2d correspondences.  Preprint arXiv:1006.3435.
  * [DFK09] Philippe Di Francesco and Rinat Kedem.  $Q$-systems as cluster algebras. II. Cartan matrix of finite type and the polynomial property.  Lett. Math. Phys., 89(3):183–216, 2009.
  * [DFK10] Philippe Di Francesco and Rinat Kedem.  $Q$-systems, heaps, paths and cluster positivity.  Comm. Math. Phys., 293(3):727–802, 2010.
  * [DFK11] Philippe Di Francesco and Rinat Kedem.  Non-commutative integrability, paths and quasi-determinants.  Adv. Math., 228(1):97–152, 2011.
  * [DFK14] Philippe Di Francesco and Rinat Kedem.  Quantum cluster algebras and fusion products.  Int. Math. Res. Not. IMRN, (10):2593–2642, 2014.
  * [DFK17a] Philippe Di Francesco and Rinat Kedem.  Quantum Q systems: from cluster algebras to quantum current algebras.  Lett. Math. Phys., 107(2):301–341, 2017.
  * [DFK17b] Philippe Di Francesco and Rinat Kedem.  Quantum Q systems: from cluster algebras to quantum current algebras.  Lett. Math. Phys., 107(2):301–341, 2017.
  * [DFK18] P. Di Francesco and R. Kedem.  Difference equations for graded characters from quantum cluster algebra.  Transform. Groups, 23(2):391–424, 2018.
  * [DFK19] Philippe Di Francesco and Rinat Kedem.  ($q,t$)-Deformed Q-Systems, DAHA and Quantum Toroidal Algebras via Generalized Macdonald Operators.  Comm. Math. Phys., 369(3):867–928, 2019.
  * [DFK21] Philippe Di Francesco and Rinat Kedem.  Macdonald Operators and Quantum Q-Systems for classical types.  2021\.  Preprint arXiv:1908.00806 [math-ph].
  * [DFKT17] Philippe Di Francesco, Rinat Kedem, and Bolor Turmunkh.  A path model for Whittaker vectors.  J. Phys. A, 50(25):255201, 35, 2017.
  * [Eti99] Pavel Etingof.  Whittaker functions on quantum groups and $q$-deformed Toda operators.  In Differential topology, infinite-dimensional Lie algebras, and applications, volume 194 of Amer. Math. Soc. Transl. Ser. 2, pages 9–25. Amer. Math. Soc., Providence, RI, 1999.
  * [FJMM09] Boris Feigin, Michio Jimbo, Tetsuji Miwa and Evgeny Mukhin.  Fermionic formulas for eigenfunctions of the difference Toda Hamiltonian.  Lett. Math. Phys., 88(1-3), 39, 2009.
  * [FL99] Boris Feigin and Sergey Loktev.  On generalized Kostka polynomials and the quantum Verlinde rule.  In Differential topology, infinite-dimensional Lie algebras, and applications, volume 194 of Amer. Math. Soc. Transl. Ser. 2, pages 61–79. Amer. Math. Soc., Providence, RI, 1999.
  * [FT19] Michael Finkelberg and Alexander Tsymbaliuk.  Multiplicative slices, relativistic Toda and shifted quantum affine algebras.  In Representations and nilpotent orbits of Lie algebraic systems, volume 330 of Progr. Math., pages 133–304. Birkhäuser/Springer, Cham, 2019.
  * [GK13] Alexander B. Goncharov and Richard Kenyon.  Dimers and cluster integrable systems.  Ann. Sci. Éc. Norm. Supér. (4), 46(5):747–813, 2013.
  * [GSV11] Michael Gekhtman, Michael Shapiro, and Alek Vainshtein.  Generalized Bäcklund-Darboux transformations for Coxeter-Toda flows from a cluster algebra perspective.  Acta Math., 206(2):245–310, 2011.
  * [GT21] Roman Gonin and Alexander Tsymbaliuk.  On Sevostyanov’s construction of quantum difference Toda lattices.  Int. Math. Res. Not. IMRN, (12):8885–8945, 2021.
  * [HKKR00] Tim Hoffmann, Johannes Kellendonk, Nadja Kutz, and Nicolai Reshetikhin.  Factorization dynamics and Coxeter-Toda lattices.  Comm. Math. Phys., 212(2):297–321, 2000.
  * [HKO+99] G. Hatayama, A. Kuniba, M. Okado, T. Takagi, and Y. Yamada.  Remarks on fermionic formula.  In Recent developments in quantum affine algebras and related topics (Raleigh, NC, 1998), volume 248 of Contemp. Math., pages 243–291. Amer. Math. Soc., Providence, RI, 1999.
  * [HKO+02] Goro Hatayama, Atsuo Kuniba, Masato Okado, Taichiro Takagi, and Zengo Tsuboi.  Paths, crystals and fermionic formulae.  In MathPhys odyssey, 2001, volume 23 of Prog. Math. Phys., pages 205–272. Birkhäuser Boston, Boston, MA, 2002.
  * [Kac90] Victor G. Kac.  Infinite-dimensional Lie algebras.  Cambridge University Press, Cambridge, third edition, 1990.
  * [Ked08] Rinat Kedem.  $Q$-systems as cluster algebras.  J. Phys. A, 41(19):194011, 14, 2008.
  * [KN99] A. N. Kirillov and M. Noumi.  $q$-difference raising operators for Macdonald polynomials and the integrality of transition coefficients.  In Algebraic methods and $q$-special functions (Montréal, QC, 1996), volume 22 of CRM Proc. Lecture Notes, pages 227–243. Amer. Math. Soc., Providence, RI, 1999.
  * [Koo92] Tom H. Koornwinder.  Askey-Wilson polynomials for root systems of type $BC$.  In Hypergeometric functions on domains of positivity, Jack polynomials, and applications (Tampa, FL, 1991), volume 138 of Contemp. Math., pages 189–204. Amer. Math. Soc., Providence, RI, 1992.
  * [KR87] A. N. Kirillov and N. Yu. Reshetikhin.  Representations of Yangians and multiplicities of the inclusion of the irreducible components of the tensor product of representations of simple Lie algebras.  Zap. Nauchn. Sem. Leningrad. Otdel. Mat. Inst. Steklov. (LOMI), 160(Anal. Teor. Chisel i Teor. Funktsii. 8):211–221, 301, 1987.
  * [LNS20] Edwin Langmann, Masatoshi Noumi and Junichi Shiraishi.  Basic properties of non- stationary Ruijsenaars functions.  SIGMA 16, 105, 2020.
  * [Lin] Mingyan Simon Lin.  Quantum twisted q-systems and graded fermionic sums.  Preprint 2021.
  * [Lin21] Mingyan Simon Lin.  Quantum $Q$-systems and fermionic sums—the non-simply laced case.  Int. Math. Res. Not. IMRN, (2):805–854, 2021.
  * [Mac01] I. G. Macdonald.  Orthogonal polynomials associated with root systems.  Sém. Lothar. Combin., 45:Art. B45a, 40, 2000/01.
  * [Mac95] I. G. Macdonald.  Symmetric functions and Hall polynomials.  Oxford Mathematical Monographs. The Clarendon Press, Oxford University Press, New York, second edition, 1995.  With contributions by A. Zelevinsky, Oxford Science Publications.
  * [Nou95] Masatoshi Noumi.  Macdonald-Koornwinder polynomials and affine Hecke rings.  Number 919, pages 44–55. 1995.  Various aspects of hypergeometric functions (Japanese) (Kyoto, 1994).
  * [NS12] Masatoshi Noumi and Jun’ichi Shiraishi.  A direct approach to the bispectral problem for the ruijsenaars-macdonald $q$-difference operators, 2012.  arXiv:1206.5364 [math.QA].
  * [Rai05] Eric M. Rains.  ${\rm BC}_{n}$-symmetric polynomials.  Transform. Groups, 10(1):63–132, 2005.
  * [Res03] N. Reshetikhin.  Integrability of characteristic Hamiltonian systems on simple Lie groups with standard Poisson Lie structure.  Comm. Math. Phys., 242(1-2):1–29, 2003.
  * [Sah99] Siddhartha Sahi.  Nonsymmetric Koornwinder polynomials and duality.  Ann. of Math. (2), 150(1):267–282, 1999.
  * [Sh19] Junichi Shiraishi.  Affine screening operators, affine Laumon spaces and conjectures concerning non-stationary Ruijsenaars functions.  Journal of Integrable Systems, 3, 1–30, 2019.
  * [SS18] Gus Schrader and Alexander Shapiro.  On $b$-whittaker functions.  2018\.  Preprint arXiv:1806.00747 [math-ph].
  * [SS19] Gus Schrader and Alexander Shapiro.  K-theoretic coulomb branches of quiver auge theories.  2019\.  Preprint arXiv:1910.03186.
  * [S14] Jasper Stokman.  Connection coefficients for basic Harish-Chandra series.  Adv. Math., 250:351–386, 2014.
  * [vD95] J. F. van Diejen.  Commuting difference operators with polynomial eigenfunctions.  Compositio Math., 95(2):183–233, 1995.
  * [vD96] J. F. van Diejen.  Self-dual Koornwinder-Macdonald polynomials.  Invent. Math., 126(2):319–339, 1996.
  * [vDE15] J. F. van Diejen and E. Emsiz.  Integrable boundary interactions for Ruijsenaars’ difference Toda chain.  Comm. Math. Phys., 337(1):171–189, 2015.
  * [Vic18] Panupong Vichitkunakorn.  Conserved quantities of Q-systems from dimer integrable systems.  Electron. J. Combin., 25(1):Paper No. 1.36, 43, 2018.
  * [Wil15] Harold Williams.  $Q$-systems, factorization dynamics, and the twist automorphism.  Int. Math. Res. Not. IMRN, (22):12042–12069, 2015.
  * [YY21] Kohei Yamaguchi and Shintarou Yanagida.  Specializing koornwinder polynomials to macdonald polynomials of type $b,c,d$ and $bc$.  2021\.  Preprint arXiv:2105.00936 [math.QA].


INTEGER MATRIX FACTORISATIONS, SUPERALGEBRAS AND THE QUADRATIC FORM
OBSTRUCTION
Nicholas J. Highama, Matthew C. Lettingtonb and Karl Michael Schmidtb
aDepartment of Mathematics, University of Manchester, Manchester, M13 9PL, UK
bSchool of Mathematics, Cardiff University, 23 Senghennydd Road, Cardiff, UK
CF24 4AG
<EMAIL_ADDRESS><EMAIL_ADDRESS><EMAIL_ADDRESS>

###### Abstract

We identify and analyse obstructions to factorisation of integer matrices into
products $N^{T}N$ or $N^{2}$ of matrices with rational or integer entries. The
obstructions arise as quadratic forms with integer coefficients and raise the
question of the discrete range of such forms. They are obtained by considering
matrix decompositions over a superalgebra. We further obtain a formula for the
determinant of a square matrix in terms of adjugates of these matrix
decompositions, as well as identifying a _co-Latin_ symmetry space.

##  1 Introduction

The question whether a given square integer matrix $M$ can be factorised into
a product of two integer matrices, either in the form of a square $M=N^{2}$ or
(in case of a symmetric positive definite matrix) in the form $M=N^{T}N$, has
a long history in number theory. It is known that if $n\leq 7$ and $M$ is an
$n\times n$ symmetric positive definite matrix with integer entries and
determinant 1, then a factorisation $M=N^{T}N$ with $N$ an $n\times n$ matrix
with integer entries exists. However, there are examples of such matrices with
dimension $n=8$ which cannot be factorised in this way. This result is
mentioned by Taussky (see [14, p. 812]) and goes back to Hermite, Minkowski,
and Mordell Mordell [11].

000 2010 _Mathematics Subject Classification_ : 15A23, 15CA30, 11H55, 05B15.
_Key words and phrases_ : integer matrix factorisation, matrix superalgebra,
quadratic forms, latin squares, adjugate matrix.

The number theoretic properties relating to the factorisation of symmetric
positive definite $n\times n$ integer matrices $M$ with fixed determinant have
classical connections to the theory of positive definite quadratic forms in
$n$ variables, see e.g. [12] and the above references.

In particular, Mordell considered the similarity classes of $n\times n$
matrices with determinant 1, where two such matrices $L$, $M$ are in the same
class if there exists a unimodular integral matrix $N$ such that $M=N^{T}LN$.
The number of such similarity classes is denoted by $h_{n}$. Then a matrix $M$
is in the similarity class of $I_{n}$ (the identity matrix) if and only if
there exists a factorisation $M=N^{T}N$ with an integer matrix $N$. This
implies that the quadratic form classically associated with the symmetric
matrix $M$, $q(x)=x^{T}Mx$, can be written as

$q(x)=x^{T}Mx=x^{T}N^{T}Nx=y^{T}y=\sum_{j=1}^{n}y_{j}^{2},$ (1)

where $y=Nx$, and $N$ has determinant 1. Thus, the factorisation can be used
to write the quadratic form $q(x)$ as a sum of squares of $n$ linear factors.
When $n=8$, such a factorisation may not exist, as Minkowski proved in 1911
that $h_{n}\geq[1+n/8]$, so $h_{n}\geq 2$ if $n=8$. Mordell showed that
$h_{8}=2$ [12], and Ko showed that $h_{9}=2$ as well [9].

In the present paper, we revisit the question of integer matrix factorisation
in the light of recent general results on matrix decompositions [7], [8]. We
establish in Corollary 3.1 that the existence of integer solutions to a
certain quadratic equation is a necessary condition for a matrix factorisation
of the type $M=N^{2}$ or $M=N^{T}N$ (for symmetric positive definite $M$) to
exist. It is interesting to note that solutions to this new type of quadratic
equation associated with a given integer matrix $M$ can also lead to rational
matrix factors $N$ with entries in $\frac{1}{n^{2}}\mathbb{Z}$.

Throughout the paper, we use the classical example of the Wilson matrix [4],
[6], [10], [13]

$W=\left(\begin{array}[]{cccc}5&7&6&5\\\ 7&10&8&7\\\ 6&8&10&9\\\ 5&7&9&10\\\
\end{array}\right)$ (2)

to demonstrate the methodologies under consideration. This integer matrix has
determinant 1 and hence an integer inverse matrix, but is moderately ill
conditioned, despite its small size and entries. It has the integer
factorisation $W=Z^{T}Z$ discovered in [6] with

$Z=\left(\begin{array}[]{cccc}2&3&2&2\\\ 1&1&2&1\\\ 0&0&1&2\\\
0&0&1&1\end{array}\right).$ (3)

We note that the entries of $Z$ are nonnegative and, although the matrix is
not triangular, it has a block upper triangular structure and can be thought
of as a block Cholesky factor of $W$.

The quadratic form associated with the Wilson matrix can be written, using (1)
and (3), as a sum of four squares:

$\displaystyle q(x)=x^{T}Wx$
$\displaystyle=(x_{3}+x_{4})^{2}+(x_{3}+2x_{4})^{2}+(x_{1}+x_{2}+2x_{3}+x_{4})^{2}$
$\displaystyle\quad{}+(2x_{1}+3x_{2}+2x_{3}+2x_{4})^{2}.$ (4)

As $Z$ is a unimodular integer matrix and hence has an integer inverse, it
follows by Lagrange’s four-square theorem that the quadratic form $q$
generated by the Wilson matrix is universal [1] in the sense that it generates
all positive integers as $x$ ranges over $\mathbb{Z}^{4}$.

This shows that integer matrix factorisation is a valuable tool in studying
the quadratic form generated by an integer matrix. As a result of our
considerations in Section 3 below, the question of factorising the Wilson
matrix in the form $W=Z^{T}Z$ is associated with the solutions of the
quadratic equation

$2w^{2}+x_{1}^{2}+x_{1}x_{2}+x_{1}x_{3}+x_{2}^{2}+x_{2}x_{3}+x_{3}^{2}=952.$
(5)

Indeed, a necessary (but not sufficient) condition for $W$ to factorise is
that integer solutions $(w,x_{1},x_{2},x_{3})$ to the quadratic equation (5)
exist. Thus this equation can be considered a quadratic form obstruction to
integer factorisability. This approach to integer matrix factorisation was
briefly alluded to, but not fully worked out in [6].

In Section 2 we derive a useful $\mathrm{S}+\mathrm{V}$ decomposition of
square matrices, first identified in [8], and give explicit formulas for
constructing it that were not given in [8]. In Section 3, we use the concept
of matrix weight, associated with the S factor, to establish the quadratic
form obstruction to integer matrix factorisation and show how a solution of
the corresponding quadratic equation can be used to calculate the matrix
factors. In Section 4, we discuss adjugate matrices in view of the matrix
decomposition and in particular show that the type S part of a matrix is
characterised by having an adjugate with all equal entries. Finally, in
Section 5 we identify the type V part of a matrix as belonging to a space of
co-Latin squares, defined as square matrices with constant sum over all
entries carrying the same symbol in any Latin square.

##  2 The $\mathrm{S}+\mathrm{V}$ decomposition of matrices

The following symmetries of ${n\times n}$ matrices were considered in [8].

(S) A matrix $M=(m_{i,j})_{i,j=1}^{n}\in{\mathbb{R}}^{n\times n}$ has the
constant sum property (or is of type S) if there is a number
$w\in{\mathbb{R}}$, called the weight of the matrix, such that

$\sum_{j=1}^{n}m_{i,j}=\sum_{j=1}^{n}m_{j,i}=nw\qquad(i\in\\{1,\dots,n\\}).$

The vector subspace of ${\mathbb{R}}^{n\times n}$ of matrices having the
constant sum property with some weight is denoted by $S_{n}$ and can be
characterised as

$S_{n}=\bigl{\\{}M\in{\,\mathbb{R}}^{n\times
n}:1_{n}^{T}Mu=0=u^{T}M1_{n}\;(u\in\\{1_{n}\\}^{\bot})\,\bigr{\\}},$

where $1_{n}\in{\mathbb{R}}^{n}$ is the column vector with all entries equal
to 1 and orthogonality is with respect to the standard inner product,
$\\{1_{n}\\}^{\bot}=\\{u\in{\mathbb{R}}^{n}:u^{T}1_{n}=0\\}$ (cf. [8, Thm. 2.6
(a)]).

(V) A matrix $M=(m_{i,j})_{i,j=1}^{n}\in{\mathbb{R}}^{n\times n}$ has the
vertex cross sum property (or is of type V) if

$m_{i,j}+m_{k,l}=m_{i,l}+m_{k,j}\quad(i,j,k,l\in\\{1,\dots,n\\})$

and the matrix entries sum to zero, $\sum_{i,j=1}^{n}m_{i,j}=0$. The vector
subspace of ${\mathbb{R}}^{n\times n}$ of matrices having the vertex cross sum
property is denoted by $V_{n}$ and can be characterised as

$V_{n}=\bigl{\\{}\,M\in{\,\mathbb{R}}^{n\times
n}:u^{T}Mv=0\;(u,v\in\\{1_{n}\\}^{\bot}),\ 1_{n}^{T}M1_{n}=0\,\bigr{\\}}$ (6)

(cf. [8, Thm. 2.6 (e)]). (We derive a surprising alternative characterisation
of this space in Section 5.)

The spaces $S_{n}$ and $V_{n}$ only have the null matrix in common; in fact,
they complement each other and give ${\mathbb{R}}^{n\times n}$ a superalgebra
structure in the following way (cf. [8, Thm. 2.5 (a)]).

######  Theorem 2.1.

${\mathbb{R}}^{n\times n}=S_{n}\oplus V_{n}$. Moreover, $S_{n}$ is a
subalgebra of ${\mathbb{R}}^{n\times n}$, and

$S_{n}S_{n}\subset S_{n},\quad S_{n}V_{n}\subset V_{n},\quad V_{n}S_{n}\subset
V_{n},\quad V_{n}V_{n}\subset S_{n}.$

We show below in Corollary 2.1 that this decomposition of the matrix algebra
is orthogonal.

We begin by showing in Theorem 2.2 that every element $M$ of $V_{n}$ can be
written in the form $M=a1_{n}^{T}+1_{n}b^{T}$, before deriving a formula for
obtaining these vectors $a$ and $b$ in Theorem 2.3.

######  Theorem 2.2.

A matrix $M\in{\mathbb{R}}^{n\times n}$ is an element of $V_{n}$ if and only
if there exist vectors $a,b\in\\{1_{n}\\}^{\bot}$ such that
$M=a1_{n}^{T}+1_{n}b^{T}$.

###### Proof.

Let $M\in V_{n}.$ Then for any $v\in\\{1_{n}\\}^{\bot},$ we have
$Mv\in\\{1_{n}\\}^{\bot\bot}={\mathbb{R}}\,1_{n};$ hence there is a linear
form $f:\\{1_{n}\\}^{\bot}\rightarrow{\mathbb{R}}$ such that $Mv=1_{n}\,f(v)$
$(v\in\\{1_{n}\\}^{\bot}).$ By the Riesz Representation Theorem, there is a
vector $b\in\\{1_{n}\\}^{\bot}$ such that $f(v)=b^{T}v$
$(v\in\\{1_{n}\\}^{\bot}),$ and consequently $Mv=1_{n}b^{T}v$
$(v\in\\{1,n\\}^{\bot}).$ Now every $x\in{\mathbb{R}}^{n}$ can be written in
the form $x=\alpha 1_{n}+v$ with suitable $v\in\\{1,n\\}^{\bot}$ and
$\alpha\in{\mathbb{R}},$ so

$\displaystyle Mx$ $\displaystyle=\alpha M1_{n}+Mv=\alpha
b1_{n}^{T}1_{n}+1_{n}b^{T}v=(a1_{n}^{T}+1_{n}b^{T})(\alpha 1_{n}+v)$
$\displaystyle=(a1_{n}^{T}+1_{n}b^{T})x,$

where $a:=\frac{1}{n}\,M1_{n},$ bearing in mind that $1_{n}^{T}1_{n}=n.$ It
follows from the last equation in (6) that $1_{n}^{T}a=0.$ Conversely, for any
$a,b\in\\{1_{n}\\}^{\bot}$ the matrix $a1_{n}^{T}+1_{n}b^{T}\in V_{n}$. ∎

Clearly $M\in V_{n}$ is symmetric if and only if $a=b$ in the above
representation. Theorem 2.2 shows that $\dim V_{n}=2n-2$, so $\dim
S_{n}=n^{2}-2n+2$ (see also [8, p. 14]), and that the rank of elements of
$V_{n}$ cannot exceed 2. Moreover, this theorem makes the spectral
decomposition of any matrix $M\in V_{n}$ very transparent. Indeed, the range
of $M$ is spanned by the orthogonal vectors $a$ and $1_{n},$ so any
eigenvector for non-zero eigenvalue must be of the form $u=\alpha a+\beta
1_{n}$ with numbers $\alpha,\beta\in{\mathbb{C}}$. Then, bearing in mind that
$1_{n}^{T}a=0$ and $b^{T}1_{n}=0$, and denoting the eigenvalue by $\lambda$,
we find

$\displaystyle\lambda u=Mu$ $\displaystyle=(a1_{n}^{T}+1_{n}b^{T})(\alpha
a+\beta 1_{n})=n\beta a+\alpha b^{T}a1_{n},$

giving $n\beta=\lambda\alpha$ and $b^{T}a\,\alpha=\lambda\beta$. Hence any
nonzero eigenvalue of $M$ is an eigenvalue of the $2\times 2$ matrix
$\left(\begin{array}[]{cc}0&n\\\ b^{T}a&0\end{array}\right)$ and vice versa.
This gives the characteristic polynomial for $M,$

$\chi_{M}(\lambda)=\lambda^{n-2}(\lambda^{2}-nb^{T}a)$

and furthermore the eigenvectors for $M$ as
$v_{\pm}=\sqrt{n}a\pm\sqrt{b^{T}a}1_{n}$ for eigenvalues
$\lambda_{\pm}=\pm\sqrt{nb^{T}a}$ if $b^{T}a>0,$ $v_{\pm}=\sqrt{n}a\pm
i\sqrt{-b^{T}a}1_{n}$ for eigenvalues $\lambda_{\pm}=\pm i\sqrt{-nb^{T}a}$ if
$b^{T}a<0$ and the single eigenvector $a$ for eigenvalue 0 if $b^{T}a=0$ (and
$a\neq 0$; otherwise any non-null vector orthogonal to $b$ will do).

Finally, Theorem 2.2 easily yields the orthogonality with respect to the
Frobenius inner product of the direct sum decomposition of the space of
$n\times n$ matrices into $S_{n}$ and $V_{n}$. The Frobenius inner product of
two matrices $A,B\in{\mathbb{R}}^{n\times n}$ is defined as
$\left<A,B\right>=\mathop{\rm tr}A^{T}B$.

######  Corollary 2.1.

The decomposition ${\mathbb{R}}^{n\times n}=S_{n}\oplus V_{n}$ is orthogonal
with respect to the Frobenius inner product.

###### Proof.

Let $S\in S_{n},$ $V\in V_{n};$ then by Theorem 2.2 there are vectors
$a,b\in\\{1_{n}\\}^{\bot}$ such that $V=a1_{n}^{T}+1_{n}b^{T}$. Hence

$\displaystyle\left<S,V\right>$ $\displaystyle=\mathop{\rm
tr}S^{T}V=\mathop{\rm
tr}(S^{T}a1_{n}^{T}+S^{T}1_{n}b^{T})=1_{n}^{T}S^{T}a+b^{T}S^{T}1_{n}$
$\displaystyle=a^{T}S1_{n}+b^{T}S^{T}1_{n}=0,$

since $S1_{n}=S^{T}1_{n}=nw_{S}1_{n}$, where $w_{S}$ is the weight of the
matrix $S$. ∎

We can calculate the weight of any matrix in $S_{n}$ as the mean of all matrix
entries. In fact, it is meaningful to define the linear form

$\mathop{\rm wt}:{\mathbb{R}}^{n\times
n}\rightarrow{\mathbb{R}},\quad\mathop{\rm
wt}M=\frac{1}{n^{2}}\sum_{i,j=1}^{n}m_{i,j}$

giving the weight of any $n\times n$ matrix $M=(m_{i,j})_{i,j=1}^{n}$. Theorem
2.2 immediately shows that $\mathop{\rm wt}M=0$ for all $M\in V_{n}$.

Any $M\in S_{n}$ can be uniquely decomposed into

$M=M_{0}+(\mathop{\rm wt}S)\,{\cal E}_{n},$

where $M_{0}\in S_{n}$ has weight 0 and ${\cal E}_{n}=1_{n}1_{n}^{T}\in S_{n}$
is the $n\times n$ matrix with all entries equal to 1. In conjunction with the
previous observations, this gives rise to the following unique decomposition
of general $n\times n$ matrices, including an explicit formula for the
calculation of the parts.

######  Theorem 2.3.

Let $M=(m_{i,j})\in\mathbb{R}^{n\times n}$. Then there is a unique
decomposition

$M=M_{V}+M_{0}+(\mathop{\rm wt}M)\,{\cal E}_{n},$

where $M_{0}\in S_{n}$ with weight 0 and $M_{V}=a1_{n}^{T}+1_{n}b^{T}\in
V_{n}$, and the entries of the vectors $a$ and $b$ are given by

$\displaystyle a_{i}$
$\displaystyle=\frac{1}{n}\sum_{j=1}^{n}m_{i,j}-\mathop{\rm
wt}M\quad(i\in\\{1,\dots,n\\}),$ $\displaystyle b_{j}$
$\displaystyle=\frac{1}{n}\sum_{i=1}^{n}m_{i,j}-\mathop{\rm
wt}M\quad(j\in\\{1,\dots,n\\}).$

In particular, if $M$ is an integer matrix then the vectors $n^{2}a$ and
$n^{2}b$ have integer entries.

###### Proof.

Let $\\{\varepsilon_{j}:j\in\\{1,\dots,n-1\\}\\}$ be an orthonormal basis of
$\\{1_{n}\\}^{\bot}.$ Then

$\left\\{\frac{1}{\sqrt{n}}\varepsilon_{j}1_{n}^{T},\frac{1}{\sqrt{n}}1_{n}\varepsilon_{j}^{T}:j\in\\{1,\dots,n-1\\}\right\\}$

is an orthonormal basis, with respect to the Frobenius inner product, of the
($2n-2$-dimensional) space $V_{n}.$ Indeed, for
$j,k\in\\{1,\dots,~{}n~{}-~{}1\\}$ we have

$\left<\varepsilon_{j}1_{n}^{T},1_{n}\varepsilon_{k}^{T}\right>=\mathop{\rm
tr}1_{n}\varepsilon_{j}^{T}1_{n}\varepsilon_{k}^{T}=0$

since $\varepsilon_{j}^{T}1_{n}=0,$ and

$\left<\varepsilon_{j}1_{n}^{T},\varepsilon_{k}1_{n}^{T}\right>=\mathop{\rm
tr}1_{n}\varepsilon_{j}^{T}\varepsilon_{k}1_{n}^{T}=\delta_{j\,k}\,\mathop{\rm
tr}1_{n}1_{n}^{T}=n\delta_{j\,k}$

and analogously for
$\left<1_{n}\varepsilon_{j}^{T},1_{n}\varepsilon_{k}^{T}\right>.$

These observations enable us to find $M_{V}$, the $V_{n}$ part of
$M\in{\mathbb{R}}^{n\times n}$ by orthogonal projection using the above
orthonormal basis. We have

$\displaystyle M_{V}$
$\displaystyle=\sum_{j=1}^{n-1}\left(\left<\frac{1}{\sqrt{n}}\varepsilon_{j}1_{n}^{T},M\right>\frac{1}{\sqrt{n}}\varepsilon_{j}1_{n}^{T}+\left<\frac{1}{\sqrt{n}}1_{n}\varepsilon_{j}^{T},M\right>\frac{1}{\sqrt{n}}1_{n}\varepsilon_{j}^{T}\right)$
$\displaystyle=\frac{1}{n}\sum_{j=1}^{n-1}\left(\mathop{\rm
tr}(1_{n}\varepsilon_{j}^{T}M)\varepsilon_{j}1_{n}^{T}+\mathop{\rm
tr}(\varepsilon_{j}1_{n}^{T}M)1_{n}\varepsilon_{j}^{T}\right)$
$\displaystyle=\frac{1}{n}\sum_{j=1}^{n-1}(\varepsilon_{j}^{T}M1_{n})\varepsilon_{j}1_{n}^{T}+\frac{1}{n}\sum_{j=1}^{n-1}1_{n}((1_{n}^{T}M\varepsilon_{j})\varepsilon_{j})^{T}=a_{M}1_{n}^{T}+1_{m}b_{M}^{T}$

with

$a_{M}=\frac{1}{n}\sum_{j=1}^{n-1}(\varepsilon_{j}^{T}M1_{n})\varepsilon_{j},\qquad
b_{M}=\frac{1}{n}\sum_{j=1}^{n-1}(1_{n}^{T}M\varepsilon_{j})\varepsilon_{j}.$

Observing that
$\\{\varepsilon_{j}:j\in\\{1,\dots,n-1\\}\\}\cup\\{\frac{1}{\sqrt{n}}1_{n}\\}$
is an orthonormal basis of ${\mathbb{R}}^{n},$ we find

$\sum_{j=1}^{n-1}(\varepsilon_{j}^{T}v)\varepsilon_{j}=v-\frac{1}{n}(1_{n}^{T}v)1_{n}\qquad(v\in{\mathbb{R}}^{n}).$

Hence

$a_{M}=\frac{1}{n}\left(M1_{n}-\frac{1}{n}n(1_{n}^{T}M1_{n})1_{n}\right),$

and

$b_{M}=\frac{1}{n}n\left(M^{T}1_{n}-\frac{1}{n}n(1_{n}^{T}M^{T}1_{n})1_{n}\right).\qed$

######  Example: the Wilson matrix.

For the Wilson matrix $W$ this decomposition takes the form

$16W=1_{4}\left(\begin{array}[]{c}-27\\\ 9\\\ 13\\\
5\end{array}\right)^{T}+\left(\begin{array}[]{c}-27\\\ 9\\\ 13\\\
5\end{array}\right)1_{4}^{T}+\left(\begin{array}[]{cccc}15&11&-9&-17\\\
11&23&-13&-21\\\ -9&-13&15&7\\\ -17&-21&7&31\\\ \end{array}\right)+119{\cal
E}_{4}$ (7)

and for the integer matrix factor $Z$ of eq. (3)

$16Z=1_{4}\left(\begin{array}[]{c}17\\\ 1\\\ -7\\\
-11\end{array}\right)^{T}+\left(\begin{array}[]{c}-7\\\ -3\\\ 5\\\
5\end{array}\right)1_{4}^{T}+\left(\begin{array}[]{cccc}3&15&-9&-9\\\
3&-1&7&-9\\\ -5&-9&-1&15\\\ -1&-5&3&3\\\ \end{array}\right)+19{\cal E}_{4}.$
(8)

The last terms on the right-hand side of these formulae correctly yield
$\mathop{\rm wt}W=\frac{119}{16}$ and $\mathop{\rm wt}Z=\frac{19}{16}$.

##  3 Integer factorisation of matrices and the quadratic form obstruction

On the basis of the matrix decomposition established in the preceding section,
we now derive the quadratic equation arising from balancing the weights in a
matrix factorisation.

######  Theorem 3.1.

Let $M,N\in\mathbb{R}^{n\times n}$, and let $M=M_{V}+M_{0}+(\mathop{\rm
wt}M)\,{\cal E}_{n}$ and $N=N_{V}+N_{0}+(\mathop{\rm wt}N)\,{\cal E}_{n}$ be
their decompositions as in Theorem 2.3, with $N_{V}=a1_{n}^{T}+1_{n}b^{T}$.

(i) If $M=N^{T}N$ then

$\mathop{\rm wt}M=|a|^{2}+n\,(\mathop{\rm wt}N)^{2},$ (9)

$M_{V}=y1_{n}^{T}+1_{n}y^{T}$ with $y=N_{0}^{T}a+n\,(\mathop{\rm wt}N)b$, and
$M_{0}=N_{0}^{T}N_{0}+n\,bb^{T}$.

(ii) If $M=N^{2}$ then

$\mathop{\rm wt}M=b^{T}a+n\,(\mathop{\rm wt}N)^{2},$ (10)

$M_{V}=y1_{n}^{T}+1_{n}z^{T}$ with $y=N_{0}^{T}a+n(\mathop{\rm wt}N)a$ and
$z=N_{0}^{T}b+(\mathop{\rm wt}N)b$, and $M_{0}=N_{0}^{2}+n\,ab^{T}$.

###### Proof.

(i) Since $N_{0}{\cal E}_{n}=N^{T}_{0}{\cal E}_{n}={\cal O}_{n}$, we have

$\displaystyle N^{T}N=$ $\displaystyle(N_{V}^{T}N_{V}+N_{0}^{T}N_{0})$
$\displaystyle+(N_{V}^{T}N_{0}+N_{0}^{T}N_{V}+(\mathop{\rm wt}N)N_{V}^{T}{\cal
E}_{n}+(\mathop{\rm wt}N){\cal E}_{n}N_{V})+n(\mathop{\rm wt}N)^{2}{\cal
E}_{n},$

where, by the superalgebra property, the matrices in the first bracket lie in
$S_{n}$, the matrices in the second bracket in $V_{n}$. Writing
$N_{V}=a1_{n}^{T}+1_{n}b^{T}$, we find $N_{V}^{T}N_{V}=a^{T}a{\cal
E}_{n}+nbb^{T}$, $N_{V}^{T}N_{0}=1_{n}a^{T}N_{0}$,
$N_{0}^{T}N_{V}=N_{0}^{T}a1_{n}^{T}$, $N_{V}^{T}{\cal E}_{n}=nb1_{n}^{T}$ and
${\cal E}_{n}N_{V}=n1_{n}b^{T}$. Hence the uniqueness of the decomposition of
$M$ by Theorem 2.3 gives the claimed identities.

(ii) Similarly, we find

$\displaystyle N^{2}=$
$\displaystyle(N_{V}^{2}+N_{0}^{2})+(N_{V}N_{0}+N_{0}N_{V}+(\mathop{\rm
wt}N)N_{V}{\cal E}_{n}+(\mathop{\rm wt}N){\cal E}_{n}N_{V})$
$\displaystyle+n(\mathop{\rm wt}N)^{2}{\cal E}_{n},$

and, setting $N_{V}=a1_{n}^{T}+1_{n}b^{T}$, note that $N_{V}^{2}=b^{T}a{\cal
E}_{n}+nab^{T}$, $N_{V}N_{0}=1_{n}b^{T}N_{0}$, $N_{0}N_{V}=N_{0}a1_{n}^{T}$,
$N_{V}{\cal E}_{n}=na1_{n}^{T}$ and ${\cal E}_{n}N_{V}=n1_{n}b^{T}$. This
gives the claimed identities by uniqueness of decomposition. ∎

The above theorem gives the following quadratic form obstructions to the
factorisation of an integer matrix $M$ into either $M=N^{T}N$ or $M=N^{2}$
with an integer matrix $N$.

######  Corollary 3.1.

Given a matrix $M\in\mathbb{Z}^{n\times n}$, it is necessary for the existence
of a factorisation $M=N^{T}N$ with $N\in\mathbb{Z}^{n\times n}$ that
$n^{2}\mathop{\rm wt}N\in\mathbb{Z}$ and that the vector components
$n^{2}a_{j},n^{2}b_{j}\in\mathbb{Z}$ $(j\in\\{1,\dots,n\\})$, where
$N=a1_{n}^{T}+1_{n}b^{T}+N_{0}+w_{N}{\cal E}_{n}$ is the decomposition of $N$
as in Theorem 2.3, form a solution of the quadratic equation

$\displaystyle n^{4}\mathop{\rm wt}M$ $\displaystyle=n\,(n^{2}\mathop{\rm
wt}N)^{2}+\sum_{j=1}^{n}(n^{2}a_{j})^{2}.$

Similarly, for the decomposition $M=N^{2}$, the quadratic equation is

$\displaystyle n^{4}\mathop{\rm wt}M$ $\displaystyle=n\,(n^{2}\mathop{\rm
wt}N)^{2}+\sum_{j=1}^{n}(n^{2}a_{j})(n^{2}b_{j}).$

Suppose we are given a symmetric integer matrix $M\in\mathbb{Z}^{n\times n}$
and have found a solution
$(w_{N},a_{1},\dots,a_{n})\in\frac{1}{n^{2}}\mathbb{Z}^{n+1}$ of the quadratic
equation associated with the factorisation $M=N^{T}N$,

$n^{4}\mathop{\rm wt}M=n^{5}w_{N}^{2}+n^{4}\sum_{j=1}^{n}a_{j}^{2}.$

To complete the factorisation, we need to identify a vector $b$ and a matrix
$N_{0}$ satisfying the equations in Theorem 3.1 (i). Using the decomposition
$M=y1_{n}^{T}+1_{n}y^{T}+M_{0}+(\mathop{\rm wt}M)\,{\cal E}_{n}$, we have
$b=\frac{1}{nw_{N}}\,(y-N_{0}^{T}a)$, so we only need to find a solution
$N_{0}$ of the quadratic matrix equation

$N_{0}^{T}(aa^{T}+nw_{N}^{2}{\cal
I}_{n})N_{0}-N_{0}^{T}ay^{T}-ya^{T}N_{0}=nw_{N}^{2}M_{0}-yy^{T}$ (11)

or, setting $N_{0}=L-\frac{1}{a^{T}a+nw_{N}^{2}}\,ay^{T}$, the simpler
quadratic

$L^{T}(aa^{T}+nw_{N}^{2}{\cal
I}_{n})L=nw_{N}^{2}\Bigl{(}M_{0}+\frac{1}{a^{T}a+nw_{N}^{2}}yy^{T}\Bigr{)}.$
(12)

The right-hand side and the middle factor on the left-hand side of these
equations are determined in terms of the $\mathrm{S}+\mathrm{V}$ decomposition
of the given matrix $M$ and the particular solution of the factorisation
quadratic form considered. Although determining $N_{0}$ or $L$ from these
equations is a factorisation problem of similar type to the original equation
$M=N^{T}N$, we found that their solution was computationally more effective.

######  Remark 1.

We remark in passing that the matrix $aa^{T}+nw_{N}^{2}$ appearing as a middle
factor on the left-hand side of eq. (12) is symmetric positive definite and
can be written as $A^{T}A$, where
$A=\sqrt{nw_{N}^{2}}+\frac{\sqrt{a^{T}a+nw_{N}^{2}}-\sqrt{nw_{N}^{2}}}{a^{T}a}\,aa^{T}$
and
$A^{-1}=\frac{1}{\sqrt{nw_{N}^{2}}}+\frac{1}{a^{T}a}\left(\frac{1}{\sqrt{a^{T}a+nw_{N})N^{2}}}-\frac{1}{\sqrt{nw_{N}^{2}}}\right)aa^{T}$.

Thus $N_{0}=A^{-1}L-\frac{1}{a^{T}a+nw_{N}^{2}}\,ay^{T}$, where

$L^{T}L=nw_{N}^{2}\Bigl{(}M_{0}+\frac{1}{a^{T}a+nw_{N}^{2}}yy^{T}\Bigr{)},$

but due to the presence of square roots this is unlikely to give a rational
solution $N_{0}$.

######  Example: the Wilson matrix.

The quadratic equation arising from balancing the weights in the assumed
factorisation of the Wilson matrix $W=Z^{T}Z$ is

$\displaystyle\mathop{\rm wt}W=\frac{119}{16}$
$\displaystyle=a_{1}^{2}+a_{2}^{2}+a_{3}^{2}+a_{4}^{2}+4w_{Z}^{2}$
$\displaystyle=2(a_{1}^{2}+a_{2}^{2}+a_{3}^{2}+a_{1}a_{2}+a_{1}a_{3}+a_{2}a_{3})+4\,w_{Z}^{2},$

as $a_{4}=-a_{1}-a_{2}-a_{3}$. Multiplying the equation by $2^{7}$ and setting
$x_{i}=16a_{i}$, $w=16w_{Z}$, we find that a necessary condition for the
integer factorisation of the Wilson matrix is that there are integer solutions
to the quadratic equation (5)

$2w^{2}+x_{1}^{2}+x_{1}x_{2}+x_{1}x_{3}+x_{2}^{2}+x_{2}x_{3}+x_{3}^{2}=952.$

Solving this equation for $w,x_{1},x_{2},x_{3}$ in Mathematica 11.0 on a PC
with a Intel Core i7 6500CPU, gave the $1728$ solutions in just under 6
seconds. Exactly one third ($576$) of these solutions lead to rational matrix
factorisations $W=Z^{T}Z$ with $Z\in\frac{1}{16}\mathbb{Z}^{4\times 4}$.

The process of converting solutions into matrix factors, i.e. of finding
suitable vectors $b$ and matrices $N_{0}$ satisfying the equations of Theorem
3.1 (i) as outlined after Corollary 3.1, took considerably longer at $34$
minutes. Our approach involved utilising (11), in which the vector $b$ is
eliminated and the right-hand-side completely determined for a given factor
weight $w_{N}$. Substituting potential solutions for the vector $a$ and weight
$w_{N}$, thus reduces the general problem of finding the matrix $N_{0}$ to
that of an $(n-1)\times(n-1)$ unknown matrix. For the Wilson matrix this is a
$9$-dimensional problem, and then for each of the $9$-dimensional solutions
for $N_{0}$, the vector $b$ can be quickly constructed. Adding together the
elements of the decomposition we then recover the rational matrix factors $N$.

To a large part, the multiplicity of solutions to the factorisation problem is
expected. Indeed, it is clear that if $M=N^{T}N$ and $U$ is an orthogonal
matrix, then $UN$ is another solution of the factorisation problem;
conversely, if $\det M=1$ and $N$ and $N^{\prime}$ are solutions with integer
entries, then $N^{\prime}=UN$, where $U=N^{\prime}N^{-1}$ is an integer
orthogonal matrix. It is also not hard to see that any integer orthogonal
matrix is a signed permutation matrix, i.e. a matrix which has exactly one
non-zero entry, either $1$ or $-1$, in each row and in each column.

It is therefore natural to classify the factorisation matrices (integer or
rational) modulo left multiplication with integer orthogonal matrices. For the
factorisations of Wilson’s matrix obtained through the above procedure, this
gives three distinct classes, represented by the matrix $Z$ of (3) and the
further two matrices

$Z^{\prime}=\left(\begin{array}[]{cccc}\frac{1}{2}&1&0&1\\\
\frac{3}{2}&2&3&3\\\ \frac{1}{2}&1&0&0\\\ \frac{3}{2}&2&1&0\\\
\end{array}\right),\quad Z^{\prime\prime}=\left(\begin{array}[]{@{\mskip
2mu}rrrr}.0mu\lx@intercol\hfil\frac{3}{2}&2&2&2\\\
.0mu\lx@intercol\hfil\frac{3}{2}&2&2&1\\\
.0mu\lx@intercol\hfil\frac{1}{2}&1&1&2\\\
.0mu\lx@intercol\hfil-\frac{1}{2}&-1&1&1\\\ \end{array}\right).$ (13)

The three factorisations $W=Z^{T}Z=Z^{\prime T}Z^{\prime}=Z^{\prime\prime
T}Z^{\prime\prime}$ correspond to the solutions
$(w,x_{1},x_{2},x_{3})=(19,17,1,-7)$, $(w,x_{1},x_{2},x_{3})=(18,-8,20,-12)$
and $(w,x_{1},x_{2},x_{3})=(19,11,7,-1)$ of eq. (5), respectively. Note that
$Z^{\prime}$ and $Z^{\prime\prime}$ are not integer matrices, so the
equivalence class of $Z$ comprises all integer factorisations of the Wilson
matrix.

##  4 Determinants, Decompositions and Adjugates

In this section we consider how the matrix decomposition of Theorem 2.3 can be
used to obtain expressions for the determinant of a matrix. We begin by
recalling that the adjugate of an $n\times n$ matrix $M$ is defined by

$\mathop{\rm adj}M=((-1)^{i+j}\det(M_{j,i}))_{i,j=1}^{n},$

where $M_{j,i}$ denotes the submatrix of $M$ obtained by deleting row $j$ and
column $i$; $\mathop{\rm adj}M$ is the transpose of the matrix of cofactors of
$M$. If $M$ is nonsingular, then the adjugate can be written as $\mathop{\rm
adj}M=(\det M)M^{-1}$.

The adjugate of a matrix appears naturally when calculating the determinant of
a rank-1 perturbation of a matrix, as shown in the following lemma (cf. [5]).

######  Lemma 4.1.

Let $M\in\mathbb{C}^{n\times n}$ and $u,v\in\mathbb{C}^{n}$. Then

$\det\left(M+uv^{T}\right)=\det\left(M\right)+v^{T}{\rm adj}\left(M\right)u.$

###### Proof.

Suppose $M$ is a nonsingular $n\times n$ matrix. Then

$\displaystyle\det(M+uv^{T})$
$\displaystyle=\det(M)\det(I+M^{-1}uv^{T})=\det(M)(1+v^{T}M^{-1}u)$
$\displaystyle=\det(M)+v^{T}(\det(M)M^{-1})u.$

This gives the stated identity for regular matrices. The general case follows
from the facts that the set of regular matrices is dense in
$\mathbb{C}^{n\times n}$ and that the determinant and adjugate are continuous
functions of the matrix. ∎

######  Corollary 4.1.

Let $M_{0}\in S_{n}$ with weight 0, and $M=M_{0}+w_{M}{\cal E}_{n}$ with some
$w_{M}\in\mathbb{C}$. Then

$\det M=w_{M}\,1_{n}^{T}({\rm adj}M_{0})1_{n}=n^{2}(\mathop{\rm
wt}M)(\mathop{\rm wt}\mathop{\rm adj}M_{0}).$

###### Proof.

Writing ${\cal E}_{n}=1_{n}1_{n}^{T}$, we have $M=M_{0}+w_{M}1_{n}1_{n}^{T}$,
and applying Lemma 4.1 gives $\det M=\det M_{0}+w_{m}1_{n}^{T}({\rm
adj}M_{0})1_{n}$. As $M_{0}1_{n}=0$, ${\det M_{0}=0}$. ∎

######  Remark 2.

If $M=a1_{n}^{T}+1_{n}b^{T}+M_{0}+(\mathop{\rm wt}M){\cal E}_{n}$ is the
decomposition of $M\in\mathbb{C}^{n\times n}$ as in Theorem 2.3, where
$a,b\in\\{1_{n}\\}^{\bot}$ and $M_{0}\in S_{n}$ with weight 0, then Lemma 4.1
and Corollary 4.1 can be used to derive the formula

$\displaystyle\det M$ $\displaystyle=(\mathop{\rm
wt}M-1)\,1_{n}^{T}\mathop{\rm adj}\left(M_{0}-ab^{T}\right)1_{n}$
$\displaystyle\qquad+(b+1_{n})^{T}\mathop{\rm
adj}\left(M_{0}-ab^{T}+(\mathop{\rm wt}M-1){\cal E}_{n}\right)(a+1_{n}).$

######  Example: the Wilson matrix.

For the Wilson matrix we have $\det(W_{S})=\frac{357}{8}$, and ${\rm
adj}(W_{0})=\frac{3}{8}{\cal E}_{4}$, where

$\displaystyle W_{S}=W_{0}+w_{W}$
$\displaystyle=\frac{1}{8}\left(\begin{array}[]{cccc}67&65&55&51\\\
65&71&53&49\\\ 55&53&67&63\\\ 51&49&63&75\\\ \end{array}\right)$
$\displaystyle=\frac{1}{16}\left(\begin{array}[]{cccc}15&11&-9&-17\\\
11&23&-13&-21\\\ -9&-13&15&7\\\ -17&-21&7&31\\\
\end{array}\right)+\frac{119}{16}{\cal E}_{4}.$

By Corollary 4.1, $\det(W_{S})=16(\mathop{\rm wt}W)(\mathop{\rm wt}\mathop{\rm
adj}W_{0})$, and as $\mathop{\rm adj}W_{0}=\frac{3}{8}{\cal E}_{4}$, this
equates to $16\,\frac{119}{16}\,\frac{3}{8}=\frac{357}{8}$.

Thus the adjugate of the weightless type S part of the Wilson matrix is a
multiple of ${\cal E}_{4}$. As we show in the following theorem, the adjugate
of a weightless type S matrix is in fact always a scalar multiple of ${\cal
E}_{n}$. Moreover, the converse holds in the sense that a matrix whose
adjugate is a non-zero multiple of ${\cal E}_{n}$ must be a weightless type S
matrix. Of course any matrix of rank $n-2$ or less has adjugate $0=0{\cal
E}_{n}$.

######  Theorem 4.1.

(a) Let $M\in S_{n}$ with weight 0. Then there is a constant $w\in\mathbb{C}$
such that $\mathop{\rm adj}M=w\,{\cal E}_{n}$.

(b) Let $M\in\mathbb{C}^{n\times n}$ such that $\mathop{\rm adj}M=w\,{\cal
E}_{n}$ with some $w\in\mathbb{C}\setminus\\{0\\}$. Then $M\in S_{n}$,
$\mathop{\rm wt}M=0$ and $\mathop{\mathrm{rank}}M=n-1$.

###### Proof.

(a) Let $j\in\\{1,\dots,n\\}$ and denote by $m_{1}$, …,
$m_{n}\in\mathbb{R}^{n-1}$ the columns of the matrix $M$ with the $j$-th row
omitted. Let $k\in\\{1,\dots,n-1\\}$. Then

$m_{k}=-\sum_{l\in\\{1,\dots,n\\},l\neq k}m_{l},$

as the row sums of $M$ vanish; hence

$\displaystyle\det(m_{1},$
$\displaystyle\dots,m_{k-1},m_{k},m_{k+2},\dots,m_{n})$
$\displaystyle=\det(m_{1},\dots,m_{k-1},-\sum_{l\neq
k}m_{l},m_{k+2},\dots,m_{n})$ $\displaystyle=-\sum_{l\neq
k}\det(m_{1},\dots,m_{k-1},m_{l},m_{k+2},\dots,m_{n})$
$\displaystyle=-\det(m_{1},\dots,m_{k-1},m_{k+1},m_{k+2},\dots,m_{n}),$

noting that all the determinants in the sum vanish except for the term
$l=k+1$. This shows that the $(j,k+1)$ entry and the $(j,k)$ entry of
$\mathop{\rm adj}M$ are equal. Since this holds for all
$k\in\\{1,\dots,n-1\\}$ and for all $j\in\\{1,\dots,n\\}$, it follows that
$\mathop{\rm adj}M$ has constant columns.

Applying the same argument to $M^{T}$ (which also has row sums $0$), we find
that $\mathop{\rm adj}M$ also has constant rows, and hence $\mathop{\rm
adj}M=w{\cal E}_{n}$ for some $w\in\mathbb{C}$.

(b) Let $j\in\\{1,\dots,n\\}$ and let $m_{1},\dots,m_{n}$ be as in (a). Since

$\det(m_{1},\dots,m_{n-1})=(-1)^{n+j}w\neq 0,$ (14)

the vectors $m_{1},\dots,m_{n-1}$ form a basis of $\mathbb{C}^{n-1}$, so

$m_{n}=\sum_{l=1}^{n-1}\alpha_{l}m_{l}$

for suitable $\alpha_{1},\dots,\alpha_{n-1}\in\mathbb{C}$. Now for any
$k\in\\{1,\dots,n-1\\}$ we have

$\displaystyle w$
$\displaystyle=(-1)^{k+j}\det(m_{1},\dots,m_{k-1},m_{k+1},\dots,m_{n-1},m_{n})$
$\displaystyle=(-1)^{k+j+n-1-k}\det(m_{1},\dots,m_{k-1},m_{n},m_{k+1},\dots,m_{n-1}).$

Combining this with eq. (14), we find that
$0=\det(m_{1},\dots,m_{k-1},m_{k}+m_{n},m_{k+1},\dots,m_{n-1})$ . Hence, for
suitable $\beta_{l}\in\mathbb{C}(l\in\\{1,\ldots,n-1\\}\setminus\\{k\\})$ we
have

$\sum_{l\neq
k}^{n-1}\beta_{l}m_{l}=m_{k}+m_{n}=m_{k}+\sum_{l=1}^{n-1}\alpha_{l}m_{l}=(1+\alpha_{k})m_{k}+\sum_{l\neq
k}^{n-1}\alpha_{l}m_{l}.$

As $m_{1},\dots,m_{n-1}$ are linearly independent, it follows that
$1+\alpha_{k}=0$, i.e. that $\alpha_{k}=-1.$

Since this holds for all $k\in\\{1,\dots,n-1\\}$, we can conclude that
${\sum\limits_{l=1}^{n}m_{l}=0}$, i.e. that the rows of $M$, except for the
jth row, add up to $0$. Since $j\in\\{1,\dots,n\\}$ was arbitrary, this holds
in fact for all rows.

Applying the above reasoning to the transpose of $M$, we find that its columns
also add up to $0$. ∎

##  5 Co-Latin Matrices

An $n\times n$ Latin square (or Latin matrix) is an $n\times n$ matrix with
entries from $\\{1,\dots,n\\}$ such that the entries of each row and of each
column are distinct, i.e. each number in $\\{1,\dots,n\\}$ appears exactly
once in each row and in each column. Evidently any Latin matrix is of type S.
If $L=(\ell_{p,q})_{p,q=1}^{n}$ is a Latin square and $k\in\\{1,\dots,n\\}$,
we define $L^{(-1)}(k)={\\{(p,q)\in\\{1,\dots,n\\}^{2}:\ell_{p,q}=k\\}}$.

######  Definition 1.

A matrix $M=(m_{i,j})_{i,j=1}^{n}\in\mathbb{R}^{n\times n}$ has the co-Latin
property if

$\sum_{(p,q)\in L^{(-1)}(k)}m_{p,q}=0$

for all $k\in\\{1,\dots,n\\}$ and all $n\times n$ Latin squares $L$.

Clearly the set of $n\times n$ co-Latin matrices forms a subspace of
$\mathbb{R}^{n\times n}$. The co-Latin property can be considered an extreme
opposite of the weightless type S property. A matrix in $S_{n}$ with weight 0
has the property that its entries in any one row or column add to 0; a co-
Latin matrix has the property that any selection of $n$ entries such that no
two selected entries lie in the same row or the same column add to 0. We show
in the following that these properties are indeed complementary in the sense
of unique decomposability of any given weightless $n\times n$ matrix into a
weightless type S matrix and a co-Latin matrix; this is a direct consequence
of the following theorem, which identifies co-Latin matrices with $V_{n}$ in
(6).

######  Theorem 5.1.

Let $n\in\mathbb{N}$. The space of all $n\times n$ co-Latin matrices is equal
to $V_{n}$.

To prepare the proof of Theorem 5.1, we first show that there exists a Latin
square which has the entries 1 and 2 pairwise on the diagonally opposite
corners of a rectangle. By suitable row and column permutations, we can assume
without loss of generality that this rectangle is the top left $2\times 2$
square. The existence of such a Latin square is not trivial; in fact there is
none in dimension 3, as the arrangement

$\left(\begin{array}[]{ccc}1&2&*\\\ 2&1&*\\\ &*&*\end{array}\right)$

enforces two 3s in the third column and the third row, violating the defining
condition of a Latin square.

######  Lemma 5.1.

Let $n\in{\mathbb{N}}\setminus\\{1,3\\}.$ Then there exists an $n\times n$
Latin square $L=(\ell_{ij})$, such that $\ell_{1\,1}=\ell_{2\,2}=1$ and
$\ell_{1\,2}=\ell_{2\,1}=2.$

###### Proof.

We consider the cases of even and odd $n$ separately.

(i) Even $n=2m$. Consider the Hankel Latin square
$\tilde{L}=(\tilde{\ell}_{i,j})_{i,j=1}^{m}$, where
$\tilde{\ell}_{i,j}=1+((i+j-2)\mathop{\rm mod}m)$ $(i,j\in\\{1,\dots,m\\}),$
e.g. for $m=4$

$\tilde{L}=\left(\begin{array}[]{cccc}1&2&3&4\\\ 2&3&4&1\\\ 3&4&1&2\\\
4&1&2&3\end{array}\right),$

and then replace each $(i,j)$ entry in this matrix with the $2\times 2$ block

$L=\left(\begin{array}[]{cc}2\tilde{\ell}_{ij}-1&2\tilde{\ell}_{ij}\\\
2\tilde{\ell}_{ij}&2\tilde{\ell}_{ij}-1\end{array}\right)$

to create the Latin square $L$.

(ii) Odd $n$. Here we start off the first 3 antidiagonals of $L$ thus:

$\left(\begin{array}[]{cccc}1&2&3&\\\ 2&1&&\\\ 3&&&\end{array}\right),$

and then fill up the antidiagonals up to and including the main antidiagonal
in the Hankel Latin manner described above. We then fill the next three
antidiagonals with

$\begin{array}[]{ccc}3\ 1\ 2\ 1\ 2&\cdots&\ 1\ 2\ 1\ 2\ 3,\\\ 2\ 3\ 3\
3&\cdots&\ 3\ 3\ 3\ 1,\\\ 1\ 2\ 1&\cdots&\ 2\ 1\ 2,\end{array}$

respectively, and complete the remaining antidiagonals in the standard Hankel
Latin manner. For example, for $n=9$ this gives the Latin square

$L=\left(\begin{array}[]{ccccccccc}1&2&3&4&5&6&7&8&9\\\ 2&1&4&5&6&7&8&9&3\\\
3&4&5&6&7&8&9&2&1\\\ 4&5&6&7&8&9&1&3&2\\\ 5&6&7&8&9&3&3&1&4\\\
6&7&8&9&1&3&2&4&5\\\ 7&8&9&2&3&1&4&5&6\\\ 8&9&1&3&2&4&5&6&7\\\
9&3&2&1&4&5&6&7&8\end{array}\right)$

This construction works from $n=5$ onwards. ∎

We can now prove the following.

######  Lemma 5.2.

Let $n\in{\mathbb{N}},$ and let $M=(m_{i,j})_{i,j=1}^{n}$ be an $n\times n$
co-Latin square. Then $M\in V_{n}$.

###### Proof.

The statement is trivial in the case $n=1$.

Now consider $n\in\mathbb{N}\setminus\\{1,3\\}$. Let
$i,j,k,l\in\\{1,\dots,n\\}$ such that $i\neq k$ and $j\neq l$. By suitable
permutation of the rows and columns of the Latin square constructed in Theorem
5.1, there exists a Latin square $L=(\ell_{p,q})_{p,q=1}^{n}$ such that
$\ell_{i,j}=\ell_{k,l}=1$ and $\ell_{i,l}=\ell_{k,j}=2.$ By the co-Latin
property, we know that

$\sum_{(p,q)\in L^{(-1)}(1)}m_{p,q}=0,$

so

$m_{i,j}+m_{k,l}=-\sum_{(p,q)\in
L^{(-1)}(1)\setminus\\{(i,j),(k,l)\\}}m_{p,q}.$ (15)

Now consider the matrix $L^{\prime}=(\ell_{p,q}^{\prime})_{p,q=1}^{n}$ which
arises from $L$ by keeping all the same entries except that
$\ell_{i,j}^{\prime}=\ell_{k,l}^{\prime}=2$ and
$\ell_{i,l}^{\prime}=\ell_{k,j}^{\prime}=1.$ Then $L^{\prime}$ is still a
Latin square, and by the co-Latin square property we find that

$\sum_{(p,q)\in{L^{\prime}}^{(-1)}(1)}m_{p,q}=0,$

so

$m_{i,l}+m_{k,j}=-\sum_{(p,q)\in{L^{\prime}}^{(-1)}(1)\setminus\\{(i,l),(k,j)\\}}m_{p,q}.$
(16)

Noting that the index sets of the sums in $(\ref{eq:lala1})$ and
$(\ref{eq:lala2}),$ and hence the values of these sums, are the same, we
conclude that $M\in V_{n}$.

Finally, to see that the statement holds true in the case $n=3$, we note that,
up to permutations of the symbols $\\{1,2,3\\}$, there are only two different
$3\times 3$ Latin squares, namely either of the form

$\left(\begin{array}[]{ccc}A&B&C\\\ B&C&A\\\
C&A&B\end{array}\right)\qquad\hbox{\rm
or}\qquad\left(\begin{array}[]{ccc}D&E&F\\\ F&D&E\\\ E&F&D\end{array}\right).$

To show that a $3\times 3$ co-Latin square $M$ has the vertex cross sum
property, without loss of generality we can consider the case $i=j=1$, $k=l=2$
(as the other cases can be reduced to this by suitable row and column
permutations). Then the first of the above Latin matrices shows that
$m_{2,1}+m_{1,2}+m_{3,3}=0,$ the second Latin matrix shows that
$m_{1,1}+m_{2,2}+m_{3,3}=0,$ and it follows that
$m_{1,1}+m_{2,2}=m_{2,1}+m_{1,2}$. ∎

We can now complete the proof of Theorem 5.1.

###### Proof of Theorem 5.1.

By Lemma 5.2, any $n\times n$ co-Latin matrix is an element of $V_{n}$.

Conversely, by Theorem 2.2 any element of $V_{n}$ is of the form
$a1_{n}^{T}+1_{n}b^{T}$, with suitable $a,b\in\\{1_{n}\\}^{\bot}$. Since a
Latin square takes each value in $\\{1,\dots,n\\}$ exactly once in each row
and each column, it is hence straightforward to see that $a1_{n}^{T}$ and
$1_{n}b^{T}$ are co-Latin squares. ∎

######  Example: the Wilson matrix.

Using the $V$ part of our integer factorisation matrix $Z_{0}$, given by
$Z_{V}$ in (8), we have

$Z_{V}=\frac{1}{8}\left(\begin{array}[]{cccc}5&7&11&11\\\ -3&-1&3&3\\\
-7&-5&-1&-1\\\ -9&-7&-3&-3\\\ \end{array}\right).$

It can be easily verified numerically that this type V matrix satisfies all
$4!=24$ Latin selections summing to 0, and so is a co-Latin matrix.

##  6 Conclusions

Matrices with integer entries play an important role in many modern
applications of mathematics. In numerical analysis, for example, they make
convenient test matrices because they are exactly representable. Families of
matrices with bounded integer entries have recently been termed Bohemian
matrices and various aspects of them have been studied [2], [3], [4]. The
factorisation of integer matrices in the form $M=N^{T}N$ with integer $N$ is
related to a classical topic in number theory, but there has not been much
work on finding such factorisations. In this work we have provided some new
ideas to determine under what circumstances integer factorisations exist and
to develop an approach to computing them. Our results are founded upon the
orthogonal decomposition of the algebra of square matrices into two parts, of
which one part is the subalgebra of constant (row and column) sum matrices,
while we identified the other part as a space of co-Latin square matrices with
symmetries determined by the properties of Latin squares.

## Acknowledgements

The first author was supported by Engineering and Physical Sciences Research
Council grant EP/P020720/1 and the Royal Society.

## References

  * [1] John H. Conway, _Universal quadratic forms and the fifteen theorem_ , In Quadratic Forms and Their Applications, Contemporary Mathematics 272, Eva Bayer-Fluckiger, David Lewis, and Andrew Ranicki, editors, (1999) 23–26. https://doi.org/10.1.1.449.7217
  * [2] E. Y. S. Chan, R. M. Corless, L. Gonzalez-Vega, J. R. Sendra, J. Sendra, S. E. Thornton, _Upper Hessenberg and Toeplitz Bohemians_ , Linear Algebra Appl. 601 (2020) 72–-100. https://doi.org/10.1016/j.laa.2020.03.037.
  * [3] R. M. Corless, S. E. Thornton, _The Bohemian Eigenvalue Project_ , ACM Comm. Comput. Alg. 50 (2016) 158–-160.
  * [4] N. J. Higham, _Rhapsodizing about Bohemian matrices_ , SIAM News 51 (10) (2018) 2.
  * [5] N. J. Higham _What is the adjugate of a matrix?_ , https://nhigham.com/2020/06/16/what-is-the-adjugate-of-a-matrix/
  * [6] N. J. Higham, M. C. Lettington, _Optimizing and Factorizing the Wilson Matrix_ , (to appear in the Amer. Math. Monthly).
  * [7] S. L. Hill, M. C. Lettington, and K. M. Schmidt, _Block representation and spectral properties of constant sum matrices_ , Electronic Journal of Linear Algebra, 34 (2018) 1170–190. https://doi.org/10.13001/1081-3810.3530.
  * [8] S. L. Hill, M. C. Lettington, and K. M. Schmidt, _On superalgebras of matrices with symmetry properties_ , Linear and Multilinear Algebra, 66 (2018) 1538–1563. https://doi.org/10.1080/03081087.2017.1363153
  * [9] C. Ko, _Determination of the class number of positive quadratic forms in nine variables with determinant unity_ , J. London Math. Soc. 13 (1938) 102–-110. https://doi.org/10.1112/jlms/s1-13.2.102
  * [10] C. B. Moler, _Reviving Wilson’s matrix_ , https://blogs.mathworks.com/cleve/2018/08/20/reviving-wilsons-matrix/, Aug. 2018.
  * [11] L. J. Mordell, _On the representation of a binary quadratic form as a sum of squares of linear forms_ , Math. Zeitschrift 35 (1932) 1–-15. https://doi.org/10.1007/BF01186544.
  * [12] L. J. Mordell, _The definite quadratic forms in eight variables with determinant unity_ , J. Math. Pures Appl. 17 (1938) 41–-46.
  * [13] J. Morris, An escalator process for the solution of linear simultaneous equations, The London, Edinburgh, and Dublin Philosophical Magazine and Journal of Science, 37 (265) (1946) 106–120.
  * [14] O. Taussky, _Sums of squares_ , Amer. Math. Monthly 77 (1970) 805–-830. https://doi.org/10.2307/2317016.