Theorem: Inverse of Cauchy Matrix

Let $C_n$ be the square Cauchy matrix of order $n$:
:<nowiki>$C_n = \begin{bmatrix}
\dfrac 1 {x_1 + y_1} & \dfrac 1 {x_1 + y_2} & \cdots & \dfrac 1 {x_1 + y_n} \\
\dfrac 1 {x_2 + y_1} & \dfrac 1 {x_2 + y_2} & \cdots & \dfrac 1 {x_2 + y_n} \\
\vdots & \vdots & \ddots & \vdots \\
\dfrac 1 {x_n + y_1} & \dfrac 1 {x_n + y_2} & \cdots & \dfrac 1 {x_n + y_n} \\
\end{bmatrix}$</nowiki>
Then its inverse $C_n^{-1} = \sqbrk b_n$ can be specified as:
:$\begin{bmatrix} b_{ij} \end{bmatrix} = \begin{bmatrix} \dfrac {\ds \prod_{k \mathop = 1}^n \paren {x_j + y_k} \paren {x_k + y_i} } {\ds \paren {x_j + y_i} \paren {\prod_{\substack {1 \mathop \le k \mathop \le n \\ k \mathop \ne j} } \paren {x_j - x_k} } \paren {\prod_{\substack {1 \mathop \le k \mathop \le n \\ k \mathop \ne i} } \paren {y_i - y_k} } } \end{bmatrix}$

Proof:
'''Preliminaries''':
Vandermonde Matrix Identity for Cauchy Matrix supplies matrix equation
:$(1): \quad -C = PV_x^{-1} V_y Q^{-1}$
:Definitions of symbols:
::<nowiki>$V_x = \paren {\begin{smallmatrix}
1           & 1           & \cdots & 1 \\
x_1         & x_2         & \cdots & x_n \\
\vdots      & \vdots      & \ddots & \vdots \\
x_1^{n - 1} & x_2^{n - 1} & \cdots & x_n^{n - 1} \\
\end{smallmatrix} },\quad 
V_y = \paren {\begin{smallmatrix}
1         & 1         & \cdots & 1 \\
y_1       & y_2       & \cdots & y_n \\
\vdots    & \vdots    & \ddots & \vdots \\
y_1^{n - 1} & y_2^{n - 1} & \cdots & y_n^{n - 1} \\
\end{smallmatrix} }$</nowiki> {{Defof|Vandermonde Matrix|subdef = Formulation 1}}
::<nowiki>$P = \paren {\begin{smallmatrix}
p_1(x_1) &  \cdots & 0 \\
\vdots   & \ddots  & \vdots \\
0        & \cdots  & p_n(x_n) \\
\end{smallmatrix} }, \quad
Q = \paren {\begin{smallmatrix}
p(y_1)  & \cdots  & 0 \\
\vdots  & \ddots  & \vdots \\
0       & \cdots  & p(y_n) \\
\end{smallmatrix} }$</nowiki> {{Defof|Diagonal Matrix}}
::<nowiki>$\ds \map p x = \prod_{i \mathop = 1}^n \paren {x - x_i},
\quad \map {p_k} x = \prod_{i \mathop = 1, i \mathop \ne k}^n \paren {x - x_i},
\quad
1 \mathop \le k \mathop \le n$</nowiki> Polynomials
'''Compute the inverse''' $C^{-1}$ for:
{{begin-eqn}}
{{eqn | l = C
      | r = <nowiki>\begin{bmatrix}
            \dfrac 1 {x_1 - y_1} & \dfrac 1 {x_1 - y_2} & \cdots & \dfrac 1 {x_1 - y_n} \\
            \dfrac 1 {x_2 - y_1} & \dfrac 1 {x_2 - y_2} & \cdots & \dfrac 1 {x_2 - y_n} \\
            \vdots               & \vdots               & \ddots & \vdots \\
            \dfrac 1 {x_n - y_1} & \dfrac 1 {x_n - y_2} & \cdots & \dfrac 1 {x_n - y_n} \\
            \end{bmatrix}</nowiki>
      | c = Assume $\set {x_1,\ldots,x_n,y_1,\ldots,y_n}$ has distinct elements.
}}
{{end-eqn}}
Replacement $y_k \to -y_k$ then gives the inverse $C_n^{-1}$ in the theorem.
Inverse of Matrix Product applied to equation (1) gives: 
:$ C^{-1} = -Q V_y^{-1} V_x P^{-1}$ 
Let ${\vec K}_1, \ldots, {\vec K}_n$ denote the columns of the $n\times n$ identity matrix.
Then $n \times n$ matrix $B = C^{-1}$ has entries $b_{ij} = {\vec K}_i^T C^{-1} {\vec K}_j$.
Hold fixed until the end of the proof the row and column index symbols $i$ and $j$.
Define column vectors $\vec A$, $\vec B$ so that $b_{ij} = {\vec A}^T \vec B$:
:$\vec A = \paren { {\vec K}_i^T \, Q \, V_y^{-1} }^T$
:$\vec B = -V_x \, P^{-1} \, {\vec K}_j$
Define $u = x_j$ and simplify:
{{begin-eqn}}
{{eqn | l = \vec A 
      | r = {\map p {y_i} } \, \paren { V_y^{-1} }^T \, {\vec K}_i 
      | c = Transpose of Matrix Product
}}
{{eqn | r = \dfrac{ \map p {y_i} } {\map \det {V_y} } \paren { \adj {V_y} }^T {\vec K}_i
      | c = Matrix Product with Adjugate Matrix
}}
{{eqn | r = <nowiki>\dfrac{ \map p {y_i} } {\map \det  {V_y} }
        \begin{bmatrix} 
          {\mathbf {Cofactor} } \paren {V_y,1,i } \\ 
          \vdots \\ 
          {\mathbf {Cofactor} } \paren {V_y,n,i } \\ 
        \end{bmatrix}</nowiki>
      | c = ${\mathbf {Cofactor} } \paren {M,r,s}$ denotes cofactor $M_{rs}$
}}
{{eqn | l = \vec B
      | r = -\dfrac 1 {\map {p_j} {x_j} } \, V_x \, {\vec K}_j
}}
{{eqn | r = <nowiki>-\dfrac 1 {\map {p_j} {x_j} } \, 
        \begin{bmatrix} 
          1 \\ 
          u \\ 
          \vdots \\ 
          u^{n - 1} \\ 
        \end{bmatrix}</nowiki>
      | c = Symbol $u = x_j$.
}}
{{end-eqn}}
Then:
{{begin-eqn}}
{{eqn | l = A^T \, B
      | r = \dfrac {-\map p {y_i} } {\map {p_j} {x_j} \map \det {V_y} } \sum_{k \mathop = 1}^n {\mathbf {Cofactor} } \paren {V_y, k,i } \, u^{k - 1}
      | c = Laplace Expansion Theorem for Determinants along column $i$.
}}
{{eqn | r = <nowiki>\dfrac {-\map p {y_i} } {\map {p_j} {x_j} \map \det {V_y} } \map \det {\begin {matrix} 
        1           & \cdots & 1                 & 1         & 1                 & \cdots & 1 \\
        y_1         & \cdots & y_{i - 1}         & u         & y_{i + 1}         & \cdots & y_n \\
        \vdots      & \cdots & \vdots            &\vdots     & \vdots            & \cdots & \vdots \\
        y_1^{n - 1} & \cdots & y_{i - 1}^{n - 1} & u^{n - 1} & y_{i + 1}^{n - 1} & \cdots & y_n^{n - 1} \\
        \end {matrix} }</nowiki>
}}
{{eqn | r = \dfrac {-\map p {y_i} } {\map {p_j} {x_j} } \ \dfrac {\bigvalueat {\paren {\map \det {V_y} } } {y_i \mathop = u} } {\map \det {V_y} }
}}
{{end-eqn}}
Simplify the fraction on the right:
{{begin-eqn}}
{{eqn | n = 2
      | l = \dfrac {\bigvalueat {\paren {\map \det {V_y} } } {y_i \mathop = u} } {\map \det {V_y} }
      | r = \dfrac {\ds \valueat {\paren {\prod_{1 \mathop \le m \mathop < k \mathop \le n} \paren {y_k - y_m} } } {y_i \mathop = u} } {\ds \prod_{1 \mathop \le m \mathop \lt k \mathop \le n} \paren {y_k - y_m} }
      | c = Value of Vandermonde Determinant
}}
{{end-eqn}}
Define sets $D,A,B,C$ with disjoint decomposition $D = A \cup B \cup C$:
{{begin-eqn}}
{{eqn | l = D 
      | r = \set {\tuple {m, k} : 1 \le m < k \le n }
}}
{{eqn | l = A 
      | r = \set {\tuple {m, k} \in D : m \ne i \mbox{ and } k \ne i}
}}
{{eqn | l = B
      | r = \set {\tuple {m, k} \in D: m = i}
}}
{{eqn | l = C
      | r = \set {\tuple {m, k} \in D: k = i}
}}
{{end-eqn}}
Use $\prod_D = \prod_A \prod_B \prod_C$ to convert the numerator and denominator in the {{RHS}} of $(2)$.
Then:
{{begin-eqn}}
{{eqn | l = \dfrac {\bigvalueat {\paren {\map \det {V_y} } } {y_i \mathop = u} } {\map \det {V_y} }
      | r = \dfrac {\ds \prod_{k \mathop = 1, \, k \mathop \ne i}^n \paren {u - y_k} } {\ds \prod_{k \mathop = 1, \, k \mathop \ne i}^n \paren {y_i - y_k} }
      | c = cancelling common factors in $(2)$
}}
{{end-eqn}}
Replacement of $\map p {y_i}$ and $\map {p_j} {x_j}$ by products gives:
{{begin-eqn}}
{{eqn | n = 3
      | l = b_{ij} 
      | r = \paren {-1} \, \dfrac {\ds \prod_{k \mathop = 1}^n \paren {y_i - x_k } } {\ds \prod_{k \mathop = 1,\, k \mathop \ne j}^n \paren { x_j - x_k } } \ \prod_{k \mathop = 1,\, k \mathop \ne i}^n \paren {\dfrac {x_j - y_k} {y_i - y_k} }
      | c = substituting $u = x_j$
}}
{{end-eqn}}
After substituting $y_k \to -y_k$ and simplifying, ''Knuth's original identity'' (1997) becomes:
{{explain|ProofWiki is self-contained so we don't refer to outside sources to justify working}}
{{begin-eqn}}
{{eqn | n = 4
      | l = b_{ij}
      | r = \dfrac {\ds \prod_{k \mathop = 1}^n \paren {x_k - y_i} } {\ds \prod_{k \mathop = 1,\, k \mathop \ne j}^n \paren {x_j - x_k} } \dfrac {\ds \prod_{k \mathop = 1, \, k \mathop \ne i}^n \paren {x_j - y_k} } {\ds \prod_{k \mathop = 1,\, k \mathop \ne i}^n \paren {-y_i + y_k} }
}}
{{end-eqn}}
In $(3)$, factor $\paren {-1}^{n + 1}$ from the numerator and $\paren {-1}^{n - 1}$ from the denominator.
Then $\paren {-1}^{n + 1} = \paren {-1}^{n - 1}$ verifies that $(3)$ matches $(4)$.
{{qed}}