Fix $k\in \mathbb{N}^+$ and let $E=(e_i,f_i)_{i=1}^I\subset \mathbb{R}^n\times \mathbb{R}^m$ be a non-empty finite set with $e_i\neq e_j$ whenever $i\neq j$.
If $n=m=1$ then it's easy to see that: $$ p(x):= \sum_{i=1}^I\, \prod_{j\neq i, j=1,\dots,I}\, (e_j-x) \frac{f_i}{ \prod_{j\neq i, j=1,\dots,I} (e_j-e_i) } $$ smoothly interpolates $E$, in the sense that $p(e_i)=p_i$ and $p\in C^{\omega}(\mathbb{R},\mathbb{R})$. However, I'm more interested in the following refined+extended question.
In general (ie, for $n,m$ potentially not equal to $1$) does there exist a function $f\in C^{k,1}(\mathbb{R}^n,\mathbb{R}^m)$ satisfying: $$ f(e_i)=f_i \qquad \mbox{ for all }i=1,\dots,I ? $$ If not, are there mild conditions under which this type of interpolator must exist?
