#include <AMA.h>

Functions
long int	AMA_MltvGrdInterp (AMA_OPTIONS options, long int nind, long int ng, double *x, double z, long int degree, AMA_SPLINE *spline)
	Interpolation of Multivariate Gridded Data More...

Function Documentation

long int AMA_MltvGrdInterp	(	AMA_OPTIONS *	options,
		long int	nind,
		long int *	ng,
		double **	x,
		double *	z,
		long int *	degree,
		AMA_SPLINE **	spline
	)

Interpolation of Multivariate Gridded Data

This function employs cnspla to compute a spline interpolant of independent variable data ${\bf X}\in{\rm R}^n$ and dependent variable data ${\bf Z}\in{\rm R}^1$ . The independent variable data is given as $x_{k,\ell_k}$ , for $\ell_k=1,\ldots,N^g_k$ and $1\le k\le n$ , and it defines the rectilinear grid

$(x_{1,1},x_{1,2},\cdots,x_{1,N^g_1-1},x_{1,N^g_1})\times\cdots\times(x_{n,1},x_{n,2},\cdots,x_{n,N^g_n-1},x_{n,N^g_n}).$

The dependent variable data lies on the grid and is given as $z_{\ell_1,\cdots,\ell_n}$ , for $\ell_1=1,\ldots,N^g_1,\cdots,\ell_n=1,\ldots,N^g_n$ ; that is, there are $N=\prod_{k=1}^nN^g_k$ dependent variable values.

For a given set of data this function employs cnspla to compute the spline

$s({\bf X}) = \sum_{j_n=1}^{m_n}\cdots\sum_{j_1=1}^{m_1}\alpha_{j_1,\ldots,j_n}B_{d_1,j_1}(x_1|{\bf\Lambda}_1)\cdots B_{d_n,j_n}(x_n|{\bf\Lambda_n})$

that minimizes

$F^{(2)}(s({\bf X})) = \sum_{k=1}^n\left[\int_R\left(\partial^2 s({\bf X})\over\partial^2 x_k\right)^2dR \right]$

subject to the interpolation constraints

$s(x_{\ell_1},\cdots,x_{\ell_n}) = z_{\ell_1,\cdots,\ell_n}$

for $\ell_1=1,\ldots,N^g_1,\cdots,\ell_n=1,\ldots,N^g_n$ . This function also imposes constraints on the spline's partial derivative ${\partial s^{(p)}({\bf X})\over\partial x^p_k}$ at its left hand boundary $x_{k,1}$ and its right hand boundary $x_{k,N^g_k}$ , for all $1\le k\le n$ . These constraints are

${\partial s^{(p)}(x_{1,\ell_1},\cdots,x_{k,1},\cdots,x_{n,\ell_n})\over\partial x^p_k} = {\partial s^{(p)}(x_{1,\ell_1},\cdots,x_{k,N^g_k},\cdots,x_{n,\ell_n})\over\partial x^p_k} = 0.0$

for

$p = \cases{ {d_k\over 2 }+1,\ldots, d_k, & {\rm if the degree is even}; \cr {d_k+1\over 2} ,\ldots, d_k-1, & {\rm if the degree is odd}. \cr}$

These constraints are imposed at $x_{k_o,\ell_{k_o}}$ , for all $1\le \ell_{k_o}\le N^g_{k_o}$ where $1\le k_o \le n$ and $k_o\ne k$ .

In the above definition of $s({\bf X})$ the $\alpha_{j_1,\ldots,j_n}$ , for $j_k=1,\ldots,m_k$ and $1\le k\le n$ , are the coefficients of the tensor product B-splines

$\Phi_{j_1,\ldots,j_n}({\bf X}) = B_{d_1,j_1}(x_1|{\bf\Lambda}_1)\cdots B_{d_n,j_n}(x_n|{\bf\Lambda}_n)$

where $B_{d_k,j_k}(x_k|{\bf\Lambda}_k)$ are the $m_k$ univariate B-splines of degree $1 \le d_k\le 5$ defined by the knot vector ${\bf\Lambda}_k\in{\rm R}^{M_k}$ . The knot vectors ${\bf\Lambda}_k$ are

${\bf\Lambda}_k = ( \lambda_{k,1},\cdots,\lambda_{k,d_k+1}, \lambda_{k,d_k+2}, \cdots, \lambda_{k,m_k}, \lambda_{k,m_k+1},\cdots,\lambda_{k,m_k+d_k+1} )^T$

and they depend on the independent variable data and the degree $d_k$ . The knot vector ${\bf\Lambda}_k$ is defined by AMA_LamdaInterp() based on the points $(x_{k,1},x_{k,2},\cdots,x_{k,N^g_k-1},x_{k,N^g_k})$ .

This function employs an efficient algorithm for computing spline interpolants of multivariate gridded data. The algorithm consists of computing the splines

$\bar s_k({\bf X}) = \sum_{q_n=1}^{\bar m_n}\cdots\sum_{q_1=1}^{\bar m_1}\bar\alpha^k_{q_1,\ldots,q_n}B_{\bar d_1,q_1}(x_1|{\bf\bar\Lambda}_1)\cdots B_{\bar d_n,q_n}(x_n|{\bf\bar\Lambda_n})$

for $1\le k\le n$ . The coefficients $\bar\alpha^k_{q_1,\ldots,q_n}$ of $\bar s_k({\bf X})$ are the coefficients of the tensor product B-splines

$\Phi_{q_1,\ldots,q_n}({\bf X}) = B_{\bar d_1,q_1}(x_1|{\bf\bar\Lambda}_1)\cdots B_{\bar d_n,q_n}(x_n|{\bf\bar\Lambda}_n)$

where $B_{\bar d_k,q_k}(x_k|{\bf\bar\Lambda}_k)$ are the $\bar m_k$ univariate B-splines of degree $\bar d_k$ defined by the knot vector ${\bf\bar\Lambda}_k\in{\rm R}^{\bar m_k + \bar d_k + 1}$ . The degree $\bar d_{k_o}$ is

$\bar d_{k_o} = \cases { d_{k_o}, & for $1\le k_o\le k$, \cr 1 , & for $k < k_o\le n $; \cr }$

the number of coefficients $\bar m_{k_o}$ is

$\bar m_{k_o} = \cases { m_{k_o} , & for $1\le k_o\le k$, \cr n^g_{k_o}, & for $k < k_o\le n $; \cr }$

and the knot vector ${\bf\bar\Lambda}_{k_o}$ is

${\bf\bar\Lambda}_{k_o} = \cases { {\bf\Lambda}_{k_o} , & for $1\le k_o\le k$, \cr ( x_{k_o,1}, x_{k_o,1}, x_{k_o,2}, \cdots, x_{k_o,n^g_{k_o}-1}, x_{k_o,n^g_{k_o}}, x_{k_o,n^g_{k_o}})^T & for $k < k_o\le n $. \cr }$

The spline $s({\bf X})$ is based on the spline $\bar s_n({\bf X})$ and is defined by setting its coefficents $\alpha_{j_1,\ldots,j_n}=\bar\alpha^n_{j_1,\ldots,j_n}$ for $j_k=1,\ldots,m_k$ and $1\le k\le n$ . The spline $\bar s_n({\bf X})$ is defined through multiple invocations of AMA_UnvInterp() in the following manner.

First, the spline $\bar s_1({\bf X})$ is defined by invoking AMA_UnvInterp() $\prod_{k=2}^nN^g_{k}$ times to compute the univariate spline

$\bar s_{q_2,\cdots,q_n}(x_1) = \sum_{q_1=1}^{m_1}\bar\alpha^1_{q_1,\ldots,q_n}B_{d_1,q_1}(x_1|{\bf\Lambda}_1)$

that minimizes

$F^{(2)}(\bar s_{q_2,\cdots,q_n}(x_1)) = \int_{x_{1,1}}^{x_{1,N^g_1}}\left( \bar s^{(2)}_{q_2,\cdots,q_n}(x_1) \right)^2 dx_1$

subject to the interpolation constraints

$\bar s_{q_2,\cdots,q_n}(x_{1,\ell_1}) = z_{\ell_1,q_2,\cdots,q_n},$

for $\ell_1=1,\ldots,N^g_1$ , and the boundary constraints

$\bar s^{(p)}_{q_2,\cdots,q_n}(x_{1,1}) = \bar s^{(p)}_{q_2,\cdots,q_n}(x_{1,N^g_1}) = 0.0$

for

$p = \cases{ {d_1\over 2 }+1,\ldots, d_1, & {\rm if the degree is even}; \cr {d_1+1\over 2} ,\ldots, d_1-1, & {\rm if the degree is odd}. \cr}$

Computing the univariate splines $\bar s_{q_2,\cdots,q_n}(x_1)$ , for $q_k=1,\ldots,N^g_k$ and $2\le k\le n$ , defines the coefficients $\bar\alpha^1_{q_1,\ldots,q_n}$ of $\bar s_1({\bf X})$ .

After computing $\bar s_1({\bf X})$ the splines $\bar s_k({\bf X})$ , for $2\le k\le n$ , are defined by invoking AMA_UnvInterp() $\prod_{k_o=1}^{k-1}m_{k_o}\prod_{k_o=k+1}^nn^g_{k_o}$ times to compute the univariate spline

$\bar s_{q_1,\cdots,q_{k-1},q_{k+1},\cdots,q_n}(x_k) = \sum_{q_k=1}^{m_k}\bar\alpha^k_{q_1,\ldots,q_n}B_{d_k,q_k}(x_k|{\bf\Lambda}_k)$

that minimizes

$F^{(2)}(\bar s_{q_1,\cdots,q_{k-1},q_{k+1},\cdots,q_n}(x_k)) = \int_{x_{k,1}}^{x_{k,n^g_k}}\left( s^{(2)}_{q_1,\cdots,q_{k-1},q_{k+1},\cdots,q_n}(x_k) \right)^2 dx_k$

subject to the interpolation constraints

$\bar s_{q_1,\cdots,q_{k-1},q_{k+1},\cdots,q_n}(x_{k,\ell_k}) = \bar\alpha^{k-1}_{q_1,\cdots,q_{k-1},\ell_k,q_{k+1},\cdots,q_n}$

for $\ell_k=1,\ldots,N^g_k$ , and the boundary constraints

$\bar s^{(p)}_{q_1,\cdots,q_{k-1},q_{k+1},\cdots,q_n}(x_{k,1}) = \bar s^{(p)}_{q_1,\cdots,q_{k-1},q_{k+1},\cdots,q_n}(x_{k,N^g_k}) = 0.0$

for

$p = \cases{ {d_k\over 2 }+1,\ldots, d_k, & {\rm if the degree is even}; \cr {d_k+1\over 2} ,\ldots, d_k-1, & {\rm if the degree is odd}. \cr}$

Computing the univariate splines $\bar s_{q_1,\cdots,q_{k-1},q_{k+1},\cdots,q_n}(x_k)$ , for $q_{k_o}=1,\ldots,m_{k_o}$ when $1\le k_o\le k-1$ and for $q_{k_o}=1,\ldots,N^g_{k_o}$ when $k+1 \le k_o \le n$ , defines the coefficients $\bar\alpha^k_{q_1,\ldots,q_n}$ of $\bar s_k({\bf X})$ .

In total this algorithm consists of $\sum_{k=1}^n\left(\prod_{k_o=1}^{k-1}m_{k_o}\prod_{k_o=k+1}^nN^g_{k_o}\right)$ invocations of AMA_UnvInterp() to compute a multivariate spline which satisfies the interpolation and boundary constraints.

This function does the following:

Checks the input parameters for validity, see AMA_MltvGrdInputCheck().
Defines the spline $s({\bf X})$ based on the knot vectors ${\bf\Lambda}_k$ , for $1\le k\le n$ .
Defines the splines $\bar s_k({\bf X})$ , for $1\le k\le n$ , through multiple invocations of AMA_UnvInterp().
Stores the interpolant $s({\bf X})=\bar s_n({\bf X})$ into spline.

Note: By default the spline coefficients are initialized to zero but a different initial value for the spline coefficients can be set with AMA_OptionsSetCoefficients().; Because this function employs the univariate spline interpolation function AMA_UnvInterp(), it does not minimizes the cross partial terms of the penalty term utilized by AMA_MltvInterp().; The boundary constraints are imposed only at the grid points and not along the entire boundary.

Parameter Note: In the parameters description given below the limits on $k$ are $1\le k\le n$ , k = $k-1$ and N $=\prod_{k=1}^nN^g_k$ .

Parameters

options	[in] Pointer to AMA_OPTIONS. Must be initialized with AMA_Options() prior to calling AMA_MltvGrdInterp().
nind	[in] The number of independent variables . Must satisfy $2\le$ nind $\le$ AMA_MXNIND.
ng	[in] Array of size nind containing the number of points where ng[k] . Must satisfy ng[k] $\ge 2$ .
x	[in] Array of size nind containing arrays of size ng[k] where x[k] contains the independent variable data $x_{k,\ell_k}$ , for $\ell_k=1,\ldots,N^g_k$ . The values of x[k] must be in strictly ascending order.
z	[in] Array of size N containing the dependent variable data $z_{\ell_1,\cdots,\ell_n}$ , for $\ell_1=1,\ldots,N^g_1,\cdots,\ell_n=1,\ldots,N^g_n$ .
degree	[in] Array of size nind containing the degree where degree[k] . Must satisfy $0\le$ degree[k] $\le 5$ .
spline	[out] Pointer to AMA_SPLINE pointer containing the spline interpolant. Must satisfy spline $\ne$ NULL.

Returns: Success/Error Code.

User Callable Function - Documented 030215

Functions

Function Documentation