\magnification=\magstep1

%\baselineskip=10.5mm

\font\kis=cmr7
\font\NAGY=cmbx10 scaled \magstep3

\long\def\DEL#1{{}}

\def\gg{{\overline{g}}}
\def\pp{{\overline{p}}}

\def\T{\bigskip\noindent{\bf T\'etel.\ \rm}}
\def\L{\bigskip\noindent\bf Lemma.\ \rm}
\def\P{\bigskip\noindent\bf Propoz\'\i ci\'o.\ \rm}
\def\K{\bigskip\noindent\bf  K\"ovetkezm\'eny.\ \rm}
\def\M{\bigskip\noindent\bf Megjegyz\'es.\ \rm}
\def\pl{\bigskip\noindent\bf P\'elda.\ \rm}
\def\b{\bigskip\noindent \bf Bizony\'{\i}t\'as.\ \rm}
\def\A{{\cal A}}
\def\e{\varepsilon}
\def\eps{\varepsilon}

\def\s{{\bigskip\noindent} }

\def\RR{{\rm I\kern-1.6pt {\rm R}}}
\def\R{{\rm I\kern-1.6pt {\rm R}}}
\def\KK{{\rm I\kern-1.6pt {\rm K}}}
\def\N{{\rm I\kern-1.6pt {\rm N}}}
\def\NN{{\rm I\kern-1.6pt {\rm N}}}
%\def\P{{\rm I\kern-1.6pt {\rm P}}}
\def\P{{\bf P}}
\def\FF{{\rm I\kern-1.6pt {\rm F}}}
\def\TT{{\rm T}\kern-4.8pt {\rm T} \kern2pt}
\def\Z{{\rm Z}\kern-3.8pt {\rm Z} \kern2pt}
\def\ZZ{{\rm Z}\kern-3.8pt {\rm Z} \kern2pt}
\def\Q{{\rm Q}\kern-3.9pt {\rm \vrule height7.1pt
        width.3pt depth-2pt} \kern7.5pt}
\def\CC{{\rm \kern 3.5pt \vrule height 6.5pt
        width.7pt depth0pt \kern-4pt C \kern.8pt}}

\def\subsetneq{\buildrel{\subset}\over{\ne}}


\def\sqr#1#2{{\,\vcenter{\vbox{\hrule height.#2pt\hbox{\vrule width.#2pt
height#1pt \kern#1pt\vrule width.#2pt}\hrule height.#2pt}}\,}}
\def\square{\mathchoice\sqr34\sqr34\sqr{2.1}3\sqr{1.5}3}

\def\bizveg{\hfill\ifmmode\sqr66\else$\sqr66$\par\fi}

\def\TR{{\rm T}}





%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%




\bf 1. Introduction

\medskip
\rm
Recently some novel interest seems to be raised toward the
symbolic LU decomposition of Vandermonde matrices.
Several explicit formulas at various levels along with matrix subfactorizations are
well-known for them and their inverses.
Our aim in this paper to extend these results to the matrices
associated with interpolation problems with Hermite type.



\bigskip

\bf 2. Preliminaries

\medskip
\noindent \bf 2.1. \
\rm
First we recall briefly the results
spread in the literature concerning the LU decomposition of
the {\it Vandermonde} matrix
$$V=V({\bf x}) := \left[ \matrix{
1 & x_0 &\ldots & x_0^N \cr
1 & x_1 &\ldots & x_1^N \cr
\vdots&\vdots&\ddots&\vdots\cr
1 & x_N &\ldots & x_N^N \cr}
\right] , \quad
\eqalign{&{\bf x}:= \big[ x_0,\ldots,x_N \big] \cr
& x_i \ne x_j \ \ {\rm for} \ i\ne j \ .\cr}$$
Henceforth $N$ is a fixed positive integer,
and we shall consider $(N+1)^2$-matrices
with indices ranging from $0$ to $N$.
In terms of the formal vectors
$$\delta ({\bf x}) := \left[
\matrix{ \delta_{x_0} \cr \delta_{x_1} \cr \vdots \cr \delta_{x_N} \cr}
\right] , \quad
{\bf e}(x) := \big[ 1, x, x^2, \ldots, x^N \big] $$
where $x$ is a variable symbol and $\delta_a$ denotes the evaluation functional
$p\mapsto p(a)$
defined for polynomials in $x$ of degree$\le N$,
we can write
$$V ({\bf x}) = \delta({\bf x}) {\bf e} (x) .$$
Consider  the Lagrange interpolation polynomial
$p:=p_{{\bf x},{\bf y}}$
defined by the requirements
$p(x_n)=y_n$ $(n=0,\ldots,N)$ (with ${\bf y}:=[y_0,\ldots,y_n]^\TR)$.
Since we can write $p$ in the form $p= \sum_n p_n x^N ={\bf e}(x) {\bf p}$,
it follows
${\bf y} = \delta({\bf x}) p_{{\bf x},{\bf y}}(x) =
\delta({\bf x}) {\bf e}(x) {\bf p}_{{\bf x},{\bf y}} =
V({\bf x}) {\bf p}_{{\bf x},{\bf y}}$ that is
$$p_{{\bf x},{\bf y}} (x) = {\bf e}(x) V({\bf x})^{-1} {\bf y} .$$
The Newtonian form
$p(x)=
\sum\limits_{n=0}^N p(x_0,\ldots,x_n) \omega_n(x)$ of this polynomial with
$$\omega_n(x):= \prod_{k: k<n} (x - x_k), \qquad
p(x_0,\ldots,x_n) = \sum_{j=0}^n \, p(x_j) \!\!\!\!
\prod_{i: j\ne i\le n} \!\! (x_j - x_i)^{-1}$$
(convention: $\omega_0 \equiv 1$) yields the relation
$p_{{\bf x},{\bf y}} (x) = {\bf e}(x) \Omega({\bf x})\Delta({\bf x}) {\bf y}$
with the upper resp. lower triangular matrices
$$\Omega({\bf x}):= \big[ \hbox{coeffs of $\omega_n$ in column $n$} \big]_{n=0}^N \ ,
\qquad
\Delta({\bf x}) := \Big[ \prod_{i: j\ne i\le n} \!\! (x_j - x_i)^{-1}\Big]_{{{n,j=0}\atop n\ge j}}^N .$$
Hence the following triangular decompositions are immediate
$$V({\bf x})^{-1} = \Omega({\bf x}) \Delta({\bf x}), \quad
V({\bf x}) = \Delta({\bf x})^{-1} \Omega({\bf x})^{-1}$$
along with the explicit closed formula
$$\Delta({\bf x})^{-1} = V({\bf x}) \Omega({\bf x}) =
\delta({\bf x}) {\bf e}(x) \Omega({\bf x}) =\big[ \omega_n(x_k) \big]_{k,n=0}^N .$$
Not stated explicitly in [Monthly] but the arguments can be continued to achieve a shortcut
to a closed (i.e. recursion free) formula for the entries of the upper triangular term
$\Omega({\bf x})^{-1}$. Indeed, we have
$$\eqalign{
\Omega({\bf x})^{-1} &= \Delta({\bf x}) V({\bf x}) = \Delta({\bf x}) \delta({\bf x}) {\bf e}(x) =
\Big[ \big[ \sum_{j=0}^n
\prod_{i: j\ne i\le n} \!\! (x_j - x_i)^{-1} \big] \delta_{x_j} \Big]_{n=0}^N {\bf e}(x) =\cr
&= \Big[ \big[ \sum_{j=0}^n
\prod_{i: j\ne i\le n} \!\! (x_j - x_i)^{-1} \big] \delta_{x_j}(x^\nu) \Big]_{{{n,\nu=0}\atop n\le \nu}}^N =
\Big[ {{\hbox{Newton difference}}\atop{\hbox{$x^\nu(x_0,\ldots,x_n)$}}} \Big]_{{{n,\nu=0}\atop n\le \nu}}^N  =\cr
&= \Big[ \sum_{{{i_0+\cdots+i_n=\nu-n}\atop{i_0,i_1,\ldots,i_n \ge 0}}}
x_0^{i_0}x_1^{i_1}\cdots x_n^{i_n} \Big]_{{{n,\nu=0}\atop n\le \nu}}^N  .\cr}$$
For any fixed degree $\nu$, the last formula can be obtained by induction on $n$
%in a straightforward manner
from the identities
$p(x_k,\ldots,x_{k+s+1})= \big[ p(x_k+1,\ldots,x_{k+s+1})-p(x_k,\ldots,x_{k+s})\big]/(x_{k+s+1}-x_k)$.

%As far as we know no closed formulas were published for the
%coefficients of the inverse matrix $\Omega^{-1}$.
%It is not hard to fill this gap with a direct approach.
%Actually, in the spirit of [Monthly], we can continue the argument as follows:

%\bigskip
%\noindent
%\bf Proposition. \it
%$V_n$ LU-felbont\'asa
%$$V_n = L_n U_n \ ,$$
%ahol $L_n$ az az als\'o-triangul\'aris $n\times n$-es m\'atrix, amelynek
%f\H o\'atl\'oja $1$-eskb\H ol \'all, az $(i,j)$-edik tagja $(i> j)$ az
%$ x_1,\ldots,x_j$ faktorokkal k\'epezhet\H o
%\"osszes $(i-1)$-edfok\'u szorzatok \"osszege,
%$U_n$ pedig azaz a fels\H o-triangul\'aris m\'atrix, ahol
%az els\H o sor csupa $1$-esekb\H ol \'all, \'es $1<i\le j$ eset\'en az
%$(i,j)$-edig tag az $(x_j-x_1)\cdots(x_j-x_{i-1})$ szorzat.
%\rm

\bigskip
\noindent
\bf 2.2. \ \rm
Next we recall the concept of Hermite (or Hermite-Vandermonde) matrices
along with their relationship to Hermite approximation.
Henceforth we fix numbers $m_1,\ldots,m_r$ such that
$(m_1+1)+\ldots+(m_r+1)=N+1$.
Given an $r$-tuple ${\bf a}:=[a_1,\ldots,a_r]$,
along with
row matrices ${\bf b}_k := \big[ b_k^{(0)}, \ldots, b_k^{(m_k)} \big]$ $(k=1,\ldots,r)$,
the Hermite interpolation polynomial
$q(x):= q_{{\bf a},{\bf b}_1,\ldots,{\bf b}_r}(x)$
is defined as the unique polynomial of degree$\le N$ satisfying
$${d^s q \over d x^s} \Big\vert_{x=a_k} = b_k^{(s)} \qquad (k=1,\ldots,r;\ \ 0\le s\le m_k) .$$
For later matrix operations, we divide the integer interval $I:=\{0,1,\ldots,N\}$
into consecutive segments
$$I_k:=\{ \nu_k^{(0)},\ldots,\nu_k^{(m_k)} \} , \quad \nu_k^{(s)} := s+\sum\limits_{\ell:\ell<k} (m_\ell+1)$$
with inverse indices
$\kappa(n) :=\big[ k \! : n\!\in\! I_k \big]$,
$\sigma(n):=\big[$position of $n$ in $I_{\kappa(n)} \big] =\big[ s \!: n \!=\!\nu_{\kappa(n)}^{(s)} \big]$.
The {\it Hermite-Vandermonde matrix} over the base point system
${\bf a}=[a_1,\ldots,a_r]$ of multiorder ${\bf m}:= [m_1,\ldots,m_r]$
is the $(N+1)^2$-matrix $H=H^{{\bf m}}({\bf a})$
of the system of linear equations of the form
${\bf q} H^{{\bf m}} ({\bf a})= [{\bf b}_1, \ldots, {\bf b}_r]$
for the coefficient vector ${\bf q}:= [q_0,\ldots,q_N]$
of the polynomial $q_{{\bf a},{\bf b}_1,\ldots,{\bf b}_r}$.
In terms of the linear functionals
$$\delta_a^{(s)} : p\mapsto {d^s p\over d x^s }\Big\vert_{x=a} , \qquad
\delta^{{\bf m}} ({\bf a}) := \big[
\delta^{(0)}_{a_1},\ldots,\delta^{(m_1)}_{a_1},\ldots,
\delta^{(0)}_{a_r},\ldots,\delta^{(m_r)}_{a_r} \big]^\TR ,$$
analogously as in the Vandermonde case ($r=N+1,m_1=\cdots=m_{N+1}=0$), we can write
$$H^{{\bf m}}({\bf a}) = \delta^{{\bf m}} ({\bf a}) {\bf e}(x), \qquad
q_{{\bf a},{\bf b}_1,\ldots,{\bf b}_r}(x) = {\bf e}(x) H^{{\bf m}}({\bf a})^{-1} {\bf b} \ \
\hbox{with} \ \ {\bf b}:=[{\bf b}_1,\ldots,{\bf b}_r]^\TR .$$
It is well known that Hermite interpolation polynomials admit also a Newtonian form
%$$q(x)= \sum_{k=1}^r \sum_{s=0}^{m_k}
%q(a_{\nu_0^{(0)}},\ldots,a_{\nu_k^{(s)}})\omega_{k,j}(x) , \quad
%\widetilde{\omega}_{k,j}(x):= \prod_{i: i<k} (x-a_i)^{m_i} (x-a_k)^j$$
$$q(x)= \sum_{n=0}^N
q( a_{ \nu_0^{(0)} } ,\ldots, a_{ \nu_k^{(s)} }) \omega_n(x) , \qquad
\omega_n(x):= \prod_{i: i<n} (x-a_{ \kappa(i) } )$$
in terms of generalized Newton differences.
\DEL{
By setting $\widetilde{{\bf x}} := {\bf a}^{{\bf m}}$ (that is
$\widetilde{x}_0=a_1,\ldots,x_{m_1}=a_1,a_{m_1+1}=a_2,\ldots,x_N=a_r)$
we let recursively
$$\eqalign{
&F(\widetilde{x}_n)= \big[ b_k^{(0)}: \ \widetilde{x}_n = a_k \big] \qquad (n=0,\ldots,N), \cr
&F(\widetilde{x}_n,\widetilde{x}_{m+1},\ldots,\widetilde{x}_{n+s}):=
\Big[ {F(\widetilde{x}_{n+1},\ldots,\widetilde{x}_{n+s}) -
F(\widetilde{x}_n,\ldots,\widetilde{x}_{n+s-1})\over \widetilde{x}_{n+s}-\widetilde{x}_n} \ \
\hbox{if} \ \widetilde{x}_n\ne \widetilde{x}_{n+s},\cr
&\hskip44mm {b_k^{(s)}\over s!} \quad \hbox{if $x_n=x_{n+s}=a_k$}\Big];\cr
&\widetilde{A}_{k,j} := F\big( \widetilde{x}_0,\ldots,\widetilde{x}_{(m_1+1)+\cdots+(m_k+1)+j} \big) .\cr}$$
}%ENDDEL
As outlined in [Monthly] hence we can get again a triangular decomposition
of the form
$$H^{{\bf m}}({\bf a})^{-1} =
\Omega^{{\bf m}}({\bf a}) \Delta^{{\bf m}}({\bf a})$$
where
$\Delta^{{\bf m}}({\bf a})$
is a lower triangular matrix whose row with index $n$ contains the coefficients of
$q_{{\bf a},{\bf b}_1,\ldots,{\bf b}_r} ( a_{ \nu_0^{(0)} } ,\ldots, a_{ \nu_k^{(s)} })$
with respect to the variables $b_k^{s)}$ with $\nu(k,s) \le n$,
while column $n$ of $\Omega^{{\bf m}}({\bf a})$
consists of the coefficients of the polynomial
$\omega_n(x)$.
Similarly as in the Vandermonde case, we can conclude that
$$\Delta^{{\bf m}}({\bf a})^{-1} = \Omega^{{\bf m}}({\bf a}) H^{{\bf m}}({\bf a}) )=
\delta^{{\bf m}}({\bf a}) {\bf e}(x) \Omega^{{\bf m}}({\bf a}) =
\Big[ {d^{\sigma(i)} \over dx^{\sigma(i)}} \Big\vert_{x=a_{\kappa(i)}} \omega_n(x) \Big]_{i,n=0}^N $$
in the standard LU decomposition
$H^{{\bf m}}({\bf a})= \Delta^{{\bf m}}({\bf a})^{-1} \Omega^{{\bf m}}({\bf a})^{-1}$.

\bigskip
\noindent
\bf 3. Combinatorial formulas of the upper triangular factors
\rm
\medskip
Closed formulas for the entries of the matrices
$\Omega^{{\bf m}}({\bf a}), \Omega^{{\bf m}}({\bf a})^{-1}$
can be obtained simply by a formal substitution of the
tuple ${\bf x}=(x_0,\ldots,x_N)$ (supposed to have pairwise different entries
in Subsection 2.1 with
${\bf x}:=
{\bf a}^{{\bf m}} =
\big(\underbrace{a_0,\ldots,a_0}_{m_0+1}, \ldots,\underbrace{a_r,\ldots,a_r}_{m_r+1} \big)$.
Namely the diagonal entries are $1$ in both cases and, for for the entries of indices $i<j$ we have
$$\eqalign{
\big[ \Omega^{{\bf m}}({\bf a}) \big]_{ij} &=
\big[ \hbox{coeff. of $x^i$ in $\prod\limits_{k=0}^{j-1} (x-x_k)$}\big]\Big\vert_{{\bf x}={\bf a}^{{\bf m}}} =\cr
&= (-1)^{j-i} \sum_{0\le \ell_1<\ell_2<\cdots<\ell_{j-i} <j}
\!\!\!\!  x_{\ell_1} x_{\ell_2}\cdots x_{\ell_{j-i}} \Big\vert_{{\bf x}={\bf a}^{{\bf m}}} =\cr
&= (-1)^{j-i} \sum_{(k_0,\ldots,k_r) \in K_{ij}} \varrho_{(k_0,\ldots,k_r)}^{ij} \ a_0^{k_0} \cdots a_r^{k_r} \cr}$$
with the index sets
$$K_{ij}:= \Big\{ (k_0,\ldots,k_r) \in
\Big[ \mathop{\times}_{\alpha<\kappa(j-1)} [0,m_\alpha] \Big]
\times \big[ 0,\sigma(\alpha) \big] \times \big\{ 0 \big\}^{r-\kappa(j)} : \
k_0+\cdots+k_r=j-i \Big\}$$
and respective weight coefficients \ \
$\varrho_{(k_0,\ldots,k_r)}^{ij} :=
\Big[ \prod\limits_{\alpha<\kappa(j-1)} \left( {m_\alpha \atop k_\alpha}\right) \Big]
\left( {\sigma(j) \atop k_{\kappa(j-1)}}\right)$.

Similarly, for $n <\nu$ we have
$$\big[ \Omega^{{\bf m}}({\bf a})^{-1} \big]_{n\nu} =
\sum_{{{i_0+\cdots+i_n=\nu-1}\atop{i_0,i_1,\ldots,i_n \ge 0}}}
x_0^{i_0}x_1^{i_1}\cdots x_n^{i_n} \Big\vert_{{\bf x}={\bf a}^{{\bf m}}} =
\sum_{(k_0,\ldots,k_r) \in \widetilde{K}_{n\nu}}
\widetilde{\varrho}_{(k_0,\ldots,k_r)}^{n\nu} \ a_0^{k_0} \cdots a_r^{k_r} $$
with
$\widetilde{K}_{n\nu}:= \Big\{ (k_0,\ldots,k_r) \!\in\!
\Big[ \! \mathop{\times}\limits_{\alpha<\kappa(j-1)} \!\! [0,m_\alpha] \Big]
\!\times\! \big[ 0,\sigma(\alpha) \big] \!\times\! \big\{ 0 \big\}^{r-\kappa(j)} \!:
k_0 \!+\!\cdots\!+\! k_r \!=\! \nu \!-\! 1 \Big\}$,
$$\eqalign{\widetilde{\varrho}_{(k_0,\ldots,k_r)}^{n\nu}  &=
\hbox{$\# \Big\{$functions $\phi:\{0,\ldots,n\}\to\ZZ_0$ with $\sum\limits_{i\in I_s} \phi(i)=k_s$
$(s=0,\ldots,\kappa(n)) \Big\} =$}\cr
&= \Big[ \prod\limits_{s=0}^{\kappa(n-1)} \mu(m_s+1,k_s) \Big] \mu\big( \sigma(n)+1,k_{\kappa(n)} \big) \cr}$$
where $\mu(\ell,k)$ denotes the number of all functions
$\psi:\{1,\ldots,\ell\}\to\ZZ_+$ with $\sum_i \psi(i)=k$.


\bigskip
\noindent
\bf 4.  Combinatorial formulas of the lower triangular factors
\rm
\medskip
\noindent
{\bf 4.1.} We know already that
$\big[ \Delta^{{\bf m}}({\bf a})^{-1} \big]_{i,n}=
{d^{\sigma(i)} \over dx^{\sigma(i)}} \Big\vert_{x=a_{\kappa(i)}} \prod_{j<n}[x - a_{\kappa(j)}]$.
Their combinatorial forms are immediate from Leibniz's rule for differentiations
of product functions:
$$\eqalign{
\big[ \Delta^{{\bf m}}({\bf a})^{-1} \big]_{i,n} &=
\sum_{0\le j_1 < \cdots < j_{n-\sigma(i)}<n}
[a_{\kappa(i)}-a_{\kappa(j_1)}] \cdots [ a_{\kappa(i)}-a_{\kappa(n-\sigma(i))}] =\cr
&= \sum_{{{\ell_1+\dots+\ell_r = n-\sigma(i)}\atop{\ell_1,\cdots,\ell_r \ge 0 = \ell_j \ (j>\kappa(n-1)) }}}
\Big[ \prod_{s=1}^{\kappa(n-1)-1} \Big( {m_s+1 \atop \ell_s}\Big) \Big]
\Big( {\min\{ m_{\kappa(n-1)},\sigma(i)\} \atop \ell_{\kappa(n-1)}}\Big) \times\cr
&\hskip35mm \times [a_{\kappa(i)}-a_1]^{\ell_1} \cdots [a_{\kappa(i)}-a_{\kappa(n-1)}]^{\ell_{\kappa(n-1)}} .\cr}$$

\medskip
\noindent
{\bf 4.2.} In accordance with the partition $I=\bigcup\limits_{k=1}^r I_k$,
we divide the vectors ${\bf z} := [z_0,\ldots,z_N]$ into subvectors
${\bf z}_k :=\big[ z_{\nu(k,0)}, z_{\nu(k,1)},\ldots, z_{\nu(k,m_k)} \big]$
$(k=1,\ldots,r)$,
consider the corresponding Newton difference matrices
$$\Delta({\bf x}_k) := \Big[
\prod_{i: j\ne i\le n} \!\! (x_{\nu(k,j)} - x_{\nu(k,i)})^{-1}
\Big]_{{{n,j=0}\atop n\ge j}}^{m_k} .$$
It is well-known from classical analysis [??] that
differentiations can be obtained as limits of Newton differences:
given any smooth function $\phi:\RR\to\RR$ along with a net
${\bf x}^{(\alpha)} \to {\bf a}^{{\bf m}}$ with pairwise different terms
$x_n^{(\alpha)} \ne x_\nu^{(\alpha)}$ $(0\le n<\nu\le N)$,
we have
$$\Delta_k \big( {\bf x}_k^{(\alpha)}\big)
\Big[ \phi\big( x_{\nu(k,0)}^{(\alpha)} \big), \ldots, \phi\big( x_{\nu(k,m_k)}^{(\alpha)}\big) \Big]^\TR
\longrightarrow
\Big[ \phi(a_k), \phi^\prime(a_k),\ldots, {d^{m_k}\over d x^{m_k}} \Big\vert_{a=a_k} \phi(x) \Big]^\TR .$$
For $n=0,\ldots,N$, applying this observation with $\phi_n(x):=x^n$
to column $n$ of the Vandermonde matrices $V({\bf x}^{(\alpha)}$, we get
${\rm diag}\big(\Delta_1({\bf x}_1^{(\alpha)}), \ldots, \Delta_1({\bf x}_1^{(\alpha)}\big)
V\big ({\bf x}^{(\alpha)}\big) \to
H^{{\bf m}} ({\bf a})$ and hence
$$V\big ({\bf x}^{(\alpha)}\big)
{\rm diag}\big(\Delta({\bf x}_1^{(\alpha)})^{-1}, \ldots, \Delta({\bf x}_r^{(\alpha)})^{-1}\big)
\to H^{{\bf m}} ({\bf a})^{-1} .$$
Taking into account that $\Delta({\bf x}^{(\alpha)}$ along with all the matrices $\Delta({\bf x}_k^{(\alpha)})$
are lower triangular, while $\Omega^{{\bf m}}({\bf a})$ and  all the $\Omega({\bf x}^{(\alpha)}$
are lower triangular with unit diagonal, we conclude that
$$\Delta^{{\bf m}} ({\bf a}) = \lim_\alpha \Delta\big ({\bf x}^{(\alpha)} \big) \,
{\rm diag}\big(\Delta({\bf x}_1^{(\alpha)})^{-1}, \ldots, \Delta({\bf x}_r^{(\alpha)})^{-1}\big) \quad
({\bf x} \to {\bf a}^{{\bf m}})$$
for the lower triangular factor of $H^{{\bf m}} ({\bf a})^{-1}$.
Let us take $\alpha:={\bf t}=(t_1,\ldots,t_r)\to 0$ and
$${\bf x}_k^{({\bf t})} := \big[ a_k, a_k+t_k,a_k+2t_k, \ldots, a_k+m_k t_k \big] \quad (k=1,\ldots,r).$$
In this case we have
$$\Delta\big( {\bf x}_k^{({\bf t})} \big)^{-1} = .................................... $$


\bigskip
\noindent
\bf 5. Hermite interpolation as limit of Lagrange interpolations

\rm
\medskip
As far as we know no explicit formulas were published for the entries of both
$\Delta^{{\bf m}}({\bf a})$ and its inverse.
We achieve them below by a limiting process from Vandermonde cases.

\bigskip
\noindent
\bf Proposition. \it
Suppose $f:\RR\to\RR$ is a function being ${\cal C}^{m_k}$-smooth in suitable neighborhoods of the
points $a_k$ with $d^j/dx^j\vert_{x=a_k} f(x) = b_k^{(j)}$ $(k=1,\ldots,r;\ 0\le j\le m_k)$.
Assume a net of $(N+1)$-tuples ${\bf x}^{[\alpha]}$ consists of points with pairwise different coordinates
$(x_i^{[\alpha]}\ne x_j^{[\alpha]}$ for $i\ne j)$ and converges to ${\bf a}^{m}$.
Then ${\cal L}f\vert\{x_0^{[\alpha]},\ldots,x_N^{[\alpha]}\} (x) \to {\cal H}{\bf F}$.

\bigskip
\noindent
\bf Proof. \rm
Let $g:= f-{\cal H}{\bf F}$. Notice that $d^j/dx^j\vert_{x=a_k} f(x) = 0$ $(k=1,\ldots,r;\ 0\le j\le m_k)$.
Since ${\cal H}{\bf F}$ is a  polynomial of degree $N$, we have
${\cal L}{\cal H}{\bf F} ={\cal H}{\bf F}$ and hence
it suffices to see that
$${\cal L}g\vert\{x_0^{[\alpha]},\ldots,x_N^{[\alpha]}\} (x) \to 0.$$
This can be done by showing that for all its Newton differences,
$$g\big( x_n^{(\alpha)},\ldots,x_{n+s}^{(\alpha)} \big) \to 0 \qquad  (s=0\ldots,N;\ 0\le n\le N-s) .
\leqno(*)$$
We verify this statement by induction on the order index $s$.
For $s=0$ and fixed $n=(m_1+1)+\cdots+(m_k+1)+j-1$ we have
$g(x_n^{(\alpha)})\to g(a_{k+1})=0$ because of the continuity of $g$ at the points $a_1,\ldots,a_r$
and since $x_{(m_1+1)+\cdots+(m_k+1)+j-1}^{(\alpha)}\to a_{k+1}$ by assumption.
Assuming $(*)$ for some $s$, we consider the behavior of
$g\big( x_n^{(\alpha)},\ldots,x_{n+s+1}^{(\alpha)} \big)$ in two cases:
(1) if $x_n^{(\alpha)}\to a_i,x_{n+s+1}^{(\alpha)}\to a_j$ with $i\ne j$ (i.e. $a_i\ne a_j)$;
(2) if $x_n^{(\alpha)},x_{n+s+1}^{(\alpha)}\to a_i$.
In this case also $x_{n+1}^{(\alpha)},\ldots,x_{n+s}^{(\alpha)}\to a_i$.
In case (1) we have

$$\eqalign{
g\big( x_{n+1}^{(\alpha)},\ldots,x_{n+s+1}^{(\alpha)} \big) &=
{ g\big( x_{n+1}^{(\alpha)},\ldots,x_{n+s+1}^{(\alpha)} \big) -
g\big( x_n^{(\alpha)},\ldots,x_{n+s}^{(\alpha)} \big) \over x_{n+s+1}^{(\alpha)} - x_n^{(\alpha)} } \to \cr
&\to {0-0 \over a_j- i a_i} = 0 .\cr}$$
In case (2) we apply the fact that a Newton difference of order $(s+1)$ can be expressed by a derivation of
or order $(s+1)$ taken at some location between the most left and right base points:
$$g\big( x_{n+1}^{(\alpha)},\ldots,x_{n+s+1}^{(\alpha)} \big) =
{d^{s+1}\over dx^{s+1}}\Big\vert_{x=\theta_\alpha} g(x) \to
{d^{s+1}\over dx^{s+1}}\Big\vert_{x=a_i} g(x) =0$$
with a suitable net $\theta_\alpha \to a_i$.

......................................................

......................................................

......................................................

\bigskip
\noindent
\bf Corollary. \it
For any $r\in\RR$, let ${\bf x}^{[t]} :=\big[ a_k+jt: \ k=1,\ldots,r; \ 0\le j\le m_k\big]$
$\big($that is ${\bf x}^{[t]}=
(a_1,a_1+t,\ldots,a_1+m_1 t,a_2,\ldots,a_r+m_r t) \big)$
and let
${\bf y}^{[t]} := \big[ \sum\limits_{i=0}^{m_k} b_k^{(i)} {[(j-1)t]^i \over i!} : \ 0\le j\le m_k\big]$.
Then, for some $\eps>0$, the components of the tuples ${\bf x}^{[t]}$ are different if $\vert t\vert<\eps$.
Fixing such a value of $\eps$, with the function germs
${\bf f}^{[t]} := \big[ x_0^{[t]} \mapsto y_0^{[t]}, \ldots, x_N^{[t]} \mapsto y_N^{[t]} \big]$ $(-\eps<t<\eps)$
we have ${\cal H}{\bf F} =\lim\limits_{t\to 0} {\cal L}{\bf f}^{[t]}$ and
$H({\bf a}^{m})^{-1} {\bf b} = \lim\limits_{t\to 0} V({\bf x}^{[t]})^{-1} {\bf y}^{[t]}$.

\bigskip
\bf Proof.\rm

......................................................

......................................................

......................................................

\bigskip
\noindent
\bf Lemma. \it With the upper triangular Newton matrices
$N_k:= \Big[ \left( {j \atop i} \right) \Big]_{0\le i \le j \le m_k}$
we have
$$H({\bf a}^{\bf m}) = \lim_{t\to 0} V({\bf x}^{[t]})
\Big[ \mathop{\oplus}\limits_{k=1}^r \big[ N_k \ {\rm diag}(1^{-1},t,\ldots,t^{-m_k})\big]\Big] .$$

\bigskip
\bf Proof.\rm

......................................................

......................................................

......................................................

\bigskip
\noindent
\bf 5. Closed formulas for the entries of the upper triangular factors.

\bigskip
\noindent
\bf Lemma. \it For the subdiagonal matrix
$S:=\big[ 1_{\{ k\}}(\ell+1) k \big]_{k,\ell=0}^N =
\sum\limits_{k=1}^N k \ {\bf e}_k {\bf e}_{k-1}^\TR$
we have
$$\exp(tS) = \Big[ t^{k-\ell} \Big( {k\atop\ell}\Big)\Big]_{k,\ell=0}^N  =
\sum_{k=0}^N \sum_{\ell=0}^k t^{k-\ell} \Big( {k\atop\ell}\Big) {\bf e}_k{\bf e}_\ell^\TR .$$
\bf Proof. \rm
The $n$-th power of a subdiagonal matrix has non-zero entries
only in the $n$-th skew row below the diagonal.
Thus we can write
$$S^n:=\Big[ 1_{\{ k\}}(\ell+n) \sigma_\ell^{(n)} \Big]_{k,\ell=0}^N =
\sum_{\ell=0}^{N-n} \sigma_\ell^{(n)} {\bf e}_{\ell+n} {\bf e}_\ell^\TR .$$
We have the recursion
$$\eqalign{
&\sum_{\ell=0}^{N-(n+1)} \sigma_\ell^{(n+1)} {\bf e}_{\ell+n+1} {\bf e}_\ell^\TR = S^{n+1} = \cr
&= S^n S =
\Big[ \sum_{\ell=0}^{N-n} \sigma_\ell^{(n)} {\bf e}_{\ell+n} {\bf e}_\ell^\TR \Big]
\Big[ \sum_{k=1}^N k \ {\bf e}_k {\bf e}_{k-1}^\TR \Big]=
\sum_{\ell=1}^{N-n}  \sigma_\ell^{(n)} \ell \ {\bf e}_{\ell+n} {\bf e}_{\ell-1}^\TR ;\cr
&\sigma_\ell^{(n+1)} = \sigma_{\ell+1}^{(n)} \qquad (\ell=0,\ldots,N-(n+1)) .\cr}$$
Taking into account the definition of $S$ implying $\sigma_\ell^{(1)} = \ell+1$,
we conclude by induction on $n$ that
$\sigma_\ell^{(n)} = (\ell+1)(\ell+2)\cdots(\ell+n)=(\ell+n)!/\ell!$ in all cases.
It follows
$\exp(t S) = \sum\limits_{n=0}^N {t^n\over n!} S^n = \sum\limits_{n=0}^N t^n
\sum\limits_{\ell=0}^{N-n} {(\ell+n)!\over n!\ell!} {\bf e}_{\ell+n} {\bf e}_\ell^\TR =
\sum\limits_{k=0}^N \sum\limits_{\ell=0}^k t^{k-\ell} \Big( {k\atop\ell}\Big) {\bf e}_k{\bf e}_\ell^\TR$.



......................................................

......................................................

......................................................

\bigskip
\noindent
\bf References















%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
\bigskip
\bigskip

\hrule
\vskip1mm
\hrule
\bigskip
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%






${\bf v}(x):=\big[ 1,x,x^2,\ldots,x^N\big]$

${\bf x}:=\big[ x_0,x_1,x_2,\ldots,x_N\big]$

${\bf y}:=\big[ y_0,y_1,y_2,\ldots,y_N\big]$

${\bf f}:=\big[ (x_0,y_0),(x_1,y_1),\ldots,(x_N,y_N)\big]$ \quad function germ

$V({\bf x}) := \left[ \matrix{
{\bf v}(x_0) \cr \vdots \cr {\bf v}(x_N) \cr} \right]$ \quad Vandermonde matrix

${\cal L}{\bf f}(x) := \big[$Lagrange polynomial of ${\bf f}$ with symbolic variable $x\big]$

${\cal L}{\bf f}(x) = c_0+c_1 x+\cdots c_N x^N = {\bf v}(x) {\bf c}$ \ where \
${\bf c} := \big[ c_0,c_1,c_2,\ldots,c_N\big]^{{\rm T}}$

$V({\bf x}) {\bf c} = {\bf y}$ \quad interpolation equations

${\cal L}{\bf f}(x) = {\bf v}(x) V({\bf x})^{-1} {\bf y}$

${\cal L}{\bf f}(x) \!=\! {\bf f}(x_0) \!+\! {\bf f}(x_0,x_1)(x \!-\! x_0) \!+\!\cdots\!+\!
{\bf f}(x_0,\!\ldots\! ,x_N) (x \!-\! x_0)\!\cdots\! (x \!-\! x_{N \!-\! 1})$  Newton form

${\bf f}(x_0,\ldots,x_k) = \sum_{j=0}^k y_j \prod_{i: j\ne i \le k} (x_k-x_i)^{-1} =
{\bf d}_k({\bf x}) {\bf y}$ \quad Newton difference quotients

$\omega_k(x) := (x-x_0)\cdots(x-x_{k-1})$ , \ \ $\omega_0(x):=1$

$\omega_k(x) = \sum_{j=0}^k w_{k,j} x^j = {\bf v}(x) {\bf w}_k$

$\Omega({\bf x}) := \big[ {\bf w}_0({\bf x}) ,\ldots, {\bf w}_N({\bf x}) \big]$ \ \ upper triangular matrix

$D{\bf f}({\bf x}) := \left[ \matrix{
{\bf d}_0({\bf x}) \cr \vdots \cr {\bf d}_N ({\bf x}) \cr} \right]$ \ \ lower  triangular matrix

${\cal L}{\bf f}(x) = {\bf v}(x) \Omega({\bf x}) D{\bf f}({\bf x}) {\bf y}$

$V({\bf x})^{-1} = \Omega({\bf x}) D{\bf f}({\bf x})$ \quad triangular decomposition


\bigskip
\hrule
\bigskip

${\bf a} := \big[ a_1,\ldots, a_r \big]$, \quad $r\le N+1$

${\bf b} := \left[ \matrix{
{\bf b}_1 \cr \vdots \cr {\bf b}_r \cr} \right]$, \quad
${\bf b}_k := \left[ \matrix{
b_k^{(0)} \cr \vdots \cr b_k^{(m_k)} \cr} \right]$

$(m_1+1)+\cdots+(m_r+1) =N+1$

${\bf F}:= \big[ \big( (a_1, {\bf b}_1), \ldots, (a_r,{\bf b}_r) \big]$ \ function germ with derivatives

${\cal H}{\bf F}(x) :=  \big[$Hermite polynomial of ${\bf F}$ with symbolic variable $x\big]$

$\left. {d^j\over dx^j} \right\vert_{x=a_k}  {\cal H}{\bf F}(x) = b_k^{(j)}$
\quad
$(k=1,\ldots,r;\ 0\le j\le m_k)$ \ \
interpolation equations

Lagrangian case: $r=N+1$, \ $m_1=\cdots=m_{N+1}=0$

${\cal H}{\bf F}(x) = {\bf v}(x) {\bf c}$ \ \ classical polynomial form, \ ${\bf c} := \big[ c_0,\ldots,c_N\big]$

${\bf m}:= \big[ m_1,\ldots,m_r \big]$

$T_m(a) := \left[ \matrix{
{\bf v}(a) \cr d/dx \vert_{x=a} {\bf v}(x) \cr \vdots \cr d^m/dx^m \vert_{x=a} {\bf v}(x) \cr} \right]$ \ \
Taylor matrices

$H({\bf a}^{{\bf m}} ) := \left[ \matrix{
T_{m_1}(a_1) \cr \vdots \cr T_{m_r}(a_r) \cr} \right]$ \quad Hermite matrix

$H({\bf a}^{{\bf m}} ) {\bf c} = {\bf b}$ \quad interpolation equations

${\cal H}{\bf F}(x) = {\bf v}(x) H({\bf a}^{{\bf m}} )^{-1} {\bf b}$

We get the Newtonian form of ${\cal H}{\bf F}(x)$ as a limit of Vandermode cases

${\bf a}^{{\bf m}} (t) := \big[ a_1, a_1+t, \ldots, a_1+m_1 t, \ldots, a_r,a_r+t,\ldots,a_r+m_r t \big]$

${\bf b}(t) := \Big[ b_1^{(0)} , b_1^{(0)} + b_1^{(1)} t, \ldots, \sum_{j=0}^{m_1} b_1^{(j)} t^j/j!, \ldots,
b_r^{(0)} , b_r^{(0)} + b_r^{(1)} t, \ldots, \sum_{j=0}^{m_r} b_r^{(j)} t^j/j! \Big]^{{\rm T}}$

${\bf f}_t \!:=\! \Big[ \big( a_1,b_1^{(0)} \big), \big(a_1, b_1^{(0)}\!+\! b_1^{(1)}(t) \big),\ldots,
\big( a_r \!+\! m_r t, \sum\limits_{j=0}^{m_1} b_1^{(j)} t^j/j!) \Big]$ \ \
pairing of ${\bf a}^{{\bf m}}(t), {\bf b}(t)$

${\cal H}{\bf F}(x) = \lim\limits_{t\to 0} {\cal L}{\bf f}_t (x)$

\bye
