Solution (source code)

= Solution

If $x$ is a node, both sides vanish and any $\xi\in(a,b)$ works. Otherwise put $\omega(t)=\prod_i(t-x_i)$ and $C=(f(x)-p_n(x))/\omega(x)$. The function $H(t)=f(t)-p_n(t)-C\omega(t)$ has $n+2$ distinct zeros, at the nodes and at $x$. Repeated application of <Rolle's theorem> gives an interior point $\xi$ where $H^{(n+1)}(\xi)=0$. As $p_n^{(n+1)}=0$ and $\omega^{(n+1)}=(n+1)!$, this gives $C=f^{(n+1)}(\xi)/(n+1)!$. Therefore
$$
\boxed{f(x)-p_n(x)=\frac{f^{(n+1)}(\xi)}{(n+1)!}\prod_{i=0}^n(x-x_i).}
$$
The distinct nodes and $n\ge1$ ensure $a<b$, so the interior point is meaningful.