Higher order interpolation

Question

Let we have some points $(x_i,y_i),~i=1,2,...n,$ that are scattered on a $2D$ domain like square. And $u_h$ is the solution of a numerical method at these points. Now, if I want to construct a higher order interpolation of $u_h$ at these points, how can I do that? Should I find $u_h^{new}$ with some more points and then use the value of $u_h^{new}$ at these points? Is there any higher order interpolation that can construct with the same points in MATLAB? I could not use interp2.

Answer 1

The griddata function is finding a surface. Here is a way to see how the fit works.

Problem statement

The data is a sequence of $m$ points of the form $$\left\{ x_{k}, y_{k}, z_{k} \right\}_{k=1}^{m}$$

The model surface is an approximation of order $d$. Consider $d=2$: $$ z(x,y) = a_{00} + a_{10}x + a_{01}y + a_{20}x^{2} + a_{11}xy + a_{02}y^{2} $$

Linear system

$$ % \begin{align} % \mathbf{A} a &= z \\ % \left[ \begin{array}{cccccc} 1 & x_{1} & y_{1} & x_{1}^{2} & x_{1} y_{1} & y_{1}^{2} \\ \vdots & \vdots & \vdots & \vdots & \vdots & \vdots \\ 1 & x_{m} & y_{m} & x_{m}^{2} & x_{m} y_{m} & y_{m}^{2} \\ \end{array} \right] % \left[ \begin{array}{c} a_{00} \\ a_{10} \\ a_{01} \\ a_{20} \\ a_{11} \\ a_{02} \end{array} \right] % &= \left[ \begin{array}{c} z_{1} \\ \vdots \\ z_{m} \end{array} \right] %% \end{align} % $$

Note: the matrix $\mathbf{A}$ is a Vandermonde matrix, the poster child for ill-conditioned matrices.

Normal equations

Forming the normal equations squares the condition number of the system. The solution here is for theoretical insight. If numerical problems manifest, investigate other methods like SVD or QR.

As $\mathbf{A}^{*}\mathbf{A}$ is symmetric, only the upper half is shown.

$$ % \begin{align} % \mathbf{A}^{*}\mathbf{A} a &= \mathbf{A}^{*}z \\ % \left[ \begin{array}{llllll} m & \sum x & \sum y & \sum x^{2} & \sum xy & \sum y^{2} \\ & \sum x^{2} & \sum yx & \sum x^{3} & \sum x^{2}y & \sum xy^{2} \\ & & \sum y^{2} & \sum x^{2}y & \sum xy^{2} & \sum y^{3} \\ & & & \sum x^{3} & \sum x^{4}y & \sum x^{2}y^{2} \\ & & & & \sum x^{2}y^{2} & \sum xy^{3} \\ & & & & & \sum y^{4} \\ \end{array} \right] % \left[ \begin{array}{c} a_{00} \\ a_{10} \\ a_{01} \\ a_{20} \\ a_{11} \\ a_{02} \end{array} \right] % &= \left[ \begin{array}{l} \sum z \\ \sum xz \\ \sum yz \\ \sum x^{2}z \\ \sum xyz \\ \sum y^{2}z \end{array} \right] %% \end{align} % $$

Normal equations solution

$$ a_{LS} = \left( \mathbf{A}^{*}\mathbf{A} \right)^{-1} \mathbf{A}^{*}z $$

Increase degree of fit

Increase $d$ from $2$ to $3$: $$ z(x,y) = a_{00} + a_{10}x + a_{01}y + a_{20}x^{2} + a_{11}xy + a_{02}y^{2} + a_{30}x^{3} + a_{21}x^{2}y + a_{12}xy^{2} + a_{03} y^{3} $$ The new system matrix is has new columns $$ \mathbf{A} = \left[ \begin{array}{ccccccc|ccc} 1 & x_{1} & y_{1} & x_{1}^{2} & x_{1} y_{1} & y_{1}^{2} & x_{1}^{3} & x_{1}^{2} y_{1} & x_{1}y_{1}^{2} & y_{1}^{3} \\ \vdots & \vdots & \vdots & \vdots & \vdots & \vdots & \vdots & \vdots & \vdots & \vdots \\ 1 & x_{m} & y_{m} & x_{m}^{2} & x_{m} y_{m} & y_{m}^{2} & x_{m}^{3} & x_{m}^{2} y_{m} & x_{m}y_{m}^{2} & y_{m}^{3} \\ \end{array} \right] $$ Algorithmically, increasing the fit order is a simple matter. Again, the numerics may present trouble.

Also, avoid the indeterminate form $0^{0}$. Force these terms to unity.