Multivariate interpolation matching a scattered dataset with known first derivatives: is there an easy way? can someone indicate me a reasonably simple approach to the following problem?
I would like to approximate a function of 2 (or, if possible, more) variables $f(x,y)$.
I have an irregular grid $(x_i,y_i), \ i=1,2,...,n$ at which I know the value of the function $f(x_i,y_i)$ and its first derivatives $\partial f/\partial x | _{(x,y)=(x_i,y_i)}$, $\partial f/\partial y | _{(x,y)=(x_i,y_i)}$. I do not need a smooth surface: I use the result of this interpolation only to find a good estimatimate of $f$ for points outside the grid.
For all I know, approaches like nearest-neighbor interpolation cannot be modified to take into account the knowledge of the derivatives. On the other hand, multivariate Hermite interpolation on irregular grids is quite problematic.
Any suggestion? Thanks a lot!
 A: Here is one way. The basic idea is SSE(sum of squared terms error) minimisation. The only difference is that you add the derivatives information as new squared terms.
With this way, you need to have a general function (a model), with a few parameters to be determined. For example 
$f(x,y)=ax^2+by^2+cxy+d$.
You need to find a suitable model for your problem. Furthermore, you might wanna put a restriction on the model to make the computation more easy. choose models that are linear, with respect to the model parameters. In the example above, $a$, $b$ and $c$ are parameters and if there was no $d$, the model would have been linear. However, $d$ is necessary for adjustment.
Having a model $f(x,y,\vec p)$ of $m$ parameters $\vec p=\{p_i|i=1,...,m\}$, make a function as below.
$g(\vec p)=\sum_{i=1}^{n}(f(x_i,y_i,\vec p)-z_i)^2+\sum_{i=1}^{n}(f_x(x_i,y_i,\vec p)-v_i)^2+\sum_{i=1}^{n}(f_y(x_i,y_i,\vec p)-w_i)^2$.
$z_i$, $v_i$ and $w_i $ are expected values of the function, the derivative of the function, with respect to $x$, and the derivative of the function, with respect to $y$, for the $i$th sample, respectively. pairs $(x_i,y_i)$ give the $i$th sample points.
We would like to minimise the function $g(\vec p)$.Taking the derivative of $g(\vec p)$, with respect to its $m$ parameters, and putting them equal to $0$, would give $m$ linear functions. Then, it is left to solve a system of linear equations.
