Functions
AMA_Mltv.h File Reference

Go to the source code of this file.

Functions

long int AMA_MltvApprox (AMA_OPTIONS *options, long int nind, long int n, double **x, double *z, double *epsilon, long int *degree, AMA_SPLINE **spline)
 Approximation of Multivariate Data More...
 
long int AMA_MltvGrdApprox (AMA_OPTIONS *options, long int nind, long int *ng, double **x, double *z, double *wht, long int *degree, AMA_SPLINE **spline)
 Approximation of Multivariate Gridded Data More...
 
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...
 
long int AMA_MltvGrdLstsqr (AMA_OPTIONS *options, long int nind, long int *ng, double **x, double *z, double *wht, long int *degree, long int *mlamda, double **lamda, double theta, AMA_SPLINE **spline)
 Least Squares Approximation of Multivariate Gridded Data More...
 
long int AMA_MltvGrdMonoApprox (AMA_OPTIONS *options, long int nind, long int *ng, double **x, double *z, double *epsilon, long int *degree, AMA_SPLINE **spline)
 Monotonic Approximation of Multivariate Gridded Data More...
 
long int AMA_MltvGrdMonoInterp (AMA_OPTIONS *options, long int nind, long int *ng, double **x, double *z, long int *degree, AMA_SPLINE **spline)
 Monotonic Interpolation of Multivariate Gridded Data More...
 
long int AMA_MltvMonoApprox (AMA_OPTIONS *options, long int nind, long int n, double **x, double *z, double *epsilon, long int *degree, AMA_SPLINE **spline)
 Monotonic Approximation of Multivariate Data More...
 
long int AMA_MltvMonoInterp (AMA_OPTIONS *options, long int nind, long int n, double **x, double *z, long int *degree, AMA_SPLINE **spline)
 Monotonic Interpolation of Multivariate Data More...
 
long int AMA_MltvInterp (AMA_OPTIONS *options, long int nind, long int n, double **x, double *z, long int *degree, AMA_SPLINE **spline)
 Interpolation of Multivariate Data More...
 
long int AMA_MltvLstsqr (AMA_OPTIONS *options, long int nind, long int n, double **x, double *z, double *wht, long int *degree, long int *mlamda, double **lamda, double theta, AMA_SPLINE **spline)
 Least Squares Approximation of Multivariate Data More...
 
long int AMA_MltvBoundaryConditions (AMA_OPTIONS *options, enum AMA_Boolean interp, CNSPLA_SPLFUN *splfun)
 Define boundary conditions for multivariate approximation and interpolation functions. More...
 
long int AMA_MltvData (AMA_OPTIONS *options, const char *datfile, long int *nind, long int *degree, long int *n, double ***x, double **z, double **epsilon, double **wht, double *theta, long int *mlamda, double ***lamda)
 Read data and approximation options for AMA Spline Library Multivariate Random Data Functions. More...
 
long int AMA_MltvDataFree (AMA_OPTIONS *options, long int nind, double ***x, double **z, double **epsilon, double **wht, double ***lamda)
 Free Multivariate Data Functions data arrays. More...
 
long int AMA_MltvGrdData (AMA_OPTIONS *options, const char *datfile, long int *nind, long int *degree, long int *ng, double ***x, double **z, double **epsilon, double **wht, double *theta, long int *mlamda, double ***lamda)
 Read data and approximation options for AMA Spline Library Multivariate Gridded Data Functions. More...
 
long int AMA_MltvGrdInputCheck (AMA_OPTIONS *options, long int nind, long int *ng, double **x, double *z, double *wht, double *epsilon, long int *degree, long int *mlamda, double **lamda, enum AMA_Boolean lstsqr, enum AMA_Boolean *mltknt)
 Perform input check for multivariate gridded data. More...
 
long int AMA_MltvInputCheck (AMA_OPTIONS *options, long int nind, long int n, double **x, double *z, double *wht, double *epsilon, long int *degree, long int *mlamda, double **lamda, enum AMA_Boolean lstsqr)
 Perform input check for multivariate data approximation functions. More...
 
long int AMA_MltvPnltrm (AMA_OPTIONS *options, double theta, long int porder, CNSPLA_SPLFUN *splfun)
 Define penalty term on a multivariate spline. More...
 
long int AMA_MltvConpnt (AMA_OPTIONS *options, long int nind, long int n, double **x, double *z, double *epsilon, CNSPLA_SPLFUN *splfun)
 Define approximation or interpolation constraints on a multivariate spline. More...
 
long int AMA_MltvConreg (AMA_OPTIONS *options, long int nind, long int n, double **x, double *z, double *epsilon, long int *degree, CNSPLA_SPLFUN *splfun)
 Define monotonicity constraints on a multivariate spline. More...
 
long int AMA_MltvRglrze (AMA_OPTIONS *options, CNSPLA_SPLFUN *splfun)
 

Function Documentation

long int AMA_MltvApprox ( AMA_OPTIONS *  options,
long int  nind,
long int  n,
double **  x,
double *  z,
double *  epsilon,
long int *  degree,
AMA_SPLINE **  spline 
)

Approximation of Multivariate Data

This function employs cnspla to compute a spline approximation of independent variable data ${\bf X}\in{\rm R}^n$ and dependent variable data ${\bf Z}\in{\rm R}^1$. For a given set of independent variable data ${\bf X}_\ell=(x_{1,\ell},\cdots,x_{n,\ell})$, dependent variable data $z_\ell$ and approximation tolerances $\epsilon_\ell$, for $\ell=1,\ldots,N$, this function computes 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 + 2.0\sum_{l=k+1}^n\int_R\left(\partial^2 s({\bf X})\over\partial x_k \partial x_l\right)^2dR\right] \]

subject to the approximation constraints

\[ z_\ell - \epsilon_\ell \le s({\bf X}_\ell) \le z_\ell + \epsilon_\ell \]

for $\ell=1,\ldots,N$. The integration region $R = [x^L_1,x^U_1] \times [x^L_2,x^U_2] \times\cdots\times [x^L_n,x^U_n]$ where $x^L_k=\min_\ell\lbrace{x_{k,\ell}\rbrace}$ and $x^U_k=\max_\ell\lbrace{x_{k,\ell}\rbrace}$, for $1\le k\le n$.

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$. They are based on the rectilinear grid $(x^g_{1,1},x^g_{1,2},\cdots,x^g_{1,N^g_1-1},x^g_{1,N^g_1})\times\cdots\times(x^g_{n,1},x^g_{n,2},\cdots,x^g_{n,N^g_n-1},x^g_{n,N^g_n})$ upon which the independent variable data lies and they are defined by AMA_LamdaInterp(). That is, the knot vector ${\bf\Lambda}_k$ is defined by AMA_LamdaInterp() based on the points $(x^g_{k,1},x^g_{k,2},\cdots,x^g_{k,N^g_k-1},x^g_{k,N^g_k})$.

Additionally, the spline $s({\bf X})$ is subject to the bounds

\[ \alpha_l \le s({\bf X}) \le \alpha_u \]

where $\alpha_l=-\alpha_\infty$, $\alpha_u=\alpha_\infty$ and $\alpha_\infty$ equals AMA_SplineInfbnd(). By default the spline is unbounded but finite bounds can be set with AMA_OptionsSetBounds(). If the bounds specified with AMA_OptionsSetBounds() can not be satisfied in conjunction with the aforementioned approximation constraints, then AMA_MltvApprox() imposes the bounds

\[ \tilde\alpha_l \le s({\bf X}) \le \tilde\alpha_u \]

where $\tilde\alpha_l = \alpha_l - \alpha^s_l$ and $\tilde\alpha_u = \alpha_u + \alpha^s_u$. The values of $\alpha^s_l$ and $\alpha^s_u$ are defined by minimizing

\[ f(\alpha^s_l,\alpha^s_u) = ( \alpha^s_l )^2 + ( \alpha^s_u )^2 \]

subject to the constraints

\[ \begin{array}{ccc} \alpha^s_l &\ge &0.0, \cr \alpha^s_u &\ge &0.0, \cr \alpha_l - \alpha^s_l &\le &s({\bf X}), \cr \alpha_u + \alpha^s_u &\ge &s({\bf X}) \cr \end{array} \]

and the aforementioned approximation constraints. If either $\alpha^s_l > 0.0$ or $\alpha^s_u > 0.0$, then AMA_MltvApprox() reports AMA_WARNING_ERROR messages which specify the bounds $\tilde\alpha_l$ and $\tilde\alpha_u$ imposed on the spline. Additionally, it sets AMA_OPTIONS::lwrbnd $=\alpha^s_l$ and AMA_OPTIONS::uprbnd $=\alpha^s_u$.

This function does the following:

  • Checks the input parameters for validity, see AMA_MltvInputCheck().
  • Defines the spline $s({\bf X})$ based on the knot vectors ${\bf\Lambda}_k$, for $1\le k\le n$.
  • Defines CNSPLA_CONPNT constraints to represent the approximation constraints, see AMA_MltvConpnt().
  • Defines CNSPLA_PNLTRM condition to represent the function $F^{(2)}(s({\bf X}))$, see AMA_MltvPnltrm().
  • Invokes cnspla to compute $\alpha^s_l$ and/or $\alpha^s_u$ if lower and/or upper bounds are specified with AMA_OptionsSetBounds().
  • Invokes cnspla to compute a spline which minimizes $F^{(2)}(s({\bf X}))$ subject to the bounds and approximation constraints.
  • Stores approximation 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().
By default the spline is unbounded but finite bounds can be set with AMA_OptionsSetBounds().
By default the cross partial terms are included in the penalty term but they can be excluded with AMA_OptionsSetPenaltyTerm().
Because this function defines knot vectors based on the rectilinear grid upon which the independent variable data lies, it should only be used when the independent variable data either defines a gridded data with holes distribution, defines a gridded data distribution or consists of a small randomly distributed data set. See Multivariate Data Functions for a description of these data distributions.

Parameter Note: In the parameter definitions given below the limits on $k$ are $1\le k\le n$ and k = $k-1$.

Parameters
options[in] Pointer to AMA_OPTIONS. Must be initialized with AMA_Options() prior to calling AMA_MltvApprox().
nind[in] The number of independent variables $n$. Must satisfy $2\le$ nind $\le$ AMA_MXNIND.
n[in] The number of data points $N$. Must satisfy n $\ge 2$.
x[in] Array of size nind containing arrays of size n where x[k] contains the independent variable data $x_{k,\ell}$, for $\ell=1,\ldots,N$.
z[in] Array of size n containing the dependent variable data $z_\ell$, for $\ell=1,\ldots,N$.
epsilon[in] Array of size n containing the approximation tolerances $\epsilon_\ell$, for $\ell=1,\ldots,N$. Must satisfy $\epsilon_\ell\ge 0.0$ for all $\ell=1,\ldots,N$. If epsilon = NULL, then this function uses $\epsilon_\ell = 0.0$ for all $\ell=1,\ldots,N$.
degree[in] Array of size nind containing the degree $d_k$ where degree[k] $= d_k$. Must satisfy $1\le$ degree[k] $\le 5$.
spline[out] Pointer to AMA_SPLINE pointer containing the spline approximation. Must satisfy spline $\ne$ NULL.
Returns
Success/Error Code.

User Callable Function - Documented 022416

long int AMA_MltvBoundaryConditions ( AMA_OPTIONS *  options,
enum AMA_Boolean  interp,
CNSPLA_SPLFUN *  splfun 
)

Define boundary conditions for multivariate approximation and interpolation functions.

Parameters
options[in] Pointer to AMA_OPTIONS. Should be initialized with AMA_Options() prior to calling AMA_MltvBoundaryConditions().
interp[in] Interpolation function flag. It has one of the following two values:
splfun[in] Pointer to CNSPLA_SPLFUN containing the spline $s(x)$ upon which the penalty term is imposed. Should satisfy splfun $\ne$ NULL.
Returns
Success/Error Code.

User Support Function - Documented nnnnnn - !!!THIS IS NOT A USER CALLABLE FUNCTION - DOCUMENT IS INCLUDED FOR COMPLETENESS!!!

long int AMA_MltvConpnt ( AMA_OPTIONS *  options,
long int  nind,
long int  n,
double **  x,
double *  z,
double *  epsilon,
CNSPLA_SPLFUN *  splfun 
)

Define approximation or interpolation constraints on a multivariate spline.

For a given set of independent variable data ${\bf X}_\ell=(x_{1,\ell},\cdots,x_{n,\ell})$, dependent variable data $z_\ell$ and approximation tolerances $\epsilon_\ell$, for $\ell=1,\ldots,N$, this function defines either the approximation constraints

\[ z_\ell - \epsilon_\ell \le s({\bf X}_\ell) \le z_\ell + \epsilon_\ell \]

or the interpolation constraints

\[ s({\bf X}_\ell) = z_\ell, \]

for $\ell=1,\ldots,N$, on a multivariate spline $s({\bf X}):{\rm R}^n\rightarrow {\rm R}^1$ given as

\[ 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}). \]

Parameter Note: In the parameter definitions given below the limits on $k$ are $1\le k\le n$ and k = $k-1$.

Parameters
options[in] Pointer to AMA_OPTIONS. Should be initialized with AMA_Options() prior to calling AMA_MltvConpnt().
nind[in] The number of independent variables $n$. Should satisfy nind $\ge 2$.
n[in] The number of data points $N$. Should satisfy n $\ge 2$.
x[in] Array of size nind containing arrays of size n where x[k] contains the independent variable data $x_{k,\ell}$, for $\ell=1,\ldots,N$.
z[in] Array of size n containing the dependent variable data $z_\ell$, for $\ell=1,\ldots,N$.
epsilon[in] Array of size n containing the approximation tolerances $\epsilon_\ell$, for $\ell=1,\ldots,N$. Should satisfy $\epsilon_\ell\ge 0.0$ for all $\ell=1,\ldots,N$. If epsilon = NULL, then the interpolation constraints are defined. Otherwise, the approximation constraints are defined.
splfun[in] Pointer to CNSPLA_SPLFUN containing spline $s({\bf X})$ upon which approximation or interpolation constraints are imposed. Should satisfy splfun $\ne$ NULL.
Returns
Success/Error Code.

User Support Function - Documented 110914 - !!!THIS IS NOT A USER CALLABLE FUNCTION - DOCUMENT IS INCLUDED FOR COMPLETENESS!!!

long int AMA_MltvConreg ( AMA_OPTIONS *  options,
long int  nind,
long int  n,
double **  x,
double *  z,
double *  epsilon,
long int *  degree,
CNSPLA_SPLFUN *  splfun 
)

Define monotonicity constraints on a multivariate spline.

For a given set of independent variable data ${\bf X}_\ell=(x_{1,\ell},\cdots,x_{n,\ell})$, dependent variable data $z_\ell$ and approximation tolerances $\epsilon_\ell$, for $\ell=1,\ldots,N$, this function defines local monotonicity constraints on a multivariate spline $s({\bf X}):{\rm R}^n\rightarrow {\rm R}^1$ given as

\[ 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}). \]

This function first defines the rectilinear grid $(x^g_{1,1},x^g_{1,2},\cdots,x^g_{1,N^g_1-1},x^g_{1,N^g_1})\times\cdots\times(x^g_{n,1},x^g_{n,2},\cdots,x^g_{n,N^g_n-1},x^g_{n,N^g_n})$ upon which the independent variable data lies and then it defines the gridded dependent variable data $z^g_{l_1,\cdots,l_n}$ and the gridded approximation tolerances $\epsilon^g_{l_1,\cdots,l_n}$ for $1\le l_k\le N^g_k$ and $1\le k\le n$. If the independent variable data defines a gridded data distribution, then $N=\prod_{k=1}^n N^g_k$. If it defines a gridded data with holes distribution, then $N<\prod_{k=1}^n N^g_k$; that is, dependent variable data does not lie everywhere on the grid. In the latter case the local monotonicity constraints are not defined in regions of missing data.

After determining the data's underlying rectilinear grid this function imposes constraints on the first order partials ${\partial s({\bf X})\over \partial x_k}$, for all $1\le k\le n$, which insure the spline satisfies the local monotonicity constraints. The form of the constraints depends on the approximation tolerances.

If $\max_\ell\lbrace{\epsilon_\ell\rbrace}=0.0$, then for some $1\le k_o\le n$ a local monotonicity constraint is defined for all $1\le l_k\le N^g_k$ with $k \ne k_o$ and it is of the form

\[ {\partial s(x^g_{l_1},\cdots,x^g_{l_{k_o-1}},x_{k_o},x^g_{l_{k_o}+1},\cdots,x^g_{l_n})\over\partial x_{k_o}} \hspace{5pt} \cases { \ge 0.0 {\rm \hspace{5pt}for\hspace{5pt}all\hspace{5pt}} x_{k_o}\in[x^g_{l_{k_o}},x^g_{l_{k_o}+1}], & if $z^g_{l_{k_o}} < z^g_{l_{k_o}+1}$; \cr \le 0.0 {\rm \hspace{5pt}for\hspace{5pt}all\hspace{5pt}} x_{k_o}\in[x^g_{l_{k_o}},x^g_{l_{k_o}+1}], & if $z^g_{l_{k_o}} > z^g_{l_{k_o}+1}$; \cr = 0.0 {\rm \hspace{5pt}for\hspace{5pt}all\hspace{5pt}} x_{k_o}\in[x^g_{l_{k_o}},x^g_{l_{k_o}+1}], & otherwise.\cr } \]

The above constraint is defined for all $l_{k_o} = 1,\ldots, N^g_{k_o}-1$. The notation $z^g_{l_{k_o}}$ and $z^g_{l_{k_o}+1}$ refers to the dependent variable data associated with the grid points $(x^g_{l_1},\cdots,x^g_{l_{k_o-1}},x^g_{l_{k_o}},x^g_{l_{k_o}+1},\cdots,x^g_{l_n})$ and $(x^g_{l_1},\cdots,x^g_{l_{k_o-1}},x^g_{l_{k_o}},x^g_{l_{k_o}+1},\cdots,x^g_{l_n})$, respectively. If data does not lie at both grid points, then the constraint is not defined.

Alternatively, if there exists an $1\le \ell\le N$ for which $\epsilon_\ell\ne 0.0$, then for some $1\le k_o\le n$ a local monotonicity constraint is defined for all $1\le l_k\le N^g_k$ with $k \ne k_o$ and it is of the form

\[ {\partial s(x^g_{l_1},\cdots,x^g_{l_{k_o}-1},x_{k_o},x^g_{l_{k_o}+1},\cdots,x^g_{l_n})\over\partial x_{k_o}} \hspace{5pt} \cases { \ge 0.0 {\rm \hspace{5pt}for\hspace{5pt}all\hspace{5pt}} x_{k_o}\in[x^g_{l_{k_o}},x^g_{l_{k_o}+1}], & if $z^g_{l_{k_o}} + \epsilon^g_{l_{k_o}} < z^g_{l_{k_o}+1} - \epsilon^g_{l_{k_o}+1}$; \cr \le 0.0 {\rm \hspace{5pt}for\hspace{5pt}all\hspace{5pt}} x_{k_o}\in[x^g_{l_{k_o}},x^g_{l_{k_o}+1}], & if $z^g_{l_{k_o}} - \epsilon^g_{l_{k_o}} > z^g_{l_{k_o}+1} + \epsilon^g_{l_{k_o}+1}$. \cr } \]

The above constraint is defined for all $l_{k_o} = 1,\ldots, N^g_{k_o}-1$. The notation $\epsilon^g_{l_{k_o}}$ and $\epsilon^g_{l_{k_o}+1}$ refers to the approximation tolerances associated with the grid points $(x^g_{l_1},\cdots,x^g_{l_{k_o}-1},x^g_{l_{k_o}},x^g_{l_{k_o}+1},\cdots,x^g_{l_n})$ and $(x^g_{l_1},\cdots,x^g_{l_{k_o}-1},x^g_{l_{k_o}+1},x^g_{l_{k_o}+1},\cdots,x^g_{l_n})$, respectively. If data does not lie at both grid points, then the constraint is not defined.

If for some $1\le l_{k_o}\le N^g_{k_o}-1$ the above conditions are not met, then the intervals $[z^g_{l_{k_o}} - \epsilon^g_{l_{k_o}}, z^g_{l_{k_o}} + \epsilon^g_{l_{k_o}}]$ and $[z^g_{l_{k_o}+1} - \epsilon^g_{l_{k_o}+1}, z^g_{l_{k_o}+1} + \epsilon^g_{l_{k_o}+1}]$ intersect and the equality constraint

\[ {\partial s(x^g_{l_1},\cdots,x^g_{l_{k_o}-1},x_{k_o},x^g_{l_{k_o}+1},\cdots,x^g_{l_n})\over\partial x_{k_o}} = 0.0 \]

can be imposed over the interval $[x^g_{l_{k_o}},x^g_{l_{k_o}+1}]$. However, it is possible for the intervals $[z^g_q - \epsilon^g_q, z^g_q + \epsilon^g_q]$ and $[z^g_{q+1} - \epsilon^g_{q+1}, z^g_{q+1} + \epsilon^g_{q+1}]$ to intersect, for each $l_{k_o}\le q\le l_{k_o}+r$ and $r\ge 1$; without an intersection existing between all the intervals. In this case the aforementioned equality constraint is inconsistent with the approximation constraints

\[ z^g_q - \epsilon^g_q \le s(x^g_{l_1},\cdots,x^g_{l_{k_o}-1},x_q,x^g_{l_{k_o}+1},\cdots,x^g_{l_n}) \le z^g_q + \epsilon^g_q, \]

for $l_{k_o}\le q\le l_{k_o}+r+1$; and, therefore, can not be imposed over the interval $[x^g_{l_{k_o}},x^g_{l_{k_o}+r+1}]$. Hence, instead of imposing an equality constraint this function minimizes

\[ \Biggr\Vert {\partial s(x^g_{l_1},\cdots,x^g_{l_{k_o}-1},x_{k_o},x^g_{l_{k_o}+1},\cdots,x^g_{l_n})\over\partial x_{k_o}} \Biggr\Vert_\infty \]

over all intervals $[x^g_{l_{k_o}},x^g_{l_{k_o}+1}]$ for which the intervals $[z^g_{l_{k_o}} - \epsilon^g_{l_{k_o}}, z^g_{l_{k_o}}+ \epsilon^g_{l_{k_o}}]$ and $[z^g_{l_{k_o}+1} - \epsilon^g_{l_{k_o}+1}, z^g_{l_{k_o}+1} + \epsilon^g_{l_{k_o}+1}]$ intersect.

Note
By default the local monotonicity constraints are imposed on the spline in all of its independent variables but the constraints can be disabled in one or more of its independent variables with AMA_OptionsSetMonotonicity().

Parameter Note: In the parameter definitions given below the limits on $k$ are $1\le k\le n$ and k = $k-1$.

Parameters
options[in] Pointer to AMA_OPTIONS. Should be initialized with AMA_Options() prior to calling AMA_MltvConreg().
nind[in] The number of independent variables $n$. Should satisfy nind $\ge 2$.
n[in] The number of data points $N$. Should satisfy n $\ge 2$.
x[in] Array of size nind containing arrays of size n where x[k] contains the independent variable data $x_{k,\ell}$, for $\ell=1,\ldots,N$.
z[in] Array of size n containing the dependent variable data $z_\ell$, for $\ell=1,\ldots,N$.
epsilon[in] Array of size n containing the approximation tolerances $\epsilon_\ell$, for $\ell=1,\ldots,N$. Should satisfy $\epsilon_\ell\ge 0.0$ for all $\ell=1,\ldots,N$. If epsilon = NULL, then this function uses $\epsilon_\ell = 0.0$ for all $\ell=1,\ldots,N$.
degree[in] Array of size nind containing the degree $d_k$ where degree[k] $= d_k$. Should satisfy $1\le$ degree[k] $\le 5$.
splfun[in] Pointer to CNSPLA_SPLFUN containing spline $s({\bf X})$ upon which monotonicity constraints are imposed. Should satisfy splfun $\ne$ NULL.
Returns
Success/Error Code.

User Support Function - Documented 121114 - !!!THIS IS NOT A USER CALLABLE FUNCTION - DOCUMENT IS INCLUDED FOR COMPLETENESS!!!

long int AMA_MltvData ( AMA_OPTIONS *  options,
const char *  datfile,
long int *  nind,
long int *  degree,
long int *  n,
double ***  x,
double **  z,
double **  epsilon,
double **  wht,
double *  theta,
long int *  mlamda,
double ***  lamda 
)

Read data and approximation options for AMA Spline Library Multivariate Random Data Functions.

This function reads data and approximation options for AMA Spline Library functions which compute spline approximations of independent variable data ${\bf X}\in{\rm R}^n$ and dependent variable data ${\bf Z}\in{\rm R}^1$.

The argument datfile should reference a readable file which consists of a data section and several, optional, approximation options sections. The Data section must preceed the approximation options sections and has the following structure:

\[ \begin{array}{llllll} {\bf Data} & & & & & \cr {\bf Nind:}n &{\bf N:}N & & & & \cr {\bf Degree:}d_1 & & & & & \cr \vdots & & & & & \cr {\bf Degree:}d_n & & & & & \cr {\bf X_1} &\cdots &{\bf X}_n &{\bf Z} &{\bf Epsilon} &{\bf Wht} \cr x_{1,1} &\cdots &x_{n,1} &z_1 &\epsilon_1 &w_1 \cr x_{1,2} &\cdots &x_{n,2} &z_2 &\epsilon_2 &w_2 \cr \vdots &\vdots &\vdots &\vdots &\vdots &\vdots \cr x_{1,N} &\cdots &x_{n,N} &z_N &\epsilon_N &w_N \cr {\bf End\_Data} & & & & & \cr \end{array} \]

where the Epsilon and Wht columns are optional. The data must satisfy the following conditions:

  • $2\le n\le$ AMA_MXNIND;
  • $1\le d_k \le 5$, for $1\le k\le n$;
  • $N\ge 2$;
  • the approximation tolerances must satisfy $\epsilon_\ell\ge 0.0$, for $\ell=1,\cdots,N$;
  • the weights must satisfy $w_\ell\ge 0.0$, for $\ell=1,\cdots,N$.

Following the Data section may be one or more of the Bounds, Least_Squares, Monotonicity or Penalty_Term sections. If an approximation options section is not defined, then this function sets the options to their default values. See Table Approximation Options Defaults for a list of Multivariate Random Data Functions approximation options default values.

The Bounds section specifies the lower and upper bounds. It has the following structure:

\[ \begin{array}{ll} {\bf Bounds} & \cr {\bf Lwrbnd:}\alpha_l &{\bf Uprbnd:}\alpha_u \cr {\bf End\_Bounds} & \cr \end{array} \]

and its options must satisfy the following conditions:

The lower bound is read as a string lwrstr and based on the value of lwrstr the value of $\alpha_l$ is set as follows:

  • If lwrstr equals infbnd, then $\alpha_l = -\alpha_\infty$ where $\alpha_\infty$ = AMA_SplineInfbnd(); or,
  • If lwrstr equals zmin , then $\alpha_l = \min_\ell\lbrace{z_\ell\rbrace}$.
  • Otherwise, $\alpha_l =$ atof(lwrstr).

Similarly, the upper bound is read as a string uprstr and based on the value of uprstr the value of $\alpha_u$ is set as follows:

  • If uprstr equals infbnd, then $\alpha_u = \alpha_\infty$ where $\alpha_\infty$ = AMA_SplineInfbnd(); or,
  • If uprstr equals zmax , then $\alpha_u = \max_\ell\lbrace{z_\ell\rbrace}$.
  • Otherwise, $\alpha_u =$ atof(uprstr).

The Least_Squares section specifies the approximation options and knots employed by AMA_MltvLstsqr(). It has the following structure:

\[ \begin{array}{ll} {\bf Least\_Squares} & \cr {\bf Theta:}\theta &{\bf Minorm:}\rm minorm \cr {\bf Mlamda:}m_1 & \cr \lambda^o_{1,1} & \cr \lambda^o_{1,2} & \cr \vdots & \cr \lambda^o_{1,m_1} & \cr \vdots & \cr {\bf Mlamda:}m_n & \cr \lambda^o_{n,1} & \cr \lambda^o_{n,2} & \cr \vdots & \cr \lambda^o_{n,m_n} & \cr {\bf End\_Least\_Squares} & \cr \end{array} \]

and its options must satisfy the following conditions:

  • $0.0\le\theta< 1.0$;
  • ${\rm minorm}$ must be Enabled or Disabled;
  • $m_k\ge 2$, for $1\le k\le n$;
  • the knots $\lambda^o_{k,i}$, for $i=1,\cdots,m_k$ and $1\le k\le n$, must be in increasing order;
  • $\lambda^o_{k,1} \le x_{k,1}$ and $\lambda^o_{k,m_k} \ge x_{k,N^g_k}$, for all $1\le k\le n$.

The Monotonicity section specifies the monotonicity constraints and continuity conditions employed by AMA_MltvMonoApprox() and AMA_MltvMonoInterp(). It has the the following structure:

\[ \begin{array}{llll} {\bf Monotonicity} & & & \cr {\bf Monpos:\rm monpos[1]} &{\bf Monneg:\rm monneg[1]} &{\bf Monzer:\rm monzer[1]} &{\bf Concnd:\rm concnd[1]} \cr \vdots &\vdots &\vdots &\vdots \cr {\bf Monpos:\rm monpos[n]} &{\bf Monneg:\rm monneg[n]} &{\bf Monzer:\rm monzer[n]} &{\bf Concnd:\rm concnd[n]} \cr {\bf End\_Monotonicity} & & & \cr \end{array} \]

and its options must satisfy the following conditions:

  • the positive monotonicity constraint flag ${\rm monpos}[k]$, negative monotonicity constraint flag ${\rm monneg}[k]$ and zero monotonicity constraint flag ${\rm monzer}[k]$ must be Enabled or Disabled, for $1\le k\le n$;
  • the continuity condition ${\rm concnd}[k]$ must be Full or Reduced, for $1\le k\le n$;

The Penalty_Term section specifies which cross partial terms are used in the penalty term. It has the following structure:

\[ \begin{array}{l} {\bf Penalty\_Term} \cr {\bf Partial\_1\_2:\rm pnltrm[1][2]} \cr \vdots \cr {\bf Partial\_1\_n:\rm pnltrm[1][n]} \cr \vdots \cr {\bf Partial\_n-1\_n:\rm pnltrm[n-1][n]} \cr {\bf End\_Penalty\_Term} \cr \end{array} \]

and its options must satisfy the following condition:

  • ${\rm pnltrm[k][l]}$ must be Include or Exclude; for $1\le k\le n-1$ and $k+1\le l\le n$.

The bold keywords are case sensitive and the string values for the approximation options are case insensitive.

This function performs the following tasks:

Parameter Note: In the parameter definitions given below the limits on $k$ are $1\le k\le n$ and k = $k-1$.

Parameters
options[in] Pointer to AMA_OPTIONS. Must be initialized with AMA_Options() prior to calling AMA_MltvData().
datfile[in] The data file name. Must satisfy datfile $\ne$ NULL.
nind[out] The number of independent variables $n$. Must satisfy nind $\ne$ NULL.
degree[out] Array of size nind containing the degree $d_k$ where degree[k] $= d_k$. Must satisfy degree $\ne$ NULL.
n[out] The number of data points $N$. Must satisfy n $\ne$ NULL.
x[out] Array of size nind containing arrays of size n where x[k] contains the independent variable data $x_{k,\ell}$, for $\ell=1,\ldots,N$. Must satisfy x $\ne$ NULL.
z[out] Array of size n containing the dependent variable data $z_\ell$, for $\ell=1,\ldots,N$. Must satisfy z $\ne$ NULL.
epsilon[out] Array of size n containing the approximation tolerances $\epsilon_\ell$, for $\ell=1,\ldots,N$. Must satisfy epsilon $\ne$ NULL.
wht[out] Array of size n containing the weights $w_\ell$, for $\ell=1,\ldots,N$. Must satisfy wht $\ne$ NULL.
theta[out] The penalty term weight $\theta$. Must satisfy theta $\ne$ NULL.
mlamda[out] Array of size nind containing the number of knots $m_k$ where mlamda[k] $= m_k$. Must satisfy mlamda $\ne$ NULL.
lamda[out] Array of size nind containing arrays of size mlamda[k] where lamda[k] contains the knots $\lambda^o_{i,k}$, for $i=1,\ldots,m_k$. Must satisfy lambda $\ne$ NULL.
Returns
Success/Error Code.

User Support Function - Documented 110315 - !!!THIS IS NOT A USER CALLABLE FUNCTION - DOCUMENT IS INCLUDED FOR COMPLETENESS!!!

long int AMA_MltvDataFree ( AMA_OPTIONS *  options,
long int  nind,
double ***  x,
double **  z,
double **  epsilon,
double **  wht,
double ***  lamda 
)

Free Multivariate Data Functions data arrays.

This function frees the arrays allocated by AMA_MltvData() and AMA_MltvGrdData().

Parameters
options[in] Pointer to AMA_OPTIONS. Must be initialized with AMA_Options() prior to calling AMA_MltvDataFree().
nind[out] The number of independent variables $n$. Must satisfy nind $\ge 2$.
x[out] Pointer to array of size nind containing the independent variable data arrays.
z[out] Pointer to array containing the dependent variable data.
epsilon[out] Pointer to array containing the approximation tolerances.
wht[out] Pointer to array containing the weights.
lamda[out] Pointer to array of size nind containing the knot arrays.
Returns
Success/Error Code.

User Support Function - Documented 110515 - !!!THIS IS NOT A USER CALLABLE FUNCTION - DOCUMENT IS INCLUDED FOR COMPLETENESS!!!

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

Approximation of Multivariate Gridded Data

This function employs cnspla to compute a spline approximation 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. Similarly, the approximation tolerances are given as $\epsilon_{\ell_1,\cdots,\ell_n}$ for $\ell_1=1,\ldots,N^g_1,\cdots,\ell_n=1,\ldots,N^g_n$.

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 approximation constraints

\[ z_{\ell_1,\cdots,\ell_n} - \epsilon_{\ell_1,\cdots,\ell_n} \le s(x_{1,\ell_1},\cdots,x_{n,\ell_n}) \le z_{\ell_1,\cdots,\ell_n} + \epsilon_{\ell_1,\cdots,\ell_n}, \]

for $\ell_1=1,\ldots,N^g_1,\cdots,\ell_n=1,\ldots,N^g_n$.

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 approximations 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_UnvApprox() in the following manner.

First, the spline $\bar s_1({\bf X})$ is defined by invoking AMA_UnvApprox() $\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 approximation constraints

\[ z_{\ell_1,q_2,\cdots,q_n} - \epsilon_{\ell_1,q_2,\cdots,q_n} \le \bar s_{q_2,\cdots,q_n}(x_{1,\ell_1}) \le z_{\ell_1,q_2,\cdots,q_n} + \epsilon_{\ell_1,q_2,\cdots,q_n} \]

for $\ell_1=1,\ldots,N^g_1$. 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_UnvApprox() $\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$. 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_UnvApprox() to compute a multivariate spline which satisfies the approximation 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_UnvApprox().
  • Stores the approximation $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 approximation function AMA_UnvApprox(), it does not minimizes the cross partial terms of the penalty term utilized by AMA_MltvApprox().
Since the approximation constraints are only imposed during the computation of $\bar s_1({\bf X})$, the approximation depends on the order of the independent variables.

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_MltvGrdApprox().
nind[in] The number of independent variables $n$. Must satisfy $2\le$ nind $\le$ AMA_MXNIND.
ng[in] Array of size nind containing the number of points $N^g_k$ where ng[k] $= N^g_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$.
epsilon[in] Array of size N containing the approximation toleraneces $\epsilon_{\ell_1,\cdots,\ell_n}$, for $\ell_1=1,\ldots,N^g_1,\cdots,\ell_n=1,\ldots,N^g_n$. Must satisfy $\epsilon_{\ell_1,\cdots,\ell_n}\ge 0.0$ for all $\ell_1=1,\ldots,N^g_1,\cdots,\ell_n=1,\ldots,N^g_n$. If epsilon = NULL, then this function uses $\epsilon_{\ell_1,\cdots,\ell_n}=0.0$.
degree[in] Array of size nind containing the degree $d_k$ where degree[k] $= d_k$. Must satisfy $0\le$ degree[k] $\le 5$.
spline[out] Pointer to AMA_SPLINE pointer containing the spline approximation. Must satisfy spline $\ne$ NULL.
Returns
Success/Error Code.

User Callable Function - Documented 103015

long int AMA_MltvGrdData ( AMA_OPTIONS *  options,
const char *  datfile,
long int *  nind,
long int *  degree,
long int *  ng,
double ***  x,
double **  z,
double **  epsilon,
double **  wht,
double *  theta,
long int *  mlamda,
double ***  lamda 
)

Read data and approximation options for AMA Spline Library Multivariate Gridded Data Functions.

This function reads data for AMA Spline Library functions which compute spline approximations of independent variable data ${\bf X}\in{\rm R}^n$ and dependent variable data ${\bf Z}\in{\rm R}^1$. The independent variable data 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}) \]

where $N^g_k\ge 2$ and $1\le k\le 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. Similarly, the approximation tolerances $\epsilon_{\ell_1,\cdots,\ell_n}$ for $\ell_1=1,\ldots,N^g_1,\cdots,\ell_n=1,\ldots,N^g_n$ lie on the grid.

The argument datfile should reference a readable file which consists of a data section and several, optional, approximation options sections. The Data section must preceed the approximation options sections and has the following structure:

\[ \begin{array}{lll} {\bf Data} & & \cr {\bf Nind:}n & & \cr {\bf Degree:}d_1 & & \cr \vdots & & \cr {\bf Degree:}d_n & & \cr {\bf Ng:}N^g_1 & & \cr x_{1,1} & & \cr \vdots & & \cr x_{1,N^g_1} & & \cr {\bf Ng:}N^g_n & & \cr x_{n,1} & & \cr \vdots & & \cr x_{n,N^g_n} & & \cr {\bf Z} &{\bf Epsilon} &{\bf Wht} \cr z_{1,\cdots,1} &\epsilon_{1,\cdots,1} &w_{1,\cdots,1} \cr z_{2,\cdots,1} &\epsilon_{2,\cdots,1} &w_{2,\cdots,1} \cr \vdots &\vdots &\vdots \cr z_{N^g_1,\cdots,1} &\epsilon_{N^g_1,\cdots,1} &w_{N^g_1,\cdots,1} \cr \vdots &\vdots &\vdots \cr z_{1,\cdots,N^g_n} &\epsilon_{1,\cdots,N^g_n} &w_{1,\cdots,N^g_n} \cr z_{2,\cdots,N^g_n} &\epsilon_{2,\cdots,N^g_n} &w_{2,\cdots,N^g_n} \cr \vdots &\vdots &\vdots \cr z_{N^g_1,\cdots,N^g_n} &\epsilon_{N^g_1,\cdots,N^g_n} &w_{N^g_1,\cdots,N^g_n} \cr {\bf End\_Data} & & \cr \end{array} \]

where the Epsilon and Wht columns are optional. The data must satisfy the following conditions:

  • $2\le n\le$ AMA_MXNIND;
  • $1\le d_k \le 5$, for $1\le k\le n$;
  • $N^g_k\ge 2$, for $1\le k\le n$;
  • the independent variable data $x_{k,\ell_k}$, for $\ell_k=1,\cdots,N^g_k$ and $1\le k\le n$, must be in increasing order;
  • the approximation tolerances must satisfy $\epsilon_{\ell_1,\cdots,\ell_n}\ge 0.0$, for $\ell_1=1,\ldots,N^g_1,\cdots,\ell_n=1,\ldots,N^g_n$;
  • the weights must satisfy $w_{\ell_1,\cdots,\ell_n}\ge 0.0$, for $\ell_1=1,\ldots,N^g_1,\cdots,\ell_n=1,\ldots,N^g_n$.

Following the Data section may be one or more of the Bounds, Least_Squares or Monotonicity sections. If an approximation options section is not defined, then this function sets the options to their default values. See Table Approximation Options Defaults for a list of Multivariate Gridded Data Functions approximation options default values.

The Bounds section specifies the lower and upper bounds employed by AMA_MltvGrdLstsqr(). It has the following structure:

\[ \begin{array}{ll} {\bf Bounds} & \cr {\bf Lwrbnd:}\alpha_l &{\bf Uprbnd:}\alpha_u \cr {\bf End\_Bounds} & \cr \end{array} \]

and its options must satisfy the following conditions:

The lower bound is read as a string lwrstr and based on the value of lwrstr the value of $\alpha_l$ is set as follows:

  • If lwrstr equals infbnd, then $\alpha_l = -\alpha_\infty$ where $\alpha_\infty$ = AMA_SplineInfbnd(); or,
  • If lwrstr equals zmin , then $\alpha_l = \min_{\ell_1,\cdots,\ell_n}\lbrace{z_{\ell_1,\cdots,\ell_n}\rbrace}$.
  • Otherwise, $\alpha_l =$ atof(lwrstr).

Similarly, the upper bound is read as a string uprstr and based on the value of uprstr the value of $\alpha_u$ is set as follows:

  • If uprstr equals infbnd, then $\alpha_u = \alpha_\infty$ where $\alpha_\infty$ = AMA_SplineInfbnd(); or,
  • If uprstr equals zmax , then $\alpha_u = \max_{\ell_1,\cdots,\ell_n}\lbrace{z_{\ell_1,\cdots,\ell_n}\rbrace}$.
  • Otherwise, $\alpha_u =$ atof(uprstr).

The Least_Squares section specifies the approximation options and knots employed by AMA_MltvGrdLstsqr(). It has the following structure:

\[ \begin{array}{ll} {\bf Least\_Squares} & \cr {\bf Theta:}\theta &{\bf Minorm:}\rm minorm \cr {\bf Mlamda:}m_1 & \cr \lambda^o_{1,1} & \cr \lambda^o_{1,2} & \cr \vdots & \cr \lambda^o_{1,m_1} & \cr \vdots & \cr {\bf Mlamda:}m_n & \cr \lambda^o_{n,1} & \cr \lambda^o_{n,2} & \cr \vdots & \cr \lambda^o_{n,m_n} & \cr {\bf End\_Least\_Squares} & \cr \end{array} \]

and its options must satisfy the following conditions:

  • $0.0\le\theta< 1.0$;
  • ${\rm minorm}$ must be Enabled or Disabled;
  • $m_k\ge 2$, for $1\le k\le n$;
  • the knots $\lambda^o_{k,i}$, for $i=1,\cdots,m_k$ and $1\le k\le n$, must be in increasing order;
  • $\lambda^o_{k,1} \le x_{k,1}$ and $\lambda^o_{k,m_k} \ge x_{k,N^g_k}$, for all $1\le k\le n$.

The Monotonicity section specifies the monotonicity constraints and continuity conditions employed by AMA_MltvGrdMonoApprox() and AMA_MltvGrdMonoInterp(). It has the the following structure:

\[ \begin{array}{llll} {\bf Monotonicity} & & & \cr {\bf Monpos:\rm monpos[1]} &{\bf Monneg:\rm monneg[1]} &{\bf Monzer:\rm monzer[1]} &{\bf Concnd:\rm concnd[1]} \cr \vdots &\vdots &\vdots &\vdots \cr {\bf Monpos:\rm monpos[n]} &{\bf Monneg:\rm monneg[n]} &{\bf Monzer:\rm monzer[n]} &{\bf Concnd:\rm concnd[n]} \cr {\bf End\_Monotonicity} & & & \cr \end{array} \]

and its options must satisfy the following conditions:

  • the positive monotonicity constraint flag ${\rm monpos}[k]$, negative monotonicity constraint flag ${\rm monneg}[k]$ and zero monotonicity constraint flag ${\rm monzer}[k]$ must be Enabled or Disabled, for $1\le k\le n$;
  • the continuity condition ${\rm concnd}[k]$ must be Full or Reduced, for $1\le k\le n$;

The bold keywords are case sensitive and the string values for the approximation options are case insensitive.

This function performs the following tasks:

Parameter Note: In the parameter definitions given below the limits on $k$ are $1\le k\le n$ and k = $k-1$.

Parameters
options[in] Pointer to AMA_OPTIONS. Must be initialized with AMA_Options() prior to calling AMA_MltvGrdData().
datfile[in] The data file name. Must satisfy datfile $\ne$ NULL.
nind[out] The number of independent variables $n$. Must satisfy nind $\ne$ NULL.
degree[out] Array of size nind containing the degree $d_k$ where degree[k] $= d_k$. Must satisfy degree $\ne$ NULL.
ng[out] Array of size nind containing the number of points $N^g_k$ where ng[k] $= N^g_k$. Must satisfy ng $\ne$ NULL.
x[out] Array of size nind containing arrays of size ng[k] where x[k] contains the independent variable data $x_{k,\ell_k}$, $\ell_k=1,\ldots,N^g_k$. Must satisfy x $\ne$ NULL.
z[out] 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$. Must satisfy z $\ne$ NULL.
epsilon[out] Array of size $N$ containing the approximation tolerances $\epsilon_{\ell_1,\cdots,\ell_n}$, for $\ell_1=1,\ldots,N^g_1,\cdots,\ell_n=1,\ldots,N^g_n$. Must satisfy epsilon $\ne$ NULL.
wht[out] Array of size $N$ containing the weights $w_{\ell_1,\cdots,\ell_n}$, for $\ell_1=1,\ldots,N^g_1,\cdots,\ell_n=1,\ldots,N^g_n$. Must satisfy wht $\ne$ NULL.
theta[out] The penalty term weight $\theta$. Must satisfy theta $\ne$ NULL.
mlamda[out] Array of size nind containing the number of knots $m_k$ where mlamda[k] $= m_k$. Must satisfy mlamda $\ne$ NULL.
lamda[out] Array of size nind containing arrays of size mlamda[k] where lamda[k] contains the knots $\lambda^o_{i,k}$, for $i=1,\ldots,m_k$. Must satisfy lamda $\ne$ NULL.
Returns
Success/Error Code.

User Support Function - Documented 110315 - !!!THIS IS NOT A USER CALLABLE FUNCTION - DOCUMENT IS INCLUDED FOR COMPLETENESS!!!

long int AMA_MltvGrdInputCheck ( AMA_OPTIONS *  options,
long int  nind,
long int *  ng,
double **  x,
double *  z,
double *  wht,
double *  epsilon,
long int *  degree,
long int *  mlamda,
double **  lamda,
enum AMA_Boolean  lstsqr,
enum AMA_Boolean *  mltknt 
)

Perform input check for multivariate gridded data.

This function checks the input parameters for the Multivariate Gridded Data Functions which 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}) \]

of degree $d_k$ which is based on the knot vectors ${\bf\Lambda_k}$, for $1\le k\le n$.

Along with the degree the Multivariate Gridded Data Functions require the independent variable data $x_{k,\ell}$, for $\ell=1,\ldots,N^g_k$ and $1\le k\le n$, which 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}) \]

and either the weights $w_{\ell_1,\cdots,\ell_n}$ or the approximation tolerances $\epsilon_{\ell_1,\cdots,\ell_n}$, for $\ell_1=1,\ldots,N^g_1,\cdots,\ell_n=1,\ldots,N^g_n$. There are a total of N $=\prod_{k=1}^nN^g_k$ points. Also, for least squares approximation functions the knot vectors ${\bf\Lambda}_k\in{\rm R}^{M_k}$, for $1\le k\le n$, are required. This function insures the following conditions are satisfied:

  • $1\le n\le$ AMA_MXNIND;
  • $N^g_k\ge 2$, for all $1\le k\le n$;
  • $1\le d_k\le 5$, for all $1\le k\le n$;
  • $x_{k,\ell}$, for all $1\le k\le n$, is in ascending order and defines at least two distinct points;
  • $w_{\ell_1,\cdots,\ell_n}\ge 0.0$ for all $\ell_1=1,\ldots,N^g_1,\cdots,\ell_n=1,\ldots,N^g_n$ and $w_{\ell_1,\cdots,\ell_n} >0.0$ for at least one $\ell_1=1,\ldots,N^g_1,\cdots,\ell_n=1,\ldots,N^g_n$;
  • $\epsilon_{\ell_1,\cdots,\ell_n}\ge 0.0$ for all $\ell_1=1,\ldots,N^g_1,\cdots,\ell_n=1,\ldots,N^g_n$;
  • $M_k\ge 2$, for all $1\le k\le n$.

Also, if lstsqr equals AMA_Boolean_True and the knot vector ${\bf\Lambda}_k$ is represented by its distinct knot vector

\[ {\bf\rm K}^{{\bf\Lambda}_k} = ( \kappa_{k,1}, \kappa_{k,2}, \cdots, \kappa_{k,m^d_k-1}, \kappa_{k,m^d_k} )^T \]

and its knot multiplicity vector

\[ {\bf\rm H}^{{\bf\Lambda}_k} = ( \eta_{k,1}, \eta_{k,2}, \cdots, \eta_{k,m^d_k-1}, \eta_{k,m^d_k} )^T, \]

then, for all $1\le k\le n$, this function insures either the conditions

  • $\kappa_{k,i} < \kappa_{k,i+1}$, for $1\le i\le m^d_k - 1$;
  • $\eta_{k,i} \le d_k + 1$, for $2\le i\le m^d_k-1$;
  • $\eta_{k,1} = \eta_{k,m^d_k} = d_k+1$;

or

  • $\kappa_{k,i} < \kappa_{k,i+1}$, for $1\le i\le m^d_k - 1$;
  • $\eta_{k,i} \le d_k + 1$, for $2\le i\le m^d_k-1$;
  • $\eta_{k,1} = \eta_{k,m^d_k} = 1$;

are satisfied. It also insures the first and last distinct knots satisfy the conditions $\kappa_{k,1}\le x_{k,1}$ and $\kappa_{k,m^d_k}\ge x_{k,N^g_k}$, respectively.

Parameter Note: In the parameter definitions 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_MltvGrdInputCheck().
nind[in] The number of independent variables $n$. Must satisfy nind $\ge 2$.
ng[in] Array of size nind containing the number of points $N^g_k$ where ng[k] $= N^g_k$. Must satisfy ng $\ne$ NULL and 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$. Must satisfy x $\ne$ NULL and x[k] $\ne$ NULL. Also, x[k] must be in 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$. Must satisfy z $\ne$ NULL.
wht[in] Array of size N containing the weights $w_{\ell_1,\cdots,\ell_n}$, for $\ell_1=1,\ldots,N^g_1,\cdots,\ell_n=1,\ldots,N^g_n$. If wht $\ne$ NULL, then must satisfy $w_\ell\ge 0.0$ for all $\ell=1,\ldots,N$.
epsilon[in] Array of size N containing the approximation tolerances $\epsilon_{\ell_1,\cdots,\ell_n}$, for $\ell_1=1,\ldots,N^g_1,\cdots,\ell_n=1,\ldots,N^g_n$. If epsilon $\ne$ NULL, then must satisfy $\epsilon_\ell\ge 0.0$ for all $\ell=1,\ldots,N$.
degree[in] Array of size nind containing the degree $d_k$ where degree[k] $= d_k$. Must satisfy degree $\ne$ NULL and $1\le$ degree[k] $\le 5$.
mlamda[in] Array of size nind containing the number of knots $M_k$ where mlamda[k] $= M_k$. If lstsqr AMA_Boolean_True, then must satisfy mlamda $\ne$ NULL and mlamda[k] $\ge 2$.
lamda[in] Array of size nind containing arrays of size mlamda[k] where lamda[k] contains the knot vector ${\bf\Lambda}_k$. If lstsqr equals AMA_Boolean_True, then must satisfy lamda $\ne$ NULL and lamda[k] $\ne$ NULL. Also, lamda[k] must be a valid knot vector.
lstsqr[in] Least squares approximation flag. If it equals AMA_Boolean_True, then this function is being called by a least squares approximation function.
mltknt[out] Array of size nind containing the multiple knot flag. Defined only if lstsqr equals AMA_Boolean_True. It has one of the following two values:
  • If mltknt[k] equals AMA_Boolean_True, then ${\bf\Lambda}_k$ has $d_k+1$ order knots at its left and right hand boundaries.
  • If mltknt[k] equals AMA_Boolean_False, then ${\bf\Lambda}_k$ has single knots at its left and right hand boundaries.
Returns
Success/Error Code.

User Support Function - Documented 021116 - !!!THIS IS NOT A USER CALLABLE FUNCTION - DOCUMENT IS INCLUDED FOR COMPLETENESS!!!

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 $n$. Must satisfy $2\le$ nind $\le$ AMA_MXNIND.
ng[in] Array of size nind containing the number of points $N^g_k$ where ng[k] $= N^g_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 $d_k$ where degree[k] $= d_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

long int AMA_MltvGrdLstsqr ( AMA_OPTIONS *  options,
long int  nind,
long int *  ng,
double **  x,
double *  z,
double *  wht,
long int *  degree,
long int *  mlamda,
double **  lamda,
double  theta,
AMA_SPLINE **  spline 
)

Least Squares Approximation of Multivariate Gridded Data

This function employs cnspla to compute a least squares spline approximation 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

\[ (1.0-\theta)\sum_{\ell_n=1}^{N^g_n}\cdots\sum_{\ell_1=1}^{N^g_1} w_{\ell_1,\cdots,\ell_n}( s(x_{1,\ell_1},\cdots,x_{n,\ell_n}) - z_{\ell_1,\cdots,\ell_n} )^2 + \sum_{k=1}^n\left[\int_{\lambda_{k,1}}^{\lambda_{k,m_k+d_k+1}}\left(\partial^2 s({\bf X})\over\partial^2 x_k\right)^2dx_k\right] \]

for $0.0\le\theta< 1.0$.

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 vector ${\bf\Lambda}_k$ is given as

\[ {\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 \]

where the knots $\lambda_{k,i}$ must satisfy the conditions

  • $\lambda_{k,1} = \ldots = \lambda_{k,d_k+1}$;
  • $\lambda_{k,i} \le \lambda_{k,i+1}$, for $1\le i\le m_k+d_k$;
  • $\lambda_{k,i} < \lambda_{i+d_k+1}$, for $1\le i\le m_k$;
  • $\lambda_{k,m_k+1} = \ldots = \lambda_{k,m_k+d_k+1}$.

If the knot vector ${\bf\Lambda}_k$ is represented by its distinct knot vector

\[ {\rm\bf K}^{{\bf\Lambda}_k} = ( \kappa_{k,1}, \kappa_{k,2}, \cdots, \kappa_{k,m^d_k-1}, \kappa_{k,m^d_k} )^T \]

and its knot multiplicity vector

\[ {\rm\bf H}^{{\bf\Lambda}_k} = ( \eta_{k,1}, \eta_{k,2}, \cdots, \eta_{k,m^d_k-1}, \eta_{k,m^d_k} )^T \]

then the aforementioned conditions are

  • $\kappa_{k,i} < \kappa_{k,i+1}$, for $1\le i\le m^d_k - 1$;
  • $\eta_{k,i} \le d_k + 1$, for $2\le i\le m^d_k-1$;
  • $\eta_{k,1} = \eta_{k,m^d_k} = d_k+1$.

Also the first and last distinct knots must satisfy the conditions $\kappa_{k,1}\le x_{k,1}$ and $\kappa_{k,m^d_k}\ge x_{k,n^g_k}$, respectively. In this case the number of knots is $M_k = m_k + d_k + 1$.

For convenience this function does not require the definition of the $d_k+1$ order knots at the left and right hand endpoints and also accepts the knot vector ${\bf\Lambda}^o_k\in{\rm R}^{M_k}$ given as

\[ {\bf\Lambda}^o_k = ( \lambda_{k,d_k+1}, \lambda_{k,d_k+2}, \cdots, \lambda_{k,m_k}, \lambda_{k,m_k+1} )^T \]

where the knots $\lambda_{k,i}$ satisfy the conditions

  • $\lambda_{k,i} \le \lambda_{k,i+1}$, for $d_k+2\le i\le m_k-1$;
  • $\lambda_{k,i} < \lambda_{i+d_k+1}$, for $d_k+1\le i\le m_k-d_k$;
  • $\lambda_{k,d_k+1} < \lambda_{k,d_k+2}$ and $\lambda_{k,m_k} < \lambda_{k,m_k+1}$.

If the knot vector ${\bf\Lambda}^o$ is represented by its distinct knot vector

\[ {\rm\bf K}^{{\bf\Lambda}^o} = ( \kappa_{k,1}, \kappa_{k,2}, \cdots, \kappa_{k,m^d_k-1}, \kappa_{k,m^d_k} )^T \]

and its knot multiplicity vector

\[ {\rm\bf H}^{{\bf\Lambda}^o} = ( \eta_{k,1}, \eta_{k,2}, \cdots, \eta_{k,m^d_k-1}, \eta_{k,m^d_k} )^T \]

then the aforementioned conditions are

  • $\kappa_{k,i} < \kappa_{k,i+1}$, for $1\le i\le m^d_k - 1$;
  • $\eta_{k,i} \le d_k + 1$, for $2\le i\le m^d_k-1$;
  • $\eta_{k,1} = \eta_{k,m^d_k} = 1$.

As before, the first and last distinct knots must satisfy the conditions $\kappa_{k,1}\le x_{k,1}$ and $\kappa_{k,m^d_k}\ge x_{k,n^g_k}$, respectively. In this case the number of knots is $M_k = m_k - d_k + 1$.

This function employs an efficient algorithm for computing least squares spline approximations 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_UnvLstsqr() in the following manner.

First, the spline $\bar s_1({\bf X})$ is defined by invoking AMA_UnvLstsqr() $\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

\[ (1.0-\theta) \sum_{\ell_1=1}^{N^g_1} w_{\ell_1,q_2,\cdots,q_n} ( \bar s_{q_2,\cdots,q_n}(x_{1,\ell_1}) - z_{\ell_1,q_2,\cdots,q_n} )^2 + \theta \int_{\lambda_{1,1}}^{\lambda_{1,m_1+d_1+1}}\left( \bar s^{(2)}_{q_2,\cdots,q_n}(x_1) \right)^2 dx_1. \]

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_UnvLstsqr() $\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

\[ (1.0-\theta) \sum_{\ell_k=1}^{N^g_k} ( \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} )^2 + \theta \int_{\lambda_{k,1}}^{\lambda_{k,m_k+d_k+1}}\left( \bar s^{(2)}_{q_1,\cdots,q_{k-1},q_{k+1},\cdots,q_n}(x_k) \right)^2 dx_k. \]

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_UnvLstsqr() to compute a least squares spline approximation.

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_UnvLstsqr().
  • Stores the least squares approximation $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().
By default the minimum norm optimization is not performed but it can be enabled with AMA_OptionsSetMinimumNorm().
Because this function employs the univariate spline approximation function AMA_UnvLstsqr(), it does not minimizes the cross partial terms of the penalty term utilized by AMA_MltvLstsqr().
If the minimum norm optimization option is enabled with AMA_OptionsSetMinimumNorm(), then it is performed during each invocation of AMA_UnvLstsqr().
The penalty term and minimum norm optimization options are mutually exclusive; that is, the minimum norm optimization is not performed if $\theta > 0.0$.
Since the weights are only used during the computation of $\bar s_1({\bf X})$, the least squares approximation depends on the order of the independent variables.

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_MltvGrdLstsqr().
nind[in] The number of independent variables $n$. Must satisfy $2\le$ nind $\le$ AMA_MXNIND.
ng[in] Array of size nind containing the number of points $N^g_k$ where ng[k] $= N^g_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 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$.
wht[in] Array of size N containing the weights $w_{\ell_1,\cdots,\ell_n}$, for $\ell_1=1,\ldots,N^g_1,\cdots,\ell_n=1,\ldots,N^g_n$. Must satisfy $w_{\ell_1,\cdots,\ell_n}\ge 0.0$ for all $\ell_1=1,\ldots,N^g_1,\cdots,\ell_n=1,\ldots,N^g_n$. If wht = NULL, then this function uses $w_{\ell_1,\cdots,\ell_n}=1.0$.
degree[in] Array of size nind containing the degree $d_k$ where degree[k] $= d_k$. Must satisfy $1\le$ degree[k] $\le 5$.
mlamda[in] Array of size nind containing the number of knots $M_k$ where mlamda[k] $= M_k$. Must satisfy mlamda[k] $\ge 2$.
lamda[in] Array of size nind containing arrays of size mlamda[k] where lamda[k] contains the knot vector ${\bf\Lambda}_k$ or ${\bf\Lambda}^o_k$.
theta[in] The penalty term weight $\theta$. Must satisfy $0.0\le$ theta $< 1.0$.
spline[out] Pointer to AMA_SPLINE pointer containing the least squares spline approximation. Must satisfy spline $\ne$ NULL.
Returns
Success/Error Code.

User Callable Function - Documented 103015

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

Monotonic Approximation of Multivariate Gridded Data

This function employs cnspla to compute a spline approximation 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. Similarly, the approximation tolerances are given as $\epsilon_{\ell_1,\cdots,\ell_n}$, for $\ell_1=1,\ldots,N^g_1,\cdots,\ell_n=1,\ldots,N^g_n$.

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 approximation constraints

\[ z_{\ell_1,\cdots,\ell_n} - \epsilon_{\ell_1,\cdots,\ell_n} \le s(x_{\ell_1},\cdots,x_{\ell_n}) \le z_{\ell_1,\cdots,\ell_n} + \epsilon_{\ell_1,\cdots,\ell_n}, \]

for $\ell_1=1,\ldots,N^g_1,\cdots,\ell_n=1,\ldots,N^g_n$. Additionally, for each $1\le k_o \le n$, this function provides the option to impose the local monotonicity constraints

\[ {\partial s({\bf X})\over\partial x_{k_o}} \hspace{5pt} \cases { \ge 0.0 {\rm \hspace{5pt}for\hspace{5pt}all\hspace{5pt}} x_{k_o}\in[x_{k_o,\ell_{k_o}},x_{k_o,\ell_{k_o}+1}], & if $z_{\ell_1,\cdots,\ell_{k_o},\cdots,\ell_n} + \epsilon_{\ell_1,\cdots,\ell_{k_o},\cdots,\ell_n} < z_{\ell_1,\cdots,\ell_{k_o}+1,\cdots,\ell_n} - \epsilon_{\ell_1,\cdots,\ell_{k_o}+1,\cdots,\ell_n}$; \cr \le 0.0 {\rm \hspace{5pt}for\hspace{5pt}all\hspace{5pt}} x_{k_o}\in[x_{k_o,\ell_{k_o}},x_{k_o,\ell_{k_o}+1}], & if $z_{\ell_1,\cdots,\ell_{k_o},\cdots,\ell_n} + \epsilon_{\ell_1,\cdots,\ell_{k_o},\cdots,\ell_n} > z_{\ell_1,\cdots,\ell_{k_o}+1,\cdots,\ell_n} - \epsilon_{\ell_1,\cdots,\ell_{k_o}+1,\cdots,\ell_n}$; \cr = 0.0 {\rm \hspace{5pt}for\hspace{5pt}all\hspace{5pt}} x_{k_o}\in[x_{k_o,\ell_{k_o}},x_{k_o,\ell_{k_o}+1}], & otherwise. \cr } \]

These constraints are imposed at $x_{\ell_1,\cdots,\ell_{k_o},\cdots,\ell_n}$, for all $1\le \ell_k \le N^g_k$ and $1\le k\le n$ with $k\ne k_o$. By default the local monotonicity constraints are imposed on the spline in all of its independent variables but the constraints can be disabled in one or more of its independent variables with AMA_OptionsSetMonotonicity().

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$. If the local monotonicity constraints are requested in the $k$-th independent variable, then the knot vector ${\bf\Lambda}_k$ is defined by AMA_LamdaMonoInterp(); otherwise, it is defined by AMA_LamdaInterp(). In either case the knot vector ${\bf\Lambda}_k$ is based on the points $(x_{k,1},x_{k,2},\cdots,x_{k,N^g_k-1},x_{k,N^g_k})$.

This function computes a spline which satisfies either the full or reduced continuity condition at the interior grid points, as well as, the approximation and local monotonicity constraints. The following table summarizes the continuity of $s({\bf X})$ with respect to $x_k$ at the interior grid points $x_{k,l_k}$, for $l_k=2,\ldots,n^g_k-1$, and the number of coefficients, $m_k$, for a spline that satisfies either the full or reduced continuity condition. In the table the term

\[ s^{(p)}(x_{k,l_k}) = {\partial^p s(x_1,\cdots,x_{k,l_k},\cdots,x_n)\over \partial x^p_k} \]

is the $p$-th order partial of $s({\bf X})$ with respect to $x_k$ evaluated at $x_{k,l_k}$ and defined over the region $[x_{1,1},x_{1,n^g_1}]\times\cdots\times[x_{k-1,1},x_{k-1,n^g_{k-1}}]\times[x_{k+1,1},x_{k+1,n^g_{k+1}}]\times\cdots\times[x_{n,1},x_{n,n^g_n}]$.

Degree Full Continuity $m_k$ Reduced Continuity $m_k$
Linear $ s^{(p)}(x_{k,l_k})$ is continuous for $p=0$ $n^g_k$ $s^{(p)}(x_{k,l_k})$ is continuous for $p=0$ $n^g_k$
Quadratic $s^{(p)}(x_{k,l_k})$ is continuous for $0\le p\le 1$ $2n^g_k$ $s^{(p)}(x_{k,l_k})$ is continuous for $p=0$ $2n^g_k-1$
Cubic $s^{(p)}(x_{k,l_k})$ is continuous for $0\le p\le 2$ $3n^g_k$ $s^{(p)}(x_{k,l_k})$ is continuous for $0\le p\le 1$ $2n^g_k$
Quartic $s^{(p)}(x_{k,l_k})$ is continuous for $0\le p\le 3$ $4n^g_k$ $s^{(p)}(x_{k,l_k})$ is continuous for $0\le p\le 2$ $3n^g_k$
Quintic $s^{(p)}(x_{k,l_k})$ is continuous for $0\le p\le 4$ $5n^g_k$ $s^{(p)}(x_{k,l_k})$ is continuous for $0\le p\le 3$ $4n^g_k$

This function employs an efficient algorithm for computing monotonic spline approximations 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_UnvMonoApprox() in the following manner.

First, the spline $\bar s_1({\bf X})$ is defined by invoking AMA_UnvMonoApprox() $\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 approximation constraints

\[ z_{\ell_1,q_2,\cdots,q_n} - \epsilon_{\ell_1,q_2,\cdots,q_n} \le \bar s_{q_2,\cdots,q_n}(x_{1,\ell_1}) \le z_{\ell_1,q_2,\cdots,q_n} + \epsilon_{\ell_1,q_2,\cdots,q_n} \]

and the local monotonicity constraints

\[ {\partial s_{q_2,\cdots,q_n}(x_1)\over\partial x_1} \hspace{5pt} \cases { \ge 0.0 {\rm \hspace{5pt}for\hspace{5pt}all\hspace{5pt}} x_1\in[x_{1,\ell_1},x_{1,\ell_1+1}], & if $z_{\ell_1,q_2,\cdots,q_n} + \epsilon_{\ell_1,q_2,\cdots,q_n} < z_{\ell_1+1,q_2,\cdots,q_n} - \epsilon_{\ell_1+1,q_2,\cdots,q_n}$; \cr \le 0.0 {\rm \hspace{5pt}for\hspace{5pt}all\hspace{5pt}} x_{1}\in[x_{1,\ell_1},x_{1,\ell_1+1}], & if $z_{\ell_1,q_2,\cdots,q_n} + \epsilon_{\ell_1,q_2,\cdots,q_n} > z_{\ell_1+1,q_2,\cdots,q_n} - \epsilon_{\ell_1+1,q_2,\cdots,q_n}$; \cr = 0.0 {\rm \hspace{5pt}for\hspace{5pt}all\hspace{5pt}} x_1\in[x_{1,\ell_1},x_{1,\ell_1+1}], & otherwise; \cr } \]

for $\ell_1=1,\ldots,N^g_1$. 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_UnvMonoApprox() $\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} \]

and the local monotonicity constraints

\[ {\partial s_{q_1,\cdots,q_{k-1},q_{k+1},\cdots,q_n}(x_k)\over\partial x_k} \hspace{5pt} \cases { \ge 0.0 {\rm \hspace{5pt}for\hspace{5pt}all\hspace{5pt}} x_k\in[x_{k,\ell_k},x_{k,\ell_k+1}], & if $\bar\alpha^{k-1}_{q_1,\cdots,q_{k-1},\ell_k,q_{k+1},\cdots,q_n} < \bar\alpha^{k-1}_{q_1,\cdots,q_{k-1},\ell_k+1,q_{k+1},\cdots,q_n}$; \cr \le 0.0 {\rm \hspace{5pt}for\hspace{5pt}all\hspace{5pt}} x_{k}\in[x_{k,\ell_k},x_{k,\ell_k+1}], & if $\bar\alpha^{k-1}_{q_1,\cdots,q_{k-1},\ell_k,q_{k+1},\cdots,q_n} > \bar\alpha^{k-1}_{q_1,\cdots,q_{k-1},\ell_k+1,q_{k+1},\cdots,q_n}$; \cr = 0.0 {\rm \hspace{5pt}for\hspace{5pt}all\hspace{5pt}} x_k\in[x_{k,\ell_k},x_{k,\ell_{k+1}}], & otherwise; \cr } \]

for $\ell_k=1,\ldots,N^g_k$. 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_UnvMonoApprox() to compute a multivariate spline which satisfies the approximation and local monotonicity 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_UnvMonoApprox().
  • Stores the monotonic approximation $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 approximation function AMA_UnvMonoApprox(), it does not minimizes the cross partial terms of the penalty term utilized by AMA_MltvMonoApprox().
Since the approximation constraints are only imposed during the computation of $\bar s_1({\bf X})$, the approximation depends on the order of the independent variables.

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_MltvGrdMonoApprox().
nind[in] The number of independent variables $n$. Must satisfy $2\le$ nind $\le$ AMA_MXNIND.
ng[in] Array of size nind containing the number of points $N^g_k$ where ng[k] $= N^g_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$.
epsilon[in] Array of size N containing the approximation toleraneces $\epsilon_{\ell_1,\cdots,\ell_n}$, for $\ell_1=1,\ldots,N^g_1,\cdots,\ell_n=1,\ldots,N^g_n$. Must satisfy $\epsilon_{\ell_1,\cdots,\ell_n}\ge 0.0$ for all $\ell_1=1,\ldots,N^g_1,\cdots,\ell_n=1,\ldots,N^g_n$. If epsilon = NULL, then this function uses $\epsilon_{\ell_1,\cdots,\ell_n}=0.0$.
degree[in] Array of size nind containing the degree $d_k$ where degree[k] $= d_k$. Must satisfy $0\le$ degree[k] $\le 5$.
spline[out] Pointer to AMA_SPLINE pointer containing the monotonic spline approximation. Must satisfy spline $\ne$ NULL.
Returns
Success/Error Code.

User Callable Function - Documented 030215

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

Monotonic 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$. Additionally, for each $1\le k_o \le n$, this function provides the option to impose the local monotonicity constraints

\[ {\partial s({\bf X})\over\partial x_{k_o}} \hspace{5pt} \cases { \ge 0.0 {\rm \hspace{5pt}for\hspace{5pt}all\hspace{5pt}} x_{k_o}\in[x_{k_o,\ell_{k_o}},x_{k_o,\ell_{k_o}+1}], & if $z_{\ell_1,\cdots,\ell_{k_o},\cdots,\ell_n} < z_{\ell_1,\cdots,\ell_{k_o}+1,\cdots,\ell_n}$; \cr \le 0.0 {\rm \hspace{5pt}for\hspace{5pt}all\hspace{5pt}} x_{k_o}\in[x_{k_o,\ell_{k_o}},x_{k_o,\ell_{k_o}+1}], & if $z_{\ell_1,\cdots,\ell_{k_o},\cdots,\ell_n} > z_{\ell_1,\cdots,\ell_{k_o}+1,\cdots,\ell_n}$; \cr = 0.0 {\rm \hspace{5pt}for\hspace{5pt}all\hspace{5pt}} x_{k_o}\in[x_{k_o,\ell_{k_o}},x_{k_o,\ell_{k_o}+1}], & otherwise. \cr } \]

These constraints are imposed at $x_{\ell_1,\cdots,\ell_{k_o},\cdots,\ell_n}$, for all $1\le \ell_k \le N^g_k$ and $1\le k\le n$ with $k\ne k_o$. By default the local monotonicity constraints are imposed on the spline in all of its independent variables but the constraints can be disabled in one or more of its independent variables with AMA_OptionsSetMonotonicity().

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$. If the local monotonicity constraints are requested in the $k$-th independent variable, then the knot vector ${\bf\Lambda}_k$ is defined by AMA_LamdaMonoInterp(); otherwise, it is defined by AMA_LamdaInterp(). In either case the knot vector ${\bf\Lambda}_k$ is based on the points $(x_{k,1},x_{k,2},\cdots,x_{k,N^g_k-1},x_{k,N^g_k})$.

This function computes a spline which satisfies either the full or reduced continuity condition at the interior grid points, as well as, the interpolation and local monotonicity constraints. The following table summarizes the continuity of $s({\bf X})$ with respect to $x_k$ at the interior grid points $x_{k,l_k}$, for $l_k=2,\ldots,n^g_k-1$, and the number of coefficients, $m_k$, for a spline that satisfies either the full or reduced continuity condition. In the table the term

\[ s^{(p)}(x_{k,l_k}) = {\partial^p s(x_1,\cdots,x_{k,l_k},\cdots,x_n)\over \partial x^p_k} \]

is the $p$-th order partial of $s({\bf X})$ with respect to $x_k$ evaluated at $x_{k,l_k}$ and defined over the region $[x_{1,1},x_{1,n^g_1}]\times\cdots\times[x_{k-1,1},x_{k-1,n^g_{k-1}}]\times[x_{k+1,1},x_{k+1,n^g_{k+1}}]\times\cdots\times[x_{n,1},x_{n,n^g_n}]$.

Degree Full Continuity $m_k$ Reduced Continuity $m_k$
Linear $ s^{(p)}(x_{k,l_k})$ is continuous for $p=0$ $n^g_k$ $s^{(p)}(x_{k,l_k})$ is continuous for $p=0$ $n^g_k$
Quadratic $s^{(p)}(x_{k,l_k})$ is continuous for $0\le p\le 1$ $2n^g_k$ $s^{(p)}(x_{k,l_k})$ is continuous for $p=0$ $2n^g_k-1$
Cubic $s^{(p)}(x_{k,l_k})$ is continuous for $0\le p\le 2$ $3n^g_k$ $s^{(p)}(x_{k,l_k})$ is continuous for $0\le p\le 1$ $2n^g_k$
Quartic $s^{(p)}(x_{k,l_k})$ is continuous for $0\le p\le 3$ $4n^g_k$ $s^{(p)}(x_{k,l_k})$ is continuous for $0\le p\le 2$ $3n^g_k$
Quintic $s^{(p)}(x_{k,l_k})$ is continuous for $0\le p\le 4$ $5n^g_k$ $s^{(p)}(x_{k,l_k})$ is continuous for $0\le p\le 3$ $4n^g_k$

This function employs an efficient algorithm for computing monotonic 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_UnvMonoInterp() in the following manner.

First, the spline $\bar s_1({\bf X})$ is defined by invoking AMA_UnvMonoInterp() $\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} \]

and the local monotonicity constraints

\[ {\partial s_{q_2,\cdots,q_n}(x_1)\over\partial x_1} \hspace{5pt} \cases { \ge 0.0 {\rm \hspace{5pt}for\hspace{5pt}all\hspace{5pt}} x_1\in[x_{1,\ell_1},x_{1,\ell_1+1}], & if $z_{\ell_1,q_2,\cdots,q_n} < z_{\ell_1+1,q_2,\cdots,q_n}$; \cr \le 0.0 {\rm \hspace{5pt}for\hspace{5pt}all\hspace{5pt}} x_{1}\in[x_{1,\ell_1},x_{1,\ell_1+1}], & if $z_{\ell_1,q_2,\cdots,q_n} > z_{\ell_1+1,q_2,\cdots,q_n}$; \cr = 0.0 {\rm \hspace{5pt}for\hspace{5pt}all\hspace{5pt}} x_1\in[x_{1,\ell_1},x_{1,\ell_1+1}], & otherwise; \cr } \]

for $\ell_1=1,\ldots,N^g_1$. 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_UnvMonoInterp() $\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} \]

and the local monotonicity constraints

\[ {\partial s_{q_1,\cdots,q_{k-1},q_{k+1},\cdots,q_n}(x_k)\over\partial x_k} \hspace{5pt} \cases { \ge 0.0 {\rm \hspace{5pt}for\hspace{5pt}all\hspace{5pt}} x_k\in[x_{k,\ell_k},x_{k,\ell_k+1}], & if $\bar\alpha^{k-1}_{q_1,\cdots,q_{k-1},\ell_k,q_{k+1},\cdots,q_n} < \bar\alpha^{k-1}_{q_1,\cdots,q_{k-1},\ell_k+1,q_{k+1},\cdots,q_n}$; \cr \le 0.0 {\rm \hspace{5pt}for\hspace{5pt}all\hspace{5pt}} x_{k}\in[x_{k,\ell_k},x_{k,\ell_k+1}], & if $\bar\alpha^{k-1}_{q_1,\cdots,q_{k-1},\ell_k,q_{k+1},\cdots,q_n} > \bar\alpha^{k-1}_{q_1,\cdots,q_{k-1},\ell_k+1,q_{k+1},\cdots,q_n}$; \cr = 0.0 {\rm \hspace{5pt}for\hspace{5pt}all\hspace{5pt}} x_k\in[x_{k,\ell_k},x_{k,\ell_{k+1}}], & otherwise; \cr } \]

for $\ell_k=1,\ldots,N^g_k$. 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_UnvMonoInterp() to compute a multivariate spline which satisfies the interpolation and local monotonicity 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_UnvMonoInterp().
  • Stores the monotonic 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_UnvMonoInterp(), it does not minimizes the cross partial terms of the penalty term utilized by AMA_MltvInterp().

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_MltvGrdMonoInterp().
nind[in] The number of independent variables $n$. Must satisfy $2\le$ nind $\le$ AMA_MXNIND.
ng[in] Array of size nind containing the number of points $N^g_k$ where ng[k] $= N^g_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 $d_k$ where degree[k] $= d_k$. Must satisfy $0\le$ degree[k] $\le 5$.
spline[out] Pointer to AMA_SPLINE pointer containing the monotonic spline interpolant. Must satisfy spline $\ne$ NULL.
Returns
Success/Error Code.

User Callable Function - Documented 030215

long int AMA_MltvInputCheck ( AMA_OPTIONS *  options,
long int  nind,
long int  n,
double **  x,
double *  z,
double *  wht,
double *  epsilon,
long int *  degree,
long int *  mlamda,
double **  lamda,
enum AMA_Boolean  lstsqr 
)

Perform input check for multivariate data approximation functions.

This function checks the input parameters for the Multivariate Random Data Functions which 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}) \]

of degree $d_k$ which is based on the knot vectors ${\bf\Lambda_k}$, for $1\le k\le n$. Along with the degree the Multivariate Random Data Functions require the independent variable data data ${\bf X}_\ell=(x_{1,\ell},\cdots,x_{n,\ell})$, for $\ell=1,\ldots,N$, and either the weights $w_\ell$ or the approximation tolerances $\epsilon_\ell$, for $\ell=1,\ldots,N$. Also, for least squares approximation functions the knot vectors ${\bf\Lambda}_k\in{\rm R}^{M_k}$, for $1\le k\le n$, are required. This function insures the following conditions are satisfied:

Also, if lstsqr equals AMA_Boolean_True and the knot vector ${\bf\Lambda}_k$ is represented by its distinct knot vector

\[ {\bf\rm K}^{{\bf\Lambda}_k} = ( \kappa_{k,1}, \kappa_{k,2}, \cdots, \kappa_{k,m^d_k-1}, \kappa_{k,m^d_k} )^T \]

and its knot multiplicity vector

\[ {\bf\rm H}^{{\bf\Lambda}_k} = ( \eta_{k,1}, \eta_{k,2}, \cdots, \eta_{k,m^d_k-1}, \eta_{k,m^d_k} )^T, \]

then, for all $1\le k\le n$, this function insures either the conditions

  • $\kappa_{k,i} < \kappa_{k,i+1}$, for $1\le i\le m^d_k - 1$;
  • $\eta_{k,i} \le d_k + 1$, for $2\le i\le m^d_k-1$;
  • $\eta_{k,1} = \eta_{k,m^d_k} = d_k+1$;

or

  • $\kappa_{k,i} < \kappa_{k,i+1}$, for $1\le i\le m^d_k - 1$;
  • $\eta_{k,i} \le d_k + 1$, for $2\le i\le m^d_k-1$;
  • $\eta_{k,1} = \eta_{k,m^d_k} = 1$;

are satisfied. It also insures the first and last distinct knots satisfy the conditions $\kappa_{k,1}\le \min_\ell\lbrace{x_{k,\ell}\rbrace}$ and $\kappa_{k,m^d_k}\ge \max_\ell\lbrace{x_{k,\ell}\rbrace}$, respectively.

Parameter Note: In the parameter definitions given below the limits on $k$ are $1\le k\le n$ and k = $k-1$.

Parameters
options[in] Pointer to AMA_OPTIONS. Must be initialized with AMA_Options() prior to calling AMA_MltvInputCheck().
nind[in] The number of independent variables $n$. Must satisfy nind $\ge 2$.
n[in] The number of data points $N$. Must satisfy n $\ge 2$.
x[in] Array of size nind containing arrays of size n where x[k] contains the independent variable data $x_{k,\ell}$, for $\ell=1,\ldots,N$. Must satisfy x $\ne$ NULL and x[k] $\ne$ NULL. Also, x[k] must be in ascending order.
z[in] Array of size n containing the dependent variable data $z_\ell$, for $\ell=1,\ldots,N$. Must satisfy z $\ne$ NULL.
wht[in] Array of size n containing the weights $w_\ell$, for $\ell=1,\ldots,N$. If wht $\ne$ NULL, then must satisfy $w_\ell\ge 0.0$ for all $\ell=1,\ldots,N$.
epsilon[in] Array of size n containing the approximation tolerances $\epsilon_\ell$, for $\ell=1,\ldots,N$. If epsilon $\ne$ NULL, then must satisfy $\epsilon_\ell\ge 0.0$ for all $\ell=1,\ldots,N$.
degree[in] Array of size nind containing the degree $d_k$ where degree[k] $= d_k$. Must satisfy degree $\ne$ NULL and $1\le$ degree[k] $\le 5$.
mlamda[in] Array of size nind containing the number of knots $M_k$ where mlamda[k] $= M_k$. If lstsqr AMA_Boolean_True, then must satisfy mlamda $\ne$ NULL and mlamda[k] $\ge 2$.
lamda[in] Array of size nind containing arrays of size mlamda[k] where lamda[k] contains the knot vector ${\bf\Lambda}_k$. If lstsqr equals AMA_Boolean_True, then must satisfy lamda $\ne$ NULL and lamda[k] $\ne$ NULL. Also, lamda[k] must be a valid knot vector.
lstsqr[in] Least squares approximation flag. If it equals AMA_Boolean_True, then this function is being called by a least squares approximation function.
Returns
Success/Error Code.

User Support Function - Documented 021116 - !!!THIS IS NOT A USER CALLABLE FUNCTION - DOCUMENT IS INCLUDED FOR COMPLETENESS!!!

long int AMA_MltvInterp ( AMA_OPTIONS *  options,
long int  nind,
long int  n,
double **  x,
double *  z,
long int *  degree,
AMA_SPLINE **  spline 
)

Interpolation of Multivariate 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$. For a given set of independent variable data ${\bf X}_\ell=(x_{1,\ell},\cdots,x_{n,\ell})$ and dependent variable data $z_\ell$, for $\ell=1,\ldots,N$, this function computes 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 + 2.0\sum_{l=k+1}^n\int_R\left(\partial^2 s({\bf X})\over\partial x_k \partial x_l\right)^2dR\right] \]

subject to the interpolation constraints

\[ s({\bf X}_\ell) = z_\ell \]

for $\ell=1,\ldots,N$. The integration region $R = [x^L_1,x^U_1] \times [x^L_2,x^U_2] \times\cdots\times [x^L_n,x^U_n]$ where $x^L_k=\min_\ell\lbrace{x_{k,\ell}\rbrace}$ and $x^U_k=\max_\ell\lbrace{x_{k,\ell}\rbrace}$, for $1\le k\le n$.

This function also imposes constraints on the spline's partial derivative ${\partial s^{(p)}(x_1,\cdots,x_n)\over\partial x^p_k}$ at its left hand boundary $x^L_k$ and its right hand boundary $x^U_k$, for all $1\le k\le n$. These constraints are

\[ {\partial s^{(p)}(x_1,\cdots,x^L_k,\cdots,x_n)\over\partial x^p_k} = {\partial s^{(p)}(x_1,\cdots,x^U_k,\cdots,x_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 are the "so-called" natural boundary conditions and they are imposed at all $x_{k_o}\in [x^L_{k_o},x^U_{k_o}]$ for $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$. They are based on the rectilinear grid $(x^g_{1,1},x^g_{1,2},\cdots,x^g_{1,n^g_1-1},x^g_{1,n^g_1})\times\cdots\times(x^g_{n,1},x^g_{n,2},\cdots,x^g_{n,n^g_n-1},x^g_{n,n^g_n})$ upon which the independent variable data lies and they are defined by AMA_LamdaInterp(). That is, the knot vector ${\bf\Lambda}_k$ is defined by AMA_LamdaInterp() based on the points $(x^g_{k,1},x^g_{k,2},\cdots,x^g_{k,n^g_k-1},x^g_{k,n^g_k})$.

Additionally, the spline $s({\bf X})$ is subject to the bounds

\[ \alpha_l \le s({\bf X}) \le \alpha_u \]

where $\alpha_l=-\alpha_\infty$, $\alpha_u=\alpha_\infty$ and $\alpha_\infty$ equals AMA_SplineInfbnd(). By default the spline is unbounded but finite bounds can be set with AMA_OptionsSetBounds(). If the bounds specified with AMA_OptionsSetBounds() can not be satisfied in conjunction with the aforementioned interpolation and boundary constraints, then AMA_MltvInterp() imposes the bounds

\[ \tilde\alpha_l \le s({\bf X}) \le \tilde\alpha_u \]

where $\tilde\alpha_l = \alpha_l - \alpha^s_l$ and $\tilde\alpha_u = \alpha_u + \alpha^s_u$. The values of $\alpha^s_l$ and $\alpha^s_u$ are defined by minimizing

\[ f(\alpha^s_l,\alpha^s_u) = ( \alpha^s_l )^2 + ( \alpha^s_u )^2 \]

subject to the constraints

\[ \begin{array}{ccc} \alpha^s_l &\ge &0.0, \cr \alpha^s_u &\ge &0.0, \cr \alpha_l - \alpha^s_l &\le &s({\bf X}), \cr \alpha_u + \alpha^s_u &\ge &s({\bf X}) \cr \end{array} \]

and the aforementioned interpolation and boundary constraints. If either $\alpha^s_l > 0.0$ or $\alpha^s_u > 0.0$, then AMA_MltvInterp() reports AMA_WARNING_ERROR messages which specify the bounds $\tilde\alpha_l$ and $\tilde\alpha_u$ imposed on the spline. Additionally, it sets AMA_OPTIONS::lwrbnd $=\alpha^s_l$ and AMA_OPTIONS::uprbnd $=\alpha^s_u$.

This function does the following:

  • Checks the input parameters for validity, see AMA_MltvInputCheck().
  • Defines the spline $s({\bf X})$ based on the knot vectors ${\bf\Lambda}_k$, for $1\le k\le n$.
  • Defines CNSPLA_CONPNT constraints to represent the interpolation constraints, see AMA_MltvConpnt().
  • Defines CNSPLA_CONREG constraints to represent the boundary constraints.
  • Defines CNSPLA_PNLTRM condition to represent the function $F^{(2)}(s({\bf X}))$, see AMA_MltvPnltrm().
  • Invokes cnspla to compute a spline which minimizes $F^{(2)}(s({\bf X}))$ subject to the interpolation and boundary constraints.
  • Stores interpolant 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().
By default the spline is unbounded but finite bounds can be set with AMA_OptionsSetBounds().
By default the cross partial terms are included in the penalty term but they can be excluded with AMA_OptionsSetPenaltyTerm().
Because this function defines knot vectors based on the rectilinear grid upon which the independent variable data lies, it should only be used when the independent variable data either defines a gridded data with holes distribution, defines a gridded data distribution or consists of a small randomly distributed data set. See Multivariate Data Functions for a description of these data distributions.

Parameter Note: In the parameter definitions given below the limits on $k$ are $1\le k\le n$ and k = $k-1$.

Parameters
options[in] Pointer to AMA_OPTIONS. Must be initialized with AMA_Options() prior to calling AMA_MltvInterp().
nind[in] The number of independent variables $n$. Must satisfy $2\le$ nind $\le$ AMA_MXNIND.
n[in] The number of data points $N$. Must satisfy n $\ge 2$.
x[in] Array of size nind containing arrays of size n where x[k] contains the independent variable data $x_{k,\ell}$, for $\ell=1,\ldots,N$.
z[in] Array of size n containing the dependent variable data $z_\ell$, for $\ell=1,\ldots,N$.
degree[in] Array of size nind containing the degree $d_k$ where degree[k] $= d_k$. Must satisfy $1\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 022416

long int AMA_MltvLstsqr ( AMA_OPTIONS *  options,
long int  nind,
long int  n,
double **  x,
double *  z,
double *  wht,
long int *  degree,
long int *  mlamda,
double **  lamda,
double  theta,
AMA_SPLINE **  spline 
)

Least Squares Approximation of Multivariate Data

This function employs cnspla to compute a least squares spline approximation of independent variable data ${\bf X}\in{\rm R}^n$ and dependent variable data ${\bf Z}\in{\rm R}^1$. For a given set of independent variable data ${\bf X}_\ell=(x_{1,\ell},\cdots,x_{n,\ell})$, dependent variable data $z_\ell$ and weights $w_\ell$, for $\ell=1,\ldots,N$, this function computes 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

\[ (1.0-\theta)\sum_\ell^N w_\ell ( s({\bf X}_\ell) - z_\ell )^2 + \theta F^{(2)}\left(s({\bf X})\right) \]

for $0.0\le\theta< 1.0$. The function $F^{(2)}\left(s({\bf X})\right)$ imposes a penalty term on the least squares approximation.

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 vector ${\bf\Lambda}_k$ is given as

\[ {\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 \]

where the knots $\lambda_{k,i}$ must satisfy the conditions

  • $\lambda_{k,1} = \ldots = \lambda_{k,d_k+1}$;
  • $\lambda_{k,i} \le \lambda_{k,i+1}$, for $1\le i\le m_k+d_k$;
  • $\lambda_{k,i} < \lambda_{i+d_k+1}$, for $1\le i\le m_k$;
  • $\lambda_{k,m_k+1} = \ldots = \lambda_{k,m_k+d_k+1}$.

If the knot vector ${\bf\Lambda}_k$ is represented by its distinct knot vector

\[ {\bf\rm K}^{{\bf\Lambda}_k} = ( \kappa_{k,1}, \kappa_{k,2}, \cdots, \kappa_{k,m^d_k-1}, \kappa_{k,m^d_k} )^T \]

and its knot multiplicity vector

\[ {\bf\rm H}^{{\bf\Lambda}_k} = ( \eta_{k,1}, \eta_{k,2}, \cdots, \eta_{k,m^d_k-1}, \eta_{k,m^d_k} )^T \]

then the aforementioned conditions are

  • $\kappa_{k,i} < \kappa_{k,i+1}$, for $1\le i\le m^d_k - 1$;
  • $\eta_{k,i} \le d_k + 1$, for $2\le i\le m^d_k-1$;
  • $\eta_{k,1} = \eta_{k,m^d_k} = d_k+1$.

Also the first and last distinct knots must satisfy the conditions $\kappa_{k,1}\le \min_\ell\lbrace x_{k,\ell}\rbrace$ and $\kappa_{k,m^d_k}\ge \max_\ell\lbrace x_{k,\ell}\rbrace$, respectively. In this case the number of knots is $M_k = m_k + d_k + 1$.

For convenience this function does not require the definition of the $d_k+1$ order knots at the left and right hand endpoints and also accepts the knot vector ${\bf\Lambda}^o_k\in{\rm R}^{M_k}$ given as

\[ {\bf\Lambda}^o_k = ( \lambda_{k,d_k+1}, \lambda_{k,d_k+2}, \cdots, \lambda_{k,m_k}, \lambda_{k,m_k+1} )^T \]

where the knots $\lambda_{k,i}$ satisfy the conditions

  • $\lambda_{k,i} \le \lambda_{k,i+1}$, for $d_k+2\le i\le m_k-1$;
  • $\lambda_{k,i} < \lambda_{i+d_k+1}$, for $d_k+1\le i\le m_k-d_k$;
  • $\lambda_{k,d_k+1} < \lambda_{k,d_k+2}$ and $\lambda_{k,m_k} < \lambda_{k,m_k+1}$.

If the knot vector ${\bf\Lambda}^o$ is represented by its distinct knot vector

\[ {\bf\rm K}^{{\bf\Lambda}^o} = ( \kappa_{k,1}, \kappa_{k,2}, \cdots, \kappa_{k,m^d_k-1}, \kappa_{k,m^d_k} )^T \]

and its knot multiplicity vector

\[ {\bf\rm H}^{{\bf\Lambda}^o} = ( \eta_{k,1}, \eta_{k,2}, \cdots, \eta_{k,m^d_k-1}, \eta_{k,m^d_k} )^T \]

then the aforementioned conditions are

  • $\kappa_{k,i} < \kappa_{k,i+1}$, for $1\le i\le m^d_k - 1$;
  • $\eta_{k,i} \le d_k + 1$, for $2\le i\le m^d_k-1$;
  • $\eta_{k,1} = \eta_{k,m^d_k} = 1$.

As before, the first and last distinct knots must satisfy the conditions $\kappa_{k,1}\le \min_\ell\lbrace x_{k,\ell}\rbrace$ and $\kappa_{k,m^d_k}\ge \max_\ell\lbrace x_{k,\ell}\rbrace$, respectively. In this case the number of knots is $M_k = m_k - d_k + 1$.

The function $F^{(2)}(s({\bf X}))$ is used to impose a penalty term on the least squares approximation and is defined as

\[ 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 + 2.0\sum_{l=k+1}^n\int_R\left(\partial^2 s({\bf X})\over\partial x_k \partial x_l\right)^2dR\right] \]

where $R = [\lambda_{1,1},\lambda_{1,m_1+d_1+1}] \times [\lambda_{2,1},\lambda_{2,m_2+d_2+1}] \times\cdots\times [\lambda_{n,1},\lambda_{n,m_n+d_n+1}]$. The cross partial terms can be excluded from the penalty term $F^{(2)}(s({\bf X}))$ with AMA_OptionsSetPenaltyTerm().

Additionally, the spline $s({\bf X})$ is subject to the bounds

\[ \alpha_l \le s({\bf X}) \le \alpha_u \]

where $\alpha_l=-\alpha_\infty$, $\alpha_u=\alpha_\infty$ and $\alpha_\infty$ equals AMA_SplineInfbnd(). By default the spline is unbounded but finite bounds can be set with AMA_OptionsSetBounds().

Depending on the number and distribution of the knots relative to the independent variable data it is possible for the least squares approximation problem to have a non-unique solution. Hence, this function provides the capability of performing a minimum norm optimization subsequent to computing the least squares approximation. The minimum norm optimization computes an alternate spline $\bar s({\bf X})$ by minimizing $F^{(2)}(\bar s({\bf X}))$ subject to the approximation constraints

\[ z_l - \epsilon_\ell \le \bar s({\bf X}_\ell) \le z_l + \epsilon_\ell, \]

for $\ell=1,\ldots,N$, where $\epsilon_l = s({\bf X}_\ell)-z_\ell$ and $s({\bf X}_\ell)$ is the least squares approximation. These constraints insure the condition

\[ \sum_\ell^N w_\ell ( \bar s({\bf X}_\ell) - z_\ell )^2 \le \sum_\ell^N w_\ell ( s({\bf X}_\ell) - z_\ell )^2, \]

that is, the minimum norm optimization does not increase the least squares error.

This function does the following:

  • Checks input parameters for validity, see AMA_MltvInputCheck().
  • Defines the spline $s({\bf X})$ based on the knot vectors ${\bf\Lambda}_k$ or ${\bf\Lambda}^o_k$, for $1\le k\le n$.
  • Invokes cnspla to compute a least squares spline approximation.
  • Invokes cnspla to perform a minimum norm optimization.
  • Stores approximation 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().
By default the spline is unbounded but finite bounds can be set with AMA_OptionsSetBounds().
By default the cross partial terms are included in the penalty term but they can be excluded with AMA_OptionsSetPenaltyTerm().
By default the minimum norm optimization is not performed but it can be enabled with AMA_OptionsSetMinimumNorm().
The penalty term and minimum norm optimization options are mutually exclusive; that is, the minimum norm optimization is not performed if $\theta > 0.0$.

Parameter Note: In the parameter definitions given below the limits on $k$ are $1\le k\le n$ and k = $k-1$.

Parameters
options[in] Pointer to AMA_OPTIONS. Must be initialized with AMA_Options() prior to calling AMA_MltvLstsqr().
nind[in] The number of independent variables $n$. Must satisfy $2\le$ nind $\le$ AMA_MXNIND.
n[in] The number of data points $N$. Must satisfy n $\ge 2$.
x[in] Array of size nind containing arrays of size n where x[k] contains the independent variable data $x_{k,\ell}$, for $\ell=1,\ldots,N$.
z[in] Array of size n containing the dependent variable data $z_\ell$, for $\ell=1,\ldots,N$.
wht[in] Array of size n containing the weights $w_\ell$, for $\ell=1,\ldots,N$. Must satisfy $w_\ell\ge 0.0$ for all $\ell=1,\ldots,N$. If wht = NULL, then this function uses $w_\ell = 1.0$ for all $\ell=1,\ldots,N$.
degree[in] Array of size nind containing the degree $d_k$ where degree[k] $= d_k$. Must satisfy $1\le$ degree[k] $\le 5$.
mlamda[in] Array of size nind containing the number of knots $M_k$ where mlamda[k] $= M_k$. Must satisfy mlamda[k] $\ge 2$.
lamda[in] Array of size nind containing arrays of size mlamda[k] where lamda[k] contains the knot vector ${\bf\Lambda}_k$ or ${\bf\Lambda}^o_k$.
theta[in] The penalty term weight $\theta$. Must satisfy $0.0\le$ theta $< 1.0$.
spline[out] Pointer to AMA_SPLINE pointer containing the least squares spline approximation. Must satisfy spline $\ne$ NULL.
Returns
Success/Error Code.

User Callable Function - Documented 101715

long int AMA_MltvMonoApprox ( AMA_OPTIONS *  options,
long int  nind,
long int  n,
double **  x,
double *  z,
double *  epsilon,
long int *  degree,
AMA_SPLINE **  spline 
)

Monotonic Approximation of Multivariate Data

This function employs cnspla to compute a spline approximation of independent variable data ${\bf X}\in{\rm R}^n$ and dependent variable data ${\bf Z}\in{\rm R}^1$. For a given set of independent variable data ${\bf X}_\ell=(x_{1,\ell},\cdots,x_{n,\ell})$, dependent variable data $z_\ell$ and approximation tolerances $\epsilon_\ell$, for $\ell=1,\ldots,N$, this function computes 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 + 2.0\sum_{l=k+1}^n\int_R\left(\partial^2 s({\bf X})\over\partial x_k \partial x_l\right)^2dR\right] \]

subject to the approximation constraints

\[ z_\ell - \epsilon_\ell \le s({\bf X}_\ell) \le z_\ell + \epsilon_\ell \]

for $\ell=1,\ldots,N$. The integration region $R = [x^L_1,x^U_1] \times [x^L_2,x^U_2] \times\cdots\times [x^L_n,x^U_n]$ where $x^L_k=\min_\ell\lbrace{x_{k,\ell}\rbrace}$ and $x^U_k=\max_\ell\lbrace{x_{k,\ell}\rbrace}$, for $1\le k\le n$. This function also imposes constraints on the first order partials ${\partial s({\bf X})\over \partial x_k}$ for all $1\le k\le n$. These constraints are known as the local monotonicity constraints and they are defined by AMA_MltvConreg().

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, the degree $d_k$ and the continuity of the spline at the interior data points. They are based on the rectilinear grid $(x^g_{1,1},x^g_{1,2},\cdots,x^g_{1,n^g_1-1},x^g_{1,n^g_1})\times\cdots\times(x^g_{n,1},x^g_{n,2},\cdots,x^g_{n,n^g_n-1},x^g_{n,n^g_n})$ upon which the independent variable data lies and they are defined by AMA_LamdaMonoInterp(). That is, the knot vector ${\bf\Lambda}_k$ is defined by AMA_LamdaMonoInterp() based on the points $(x^g_{k,1},x^g_{k,2},\cdots,x^g_{k,n^g_k-1},x^g_{k,n^g_k})$.

Additionally, the spline $s({\bf X})$ is subject to the bounds

\[ \alpha_l \le s({\bf X}) \le \alpha_u \]

where $\alpha_l=-\alpha_\infty$, $\alpha_u=\alpha_\infty$ and $\alpha_\infty$ equals AMA_SplineInfbnd(). By default the spline is unbounded but finite bounds can be set with AMA_OptionsSetBounds(). If the bounds specified with AMA_OptionsSetBounds() can not be satisfied in conjunction with the aforementioned approximation and local monotonicity constraints, then AMA_MltvMonoApprox() imposes the bounds

\[ \tilde\alpha_l \le s({\bf X}) \le \tilde\alpha_u \]

where $\tilde\alpha_l = \alpha_l - \alpha^s_l$ and $\tilde\alpha_u = \alpha_u + \alpha^s_u$. The values of $\alpha^s_l$ and $\alpha^s_u$ are defined by minimizing

\[ f(\alpha^s_l,\alpha^s_u) = ( \alpha^s_l )^2 + ( \alpha^s_u )^2 \]

subject to the constraints

\[ \begin{array}{ccc} \alpha^s_l &\ge &0.0, \cr \alpha^s_u &\ge &0.0, \cr \alpha_l - \alpha^s_l &\le &s({\bf X}), \cr \alpha_u + \alpha^s_u &\ge &s({\bf X}) \cr \end{array} \]

and the aforementioned approximation and local monotonicity constraints. If either $\alpha^s_l > 0.0$ or $\alpha^s_u > 0.0$, then AMA_MltvMonoApprox() reports AMA_WARNING_ERROR messages which specify the bounds $\tilde\alpha_l$ and $\tilde\alpha_u$ imposed on the spline. Additionally, it sets AMA_OPTIONS::lwrbnd $=\alpha^s_l$ and AMA_OPTIONS::uprbnd $=\alpha^s_u$.

This function computes a spline which satisfies either the full or reduced continuity condition at the interior grid points, as well as, the approximation and local monotonicity constraints. The following table summarizes the continuity of $s({\bf X})$ with respect to $x_k$ at the interior grid points $x^g_{k,l_k}$, for $l_k=2,\ldots,N^g_k-1$, and the number of coefficients, $m_k$, for a spline that satisfies either the full or reduced continuity condition. In the table the term

\[ s^{(p)}(x^g_{k,l_k}) = {\partial^p s(x_1,\cdots,x^g_{k,l_k},\cdots,x_n)\over \partial x^p_k} \]

is the $p$-th order partial of $s({\bf X})$ with respect to $x_k$ evaluated at $x^g_{k,l_k}$ and defined over the region $[x^L_1,x^U_1]\times\cdots\times[x^L_{k-1},x^U_{k-1}]\times[x^L_{k+1},x^U_{k+1}]\times\cdots\times[x^L_n,x^U_n]$.

Degree Full Continuity $m_k$ Reduced Continuity $m_k$
Linear $ s^{(p)}(x^g_{k,l_k})$ is continuous for $p=0$ $N^g_k$ $s^{(p)}(x^g_{k,l_k})$ is continuous for $p=0$ $N^g_k$
Quadratic $s^{(p)}(x^g_{k,l_k})$ is continuous for $0\le p\le 1$ $2N^g_k$ $s^{(p)}(x^g_{k,l_k})$ is continuous for $p=0$ $2N^g_k-1$
Cubic $s^{(p)}(x^g_{k,l_k})$ is continuous for $0\le p\le 2$ $3N^g_k$ $s^{(p)}(x^g_{k,l_k})$ is continuous for $0\le p\le 1$ $2N^g_k$
Quartic $s^{(p)}(x^g_{k,l_k})$ is continuous for $0\le p\le 3$ $4N^g_k$ $s^{(p)}(x^g_{k,l_k})$ is continuous for $0\le p\le 2$ $3N^g_k$
Quintic $s^{(p)}(x^g_{k,l_k})$ is continuous for $0\le p\le 4$ $5N^g_k$ $s^{(p)}(x^g_{k,l_k})$ is continuous for $0\le p\le 3$ $4N^g_k$

This function does the following:

  • Checks the input parameters for validity, see AMA_MltvInputCheck().
  • Defines the spline $s({\bf X})$ based on the knot vectors ${\bf\Lambda}_k$, for $1\le k\le n$.
  • Defines CNSPLA_CONPNT constraints to represent the approximation constraints, see AMA_MltvConpnt().
  • Defines CNSPLA_PNLTRM condition to represent the function $F^{(2)}(s({\bf X}))$, see AMA_MltvPnltrm().
  • Defines CNSPLA_CONREG constraints to represent the local monotonicity constraints, see AMA_MltvConreg().
  • Invokes cnspla to compute a spline which minimizes $F^{(2)}(s({\bf X}))$ subject to the approximation and local monotonicity constraints.
  • Stores approximation 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().
By default the spline is unbounded but finite bounds can be set with AMA_OptionsSetBounds().
By default the spline satisfies the full continuity condition in all of its independent variables but the reduced continuity condition can be requested with AMA_OptionsSetContinuity().
By default the cross partial terms are included in the penalty term but they can be excluded with AMA_OptionsSetPenaltyTerm().
Because this function defines knot vectors based on the rectilinear grid upon which the independent variable data lies, it should only be used when the independent variable data either defines a gridded data with holes distribution, defines a gridded data distribution or consists of a small randomly distributed data set. See Multivariate Data Functions for a description of these data distributions.

Parameter Note: In the parameter definitions given below the limits on $k$ are $1\le k\le n$ and k = $k-1$.

Parameters
options[in] Pointer to AMA_OPTIONS. Must be initialized with AMA_Options() prior to calling AMA_MltvMonoApprox().
nind[in] The number of independent variables $n$. Must satisfy $2\le$ nind $\le$ AMA_MXNIND.
n[in] The number of data points $N$. Must satisfy n $\ge 2$.
x[in] Array of size nind containing arrays of size n where x[k] contains the independent variable data $x_{k,\ell}$, for $\ell=1,\ldots,N$.
z[in] Array of size n containing the dependent variable data $z_\ell$, for $\ell=1,\ldots,N$.
epsilon[in] Array of size n containing the approximation tolerances $\epsilon_\ell$, for $\ell=1,\ldots,N$. Must satisfy $\epsilon_\ell\ge 0.0$ for all $\ell=1,\ldots,N$. If epsilon = NULL, then this function uses $\epsilon_\ell = 0.0$ for all $\ell=1,\ldots,N$.
degree[in] Array of size nind containing the degree $d_k$ where degree[k] $= d_k$. Must satisfy $1\le$ degree[k] $\le 5$.
spline[out] Pointer to AMA_SPLINE pointer containing the monotonic spline approximation. Must satisfy spline $\ne$ NULL.
Returns
Success/Error Code.

User Callable Function - Documented 022416

long int AMA_MltvMonoInterp ( AMA_OPTIONS *  options,
long int  nind,
long int  n,
double **  x,
double *  z,
long int *  degree,
AMA_SPLINE **  spline 
)

Monotonic Interpolation of Multivariate 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$. For a given set of independent variable data ${\bf X}_\ell=(x_{1,\ell},\cdots,x_{n,\ell})$ and dependent variable data $z_\ell$, for $\ell=1,\ldots,N$, this function computes 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 + 2.0\sum_{l=k+1}^n\int_R\left(\partial^2 s({\bf X})\over\partial x_k \partial x_l\right)^2dR\right] \]

subject to the interpolation constraints

\[ s({\bf X}_\ell) = z_\ell \]

for $\ell=1,\ldots,N$. The integration region $R = [x^L_1,x^U_1] \times [x^L_2,x^U_2] \times\cdots\times [x^L_n,x^U_n]$ where $x^L_k=\min_\ell\lbrace{x_{k,\ell}\rbrace}$ and $x^U_k=\max_\ell\lbrace{x_{k,\ell}\rbrace}$, for $1\le k\le n$. This function also imposes constraints on the first order partials ${\partial s({\bf X})\over \partial x_k}$ for all $1\le k\le n$. These constraints are known as the local monotonicity constraints and they are defined by AMA_MltvConreg().

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, the degree $d_k$ and the continuity of the spline at the interior data points. They are based on the rectilinear grid $(x^g_{1,1},x^g_{1,2},\cdots,x^g_{1,n^g_1-1},x^g_{1,n^g_1})\times\cdots\times(x^g_{n,1},x^g_{n,2},\cdots,x^g_{n,n^g_n-1},x^g_{n,n^g_n})$ upon which the independent variable data lies and they are defined by AMA_LamdaMonoInterp(). That is, the knot vector ${\bf\Lambda}_k$ is defined by AMA_LamdaMonoInterp() based on the points $(x^g_{k,1},x^g_{k,2},\cdots,x^g_{k,n^g_k-1},x^g_{k,n^g_k})$.

Additionally, the spline $s({\bf X})$ is subject to the bounds

\[ \alpha_l \le s({\bf X}) \le \alpha_u \]

where $\alpha_l=-\alpha_\infty$, $\alpha_u=\alpha_\infty$ and $\alpha_\infty$ equals AMA_SplineInfbnd(). By default the spline is unbounded but finite bounds can be set with AMA_OptionsSetBounds(). If the bounds specified with AMA_OptionsSetBounds() can not be satisfied in conjunction with the aforementioned approximation and local monotonicity constraints, then AMA_MltvMonoApprox() imposes the bounds

\[ \tilde\alpha_l \le s({\bf X}) \le \tilde\alpha_u \]

where $\tilde\alpha_l = \alpha_l - \alpha^s_l$ and $\tilde\alpha_u = \alpha_u + \alpha^s_u$. The values of $\alpha^s_l$ and $\alpha^s_u$ are defined by minimizing

\[ f(\alpha^s_l,\alpha^s_u) = ( \alpha^s_l )^2 + ( \alpha^s_u )^2 \]

subject to the constraints

\[ \begin{array}{ccc} \alpha^s_l &\ge &0.0, \cr \alpha^s_u &\ge &0.0, \cr \alpha_l - \alpha^s_l &\le &s({\bf X}), \cr \alpha_u + \alpha^s_u &\ge &s({\bf X}) \cr \end{array} \]

and the aforementioned approximation and local monotonicity constraints. If either $\alpha^s_l > 0.0$ or $\alpha^s_u > 0.0$, then AMA_MltvMonoApprox() reports AMA_WARNING_ERROR messages which specify the bounds $\tilde\alpha_l$ and $\tilde\alpha_u$ imposed on the spline. Additionally, it sets AMA_OPTIONS::lwrbnd $=\alpha^s_l$ and AMA_OPTIONS::uprbnd $=\alpha^s_u$.

This function computes a spline which satisfies either the full or reduced continuity condition at the interior grid points, as well as, the interpolation and local monotonicity constraints. The following table summarizes the continuity of $s({\bf X})$ with respect to $x_k$ at the interior grid points $x^g_{k,l_k}$, for $l_k=2,\ldots,N^g_k-1$, and the number of coefficients, $m_k$, for a spline that satisfies either the full or reduced continuity condition. In the table the term

\[ s^{(p)}(x^g_{k,l_k}) = {\partial^p s(x_1,\cdots,x^g_{k,l_k},\cdots,x_n)\over \partial x^p_k} \]

is the $p$-th order partial of $s({\bf X})$ with respect to $x_k$ evaluated at $x^g_{k,l_k}$ and defined over the region $[x^L_1,x^U_1]\times\cdots\times[x^L_{k-1},x^U_{k-1}]\times[x^L_{k+1},x^U_{k+1}]\times\cdots\times[x^L_n,x^U_n]$.

Degree Full Continuity $m_k$ Reduced Continuity $m_k$
Linear $ s^{(p)}(x^g_{k,l_k})$ is continuous for $p=0$ $N^g_k$ $s^{(p)}(x^g_{k,l_k})$ is continuous for $p=0$ $N^g_k$
Quadratic $s^{(p)}(x^g_{k,l_k})$ is continuous for $0\le p\le 1$ $2N^g_k$ $s^{(p)}(x^g_{k,l_k})$ is continuous for $p=0$ $2N^g_k-1$
Cubic $s^{(p)}(x^g_{k,l_k})$ is continuous for $0\le p\le 2$ $3N^g_k$ $s^{(p)}(x^g_{k,l_k})$ is continuous for $0\le p\le 1$ $2N^g_k$
Quartic $s^{(p)}(x^g_{k,l_k})$ is continuous for $0\le p\le 3$ $4N^g_k$ $s^{(p)}(x^g_{k,l_k})$ is continuous for $0\le p\le 2$ $3N^g_k$
Quintic $s^{(p)}(x^g_{k,l_k})$ is continuous for $0\le p\le 4$ $5N^g_k$ $s^{(p)}(x^g_{k,l_k})$ is continuous for $0\le p\le 3$ $4N^g_k$

This function does the following:

  • Checks the input parameters for validity, see AMA_MltvInputCheck().
  • Defines the spline $s({\bf X})$ based on the knot vectors ${\bf\Lambda}_k$, for $1\le k\le n$.
  • Defines CNSPLA_CONPNT constraints to represent the interpolation constraints, see AMA_MltvConpnt().
  • Defines CNSPLA_PNLTRM condition to represent the function $F^{(2)}(s({\bf X}))$, see AMA_MltvPnltrm().
  • Defines CNSPLA_CONREG constraints to represent the local monotonicity constraints, see AMA_MltvConreg().
  • Invokes cnspla to compute a spline which minimizes $F^{(2)}(s({\bf X}))$ subject to the interpolation and local monotonicity constraints.
  • Stores interpolant 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().
By default the spline is unbounded but finite bounds can be set with AMA_OptionsSetBounds().
By default the spline satisfies the full continuity condition in all of its independent variables but the reduced continuity condition can be requested with AMA_OptionsSetContinuity().
By default the cross partial terms are included in the penalty term but they can be excluded with AMA_OptionsSetPenaltyTerm().
Because this function defines knot vectors based on the rectilinear grid upon which the independent variable data lies, it should only be used when the independent variable data either defines a gridded data with holes distribution, defines a gridded data distribution or consists of a small randomly distributed data set. See Multivariate Data Functions for a description of these data distributions.

Parameter Note: In the parameter definitions given below the limits on $k$ are $1\le k\le n$ and k = $k-1$.

Parameters
options[in] Pointer to AMA_OPTIONS. Must be initialized with AMA_Options() prior to calling AMA_MltvMonoInterp().
nind[in] The number of independent variables $n$. Must satisfy $2\le$ nind $\le$ AMA_MXNIND.
n[in] The number of data points $N$. Must satisfy n $\ge 2$.
x[in] Array of size nind containing arrays of size n where x[k] contains the independent variable data $x_{k,\ell}$, for $\ell=1,\ldots,N$.
z[in] Array of size n containing the dependent variable data $z_\ell$, for $\ell=1,\ldots,N$.
degree[in] Array of size nind containing the degree $d_k$ where degree[k] $= d_k$. Must satisfy $1\le$ degree[k] $\le 5$.
spline[out] Pointer to AMA_SPLINE pointer containing the monotonic spline interpolant. Must satisfy spline $\ne$ NULL.
Returns
Success/Error Code.

User Callable Function - Documented 022416

long int AMA_MltvPnltrm ( AMA_OPTIONS *  options,
double  theta,
long int  porder,
CNSPLA_SPLFUN *  splfun 
)

Define penalty term on a multivariate spline.

This function defines a penalty term $\theta F^{(p)}(s({\bf X}))$ for $p=0$ or $p=2$ on a multivariate spline $s({\bf X}):{\rm R}^n\rightarrow {\rm R}^1$ given as

\[ 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}) \]

where the $\alpha_{j_1,\ldots,j_n}$s, 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). \]

The $B_{d_k,j_k}(x_k|{\bf\Lambda}_k)$s are the $m_k$ univariate B-splines of degree $1 \le d_k\le 5$ defined by the knot vectors

\[ {\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. \]

For $p=0$ the penalty term $F^{(0)}(s({\bf X}))$ is

\[ F^{(0)}(s({\bf X})) = \int_R\left(s({\bf X})\right)^2dR \]

where $R = [\lambda_{1,1},\lambda_{1,m_1+1}] \times [\lambda_{2,1},\lambda_{2,m_2+1}] \times\cdots\times [\lambda_{n,1},\lambda_{n,m_n+1}]$. Alternately, for $p=2$ the penalty term $F^{(2)}(s({\bf X}))$ is

\[ 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 + 2.0\sum_{l=k+1}^n\int_R\left(\partial^2 s({\bf X})\over\partial x_k \partial x_l\right)^2dR\right]. \]

The cross partial terms can be excluded from the penalty term $F^{(2)}(s({\bf X}))$ with AMA_OptionsSetPenaltyTerm().

Parameters
options[in] Pointer to AMA_OPTIONS. Should be initialized with AMA_Options() prior to calling AMA_MltvPnltrm().
theta[in] The penalty term weight $\theta$. Should satisfy $0.0\le$ theta $\le 1.0$.
porder[in] The penalty term order $p$. Should satisfy porder $=0$ or porder $=2$.
splfun[in] Pointer to CNSPLA_SPLFUN containing spline $s({\bf X})$ upon which penalty term is imposed. Should satisfy splfun $\ne$ NULL.
Returns
Success/Error Code.

User Support Function - Documented 121114 - !!!THIS IS NOT A USER CALLABLE FUNCTION - DOCUMENT IS INCLUDED FOR COMPLETENESS!!!

long int AMA_MltvRglrze ( AMA_OPTIONS *  options,
CNSPLA_SPLFUN *  splfun 
)
AMA Spline Library Document - Version:1.0.0-a Date:02/24/16 - Generated by   doxygen 1.8.5