Source code for naginterfaces.library.interp

# -*- coding: utf-8 -*-
r"""
Module Summary
--------------
Interfaces for the NAG Mark 28.5 `interp` Chapter.

``interp`` - Interpolation

This module is concerned with the interpolation of a function of one or more variables.
When provided with the value of the function (and possibly one or more of its lowest-order derivatives) at each of a number of values of the variable(s), the NAG Library functions provide either an interpolating function or an interpolated value.
For some of the interpolating functions, there are supporting NAG Library functions to evaluate, differentiate or integrate them.

Functionality Index
-------------------

**Derivative**

  of interpolant

    from :meth:`dim1_monotonic`: :meth:`dim1_monotonic_deriv`

    from :meth:`dim2_scat_shep`: :meth:`dim2_scat_shep_eval`

    from :meth:`dim3_scat_shep`: :meth:`dim3_scat_shep_eval`

    from :meth:`dim4_scat_shep`: :meth:`dim4_scat_shep_eval`

    from :meth:`dim5_scat_shep`: :meth:`dim5_scat_shep_eval`

    from :meth:`dimn_scat_shep`: :meth:`dimn_scat_shep_eval`

**Evaluation**

  of interpolant

    from :meth:`dim1_monotonic`: :meth:`dim1_monotonic_eval`

    from :meth:`dim1_ratnl`: :meth:`dim1_ratnl_eval`

    from :meth:`dim2_scat`: :meth:`dim2_scat_eval`

    from :meth:`dim2_scat_shep`: :meth:`dim2_scat_shep_eval`

    from :meth:`dim3_scat_shep`: :meth:`dim3_scat_shep_eval`

    from :meth:`dim4_scat_shep`: :meth:`dim4_scat_shep_eval`

    from :meth:`dim5_scat_shep`: :meth:`dim5_scat_shep_eval`

    from :meth:`dimn_scat_shep`: :meth:`dimn_scat_shep_eval`

    from triangulation from :meth:`dim2_triangulate`: :meth:`dim2_triang_bary_eval`

    using variables computed by :meth:`dim1_monconv_disc`: :meth:`dim1_monconv_eval`

**Extrapolation**

  one variable

    monotonic convex: :meth:`dim1_monconv_disc`

    piecewise cubic: :meth:`dim1_monotonic`

    polynomial

      data with or without derivatives: :meth:`dim1_cheb`

      general data: :meth:`dim1_aitken`

    rational function: :meth:`dim1_ratnl`

Integration (definite) of interpolant from :meth:`dim1_monotonic`: :meth:`dim1_monotonic_intg`

**Interpolated values**

  :math:`d` variables

    from interpolant from :meth:`dimn_scat_shep`: :meth:`dimn_scat_shep_eval`

    modified Shepard method, Linear or Cubic: :meth:`dimn_grid`

  five variables

    from interpolant from :meth:`dim5_scat_shep`: :meth:`dim5_scat_shep_eval`

  four variables

    from interpolant from :meth:`dim4_scat_shep`: :meth:`dim4_scat_shep_eval`

  one variable

    from interpolant from :meth:`dim1_monotonic`: :meth:`dim1_monotonic_eval`

    from interpolant from :meth:`dim1_monotonic` (including derivative): :meth:`dim1_monotonic_deriv`

    from polynomial

      equally spaced data: :meth:`dim1_everett`

      general data: :meth:`dim1_aitken`

    from rational function: :meth:`dim1_ratnl_eval`

    using variables computed by :meth:`dim1_monconv_disc`: :meth:`dim1_monconv_eval`

  three variables

    from interpolant from :meth:`dim3_scat_shep`: :meth:`dim3_scat_shep_eval`

  two variables

    barycentric, from triangulation from :meth:`dim2_triangulate`: :meth:`dim2_triang_bary_eval`

    from interpolant from :meth:`dim2_scat`: :meth:`dim2_scat_eval`

    from interpolant from :meth:`dim2_scat_shep`: :meth:`dim2_scat_shep_eval`

**Interpolating function**

  :math:`d` variables

    modified Shepard method: :meth:`dimn_scat_shep`

  five variables

    modified Shepard method: :meth:`dim5_scat_shep`

  four variables

    modified Shepard method: :meth:`dim4_scat_shep`

  one variable

    cubic spline: :meth:`dim1_spline`

    monotonic convex piecewise polynomial: :meth:`dim1_monconv_disc`

    other piecewise polynomial: :meth:`dim1_monotonic`

    polynomial

      data with or without derivatives: :meth:`dim1_cheb`

    rational function: :meth:`dim1_ratnl`

  three variables

    modified Shepard method: :meth:`dim3_scat_shep`

  two variables

    bicubic spline: :meth:`dim2_spline_grid`

    modified Shepard method: :meth:`dim2_scat_shep`

    other piecewise polynomial: :meth:`dim2_scat`

    triangulation: :meth:`dim2_triangulate`

For full information please refer to the NAG Library document

https://www.nag.com/numeric/nl/nagdoc_28.5/flhtml/e01/e01intro.html
"""

# NAG Copyright 2017-2022.

[docs]def dim1_aitken(a, b, x): r""" ``dim1_aitken`` interpolates a function of one variable at a given point :math:`x` from a table of function values :math:`y_i` evaluated at equidistant or non-equidistant points :math:`x_i`, for :math:`\textit{i} = 1,2,\ldots,n+1`, using Aitken's technique of successive linear interpolations. .. _e01aa-py2-py-doc: For full information please refer to the NAG Library document for e01aa https://www.nag.com/numeric/nl/nagdoc_28.5/flhtml/e01/e01aaf.html .. _e01aa-py2-py-parameters: **Parameters** **a** : float, array-like, shape :math:`\left(n+1\right)` :math:`\mathrm{a}[\textit{i}-1]` must contain the :math:`x`-component of the :math:`\textit{i}`\ th data point, :math:`x_{\textit{i}}`, for :math:`\textit{i} = 1,2,\ldots,n+1`. **b** : float, array-like, shape :math:`\left(n+1\right)` :math:`\mathrm{b}[\textit{i}-1]` must contain the :math:`y`-component (function value) of the :math:`\textit{i}`\ th data point, :math:`y_{\textit{i}}`, for :math:`\textit{i} = 1,2,\ldots,n+1`. **x** : float The point :math:`x` at which the interpolation is required. Note that :math:`\mathrm{x}` may lie outside the interval defined by the minimum and maximum values in :math:`\mathrm{a}`, in which case an extrapolated value will be computed; extrapolated results should be treated with considerable caution since there is no information on the behaviour of the function outside the defined interval. **Returns** **a** : float, ndarray, shape :math:`\left(n+1\right)` :math:`\mathrm{a}[\textit{i}-1]` contains the value :math:`x_{\textit{i}}-x`, for :math:`\textit{i} = 1,2,\ldots,n+1`. **b** : float, ndarray, shape :math:`\left(n+1\right)` The contents of :math:`\mathrm{b}` are unspecified. **c** : float, ndarray, shape :math:`\left(n\times \left(n+1\right)/2\right)` :math:`\mathrm{c}[0],\ldots,\mathrm{c}[n-1]` contain the first set of linear interpolations, :math:`\mathrm{c}[n],\ldots,\mathrm{c}[2\times n-2]` contain the second set of linear interpolations, :math:`\mathrm{c}[2n-1],\ldots,\mathrm{c}[3\times n-4]` contain the third set of linear interpolations, :math:`\vdots` :math:`\mathrm{c}[{n\times \left(n+1\right)/2}-1]` contains the interpolated function value at the point :math:`x`. .. _e01aa-py2-py-errors: **Raises** **NagValueError** (`errno` :math:`6`) On entry, error in parameter :math:`\textit{n}`. Constraint: :math:`n > 0`. .. _e01aa-py2-py-notes: **Notes** ``dim1_aitken`` interpolates a function of one variable at a given point :math:`x` from a table of values :math:`x_i` and :math:`y_i`, for :math:`i = 1,2,\ldots,n+1` using Aitken's method (see Fröberg (1970)). The intermediate values of linear interpolations are stored to enable an estimate of the accuracy of the results to be made. .. _e01aa-py2-py-references: **References** Fröberg, C E, 1970, `Introduction to Numerical Analysis`, Addison--Wesley """ raise NotImplementedError
[docs]def dim1_everett(p, a): r""" ``dim1_everett`` interpolates a function of one variable at a given point :math:`x` from a table of function values evaluated at equidistant points, using Everett's formula. .. _e01ab-py2-py-doc: For full information please refer to the NAG Library document for e01ab https://www.nag.com/numeric/nl/nagdoc_28.5/flhtml/e01/e01abf.html .. _e01ab-py2-py-parameters: **Parameters** **p** : float The point :math:`p` at which the interpolated function value is required, i.e., :math:`p = \left(x-x_0\right)/h` with :math:`{-1.0} < p < 1.0`. **a** : float, array-like, shape :math:`\left(2\times n\right)` :math:`\mathrm{a}[\textit{i}-1]` must be set to the function value :math:`y_{{\textit{i}-n}}`, for :math:`\textit{i} = 1,2,\ldots,2n`. **Returns** **a** : float, ndarray, shape :math:`\left(2\times n\right)` The contents of :math:`\mathrm{a}` are unspecified. **g** : float, ndarray, shape :math:`\left(2\times n+1\right)` The array contains .. rst-class:: nag-rules-none nag-align-left +------------------------+------------------------------------------------------------------------------+ |:math:`y_0` |in :math:`\mathrm{g}[0]` | +------------------------+------------------------------------------------------------------------------+ |:math:`y_1` |in :math:`\mathrm{g}[1]` | +------------------------+------------------------------------------------------------------------------+ |:math:`\delta^{{2r}}y_0`|in :math:`\mathrm{g}[2r]` | +------------------------+------------------------------------------------------------------------------+ |:math:`\delta^{{2r}}y_1`|in :math:`\mathrm{g}[2\textit{r}+1]`, for :math:`\textit{r} = 1,2,\ldots,n-1`.| +------------------------+------------------------------------------------------------------------------+ The interpolated function value :math:`y_p` is stored in :math:`\mathrm{g}[2n]`. .. _e01ab-py2-py-errors: **Raises** **NagValueError** (`errno` :math:`1`) On entry, :math:`\mathrm{p} = \langle\mathit{\boldsymbol{value}}\rangle`. Constraint: :math:`\mathrm{p}\leq 1.0`. (`errno` :math:`1`) On entry, :math:`\mathrm{p} = \langle\mathit{\boldsymbol{value}}\rangle`. Constraint: :math:`\mathrm{p}\geq {-1.0}`. (`errno` :math:`2`) On entry, :math:`n = \langle\mathit{\boldsymbol{value}}\rangle`. Constraint: :math:`n > 0`. .. _e01ab-py2-py-notes: **Notes** ``dim1_everett`` interpolates a function of one variable at a given point .. math:: x = x_0+ph\text{,} where :math:`-1\leq p\leq 1` and :math:`h` is the interval of differencing, from a table of values :math:`x_m = x_0+mh` and :math:`y_m` where :math:`m = -\left(n-1\right),-\left(n-2\right),\ldots,-1,0,1,\ldots,n`. The formula used is that of Fröberg (1970), neglecting the remainder term: .. math:: y_p = \sum_{{r = 0}}^{{n-1}}\left(\frac{{1-p+r}}{{2r+1}}\right)\delta^{{2r}}y_0+\sum_{{r = 0}}^{{n-1}}\left(\frac{{p+r}}{{2r+1}}\right)\delta^{{2r}}y_1\text{.} The values of :math:`\delta^{{2r}}y_0` and :math:`\delta^{{2r}}y_1` are stored on exit from the function in addition to the interpolated function value :math:`y_p`. .. _e01ab-py2-py-references: **References** Fröberg, C E, 1970, `Introduction to Numerical Analysis`, Addison--Wesley """ raise NotImplementedError
[docs]def dim1_cheb(xmin, xmax, x, y, ip, lwrk, liwrk, itmin=0, itmax=0): r""" ``dim1_cheb`` constructs the Chebyshev series representation of a polynomial interpolant to a set of data which may contain derivative values. .. _e01ae-py2-py-doc: For full information please refer to the NAG Library document for e01ae https://www.nag.com/numeric/nl/nagdoc_28.5/flhtml/e01/e01aef.html .. _e01ae-py2-py-parameters: **Parameters** **xmin** : float The lower and upper end points, respectively, of the interval :math:`\left[x_{\mathrm{min}}, x_{\mathrm{max}}\right]`. If they are not determined by your problem, it is recommended that they be set respectively to the smallest and largest values among the :math:`x_i`. **xmax** : float The lower and upper end points, respectively, of the interval :math:`\left[x_{\mathrm{min}}, x_{\mathrm{max}}\right]`. If they are not determined by your problem, it is recommended that they be set respectively to the smallest and largest values among the :math:`x_i`. **x** : float, array-like, shape :math:`\left(m\right)` The value of :math:`x_{\textit{i}}`, for :math:`\textit{i} = 1,2,\ldots,m`. The :math:`\mathrm{x}[i-1]` need not be ordered. **y** : float, array-like, shape :math:`\left(n\right)` The given values of the dependent variable, and derivatives, as follows: The first :math:`p_1+1` elements contain :math:`y_1,y_1^{\left(1\right)},\ldots,y_1^{\left(p_1\right)}` in that order. The next :math:`p_2+1` elements contain :math:`y_2,y_2^{\left(1\right)},\ldots,y_2^{\left(p_2\right)}` in that order. :math:`\quad \text{ }\quad \vdots` The last :math:`p_m+1` elements contain :math:`y_m,y_m^{\left(1\right)},\ldots,y_m^{\left(p_m\right)}` in that order. **ip** : int, array-like, shape :math:`\left(m\right)` :math:`p_{\textit{i}}`, the order of the highest-order derivative whose value is given at :math:`x_{\textit{i}}`, for :math:`\textit{i} = 1,2,\ldots,m`. If the value of :math:`y` only is given for some :math:`x_i` then the corresponding value of :math:`\mathrm{ip}[i-1]` must be zero. **lwrk** : int The dimension of the array :math:`\mathrm{wrk}`. **liwrk** : int The dimension of the array :math:`\mathrm{iwrk}`. **itmin** : int, optional Respectively the minimum and maximum number of iterations to be performed by the function (for full details see `Further Comments <https://www.nag.com/numeric/nl/nagdoc_28.5/flhtml/e01/e01aef.html#fcomments4>`__). Setting :math:`\mathrm{itmin}` and/or :math:`\mathrm{itmax}` negative or zero invokes default value(s) of :math:`2` and/or :math:`10`, respectively. The default values will be satisfactory for most problems, but occasionally significant improvement will result from using higher values. **itmax** : int, optional Respectively the minimum and maximum number of iterations to be performed by the function (for full details see `Further Comments <https://www.nag.com/numeric/nl/nagdoc_28.5/flhtml/e01/e01aef.html#fcomments4>`__). Setting :math:`\mathrm{itmin}` and/or :math:`\mathrm{itmax}` negative or zero invokes default value(s) of :math:`2` and/or :math:`10`, respectively. The default values will be satisfactory for most problems, but occasionally significant improvement will result from using higher values. **Returns** **a** : float, ndarray, shape :math:`\left(n\right)` :math:`\mathrm{a}[\textit{i}]` contains the coefficient :math:`a_{\textit{i}}` in the Chebyshev series representation of :math:`q\left(x\right)`, for :math:`\textit{i} = 0,\ldots,n-1`. **wrk** : float, ndarray, shape :math:`\left(\mathrm{lwrk}\right)` Used as workspace, but see also `Further Comments <https://www.nag.com/numeric/nl/nagdoc_28.5/flhtml/e01/e01aef.html#fcomments5>`__. **iwrk** : int, ndarray, shape :math:`\left(\mathrm{liwrk}\right)` Used as workspace, but see also `Further Comments <https://www.nag.com/numeric/nl/nagdoc_28.5/flhtml/e01/e01aef.html#fcomments5>`__. .. _e01ae-py2-py-errors: **Raises** **NagValueError** (`errno` :math:`1`) On entry, :math:`\mathrm{liwrk}` is too small. :math:`\mathrm{liwrk} = \langle\mathit{\boldsymbol{value}}\rangle`. Minimum possible dimension: :math:`\langle\mathit{\boldsymbol{value}}\rangle`. (`errno` :math:`1`) On entry, :math:`\mathrm{lwrk}` is too small. :math:`\mathrm{lwrk} = \langle\mathit{\boldsymbol{value}}\rangle`. Minimum possible dimension: :math:`\langle\mathit{\boldsymbol{value}}\rangle`. (`errno` :math:`1`) On entry, :math:`n = \langle\mathit{\boldsymbol{value}}\rangle` and :math:`m+\mathrm{ip}[0]+\mathrm{ip}[1] + \cdots +\mathrm{ip}[m-1] = \langle\mathit{\boldsymbol{value}}\rangle`. Constraint: :math:`n = m+\mathrm{ip}[0]+\mathrm{ip}[1] + \cdots +\mathrm{ip}[m-1]`. (`errno` :math:`1`) On entry, :math:`m = \langle\mathit{\boldsymbol{value}}\rangle`. Constraint: :math:`m\geq 1`. (`errno` :math:`2`) On entry, :math:`\textit{I} = \langle\mathit{\boldsymbol{value}}\rangle` and :math:`\mathrm{ip}[\textit{I}-1] = \langle\mathit{\boldsymbol{value}}\rangle`. Constraint: :math:`\mathrm{ip}[\textit{I}-1]\geq 0`. (`errno` :math:`3`) On entry, :math:`\textit{I} = \langle\mathit{\boldsymbol{value}}\rangle`, :math:`\textit{J} = \langle\mathit{\boldsymbol{value}}\rangle` and :math:`\mathrm{x}[\textit{I}-1] = \langle\mathit{\boldsymbol{value}}\rangle`. Constraint: :math:`\mathrm{x}[\textit{I}-1]\neq \mathrm{x}[\textit{J}-1]`. (`errno` :math:`3`) On entry, :math:`\textit{I} = \langle\mathit{\boldsymbol{value}}\rangle`, :math:`\mathrm{x}[\textit{I}-1] = \langle\mathit{\boldsymbol{value}}\rangle`, :math:`\mathrm{xmin} = \langle\mathit{\boldsymbol{value}}\rangle` and :math:`\mathrm{xmax} = \langle\mathit{\boldsymbol{value}}\rangle`. Constraint: :math:`\mathrm{xmin}\leq \mathrm{x}[\textit{I}-1]\leq \mathrm{xmax}`. (`errno` :math:`3`) On entry, :math:`\mathrm{xmin} = \langle\mathit{\boldsymbol{value}}\rangle` and :math:`\mathrm{xmax} = \langle\mathit{\boldsymbol{value}}\rangle`. Constraint: :math:`\mathrm{xmin} < \mathrm{xmax}`. (`errno` :math:`5`) The computation has been terminated because the iterative process appears to be diverging. **Warns** **NagAlgorithmicWarning** (`errno` :math:`4`) Not all the performance indices are less than eight times the machine precision, although :math:`\mathrm{itmax}` iterations have been performed. A more accurate solution may possibly be obtained by increasing :math:`\mathrm{itmax}` and recalling the function. .. _e01ae-py2-py-notes: **Notes** Let :math:`m` distinct values :math:`x_{\textit{i}}` of an independent variable :math:`x` be given, with :math:`x_{\mathrm{min}}\leq x_{\textit{i}}\leq x_{\mathrm{max}}`, for :math:`\textit{i} = 1,2,\ldots,m`. For each value :math:`x_i`, suppose that the value :math:`y_i` of the dependent variable :math:`y` together with the first :math:`p_i` derivatives of :math:`y` with respect to :math:`x` are given. Each :math:`p_i` must, therefore, be a non-negative integer, with the total number of interpolating conditions, :math:`n`, equal to :math:`m+\sum_{{i = 1}}^mp_i`. ``dim1_cheb`` calculates the unique polynomial :math:`q\left(x\right)` of degree :math:`n-1` (or less) which is such that :math:`q^{\left(\textit{k}\right)}\left(x_{\textit{i}}\right) = y_{\textit{i}}^{\left(\textit{k}\right)}`, for :math:`\textit{k} = 0,1,\ldots,p_{\textit{i}}`, for :math:`\textit{i} = 1,2,\ldots,m`. Here :math:`q^{\left(0\right)}\left(x_i\right)` means :math:`q\left(x_i\right)`. This polynomial is represented in Chebyshev series form in the normalized variable :math:`\bar{x}`, as follows: .. math:: q\left(x\right) = \frac{1}{2}a_0T_0\left(\bar{x}\right)+a_1T_1\left(\bar{x}\right)+ \cdots +a_{{n-1}}T_{{n-1}}\left(\bar{x}\right)\text{,} where .. math:: \bar{x} = \frac{{2x-x_{\mathrm{min}}-x_{\mathrm{max}}}}{{x_{\mathrm{max}}-x_{\mathrm{min}}}} so that :math:`-1\leq \bar{x}\leq 1` for :math:`x` in the interval :math:`x_{\mathrm{min}}` to :math:`x_{\mathrm{max}}`, and where :math:`T_i\left(\bar{x}\right)` is the Chebyshev polynomial of the first kind of degree :math:`i` with argument :math:`\bar{x}`. (The polynomial interpolant can subsequently be evaluated for any value of :math:`x` in the given range by using :meth:`fit.dim1_cheb_eval2 <naginterfaces.library.fit.dim1_cheb_eval2>`. Chebyshev series representations of the derivative(s) and integral(s) of :math:`q\left(x\right)` may be obtained by (repeated) use of :meth:`fit.dim1_cheb_deriv <naginterfaces.library.fit.dim1_cheb_deriv>` and :meth:`fit.dim1_cheb_integ <naginterfaces.library.fit.dim1_cheb_integ>`.) The method used consists first of constructing a divided-difference table from the normalized :math:`\bar{x}` values and the given values of :math:`y` and its derivatives with respect to :math:`\bar{x}`. The Newton form of :math:`q\left(x\right)` is then obtained from this table, as described in Huddleston (1974) and Krogh (1970), with the modification described in `Further Comments <https://www.nag.com/numeric/nl/nagdoc_28.5/flhtml/e01/e01aef.html#fcomments2>`__. The Newton form of the polynomial is then converted to Chebyshev series form as described in `Conversion to Chebyshev Form <https://www.nag.com/numeric/nl/nagdoc_28.5/flhtml/e01/e01aef.html#fcomments3>`__. Since the errors incurred by these stages can be considerable, a form of iterative refinement is used to improve the solution. This refinement is particularly useful when derivatives of rather high order are given in the data. In reasonable examples, the refinement will usually terminate with a certain accuracy criterion satisfied by the polynomial (see `Accuracy <https://www.nag.com/numeric/nl/nagdoc_28.5/flhtml/e01/e01aef.html#accuracy>`__). In more difficult examples, the criterion may not be satisfied and refinement will continue until the maximum number of iterations (as specified by the input argument :math:`\mathrm{itmax}`) is reached. In extreme examples, the iterative process may diverge (even though the accuracy criterion is satisfied): if a certain divergence criterion is satisfied, the process terminates at once. In all cases the function returns the 'best' polynomial achieved before termination. For the definition of 'best' and details of iterative refinement and termination criteria, see `Further Comments <https://www.nag.com/numeric/nl/nagdoc_28.5/flhtml/e01/e01aef.html#fcomments4>`__. .. _e01ae-py2-py-references: **References** Huddleston, R E, 1974, `CDC 6600 routines for the interpolation of data and of data with derivatives`, SLL-74-0214, Sandia Laboratories (Reprint) Krogh, F T, 1970, `Efficient algorithms for polynomial interpolation and numerical differentiation`, Math. Comput. (24), 185--190 """ raise NotImplementedError
[docs]def dim1_spline(x, y): r""" ``dim1_spline`` determines a cubic spline interpolant to a given set of data. .. _e01ba-py2-py-doc: For full information please refer to the NAG Library document for e01ba https://www.nag.com/numeric/nl/nagdoc_28.5/flhtml/e01/e01baf.html .. _e01ba-py2-py-parameters: **Parameters** **x** : float, array-like, shape :math:`\left(m\right)` :math:`\mathrm{x}[\textit{i}-1]` must be set to :math:`x_{\textit{i}}`, the :math:`\textit{i}`\ th data value of the independent variable :math:`x`, for :math:`\textit{i} = 1,2,\ldots,m`. **y** : float, array-like, shape :math:`\left(m\right)` :math:`\mathrm{y}[\textit{i}-1]` must be set to :math:`y_{\textit{i}}`, the :math:`\textit{i}`\ th data value of the dependent variable :math:`y`, for :math:`\textit{i} = 1,2,\ldots,m`. **Returns** **lamda** : float, ndarray, shape :math:`\left(m+4\right)` The value of :math:`\lambda_{\textit{i}}`, the :math:`\textit{i}`\ th knot, for :math:`\textit{i} = 1,2,\ldots,m+4`. **c** : float, ndarray, shape :math:`\left(m\right)` The coefficient :math:`c_{\textit{i}}` of the B-spline :math:`N_{\textit{i}}\left(x\right)`, for :math:`\textit{i} = 1,2,\ldots,m`. The remaining elements of the array are not used. .. _e01ba-py2-py-errors: **Raises** **NagValueError** (`errno` :math:`1`) On entry, :math:`\textit{lck} = \langle\mathit{\boldsymbol{value}}\rangle` and :math:`m = \langle\mathit{\boldsymbol{value}}\rangle`. Constraint: :math:`\textit{lck}\geq m+4`. (`errno` :math:`1`) On entry, :math:`m = \langle\mathit{\boldsymbol{value}}\rangle`. Constraint: :math:`m\geq 4`. (`errno` :math:`2`) On entry, :math:`\textit{I} = \langle\mathit{\boldsymbol{value}}\rangle`, :math:`\mathrm{x}[\textit{I}-1] = \langle\mathit{\boldsymbol{value}}\rangle` and :math:`\mathrm{x}[\textit{I}-2] = \langle\mathit{\boldsymbol{value}}\rangle`. Constraint: :math:`\mathrm{x}[\textit{I}-1] > \mathrm{x}[\textit{I}-2]`. .. _e01ba-py2-py-notes: **Notes** `In the NAG Library the traditional C interface for this routine uses a different algorithmic base. Please contact NAG if you have any questions about compatibility.` ``dim1_spline`` determines a cubic spline :math:`s\left(x\right)`, defined in the range :math:`x_1\leq x\leq x_m`, which interpolates (passes exactly through) the set of data points :math:`\left(x_{\textit{i}}, y_{\textit{i}}\right)`, for :math:`\textit{i} = 1,2,\ldots,m`, where :math:`m\geq 4` and :math:`x_1 < x_2 < \cdots < x_m`. Unlike some other spline interpolation algorithms, derivative end conditions are not imposed. The spline interpolant chosen has :math:`m-4` interior knots :math:`\lambda_5,\lambda_6,\ldots,\lambda_m`, which are set to the values of :math:`x_3,x_4,\ldots,x_{{m-2}}` respectively. This spline is represented in its B-spline form (see Cox (1975)): .. math:: s\left(x\right) = \sum_{{i = 1}}^mc_iN_i\left(x\right)\text{,} where :math:`N_i\left(x\right)` denotes the normalized B-spline of degree :math:`3`, defined upon the knots :math:`\lambda_i,\lambda_{{i+1}},\ldots,\lambda_{{i+4}}`, and :math:`c_i` denotes its coefficient, whose value is to be determined by the function. The use of B-splines requires eight additional knots :math:`\lambda_1`, :math:`\lambda_2`, :math:`\lambda_3`, :math:`\lambda_4`, :math:`\lambda_{{m+1}}`, :math:`\lambda_{{m+2}}`, :math:`\lambda_{{m+3}}` and :math:`\lambda_{{m+4}}` to be specified; ``dim1_spline`` sets the first four of these to :math:`x_1` and the last four to :math:`x_m`. The algorithm for determining the coefficients is as described in Cox (1975) except that :math:`QR` factorization is used instead of :math:`LU` decomposition. The implementation of the algorithm involves setting up appropriate information for the related function :meth:`fit.dim1_spline_knots <naginterfaces.library.fit.dim1_spline_knots>` followed by a call of that function. (See :meth:`fit.dim1_spline_knots <naginterfaces.library.fit.dim1_spline_knots>` for further details.) Values of the spline interpolant, or of its derivatives or definite integral, can subsequently be computed as detailed in `Further Comments <https://www.nag.com/numeric/nl/nagdoc_28.5/flhtml/e01/e01baf.html#fcomments>`__. .. _e01ba-py2-py-references: **References** Cox, M G, 1975, `An algorithm for spline interpolation`, J. Inst. Math. Appl. (15), 95--108 Cox, M G, 1977, `A survey of numerical methods for data and function approximation`, The State of the Art in Numerical Analysis, (ed D A H Jacobs), 627--668, Academic Press """ raise NotImplementedError
[docs]def dim1_monotonic(x, f): r""" ``dim1_monotonic`` computes a monotonicity-preserving piecewise cubic Hermite interpolant to a set of data points. .. _e01be-py2-py-doc: For full information please refer to the NAG Library document for e01be https://www.nag.com/numeric/nl/nagdoc_28.5/flhtml/e01/e01bef.html .. _e01be-py2-py-parameters: **Parameters** **x** : float, array-like, shape :math:`\left(n\right)` :math:`\mathrm{x}[\textit{r}-1]` must be set to :math:`x_{\textit{r}}`, the :math:`\textit{r}`\ th value of the independent variable (abscissa), for :math:`\textit{r} = 1,2,\ldots,n`. **f** : float, array-like, shape :math:`\left(n\right)` :math:`\mathrm{f}[\textit{r}-1]` must be set to :math:`f_{\textit{r}}`, the :math:`\textit{r}`\ th value of the dependent variable (ordinate), for :math:`\textit{r} = 1,2,\ldots,n`. **Returns** **d** : float, ndarray, shape :math:`\left(n\right)` Estimates of derivatives at the data points. :math:`\mathrm{d}[r-1]` contains the derivative at :math:`\mathrm{x}[r-1]`. .. _e01be-py2-py-errors: **Raises** **NagValueError** (`errno` :math:`1`) On entry, :math:`n = \langle\mathit{\boldsymbol{value}}\rangle`. Constraint: :math:`n\geq 2`. (`errno` :math:`2`) On entry, :math:`r = \langle\mathit{\boldsymbol{value}}\rangle`, :math:`\mathrm{x}[r-2] = \langle\mathit{\boldsymbol{value}}\rangle` and :math:`\mathrm{x}[r-1] = \langle\mathit{\boldsymbol{value}}\rangle`. Constraint: :math:`\mathrm{x}[r-2] < \mathrm{x}[r-1]` for all :math:`r`. .. _e01be-py2-py-notes: **Notes** `In the NAG Library the traditional C interface for this routine uses a different algorithmic base. Please contact NAG if you have any questions about compatibility.` ``dim1_monotonic`` estimates first derivatives at the set of data points :math:`\left(x_{\textit{r}}, f_{\textit{r}}\right)`, for :math:`\textit{r} = 1,2,\ldots,n`, which determine a piecewise cubic Hermite interpolant to the data, that preserves monotonicity over ranges where the data points are monotonic. If the data points are only piecewise monotonic, the interpolant will have an extremum at each point where monotonicity switches direction. The estimates of the derivatives are computed by a formula due to Brodlie, which is described in Fritsch and Butland (1984), with suitable changes at the boundary points. The function is derived from function PCHIM in Fritsch (1982). Values of the computed interpolant, and of its first derivative and definite integral, can subsequently be computed by calling :meth:`dim1_monotonic_eval`, :meth:`dim1_monotonic_deriv` and :meth:`dim1_monotonic_intg`, as described in `Further Comments <https://www.nag.com/numeric/nl/nagdoc_28.5/flhtml/e01/e01bef.html#fcomments>`__. .. _e01be-py2-py-references: **References** Fritsch, F N, 1982, `PCHIP final specifications`, Report UCID-30194, Lawrence Livermore National Laboratory Fritsch, F N and Butland, J, 1984, `A method for constructing local monotone piecewise cubic interpolants`, SIAM J. Sci. Statist. Comput. (5), 300--304 """ raise NotImplementedError
[docs]def dim1_monotonic_eval(x, f, d, px): r""" ``dim1_monotonic_eval`` evaluates a piecewise cubic Hermite interpolant at a set of points. .. _e01bf-py2-py-doc: For full information please refer to the NAG Library document for e01bf https://www.nag.com/numeric/nl/nagdoc_28.5/flhtml/e01/e01bff.html .. _e01bf-py2-py-parameters: **Parameters** **x** : float, array-like, shape :math:`\left(n\right)` :math:`\textit{n}`, :math:`\mathrm{x}`, :math:`\mathrm{f}` and :math:`\mathrm{d}` must be unchanged from the previous call of :meth:`dim1_monotonic`. **f** : float, array-like, shape :math:`\left(n\right)` :math:`\textit{n}`, :math:`\mathrm{x}`, :math:`\mathrm{f}` and :math:`\mathrm{d}` must be unchanged from the previous call of :meth:`dim1_monotonic`. **d** : float, array-like, shape :math:`\left(n\right)` :math:`\textit{n}`, :math:`\mathrm{x}`, :math:`\mathrm{f}` and :math:`\mathrm{d}` must be unchanged from the previous call of :meth:`dim1_monotonic`. **px** : float, array-like, shape :math:`\left(m\right)` The :math:`m` values of :math:`x` at which the interpolant is to be evaluated. **Returns** **pf** : float, ndarray, shape :math:`\left(m\right)` :math:`\mathrm{pf}[\textit{i}-1]` contains the value of the interpolant evaluated at the point :math:`\mathrm{px}[\textit{i}-1]`, for :math:`\textit{i} = 1,2,\ldots,m`. .. _e01bf-py2-py-errors: **Raises** **NagValueError** (`errno` :math:`1`) On entry, :math:`n = \langle\mathit{\boldsymbol{value}}\rangle`. Constraint: :math:`n\geq 2`. (`errno` :math:`2`) On entry, :math:`r = \langle\mathit{\boldsymbol{value}}\rangle`, :math:`\mathrm{x}[r-2] = \langle\mathit{\boldsymbol{value}}\rangle` and :math:`\mathrm{x}[r-1] = \langle\mathit{\boldsymbol{value}}\rangle`. Constraint: :math:`\mathrm{x}[r-2] < \mathrm{x}[r-1]` for all :math:`r`. (`errno` :math:`3`) On entry, :math:`m = \langle\mathit{\boldsymbol{value}}\rangle`. Constraint: :math:`m\geq 1`. **Warns** **NagAlgorithmicWarning** (`errno` :math:`4`) Warning -- some points in array :math:`\mathrm{px}` lie outside the range :math:`\mathrm{x}[0] \cdots \mathrm{x}[n-1]`. Values at these points are unreliable because computed by extrapolation. .. _e01bf-py2-py-notes: **Notes** `In the NAG Library the traditional C interface for this routine uses a different algorithmic base. Please contact NAG if you have any questions about compatibility.` ``dim1_monotonic_eval`` evaluates a piecewise cubic Hermite interpolant, as computed by :meth:`dim1_monotonic`, at the points :math:`\mathrm{px}[\textit{i}-1]`, for :math:`\textit{i} = 1,2,\ldots,m`. If any point lies outside the interval from :math:`\mathrm{x}[0]` to :math:`\mathrm{x}[n-1]`, a value is extrapolated from the nearest extreme cubic, and a warning is returned. The function is derived from function PCHFE in Fritsch (1982). .. _e01bf-py2-py-references: **References** Fritsch, F N, 1982, `PCHIP final specifications`, Report UCID-30194, Lawrence Livermore National Laboratory """ raise NotImplementedError
[docs]def dim1_monotonic_deriv(x, f, d, px): r""" ``dim1_monotonic_deriv`` evaluates a piecewise cubic Hermite interpolant and its first derivative at a set of points. .. _e01bg-py2-py-doc: For full information please refer to the NAG Library document for e01bg https://www.nag.com/numeric/nl/nagdoc_28.5/flhtml/e01/e01bgf.html .. _e01bg-py2-py-parameters: **Parameters** **x** : float, array-like, shape :math:`\left(n\right)` :math:`\textit{n}`, :math:`\mathrm{x}`, :math:`\mathrm{f}` and :math:`\mathrm{d}` must be unchanged from the previous call of :meth:`dim1_monotonic`. **f** : float, array-like, shape :math:`\left(n\right)` :math:`\textit{n}`, :math:`\mathrm{x}`, :math:`\mathrm{f}` and :math:`\mathrm{d}` must be unchanged from the previous call of :meth:`dim1_monotonic`. **d** : float, array-like, shape :math:`\left(n\right)` :math:`\textit{n}`, :math:`\mathrm{x}`, :math:`\mathrm{f}` and :math:`\mathrm{d}` must be unchanged from the previous call of :meth:`dim1_monotonic`. **px** : float, array-like, shape :math:`\left(m\right)` The :math:`m` values of :math:`x` at which the interpolant is to be evaluated. **Returns** **pf** : float, ndarray, shape :math:`\left(m\right)` :math:`\mathrm{pf}[\textit{i}-1]` contains the value of the interpolant evaluated at the point :math:`\mathrm{px}[\textit{i}-1]`, for :math:`\textit{i} = 1,2,\ldots,m`. **pd** : float, ndarray, shape :math:`\left(m\right)` :math:`\mathrm{pd}[\textit{i}-1]` contains the first derivative of the interpolant evaluated at the point :math:`\mathrm{px}[\textit{i}-1]`, for :math:`\textit{i} = 1,2,\ldots,m`. .. _e01bg-py2-py-errors: **Raises** **NagValueError** (`errno` :math:`1`) On entry, :math:`n = \langle\mathit{\boldsymbol{value}}\rangle`. Constraint: :math:`n\geq 2`. (`errno` :math:`2`) On entry, :math:`r = \langle\mathit{\boldsymbol{value}}\rangle`, :math:`\mathrm{x}[r-2] = \langle\mathit{\boldsymbol{value}}\rangle` and :math:`\mathrm{x}[r-1] = \langle\mathit{\boldsymbol{value}}\rangle`. Constraint: :math:`\mathrm{x}[r-2] < \mathrm{x}[r-1]` for all :math:`r`. (`errno` :math:`3`) On entry, :math:`m = \langle\mathit{\boldsymbol{value}}\rangle`. Constraint: :math:`m\geq 1`. **Warns** **NagAlgorithmicWarning** (`errno` :math:`4`) Warning -- some points in array :math:`\mathrm{px}` lie outside the range :math:`\mathrm{x}[0] \cdots \mathrm{x}[n-1]`. Values at these points are unreliable because computed by extrapolation. .. _e01bg-py2-py-notes: **Notes** `In the NAG Library the traditional C interface for this routine uses a different algorithmic base. Please contact NAG if you have any questions about compatibility.` ``dim1_monotonic_deriv`` evaluates a piecewise cubic Hermite interpolant, as computed by :meth:`dim1_monotonic`, at the points :math:`\mathrm{px}[\textit{i}-1]`, for :math:`\textit{i} = 1,2,\ldots,m`. The first derivatives at the points are also computed. If any point lies outside the interval from :math:`\mathrm{x}[0]` to :math:`\mathrm{x}[n-1]`, values of the interpolant and its derivative are extrapolated from the nearest extreme cubic, and a warning is returned. If values of the interpolant only, and not of its derivative, are required, :meth:`dim1_monotonic_eval` should be used. The function is derived from function PCHFD in Fritsch (1982). .. _e01bg-py2-py-references: **References** Fritsch, F N, 1982, `PCHIP final specifications`, Report UCID-30194, Lawrence Livermore National Laboratory """ raise NotImplementedError
[docs]def dim1_monotonic_intg(x, f, d, a, b): r""" ``dim1_monotonic_intg`` evaluates the definite integral of a piecewise cubic Hermite interpolant over the interval :math:`\left[a, b\right]`. .. _e01bh-py2-py-doc: For full information please refer to the NAG Library document for e01bh https://www.nag.com/numeric/nl/nagdoc_28.5/flhtml/e01/e01bhf.html .. _e01bh-py2-py-parameters: **Parameters** **x** : float, array-like, shape :math:`\left(n\right)` :math:`\textit{n}`, :math:`\mathrm{x}`, :math:`\mathrm{f}` and :math:`\mathrm{d}` must be unchanged from the previous call of :meth:`dim1_monotonic`. **f** : float, array-like, shape :math:`\left(n\right)` :math:`\textit{n}`, :math:`\mathrm{x}`, :math:`\mathrm{f}` and :math:`\mathrm{d}` must be unchanged from the previous call of :meth:`dim1_monotonic`. **d** : float, array-like, shape :math:`\left(n\right)` :math:`\textit{n}`, :math:`\mathrm{x}`, :math:`\mathrm{f}` and :math:`\mathrm{d}` must be unchanged from the previous call of :meth:`dim1_monotonic`. **a** : float The interval :math:`\left[a, b\right]` over which integration is to be performed. **b** : float The interval :math:`\left[a, b\right]` over which integration is to be performed. **Returns** **pint** : float The value of the definite integral of the interpolant over the interval :math:`\left[a, b\right]`. .. _e01bh-py2-py-errors: **Raises** **NagValueError** (`errno` :math:`1`) On entry, :math:`n = \langle\mathit{\boldsymbol{value}}\rangle`. Constraint: :math:`n\geq 2`. (`errno` :math:`2`) On entry, :math:`r = \langle\mathit{\boldsymbol{value}}\rangle`, :math:`\mathrm{x}[r-2] = \langle\mathit{\boldsymbol{value}}\rangle` and :math:`\mathrm{x}[r-1] = \langle\mathit{\boldsymbol{value}}\rangle`. Constraint: :math:`\mathrm{x}[r-2] < \mathrm{x}[r-1]` for all :math:`r`. **Warns** **NagAlgorithmicWarning** (`errno` :math:`3`) Warning -- either :math:`\mathrm{a}` or :math:`\mathrm{b}` is outside the range :math:`\mathrm{x}[0] \cdots \mathrm{x}[n-1]`. The result has been computed by extrapolation and is unreliable. :math:`\mathrm{a} = \langle\mathit{\boldsymbol{value}}\rangle` :math:`\mathrm{b} = \langle\mathit{\boldsymbol{value}}\rangle`. .. _e01bh-py2-py-notes: **Notes** `In the NAG Library the traditional C interface for this routine uses a different algorithmic base. Please contact NAG if you have any questions about compatibility.` ``dim1_monotonic_intg`` evaluates the definite integral of a piecewise cubic Hermite interpolant, as computed by :meth:`dim1_monotonic`, over the interval :math:`\left[a, b\right]`. If either :math:`a` or :math:`b` lies outside the interval from :math:`\mathrm{x}[0]` to :math:`\mathrm{x}[n-1]` computation of the integral involves extrapolation and a warning is returned. The function is derived from function PCHIA in Fritsch (1982). .. _e01bh-py2-py-references: **References** Fritsch, F N, 1982, `PCHIP final specifications`, Report UCID-30194, Lawrence Livermore National Laboratory """ raise NotImplementedError
[docs]def dim1_monconv_disc(negfor, yfor, x, y, lam=0.2): r""" ``dim1_monconv_disc`` computes, for a given set of data points, the forward values and other values required for monotone convex interpolation as defined in Hagan and West (2008). This form of interpolation is particularly suited to the construction of yield curves in Financial Mathematics but can be applied to any data where it is desirable to preserve both monotonicity and convexity. .. _e01ce-py2-py-doc: For full information please refer to the NAG Library document for e01ce https://www.nag.com/numeric/nl/nagdoc_28.5/flhtml/e01/e01cef.html .. _e01ce-py2-py-parameters: **Parameters** **negfor** : bool Determines whether or not to allow negative forward rates. :math:`\mathrm{negfor} = \mathbf{True}` Negative forward rates are permitted. :math:`\mathrm{negfor} = \mathbf{False}` Forward rates calculated must be non-negative. **yfor** : bool Determines whether the array :math:`\mathrm{y}` contains values, :math:`y`, or discrete forward rates :math:`f^d`. :math:`\mathrm{yfor} = \mathbf{True}` :math:`\mathrm{y}` contains the discrete forward rates :math:`f_i^d`, for :math:`\textit{i} = 1,2,\ldots,n`. :math:`\mathrm{yfor} = \mathbf{False}` :math:`\mathrm{y}` contains the values :math:`y_i`, for :math:`\textit{i} = 1,2,\ldots,n`. **x** : float, array-like, shape :math:`\left(n\right)` :math:`x`, the (possibly unordered) set of data points. **y** : float, array-like, shape :math:`\left(n\right)` If :math:`\mathrm{yfor} = \mathbf{True}`, the discrete forward rates :math:`f_i^d` corresponding to the data points :math:`x_i`, for :math:`\textit{i} = 1,2,\ldots,n`. If :math:`\mathrm{yfor} = \mathbf{False}`, the data values :math:`y_i` corresponding to the data points :math:`x_i`, for :math:`\textit{i} = 1,2,\ldots,n`. **lam** : float, optional :math:`\lambda`, the amelioration (smoothing) parameter. Forward rates are first computed using `(2) <https://www.nag.com/numeric/nl/nagdoc_28.5/flhtml/e01/e01cef.html#eqn2>`__ and then, if :math:`\lambda > 0`, a limiting filter is applied which depends on neighbouring discrete forward values. This filter has a smoothing effect on the curve that increases with :math:`\lambda`. **Returns** **comm** : dict, communication object Communication structure. .. _e01ce-py2-py-errors: **Raises** **NagValueError** (`errno` :math:`1`) On entry, :math:`n = \langle\mathit{\boldsymbol{value}}\rangle`. Constraint: :math:`n\geq 2`. (`errno` :math:`2`) On entry, :math:`\mathrm{lam} = \langle\mathit{\boldsymbol{value}}\rangle`. Constraint: :math:`0.0\leq \mathrm{lam}\leq 1.0`. (`errno` :math:`3`) On entry, :math:`\mathrm{x}` contains duplicate data points. .. _e01ce-py2-py-notes: **Notes** ``dim1_monconv_disc`` computes, for a set of data points, :math:`\left(x_i, y_i\right)`, for :math:`\textit{i} = 1,2,\ldots,n`, the discrete forward rates, :math:`f_i^d`, and the instantaneous forward rates, :math:`f_i`, which are used in a monotone convex interpolation method that attempts to preserve both the monotonicity and the convexity of the original data. The monotone convex interpolation method is due to Hagan and West and is described in Hagan and West (2006), Hagan and West (2008) and West (2011). The discrete forward rates are defined simply, for ordered data, by .. math:: \begin{array}{l} f_1^d = y_1\text{;} \\ f_i^d = \frac{{x_iy_i-x_{{i-1}}y_{{i-1}}}}{{x_i-x_{{i-1}}}} \text{, for } i = 2,3,\ldots,n\text{.}\end{array} The discrete forward rates, if pre-computed, may be supplied instead of :math:`y`, in which case the original values :math:`y` are computed using the inverse of `(1) <https://www.nag.com/numeric/nl/nagdoc_28.5/flhtml/e01/e01cef.html#eqn1>`__. The data points :math:`x_i` need not be ordered on input (though :math:`y_i` must correspond to :math:`x_i`); a set of ordered and scaled values :math:`\xi_i` are calculated from :math:`x_i` and stored. In its simplest form, the instantaneous forward rates, :math:`f_i`, at the data points are computed as linear interpolations of the :math:`f_i^d`: .. math:: \begin{array}{l} f_i = \frac{{x_i-x_{{i-1}}}}{{x_{{i+1}}-x_{{i-1}}}} f_{{i+1}}^d + \frac{{x_{{i+1}}-x_i}}{{x_{{i+1}}-x_{{i-1}}}} f_i^d \text{, for } i = 2,3,\ldots,{n-1} \\f_1 = f_2^d - \frac{1}{2} \left(f_2-f_2^d\right) \\ f_n = f_n^d - \frac{1}{2} \left(f_{{n-1}}-f_n^d\right)\text{.} \end{array} If it is required, as a constraint, that these values should never be negative then a limiting filter is applied to :math:`f` as described in Section 3.6 of West (2011). An ameliorated (smoothed) form of this linear interpolation for the forward rates is implemented using the amelioration (smoothing) parameter :math:`\lambda`. For :math:`\lambda \equiv 0`, equation `(2) <https://www.nag.com/numeric/nl/nagdoc_28.5/flhtml/e01/e01cef.html#eqn2>`__ is used (with possible post-process filtering); for :math:`0 < \lambda \leq 1`, the ameliorated method described fully in West (2011) is used. The values computed by ``dim1_monconv_disc`` are used by :meth:`dim1_monconv_eval` to compute, for a given value :math:`\hat{x}`, the monotone convex interpolated (or extrapolated) value :math:`\hat{y}\left(\hat{x}\right)` and the corresponding instantaneous forward rate :math:`f`; the curve gradient at :math:`\hat{x}` can be derived as :math:`y^{\prime } = \left(f-\hat{y}\right)/\hat{x}` for :math:`\hat{x}\neq 0`. .. _e01ce-py2-py-references: **References** Hagan, P S and West, G, 2006, `Interpolation methods for curve construction`, Applied Mathematical Finance (13(2)), 89--129 Hagan, P S and West, G, 2008, `Methods for constructing a yield curve`, WILLMOTT Magazine (May), 70--81 West, G, 2011, `The monotone convex method of interpolation`, Financial Modelling Agency """ raise NotImplementedError
[docs]def dim1_monconv_eval(x, comm): r""" ``dim1_monconv_eval`` evaluates a monotonic convex interpolant at a set of points. .. _e01cf-py2-py-doc: For full information please refer to the NAG Library document for e01cf https://www.nag.com/numeric/nl/nagdoc_28.5/flhtml/e01/e01cff.html .. _e01cf-py2-py-parameters: **Parameters** **x** : float, array-like, shape :math:`\left(m\right)` :math:`x`, the points at which the interpolant is to be evaluated. **comm** : dict, communication object, modified in place Communication structure. This argument must have been initialized by a prior call to :meth:`dim1_monconv_disc`. **Returns** **val** : float, ndarray, shape :math:`\left(m\right)` The values of the interpolant at :math:`x`. **fwd** : float, ndarray, shape :math:`\left(m\right)` The values of the forward rates at :math:`x`. .. _e01cf-py2-py-errors: **Raises** **NagValueError** (`errno` :math:`1`) Either :meth:`dim1_monconv_disc` was not called first or the communication array has become corrupted. .. _e01cf-py2-py-notes: **Notes** ``dim1_monconv_eval`` evaluates a monotonic convex interpolant, as setup by :meth:`dim1_monconv_disc`, at the points :math:`x`. The function is derived from the work of Hagan and West and is described in Hagan and West (2006), Hagan and West (2008) and West (2011). .. _e01cf-py2-py-references: **References** Hagan, P S and West, G, 2006, `Interpolation methods for curve construction`, Applied Mathematical Finance (13(2)), 89--129 Hagan, P S and West, G, 2008, `Methods for constructing a yield curve`, WILLMOTT Magazine (May), 70--81 West, G, 2011, `The monotone convex method of interpolation`, Financial Modelling Agency """ raise NotImplementedError
[docs]def dim2_spline_grid(x, y, f): r""" ``dim2_spline_grid`` computes a bicubic spline interpolating surface through a set of data values, given on a rectangular grid in the :math:`x`-:math:`y` plane. .. _e01da-py2-py-doc: For full information please refer to the NAG Library document for e01da https://www.nag.com/numeric/nl/nagdoc_28.5/flhtml/e01/e01daf.html .. _e01da-py2-py-parameters: **Parameters** **x** : float, array-like, shape :math:`\left(\textit{mx}\right)` :math:`\mathrm{x}[\textit{q}-1]` and :math:`\mathrm{y}[\textit{r}-1]` must contain :math:`x_{\textit{q}}`, for :math:`\textit{q} = 1,2,\ldots,m_x`, and :math:`y_{\textit{r}}`, for :math:`\textit{r} = 1,2,\ldots,m_y`, respectively. **y** : float, array-like, shape :math:`\left(\textit{my}\right)` :math:`\mathrm{x}[\textit{q}-1]` and :math:`\mathrm{y}[\textit{r}-1]` must contain :math:`x_{\textit{q}}`, for :math:`\textit{q} = 1,2,\ldots,m_x`, and :math:`y_{\textit{r}}`, for :math:`\textit{r} = 1,2,\ldots,m_y`, respectively. **f** : float, array-like, shape :math:`\left(\textit{mx}\times \textit{my}\right)` :math:`\mathrm{f}[m_y\times \left(\textit{q}-1\right)+\textit{r}-1]` must contain :math:`f_{{\textit{q},\textit{r}}}`, for :math:`\textit{r} = 1,2,\ldots,m_y`, for :math:`\textit{q} = 1,2,\ldots,m_x`. **Returns** **px** : int :math:`\mathrm{px}` and :math:`\mathrm{py}` contain :math:`m_x+4` and :math:`m_y+4`, the total number of knots of the computed spline with respect to the :math:`x` and :math:`y` variables, respectively. **py** : int :math:`\mathrm{px}` and :math:`\mathrm{py}` contain :math:`m_x+4` and :math:`m_y+4`, the total number of knots of the computed spline with respect to the :math:`x` and :math:`y` variables, respectively. **lamda** : float, ndarray, shape :math:`\left(\textit{mx}+4\right)` :math:`\mathrm{lamda}` contains the complete set of knots :math:`\lambda_i` associated with the :math:`x` variable **mu** : float, ndarray, shape :math:`\left(\textit{my}+4\right)` :math:`\mathrm{mu}` contains the complete set of knots :math:`\mu_i` associated with the :math:`y` variable **c** : float, ndarray, shape :math:`\left(\textit{mx}\times \textit{my}\right)` The coefficients of the spline interpolant. :math:`\mathrm{c}[m_y\times \left(i-1\right)+j-1]` contains the coefficient :math:`c_{{ij}}` described in :ref:`Notes <e01da-py2-py-notes>`. .. _e01da-py2-py-errors: **Raises** **NagValueError** (`errno` :math:`1`) On entry, :math:`\textit{my} = \langle\mathit{\boldsymbol{value}}\rangle`. Constraint: :math:`\textit{my}\geq 4`. (`errno` :math:`1`) On entry, :math:`\textit{mx} = \langle\mathit{\boldsymbol{value}}\rangle`. Constraint: :math:`\textit{mx}\geq 4`. (`errno` :math:`2`) On entry, the :math:`\mathrm{x}` or the :math:`\mathrm{y}` mesh points are not in strictly ascending order. (`errno` :math:`3`) An intermediate set of linear equations is singular -- the data is too ill-conditioned to compute :math:`B`-spline coefficients. .. _e01da-py2-py-notes: **Notes** `In the NAG Library the traditional C interface for this routine uses a different algorithmic base. Please contact NAG if you have any questions about compatibility.` ``dim2_spline_grid`` determines a bicubic spline interpolant to the set of data points :math:`\left(x_{\textit{q}}, y_{\textit{r}}, f_{{\textit{q},\textit{r}}}\right)`, for :math:`\textit{r} = 1,2,\ldots,m_y`, for :math:`\textit{q} = 1,2,\ldots,m_x`. The spline is given in the B-spline representation .. math:: s\left(x, y\right) = \sum_{{i = 1}}^{m_x}\sum_{{j = 1}}^{m_y}c_{{ij}}M_i\left(x\right)N_j\left(y\right)\text{,} such that .. math:: s\left(x_q, y_r\right) = f_{{q,r}}\text{,} where :math:`M_i\left(x\right)` and :math:`N_j\left(y\right)` denote normalized cubic B-splines, the former defined on the knots :math:`\lambda_i` to :math:`\lambda_{{i+4}}` and the latter on the knots :math:`\mu_j` to :math:`\mu_{{j+4}}`, and the :math:`c_{{ij}}` are the spline coefficients. These knots, as well as the coefficients, are determined by the function, which is derived from the function B2IRE in Anthony `et al.` (1982). The method used is described in `Further Comments <https://www.nag.com/numeric/nl/nagdoc_28.5/flhtml/e01/e01daf.html#fcomments2>`__. For further information on splines, see Hayes and Halliday (1974) for bicubic splines and de Boor (1972) for normalized B-splines. Values and derivatives of the computed spline can subsequently be computed by calling :meth:`fit.dim2_spline_evalv <naginterfaces.library.fit.dim2_spline_evalv>`, :meth:`fit.dim2_spline_evalm <naginterfaces.library.fit.dim2_spline_evalm>` or :meth:`fit.dim2_spline_derivm <naginterfaces.library.fit.dim2_spline_derivm>` as described in `Evaluation of Computed Spline <https://www.nag.com/numeric/nl/nagdoc_28.5/flhtml/e01/e01daf.html#fcomments3>`__. .. _e01da-py2-py-references: **References** Anthony, G T, Cox, M G and Hayes, J G, 1982, `DASL -- Data Approximation Subroutine Library`, National Physical Laboratory Cox, M G, 1975, `An algorithm for spline interpolation`, J. Inst. Math. Appl. (15), 95--108 de Boor, C, 1972, `On calculating with B-splines`, J. Approx. Theory (6), 50--62 Hayes, J G and Halliday, J, 1974, `The least squares fitting of cubic spline surfaces to general data sets`, J. Inst. Math. Appl. (14), 89--103 """ raise NotImplementedError
[docs]def dim2_triangulate(x, y): r""" ``dim2_triangulate`` generates a triangulation for a given set of two-dimensional points using the method of Renka and Cline. .. _e01ea-py2-py-doc: For full information please refer to the NAG Library document for e01ea https://www.nag.com/numeric/nl/nagdoc_28.5/flhtml/e01/e01eaf.html .. _e01ea-py2-py-parameters: **Parameters** **x** : float, array-like, shape :math:`\left(n\right)` The :math:`x` coordinates of the :math:`n` data points. **y** : float, array-like, shape :math:`\left(n\right)` The :math:`y` coordinates of the :math:`n` data points. **Returns** **triang** : int, ndarray, shape :math:`\left(7\times n\right)` A data structure defining the computed triangulation, in a form suitable for passing to :meth:`dim2_triang_bary_eval`. Details of how the triangulation is encoded in :math:`\mathrm{triang}` are given in `Further Comments <https://www.nag.com/numeric/nl/nagdoc_28.5/flhtml/e01/e01eaf.html#fcomments>`__. These details are most likely to be of use when plotting the computed triangulation. .. _e01ea-py2-py-errors: **Raises** **NagValueError** (`errno` :math:`1`) On entry, :math:`n = \langle\mathit{\boldsymbol{value}}\rangle`. Constraint: :math:`n\geq 3`. (`errno` :math:`2`) On entry, all the :math:`\left(x, y\right)` pairs are collinear. .. _e01ea-py2-py-notes: **Notes** ``dim2_triangulate`` creates a Thiessen triangulation with a given set of two-dimensional data points as nodes. This triangulation will be as equiangular as possible (Cline and Renka (1984)). See Renka and Cline (1984) for more detailed information on the algorithm, a development of that by Lawson (1977). The code is derived from Renka (1984). The computed triangulation is returned in a form suitable for passing to :meth:`dim2_triang_bary_eval` which, for a set of nodal function values, computes interpolated values at a set of points. .. _e01ea-py2-py-references: **References** Cline, A K and Renka, R L, 1984, `A storage-efficient method for construction of a Thiessen triangulation`, Rocky Mountain J. Math. (14), 119--139 Lawson, C L, 1977, `Software for` :math:`C^1` `surface interpolation`, Mathematical Software III, (ed J R Rice), 161--194, Academic Press Renka, R L, 1984, `Algorithm 624: triangulation and interpolation of arbitrarily distributed points in the plane`, ACM Trans. Math. Software (10), 440--442 Renka, R L and Cline, A K, 1984, `A triangle-based` :math:`C^1` `interpolation method`, Rocky Mountain J. Math. (14), 223--237 """ raise NotImplementedError
[docs]def dim2_triang_bary_eval(x, y, f, triang, px, py): r""" ``dim2_triang_bary_eval`` performs barycentric interpolation, at a given set of points, using a set of function values on a scattered grid and a triangulation of that grid computed by :meth:`dim2_triangulate`. .. _e01eb-py2-py-doc: For full information please refer to the NAG Library document for e01eb https://www.nag.com/numeric/nl/nagdoc_28.5/flhtml/e01/e01ebf.html .. _e01eb-py2-py-parameters: **Parameters** **x** : float, array-like, shape :math:`\left(n\right)` The coordinates of the :math:`\textit{r}`\ th data point, :math:`\left(x_r, y_r\right)`, for :math:`\textit{r} = 1,2,\ldots,n`. :math:`\mathrm{x}` and :math:`\mathrm{y}` must be unchanged from the previous call of :meth:`dim2_triangulate`. **y** : float, array-like, shape :math:`\left(n\right)` The coordinates of the :math:`\textit{r}`\ th data point, :math:`\left(x_r, y_r\right)`, for :math:`\textit{r} = 1,2,\ldots,n`. :math:`\mathrm{y}` and :math:`\mathrm{y}` must be unchanged from the previous call of :meth:`dim2_triangulate`. **f** : float, array-like, shape :math:`\left(n\right)` The function values :math:`f_{\textit{r}}` at :math:`\left(x_{\textit{r}}, y_{\textit{r}}\right)`, for :math:`\textit{r} = 1,2,\ldots,n`. **triang** : int, array-like, shape :math:`\left(7\times n\right)` The triangulation computed by the previous call of :meth:`dim2_triangulate`. See `Further Comments <https://www.nag.com/numeric/nl/nagdoc_28.5/flhtml/e01/e01eaf.html#fcomments>`__ for details of how the triangulation used is encoded in :math:`\mathrm{triang}`. **px** : float, array-like, shape :math:`\left(m\right)` The coordinates :math:`\left(\textit{px}_{\textit{i}}, \textit{py}_{\textit{i}}\right)`, for :math:`\textit{i} = 1,2,\ldots,m`, at which interpolated function values are sought. **py** : float, array-like, shape :math:`\left(m\right)` The coordinates :math:`\left(\textit{px}_{\textit{i}}, \textit{py}_{\textit{i}}\right)`, for :math:`\textit{i} = 1,2,\ldots,m`, at which interpolated function values are sought. **Returns** **pf** : float, ndarray, shape :math:`\left(m\right)` The interpolated values :math:`F\left(\textit{px}_{\textit{i}}, \textit{py}_{\textit{i}}\right)`, for :math:`\textit{i} = 1,2,\ldots,m`. .. _e01eb-py2-py-errors: **Raises** **NagValueError** (`errno` :math:`1`) On entry, :math:`n = \langle\mathit{\boldsymbol{value}}\rangle`. Constraint: :math:`n\geq 3`. (`errno` :math:`2`) On entry, :math:`m = \langle\mathit{\boldsymbol{value}}\rangle`. Constraint: :math:`m\geq 1`. (`errno` :math:`3`) On entry, the triangulation information held in the array :math:`\mathrm{triang}` does not specify a valid triangulation of the data points. :math:`\mathrm{triang}` has been corrupted since the call to :meth:`dim2_triangulate`. **Warns** **NagAlgorithmicWarning** (`errno` :math:`4`) At least one evaluation point lies outside the nodal triangulation. For each such point the value returned in :math:`\mathrm{pf}` is that corresponding to a node on the closest boundary line segment. .. _e01eb-py2-py-notes: **Notes** ``dim2_triang_bary_eval`` takes as input a set of scattered data points :math:`\left(x_{\textit{r}}, y_{\textit{r}}, f_{\textit{r}}\right)`, for :math:`\textit{r} = 1,2,\ldots,n`, and a Thiessen triangulation of the :math:`\left(x_r, y_r\right)` computed by :meth:`dim2_triangulate`, and interpolates at a set of points :math:`\left(\textit{px}_i, \textit{py}_i\right)`, for :math:`\textit{i} = 1,2,\ldots,m`. If the :math:`i`\ th interpolation point :math:`\left(\textit{px}_i, \textit{py}_i\right)` is equal to :math:`\left(x_r, y_r\right)` for some value of :math:`r`, the returned value will be equal to :math:`f_r`; otherwise a barycentric transformation will be used to calculate the interpolant. For each point :math:`\left(\textit{px}_i, \textit{py}_i\right)`, a triangle is sought which contains the point; the vertices of the triangle and :math:`f_r` values at the vertices are then used to compute the value :math:`F\left(\textit{px}_i, \textit{py}_i\right)`. If any interpolation point lies outside the triangulation defined by the input arguments, the returned value is the value provided, :math:`f_s`, at the closest node :math:`\left(x_s, y_s\right)`. ``dim2_triang_bary_eval`` must only be called after a call to :meth:`dim2_triangulate`. .. _e01eb-py2-py-references: **References** Cline, A K and Renka, R L, 1984, `A storage-efficient method for construction of a Thiessen triangulation`, Rocky Mountain J. Math. (14), 119--139 Lawson, C L, 1977, `Software for` :math:`C^1` `surface interpolation`, Mathematical Software III, (ed J R Rice), 161--194, Academic Press Renka, R L, 1984, `Algorithm 624: triangulation and interpolation of arbitrarily distributed points in the plane`, ACM Trans. Math. Software (10), 440--442 Renka, R L and Cline, A K, 1984, `A triangle-based` :math:`C^1` `interpolation method`, Rocky Mountain J. Math. (14), 223--237 """ raise NotImplementedError
[docs]def dim1_ratnl(x, f): r""" ``dim1_ratnl`` produces, from a set of function values and corresponding abscissae, the coefficients of an interpolating rational function expressed in continued fraction form. .. _e01ra-py2-py-doc: For full information please refer to the NAG Library document for e01ra https://www.nag.com/numeric/nl/nagdoc_28.5/flhtml/e01/e01raf.html .. _e01ra-py2-py-parameters: **Parameters** **x** : float, array-like, shape :math:`\left(n\right)` :math:`\mathrm{x}[\textit{i}-1]` must be set to the value of the :math:`\textit{i}`\ th data abscissa, :math:`x_{\textit{i}}`, for :math:`\textit{i} = 1,2,\ldots,n`. **f** : float, array-like, shape :math:`\left(n\right)` :math:`\mathrm{f}[\textit{i}-1]` must be set to the value of the data ordinate, :math:`f_{\textit{i}}`, corresponding to :math:`x_{\textit{i}}`, for :math:`\textit{i} = 1,2,\ldots,n`. **Returns** **m** : int :math:`m`, the number of terms in the continued fraction representation of :math:`R\left(x\right)`. **a** : float, ndarray, shape :math:`\left(n\right)` :math:`\mathrm{a}[\textit{j}-1]` contains the value of the parameter :math:`a_{\textit{j}}` in :math:`R\left(x\right)`, for :math:`\textit{j} = 1,2,\ldots,m`. The remaining elements of :math:`\mathrm{a}`, if any, are set to zero. **u** : float, ndarray, shape :math:`\left(n\right)` :math:`\mathrm{u}[\textit{j}-1]` contains the value of the parameter :math:`u_{\textit{j}}` in :math:`R\left(x\right)`, for :math:`\textit{j} = 1,2,\ldots,m-1`. The :math:`u_j` are a permuted subset of the elements of :math:`\mathrm{x}`. The remaining :math:`n-m+1` locations contain a permutation of the remaining :math:`x_i`, which can be ignored. .. _e01ra-py2-py-errors: **Raises** **NagValueError** (`errno` :math:`1`) On entry, :math:`n = \langle\mathit{\boldsymbol{value}}\rangle`. Constraint: :math:`n > 0`. (`errno` :math:`2`) On entry, :math:`\mathrm{x}[\textit{I}-1]` is very close to :math:`\mathrm{x}[\textit{J}-1]`: :math:`\textit{I} = \langle\mathit{\boldsymbol{value}}\rangle`, :math:`\mathrm{x}[\textit{I}-1] = \langle\mathit{\boldsymbol{value}}\rangle`, :math:`\textit{J} = \langle\mathit{\boldsymbol{value}}\rangle` and :math:`\mathrm{x}[\textit{J}-1] = \langle\mathit{\boldsymbol{value}}\rangle`. (`errno` :math:`3`) A continued fraction of the required form does not exist. .. _e01ra-py2-py-notes: **Notes** ``dim1_ratnl`` produces the parameters of a rational function :math:`R\left(x\right)` which assumes prescribed values :math:`f_i` at prescribed values :math:`x_i` of the independent variable :math:`x`, for :math:`\textit{i} = 1,2,\ldots,n`. More specifically, ``dim1_ratnl`` determines the parameters :math:`a_j`, for :math:`\textit{j} = 1,2,\ldots,m` and :math:`u_j`, for :math:`\textit{j} = 1,2,\ldots,m-1`, in the continued fraction .. math:: R\left(x\right) = a_1+R_m\left(x\right) where .. math:: R_i\left(x\right) = \frac{{a_{{m-i+2}}\left(x-u_{{m-i+1}}\right)}}{{1+R_{{i-1}}\left(x\right)}}\text{, for }i = m,m-1,\ldots,2\text{,} and .. math:: R_1\left(x\right) = 0\text{,} such that :math:`R\left(x_i\right) = f_i`, for :math:`\textit{i} = 1,2,\ldots,n`. The value of :math:`m` in `(1) <https://www.nag.com/numeric/nl/nagdoc_28.5/flhtml/e01/e01raf.html#eqn1>`__ is determined by the function; normally :math:`m = n`. The values of :math:`u_j` form a reordered subset of the values of :math:`x_i` and their ordering is designed to ensure that a representation of the form `(1) <https://www.nag.com/numeric/nl/nagdoc_28.5/flhtml/e01/e01raf.html#eqn1>`__ is determined whenever one exists. The subsequent evaluation of `(1) <https://www.nag.com/numeric/nl/nagdoc_28.5/flhtml/e01/e01raf.html#eqn1>`__ for given values of :math:`x` can be carried out using :meth:`dim1_ratnl_eval`. The computational method employed in ``dim1_ratnl`` is the modification of the Thacher--Tukey algorithm described in Graves--Morris and Hopkins (1981). .. _e01ra-py2-py-references: **References** Graves--Morris, P R and Hopkins, T R, 1981, `Reliable rational interpolation`, Numer. Math. (36), 111--128 """ raise NotImplementedError
[docs]def dim1_ratnl_eval(a, u, x): r""" ``dim1_ratnl_eval`` evaluates continued fractions of the form produced by :meth:`dim1_ratnl`. .. _e01rb-py2-py-doc: For full information please refer to the NAG Library document for e01rb https://www.nag.com/numeric/nl/nagdoc_28.5/flhtml/e01/e01rbf.html .. _e01rb-py2-py-parameters: **Parameters** **a** : float, array-like, shape :math:`\left(m\right)` :math:`\mathrm{a}[\textit{j}-1]` must be set to the value of the parameter :math:`a_{\textit{j}}` in the continued fraction, for :math:`\textit{j} = 1,2,\ldots,m`. **u** : float, array-like, shape :math:`\left(m\right)` :math:`\mathrm{u}[\textit{j}-1]` must be set to the value of the parameter :math:`u_{\textit{j}}` in the continued fraction, for :math:`\textit{j} = 1,2,\ldots,m-1`. (The element :math:`\mathrm{u}[m-1]` is not used). **x** : float The value of :math:`x` at which the continued fraction is to be evaluated. **Returns** **f** : float The value of the continued fraction corresponding to the value of :math:`x`. .. _e01rb-py2-py-errors: **Raises** **NagValueError** (`errno` :math:`1`) :math:`\mathrm{x}` corresponds to a pole of :math:`R\left(x\right)`, or is very close. :math:`\mathrm{x} = \langle\mathit{\boldsymbol{value}}\rangle`. .. _e01rb-py2-py-notes: **Notes** ``dim1_ratnl_eval`` evaluates the continued fraction .. math:: R\left(x\right) = a_1+R_m\left(x\right) where .. math:: R_i\left(x\right) = \frac{{a_{{m-i+2}}\left(x-u_{{m-i+1}}\right)}}{{1+R_{{i-1}}\left(x\right)}}\text{, for }i = m,m-1,\ldots,2\text{.} and .. math:: R_1\left(x\right) = 0 for a prescribed value of :math:`x`. ``dim1_ratnl_eval`` is intended to be used to evaluate the continued fraction representation (of an interpolatory rational function) produced by :meth:`dim1_ratnl`. .. _e01rb-py2-py-references: **References** Graves--Morris, P R and Hopkins, T R, 1981, `Reliable rational interpolation`, Numer. Math. (36), 111--128 """ raise NotImplementedError
[docs]def dim2_scat(x, y, f): r""" ``dim2_scat`` generates a two-dimensional surface interpolating a set of scattered data points, using the method of Renka and Cline. .. _e01sa-py2-py-doc: For full information please refer to the NAG Library document for e01sa https://www.nag.com/numeric/nl/nagdoc_28.5/flhtml/e01/e01saf.html .. _e01sa-py2-py-parameters: **Parameters** **x** : float, array-like, shape :math:`\left(m\right)` The coordinates of the :math:`\textit{r}`\ th data point, for :math:`\textit{r} = 1,2,\ldots,m`. The data points are accepted in any order, but see `Further Comments <https://www.nag.com/numeric/nl/nagdoc_28.5/flhtml/e01/e01saf.html#fcomments>`__. **y** : float, array-like, shape :math:`\left(m\right)` The coordinates of the :math:`\textit{r}`\ th data point, for :math:`\textit{r} = 1,2,\ldots,m`. The data points are accepted in any order, but see `Further Comments <https://www.nag.com/numeric/nl/nagdoc_28.5/flhtml/e01/e01saf.html#fcomments>`__. **f** : float, array-like, shape :math:`\left(m\right)` The coordinates of the :math:`\textit{r}`\ th data point, for :math:`\textit{r} = 1,2,\ldots,m`. The data points are accepted in any order, but see `Further Comments <https://www.nag.com/numeric/nl/nagdoc_28.5/flhtml/e01/e01saf.html#fcomments>`__. **Returns** **triang** : int, ndarray, shape :math:`\left(7\times m\right)` A data structure defining the computed triangulation, in a form suitable for passing to :meth:`dim2_scat_eval`. **grads** : float, ndarray, shape :math:`\left(2, m\right)` The estimated partial derivatives at the nodes, in a form suitable for passing to :meth:`dim2_scat_eval`. The derivatives at node :math:`\textit{r}` with respect to :math:`x` and :math:`y` are contained in :math:`\mathrm{grads}[0,\textit{r}-1]` and :math:`\mathrm{grads}[1,\textit{r}-1]` respectively, for :math:`\textit{r} = 1,2,\ldots,m`. .. _e01sa-py2-py-errors: **Raises** **NagValueError** (`errno` :math:`1`) On entry, :math:`m = \langle\mathit{\boldsymbol{value}}\rangle`. Constraint: :math:`m\geq 3`. (`errno` :math:`2`) All nodes are collinear. There is no unique solution. (`errno` :math:`3`) On entry, :math:`\left({\mathrm{x}[\textit{I}-1]}, {\mathrm{y}[\textit{I}-1]}\right) = \left({\mathrm{x}[\textit{J}-1]}, {\mathrm{y}[\textit{J}-1]}\right)`, for :math:`\textit{I},\textit{J} = \langle\mathit{\boldsymbol{value}}\rangle\langle\mathit{\boldsymbol{value}}\rangle`, :math:`\mathrm{x}[\textit{I}-1]`, :math:`\mathrm{y}[\textit{I}-1] = \langle\mathit{\boldsymbol{value}}\rangle\langle\mathit{\boldsymbol{value}}\rangle`. .. _e01sa-py2-py-notes: **Notes** `In the NAG Library the traditional C interface for this routine uses a different algorithmic base. Please contact NAG if you have any questions about compatibility.` ``dim2_scat`` constructs an interpolating surface :math:`F\left(x, y\right)` through a set of :math:`m` scattered data points :math:`\left(x_{\textit{r}}, y_{\textit{r}}, f_{\textit{r}}\right)`, for :math:`\textit{r} = 1,2,\ldots,m`, using a method due to Renka and Cline. In the :math:`\left(x, y\right)` plane, the data points must be distinct. The constructed surface is continuous and has continuous first derivatives. The method involves firstly creating a triangulation with all the :math:`\left(x, y\right)` data points as nodes, the triangulation being as nearly equiangular as possible (see Cline and Renka (1984)). Then gradients in the :math:`x`- and :math:`y`-directions are estimated at node :math:`\textit{r}`, for :math:`\textit{r} = 1,2,\ldots,m`, as the partial derivatives of a quadratic function of :math:`x` and :math:`y` which interpolates the data value :math:`f_r`, and which fits the data values at nearby nodes (those within a certain distance chosen by the algorithm) in a weighted least squares sense. The weights are chosen such that closer nodes have more influence than more distant nodes on derivative estimates at node :math:`r`. The computed partial derivatives, with the :math:`f_r` values, at the three nodes of each triangle define a piecewise polynomial surface of a certain form which is the interpolant on that triangle. See Renka and Cline (1984) for more detailed information on the algorithm, a development of that by Lawson (1977). The code is derived from Renka (1984). The interpolant :math:`F\left(x, y\right)` can subsequently be evaluated at any point :math:`\left(x, y\right)` inside or outside the domain of the data by a call to :meth:`dim2_scat_eval`. Points outside the domain are evaluated by extrapolation. .. _e01sa-py2-py-references: **References** Cline, A K and Renka, R L, 1984, `A storage-efficient method for construction of a Thiessen triangulation`, Rocky Mountain J. Math. (14), 119--139 Lawson, C L, 1977, `Software for` :math:`C^1` `surface interpolation`, Mathematical Software III, (ed J R Rice), 161--194, Academic Press Renka, R L, 1984, `Algorithm 624: triangulation and interpolation of arbitrarily distributed points in the plane`, ACM Trans. Math. Software (10), 440--442 Renka, R L and Cline, A K, 1984, `A triangle-based` :math:`C^1` `interpolation method`, Rocky Mountain J. Math. (14), 223--237 """ raise NotImplementedError
[docs]def dim2_scat_eval(x, y, f, triang, grads, px, py, ist=1): r""" ``dim2_scat_eval`` evaluates at a given point the two-dimensional interpolant function computed by :meth:`dim2_scat`. .. _e01sb-py2-py-doc: For full information please refer to the NAG Library document for e01sb https://www.nag.com/numeric/nl/nagdoc_28.5/flhtml/e01/e01sbf.html .. _e01sb-py2-py-parameters: **Parameters** **x** : float, array-like, shape :math:`\left(m\right)` :math:`\mathrm{x}` must be unchanged from the previous call of :meth:`dim2_scat` **y** : float, array-like, shape :math:`\left(m\right)` :math:`\mathrm{y}` must be unchanged from the previous call of :meth:`dim2_scat` **f** : float, array-like, shape :math:`\left(m\right)` :math:`\mathrm{f}` must be unchanged from the previous call of :meth:`dim2_scat` **triang** : int, array-like, shape :math:`\left(7\times m\right)` :math:`\mathrm{triang}` must be unchanged from the previous call of :meth:`dim2_scat` **grads** : float, array-like, shape :math:`\left(2, m\right)` :math:`\mathrm{grads}` must be unchanged from the previous call of :meth:`dim2_scat` **px** : float The point :math:`\left({px}, {py}\right)` at which the interpolant is to be evaluated. **py** : float The point :math:`\left({px}, {py}\right)` at which the interpolant is to be evaluated. **ist** : int, optional The index of the starting node in the search for a triangle containing the point :math:`\left({px}, {py}\right)`. On the first call to ``dim2_scat_eval``, :math:`\mathrm{ist}` must be set to :math:`1`. For efficiency on subsequent calls to ``dim2_scat_eval`` an updated value of :math:`\mathrm{ist}` as returned by ``dim2_scat_eval`` may be supplied instead. An input value outside the range :math:`1\leq \mathrm{ist}\leq m` will be treated as :math:`1`. **Returns** **ist** : int The index of one of the vertices of the triangle containing the point :math:`\left({px}, {py}\right)`. **pf** : float The value of the interpolant evaluated at the point :math:`\left({px}, {py}\right)`. .. _e01sb-py2-py-errors: **Raises** **NagValueError** (`errno` :math:`1`) On entry, :math:`m = \langle\mathit{\boldsymbol{value}}\rangle`. Constraint: :math:`m\geq 3`. (`errno` :math:`2`) On entry, :math:`\mathrm{triang}` does not contain a valid data point triangulation; :math:`\mathrm{triang}` may have been corrupted since the call to :meth:`dim2_scat`. **Warns** **NagAlgorithmicWarning** (`errno` :math:`3`) Warning -- the evaluation point :math:`\left(\langle\mathit{\boldsymbol{value}}\rangle, \langle\mathit{\boldsymbol{value}}\rangle\right)` lies outside the triangulation boundary. The returned value was computed by extrapolation. .. _e01sb-py2-py-notes: **Notes** `In the NAG Library the traditional C interface for this routine uses a different algorithmic base. Please contact NAG if you have any questions about compatibility.` ``dim2_scat_eval`` takes as input the arguments defining the interpolant :math:`F\left(x, y\right)` of a set of scattered data points :math:`\left(x_r, y_r, f_r\right)`, for :math:`\textit{r} = 1,2,\ldots,m`, as computed by :meth:`dim2_scat`, and evaluates the interpolant at the point :math:`\left({px}, {py}\right)`. If :math:`\left({px}, {py}\right)` is equal to :math:`\left(x_r, y_r\right)` for some value of :math:`r`, the returned value will be equal to :math:`f_r`. If :math:`\left({px}, {py}\right)` is not equal to :math:`\left(x_r, y_r\right)` for any :math:`r`, the derivatives in :math:`\mathrm{grads}` will be used to compute the interpolant. A triangle is sought which contains the point :math:`\left({px}, {py}\right)`, and the vertices of the triangle along with the partial derivatives and :math:`f_r` values at the vertices are used to compute the value :math:`F\left({px}, {py}\right)`. If the point :math:`\left({px}, {py}\right)` lies outside the triangulation defined by the input arguments, the returned value is obtained by extrapolation. In this case, the interpolating function :math:`\mathrm{f}` is extended linearly beyond the triangulation boundary. The method is described in more detail in Renka and Cline (1984) and the code is derived from Renka (1984). ``dim2_scat_eval`` must only be called after a call to :meth:`dim2_scat`. .. _e01sb-py2-py-references: **References** Renka, R L, 1984, `Algorithm 624: triangulation and interpolation of arbitrarily distributed points in the plane`, ACM Trans. Math. Software (10), 440--442 Renka, R L and Cline, A K, 1984, `A triangle-based` :math:`C^1` `interpolation method`, Rocky Mountain J. Math. (14), 223--237 """ raise NotImplementedError
[docs]def dim2_scat_shep(x, y, f, nw, nq): r""" ``dim2_scat_shep`` generates a two-dimensional interpolant to a set of scattered data points, using a modified Shepard method. .. _e01sg-py2-py-doc: For full information please refer to the NAG Library document for e01sg https://www.nag.com/numeric/nl/nagdoc_28.5/flhtml/e01/e01sgf.html .. _e01sg-py2-py-parameters: **Parameters** **x** : float, array-like, shape :math:`\left(m\right)` The Cartesian coordinates of the data points :math:`\left(x_{\textit{r}}, y_{\textit{r}}\right)`, for :math:`\textit{r} = 1,2,\ldots,m`. **y** : float, array-like, shape :math:`\left(m\right)` The Cartesian coordinates of the data points :math:`\left(x_{\textit{r}}, y_{\textit{r}}\right)`, for :math:`\textit{r} = 1,2,\ldots,m`. **f** : float, array-like, shape :math:`\left(m\right)` :math:`\mathrm{f}[\textit{r}-1]` must be set to the data value :math:`f_{\textit{r}}`, for :math:`\textit{r} = 1,2,\ldots,m`. **nw** : int The number :math:`N_w` of data points that determines each radius of influence :math:`R_w`, appearing in the definition of each of the weights :math:`w_{\textit{r}}`, for :math:`\textit{r} = 1,2,\ldots,m` (see :ref:`Notes <e01sg-py2-py-notes>`). Note that :math:`R_w` is different for each weight. If :math:`\mathrm{nw}\leq 0` the default value :math:`\mathrm{nw} = \mathrm{min}\left(19, {m-1}\right)` is used instead. **nq** : int The number :math:`N_q` of data points to be used in the least squares fit for coefficients defining the nodal functions :math:`q_r\left(x, y\right)` (see :ref:`Notes <e01sg-py2-py-notes>`). If :math:`\mathrm{nq}\leq 0` the default value :math:`\mathrm{nq} = \mathrm{min}\left(13, {m-1}\right)` is used instead. **Returns** **iq** : int, ndarray, shape :math:`\left(2\times m+1\right)` Integer data defining the interpolant :math:`Q\left(x, y\right)`. **rq** : float, ndarray, shape :math:`\left(6\times m+5\right)` Real data defining the interpolant :math:`Q\left(x, y\right)`. .. _e01sg-py2-py-errors: **Raises** **NagValueError** (`errno` :math:`1`) On entry, :math:`\textit{lrq}` is too small: :math:`\textit{lrq} = \langle\mathit{\boldsymbol{value}}\rangle`. (`errno` :math:`1`) On entry, :math:`\textit{liq}` is too small: :math:`\textit{liq} = \langle\mathit{\boldsymbol{value}}\rangle`. (`errno` :math:`1`) On entry, :math:`\mathrm{nw} = \langle\mathit{\boldsymbol{value}}\rangle` and :math:`m = \langle\mathit{\boldsymbol{value}}\rangle`. Constraint: :math:`\mathrm{nw}\leq \mathrm{min}\left(40, {m-1}\right)`. (`errno` :math:`1`) On entry, :math:`\mathrm{nq} = \langle\mathit{\boldsymbol{value}}\rangle` and :math:`m = \langle\mathit{\boldsymbol{value}}\rangle`. Constraint: :math:`\mathrm{nq}\leq \mathrm{min}\left(40, {m-1}\right)`. (`errno` :math:`1`) On entry, :math:`\mathrm{nq} = \langle\mathit{\boldsymbol{value}}\rangle`. Constraint: :math:`\mathrm{nq} \leq 0` or :math:`\mathrm{nq}\geq 5`. (`errno` :math:`1`) On entry, :math:`m = \langle\mathit{\boldsymbol{value}}\rangle`. Constraint: :math:`m\geq 6`. (`errno` :math:`2`) There are duplicate nodes in the dataset. :math:`\left({\mathrm{x}[\textit{I}-1]}, {\mathrm{y}[\textit{I}-1]}\right) = \left({\mathrm{x}[\textit{J}-1]}, {\mathrm{y}[\textit{J}-1]}\right)`, for :math:`\textit{I} = \langle\mathit{\boldsymbol{value}}\rangle` and :math:`\textit{J} = \langle\mathit{\boldsymbol{value}}\rangle`. The interpolant cannot be derived. (`errno` :math:`3`) All nodes are collinear. There is no unique solution. .. _e01sg-py2-py-notes: **Notes** ``dim2_scat_shep`` constructs a smooth function :math:`Q\left(x, y\right)` which interpolates a set of :math:`m` scattered data points :math:`\left(x_r, y_r, f_r\right)`, for :math:`r = 1,2,\ldots,m`, using a modification of Shepard's method. The surface is continuous and has continuous first partial derivatives. The basic Shepard (1968) method interpolates the input data with the weighted mean .. math:: Q\left(x, y\right) = \frac{{\sum_{{r = 1}}^mw_r\left(x, y\right)q_r}}{{\sum_{{r = 1}}^mw_r\left(x, y\right)}}\text{,} where :math:`q_r = f_r`, :math:`w_r\left(x, y\right) = \frac{1}{{d_r^2}}` and :math:`d_r^2 = \left(x-x_r\right)^2+\left(y-y_r\right)^2`. The basic method is global in that the interpolated value at any point depends on all the data, but this function uses a modification (see Franke and Nielson (1980) and Renka (1988a)), whereby the method becomes local by adjusting each :math:`w_r\left(x, y\right)` to be zero outside a circle with centre :math:`\left(x_r, y_r\right)` and some radius :math:`R_w`. Also, to improve the performance of the basic method, each :math:`q_r` above is replaced by a function :math:`q_r\left(x, y\right)`, which is a quadratic fitted by weighted least squares to data local to :math:`\left(x_r, y_r\right)` and forced to interpolate :math:`\left(x_r, y_r, f_r\right)`. In this context, a point :math:`\left(x, y\right)` is defined to be local to another point if it lies within some distance :math:`R_q` of it. Computation of these quadratics constitutes the main work done by this function. The efficiency of the function is further enhanced by using a cell method for nearest neighbour searching due to Bentley and Friedman (1979). The radii :math:`R_w` and :math:`R_q` are chosen to be just large enough to include :math:`N_w` and :math:`N_q` data points, respectively, for user-supplied constants :math:`N_w` and :math:`N_q`. Default values of these arguments are provided by the function, and advice on alternatives is given in `Further Comments <https://www.nag.com/numeric/nl/nagdoc_28.5/flhtml/e01/e01sgf.html#fcomments2>`__. This function is derived from the function QSHEP2 described by Renka (1988b). Values of the interpolant :math:`Q\left(x, y\right)` generated by this function, and its first partial derivatives, can subsequently be evaluated for points in the domain of the data by a call to :meth:`dim2_scat_shep_eval`. .. _e01sg-py2-py-references: **References** Bentley, J L and Friedman, J H, 1979, `Data structures for range searching`, ACM Comput. Surv. (11), 397--409 Franke, R and Nielson, G, 1980, `Smooth interpolation of large sets of scattered data`, Internat. J. Num. Methods Engrg. (15), 1691--1704 Renka, R J, 1988, `Algorithm 660: QSHEP2D: Quadratic Shepard method for bivariate interpolation of scattered data`, ACM Trans. Math. Software (14), 149--150 Renka, R J, 1988, `Multivariate interpolation of large sets of scattered data`, ACM Trans. Math. Software (14), 139--148 Shepard, D, 1968, `A two-dimensional interpolation function for irregularly spaced data`, Proc. 23rd Nat. Conf. ACM, 517--523, Brandon/Systems Press Inc., Princeton """ raise NotImplementedError
[docs]def dim2_scat_shep_eval(x, y, f, iq, rq, u, v): r""" ``dim2_scat_shep_eval`` evaluates the two-dimensional interpolating function generated by :meth:`dim2_scat_shep` and its first partial derivatives. .. _e01sh-py2-py-doc: For full information please refer to the NAG Library document for e01sh https://www.nag.com/numeric/nl/nagdoc_28.5/flhtml/e01/e01shf.html .. _e01sh-py2-py-parameters: **Parameters** **x** : float, array-like, shape :math:`\left(m\right)` :math:`\textit{m}`, :math:`\mathrm{x}`, :math:`\mathrm{y}` and :math:`\mathrm{f}` must be the same values as were supplied in the preceding call to :meth:`dim2_scat_shep`. **y** : float, array-like, shape :math:`\left(m\right)` :math:`\textit{m}`, :math:`\mathrm{x}`, :math:`\mathrm{y}` and :math:`\mathrm{f}` must be the same values as were supplied in the preceding call to :meth:`dim2_scat_shep`. **f** : float, array-like, shape :math:`\left(m\right)` :math:`\textit{m}`, :math:`\mathrm{x}`, :math:`\mathrm{y}` and :math:`\mathrm{f}` must be the same values as were supplied in the preceding call to :meth:`dim2_scat_shep`. **iq** : int, array-like, shape :math:`\left(2\times m+1\right)` Must be unchanged from the value returned from a previous call to :meth:`dim2_scat_shep`. **rq** : float, array-like, shape :math:`\left(6\times m+5\right)` Must be unchanged from the value returned from a previous call to :meth:`dim2_scat_shep`. **u** : float, array-like, shape :math:`\left(n\right)` The evaluation points :math:`\left(u_{\textit{i}}, v_{\textit{i}}\right)`, for :math:`\textit{i} = 1,2,\ldots,n`. **v** : float, array-like, shape :math:`\left(n\right)` The evaluation points :math:`\left(u_{\textit{i}}, v_{\textit{i}}\right)`, for :math:`\textit{i} = 1,2,\ldots,n`. **Returns** **q** : float, ndarray, shape :math:`\left(n\right)` The values of the interpolant at :math:`\left(u_{\textit{i}}, v_{\textit{i}}\right)`, for :math:`\textit{i} = 1,2,\ldots,n`. If any of these evaluation points lie outside the region of definition of the interpolant the corresponding entries in :math:`\mathrm{q}` are set to an extrapolated approximation, and ``dim2_scat_shep_eval`` returns with :math:`\mathrm{errno}` = 3. **qx** : float, ndarray, shape :math:`\left(n\right)` The values of the partial derivatives of the interpolant :math:`Q\left(x, y\right)` at :math:`\left(u_{\textit{i}}, v_{\textit{i}}\right)`, for :math:`\textit{i} = 1,2,\ldots,n`. If any of these evaluation points lie outside the region of definition of the interpolant, the corresponding entries in :math:`\mathrm{qx}` and :math:`\mathrm{qy}` are set to extrapolated approximations to the partial derivatives, and ``dim2_scat_shep_eval`` returns with :math:`\mathrm{errno}` = 3. **qy** : float, ndarray, shape :math:`\left(n\right)` The values of the partial derivatives of the interpolant :math:`Q\left(x, y\right)` at :math:`\left(u_{\textit{i}}, v_{\textit{i}}\right)`, for :math:`\textit{i} = 1,2,\ldots,n`. If any of these evaluation points lie outside the region of definition of the interpolant, the corresponding entries in :math:`\mathrm{qx}` and :math:`\mathrm{qy}` are set to extrapolated approximations to the partial derivatives, and ``dim2_scat_shep_eval`` returns with :math:`\mathrm{errno}` = 3. .. _e01sh-py2-py-errors: **Raises** **NagValueError** (`errno` :math:`1`) On entry, :math:`\textit{lrq}` is too small: :math:`\textit{lrq} = \langle\mathit{\boldsymbol{value}}\rangle`. (`errno` :math:`1`) On entry, :math:`\textit{liq}` is too small: :math:`\textit{liq} = \langle\mathit{\boldsymbol{value}}\rangle`. (`errno` :math:`1`) On entry, :math:`n = \langle\mathit{\boldsymbol{value}}\rangle`. Constraint: :math:`n\geq 1`. (`errno` :math:`1`) On entry, :math:`m = \langle\mathit{\boldsymbol{value}}\rangle`. Constraint: :math:`m\geq 6`. (`errno` :math:`2`) On entry, values in :math:`\mathrm{rq}` appear to be invalid. Check that :math:`\mathrm{rq}` has not been corrupted between calls to :meth:`dim2_scat_shep` and ``dim2_scat_shep_eval``. (`errno` :math:`2`) On entry, values in :math:`\mathrm{iq}` appear to be invalid. Check that :math:`\mathrm{iq}` has not been corrupted between calls to :meth:`dim2_scat_shep` and ``dim2_scat_shep_eval``. **Warns** **NagAlgorithmicWarning** (`errno` :math:`3`) On entry, at least one evaluation point lies outside the region of definition of the interpolant. At such points the corresponding values in :math:`\mathrm{q}` and :math:`\mathrm{qx}` contain extrapolated approximations. Points should be evaluated one by one to identify extrapolated values. .. _e01sh-py2-py-notes: **Notes** ``dim2_scat_shep_eval`` takes as input the interpolant :math:`Q\left(x, y\right)` of a set of scattered data points :math:`\left(x_r, y_r, f_r\right)`, for :math:`\textit{r} = 1,2,\ldots,m`, as computed by :meth:`dim2_scat_shep`, and evaluates the interpolant and its first partial derivatives at the set of points :math:`\left(u_i, v_i\right)`, for :math:`\textit{i} = 1,2,\ldots,n`. ``dim2_scat_shep_eval`` must only be called after a call to :meth:`dim2_scat_shep`. This function is derived from the function QS2GRD described by Renka (1988). .. _e01sh-py2-py-references: **References** Renka, R J, 1988, `Algorithm 660: QSHEP2D: Quadratic Shepard method for bivariate interpolation of scattered data`, ACM Trans. Math. Software (14), 149--150 """ raise NotImplementedError
[docs]def dim3_scat_shep(x, y, z, f, nw, nq): r""" ``dim3_scat_shep`` generates a three-dimensional interpolant to a set of scattered data points, using a modified Shepard method. .. _e01tg-py2-py-doc: For full information please refer to the NAG Library document for e01tg https://www.nag.com/numeric/nl/nagdoc_28.5/flhtml/e01/e01tgf.html .. _e01tg-py2-py-parameters: **Parameters** **x** : float, array-like, shape :math:`\left(m\right)` :math:`\mathrm{x}[\textit{r}-1]`, :math:`\mathrm{y}[\textit{r}-1]`, :math:`\mathrm{z}[\textit{r}-1]` must be set to the Cartesian coordinates of the data point :math:`\left(x_{\textit{r}}, y_{\textit{r}}, z_{\textit{r}}\right)`, for :math:`\textit{r} = 1,2,\ldots,m`. **y** : float, array-like, shape :math:`\left(m\right)` :math:`\mathrm{x}[\textit{r}-1]`, :math:`\mathrm{y}[\textit{r}-1]`, :math:`\mathrm{z}[\textit{r}-1]` must be set to the Cartesian coordinates of the data point :math:`\left(x_{\textit{r}}, y_{\textit{r}}, z_{\textit{r}}\right)`, for :math:`\textit{r} = 1,2,\ldots,m`. **z** : float, array-like, shape :math:`\left(m\right)` :math:`\mathrm{x}[\textit{r}-1]`, :math:`\mathrm{y}[\textit{r}-1]`, :math:`\mathrm{z}[\textit{r}-1]` must be set to the Cartesian coordinates of the data point :math:`\left(x_{\textit{r}}, y_{\textit{r}}, z_{\textit{r}}\right)`, for :math:`\textit{r} = 1,2,\ldots,m`. **f** : float, array-like, shape :math:`\left(m\right)` :math:`\mathrm{f}[\textit{r}-1]` must be set to the data value :math:`f_{\textit{r}}`, for :math:`\textit{r} = 1,2,\ldots,m`. **nw** : int The number :math:`N_w` of data points that determines each radius of influence :math:`R_w`, appearing in the definition of each of the weights :math:`w_{\textit{r}}`, for :math:`\textit{r} = 1,2,\ldots,m` (see :ref:`Notes <e01tg-py2-py-notes>`). Note that :math:`R_w` is different for each weight. If :math:`\mathrm{nw}\leq 0` the default value :math:`\mathrm{nw} = \mathrm{min}\left(32, {m-1}\right)` is used instead. **nq** : int The number :math:`N_q` of data points to be used in the least squares fit for coefficients defining the nodal functions :math:`q_r\left(x, y, z\right)` (see :ref:`Notes <e01tg-py2-py-notes>`). If :math:`\mathrm{nq}\leq 0` the default value :math:`\mathrm{nq} = \mathrm{min}\left(17, {m-1}\right)` is used instead. **Returns** **iq** : int, ndarray, shape :math:`\left(2\times m+1\right)` Integer data defining the interpolant :math:`Q\left(x, y, z\right)`. **rq** : float, ndarray, shape :math:`\left(10\times m+7\right)` Real data defining the interpolant :math:`Q\left(x, y, z\right)`. .. _e01tg-py2-py-errors: **Raises** **NagValueError** (`errno` :math:`1`) On entry, :math:`\textit{lrq}` is too small: :math:`\textit{lrq} = \langle\mathit{\boldsymbol{value}}\rangle`. (`errno` :math:`1`) On entry, :math:`\textit{liq}` is too small: :math:`\textit{liq} = \langle\mathit{\boldsymbol{value}}\rangle`. (`errno` :math:`1`) On entry, :math:`\mathrm{nw} = \langle\mathit{\boldsymbol{value}}\rangle` and :math:`m = \langle\mathit{\boldsymbol{value}}\rangle`. Constraint: :math:`\mathrm{nw}\leq \mathrm{min}\left(40, {m-1}\right)`. (`errno` :math:`1`) On entry, :math:`\mathrm{nq} = \langle\mathit{\boldsymbol{value}}\rangle` and :math:`m = \langle\mathit{\boldsymbol{value}}\rangle`. Constraint: :math:`\mathrm{nq}\leq \mathrm{min}\left(40, {m-1}\right)`. (`errno` :math:`1`) On entry, :math:`\mathrm{nq} = \langle\mathit{\boldsymbol{value}}\rangle`. Constraint: :math:`\mathrm{nq} \leq 0` or :math:`\mathrm{nq} \geq 9`. (`errno` :math:`1`) On entry, :math:`m = \langle\mathit{\boldsymbol{value}}\rangle`. Constraint: :math:`m\geq 10`. (`errno` :math:`2`) There are duplicate nodes in the dataset. :math:`\left({\mathrm{x}[\textit{I}-1]}, {\mathrm{y}[\textit{I}-1]}, {\mathrm{z}[\textit{I}-1]}\right) = \left({\mathrm{x}[\textit{J}-1]}, {\mathrm{y}[\textit{J}-1]}, {\mathrm{z}[\textit{J}-1]}\right)` for: :math:`\textit{I} = \langle\mathit{\boldsymbol{value}}\rangle` and :math:`\textit{J} = \langle\mathit{\boldsymbol{value}}\rangle`. The interpolant cannot be derived. (`errno` :math:`3`) All nodes are coplanar. There is no unique solution. .. _e01tg-py2-py-notes: **Notes** ``dim3_scat_shep`` constructs a smooth function :math:`Q\left(x, y, z\right)` which interpolates a set of :math:`m` scattered data points :math:`\left(x_r, y_r, z_r, f_r\right)`, for :math:`r = 1,2,\ldots,m`, using a modification of Shepard's method. The surface is continuous and has continuous first partial derivatives. The basic Shepard method, which is a generalization of the two-dimensional method described in Shepard (1968), interpolates the input data with the weighted mean .. math:: Q\left(x, y, z\right) = \frac{{\sum_{{r = 1}}^mw_r\left(x, y, z\right)q_r}}{{\sum_{{r = 1}}^mw_r\left(x, y, z\right)}}\text{,} where .. math:: q_r = f_r\text{ and }w_r\left(x, y, z\right) = \frac{1}{{d_r^2}}\text{ and }d_r^2 = \left(x-x_r\right)^2+\left(y-y_r\right)^2+\left(z-z_r\right)^2\text{.} The basic method is global in that the interpolated value at any point depends on all the data, but this function uses a modification (see Franke and Nielson (1980) and Renka (1988a)), whereby the method becomes local by adjusting each :math:`w_r\left(x, y, z\right)` to be zero outside a sphere with centre :math:`\left(x_r, y_r, z_r\right)` and some radius :math:`R_w`. Also, to improve the performance of the basic method, each :math:`q_r` above is replaced by a function :math:`q_r\left(x, y, z\right)`, which is a quadratic fitted by weighted least squares to data local to :math:`\left(x_r, y_r, z_r\right)` and forced to interpolate :math:`\left(x_r, y_r, z_r, f_r\right)`. In this context, a point :math:`\left(x, y, z\right)` is defined to be local to another point if it lies within some distance :math:`R_q` of it. Computation of these quadratics constitutes the main work done by this function. The efficiency of the function is further enhanced by using a cell method for nearest neighbour searching due to Bentley and Friedman (1979). The radii :math:`R_w` and :math:`R_q` are chosen to be just large enough to include :math:`N_w` and :math:`N_q` data points, respectively, for user-supplied constants :math:`N_w` and :math:`N_q`. Default values of these arguments are provided by the function, and advice on alternatives is given in `Further Comments <https://www.nag.com/numeric/nl/nagdoc_28.5/flhtml/e01/e01tgf.html#fcomments2>`__. This function is derived from the function QSHEP3 described by Renka (1988b). Values of the interpolant :math:`Q\left(x, y, z\right)` generated by this function, and its first partial derivatives, can subsequently be evaluated for points in the domain of the data by a call to :meth:`dim3_scat_shep_eval`. .. _e01tg-py2-py-references: **References** Bentley, J L and Friedman, J H, 1979, `Data structures for range searching`, ACM Comput. Surv. (11), 397--409 Franke, R and Nielson, G, 1980, `Smooth interpolation of large sets of scattered data`, Internat. J. Num. Methods Engrg. (15), 1691--1704 Renka, R J, 1988, `Multivariate interpolation of large sets of scattered data`, ACM Trans. Math. Software (14), 139--148 Renka, R J, 1988, `Algorithm 661: QSHEP3D: Quadratic Shepard method for trivariate interpolation of scattered data`, ACM Trans. Math. Software (14), 151--152 Shepard, D, 1968, `A two-dimensional interpolation function for irregularly spaced data`, Proc. 23rd Nat. Conf. ACM, 517--523, Brandon/Systems Press Inc., Princeton """ raise NotImplementedError
[docs]def dim3_scat_shep_eval(x, y, z, f, iq, rq, u, v, w): r""" ``dim3_scat_shep_eval`` evaluates the three-dimensional interpolating function generated by :meth:`dim3_scat_shep` and its first partial derivatives. .. _e01th-py2-py-doc: For full information please refer to the NAG Library document for e01th https://www.nag.com/numeric/nl/nagdoc_28.5/flhtml/e01/e01thf.html .. _e01th-py2-py-parameters: **Parameters** **x** : float, array-like, shape :math:`\left(m\right)` :math:`\textit{m}`, :math:`\mathrm{x}`, :math:`\mathrm{y}`, :math:`\mathrm{z}` and :math:`\mathrm{f}` must be the same values as were supplied in the preceding call to :meth:`dim3_scat_shep`. **y** : float, array-like, shape :math:`\left(m\right)` :math:`\textit{m}`, :math:`\mathrm{x}`, :math:`\mathrm{y}`, :math:`\mathrm{z}` and :math:`\mathrm{f}` must be the same values as were supplied in the preceding call to :meth:`dim3_scat_shep`. **z** : float, array-like, shape :math:`\left(m\right)` :math:`\textit{m}`, :math:`\mathrm{x}`, :math:`\mathrm{y}`, :math:`\mathrm{z}` and :math:`\mathrm{f}` must be the same values as were supplied in the preceding call to :meth:`dim3_scat_shep`. **f** : float, array-like, shape :math:`\left(m\right)` :math:`\textit{m}`, :math:`\mathrm{x}`, :math:`\mathrm{y}`, :math:`\mathrm{z}` and :math:`\mathrm{f}` must be the same values as were supplied in the preceding call to :meth:`dim3_scat_shep`. **iq** : int, array-like, shape :math:`\left(2\times m+1\right)` Must be unchanged from the value returned from a previous call to :meth:`dim3_scat_shep`. **rq** : float, array-like, shape :math:`\left(10\times m+7\right)` Must be unchanged from the value returned from a previous call to :meth:`dim3_scat_shep`. **u** : float, array-like, shape :math:`\left(n\right)` :math:`\mathrm{u}[\textit{i}-1]`, :math:`\mathrm{v}[\textit{i}-1]`, :math:`\mathrm{w}[\textit{i}-1]` must be set to the evaluation point :math:`\left(u_{\textit{i}}, v_{\textit{i}}, w_{\textit{i}}\right)`, for :math:`\textit{i} = 1,2,\ldots,n`. **v** : float, array-like, shape :math:`\left(n\right)` :math:`\mathrm{u}[\textit{i}-1]`, :math:`\mathrm{v}[\textit{i}-1]`, :math:`\mathrm{w}[\textit{i}-1]` must be set to the evaluation point :math:`\left(u_{\textit{i}}, v_{\textit{i}}, w_{\textit{i}}\right)`, for :math:`\textit{i} = 1,2,\ldots,n`. **w** : float, array-like, shape :math:`\left(n\right)` :math:`\mathrm{u}[\textit{i}-1]`, :math:`\mathrm{v}[\textit{i}-1]`, :math:`\mathrm{w}[\textit{i}-1]` must be set to the evaluation point :math:`\left(u_{\textit{i}}, v_{\textit{i}}, w_{\textit{i}}\right)`, for :math:`\textit{i} = 1,2,\ldots,n`. **Returns** **q** : float, ndarray, shape :math:`\left(n\right)` :math:`\mathrm{q}[\textit{i}-1]` contains the value of the interpolant, at :math:`\left(u_{\textit{i}}, v_{\textit{i}}, w_{\textit{i}}\right)`, for :math:`\textit{i} = 1,2,\ldots,n`. If any of these evaluation points lie outside the region of definition of the interpolant the corresponding entries in :math:`\mathrm{q}` are set to an extrapolated approximation, and ``dim3_scat_shep_eval`` returns with :math:`\mathrm{errno}` = 3. **qx** : float, ndarray, shape :math:`\left(n\right)` :math:`\mathrm{qx}[\textit{i}-1]`, :math:`\mathrm{qy}[\textit{i}-1]`, :math:`\mathrm{qz}[\textit{i}-1]` contains the value of the partial derivatives of the interpolant :math:`Q\left(x, y, z\right)` at :math:`\left(u_{\textit{i}}, v_{\textit{i}}, w_{\textit{i}}\right)`, for :math:`\textit{i} = 1,2,\ldots,n`. If any of these evaluation points lie outside the region of definition of the interpolant, the corresponding entries in :math:`\mathrm{qx}`, :math:`\mathrm{qy}` and :math:`\mathrm{qz}` are set to extrapolated approximations to the partial derivatives, and ``dim3_scat_shep_eval`` returns with :math:`\mathrm{errno}` = 3. **qy** : float, ndarray, shape :math:`\left(n\right)` :math:`\mathrm{qx}[\textit{i}-1]`, :math:`\mathrm{qy}[\textit{i}-1]`, :math:`\mathrm{qz}[\textit{i}-1]` contains the value of the partial derivatives of the interpolant :math:`Q\left(x, y, z\right)` at :math:`\left(u_{\textit{i}}, v_{\textit{i}}, w_{\textit{i}}\right)`, for :math:`\textit{i} = 1,2,\ldots,n`. If any of these evaluation points lie outside the region of definition of the interpolant, the corresponding entries in :math:`\mathrm{qx}`, :math:`\mathrm{qy}` and :math:`\mathrm{qz}` are set to extrapolated approximations to the partial derivatives, and ``dim3_scat_shep_eval`` returns with :math:`\mathrm{errno}` = 3. **qz** : float, ndarray, shape :math:`\left(n\right)` :math:`\mathrm{qx}[\textit{i}-1]`, :math:`\mathrm{qy}[\textit{i}-1]`, :math:`\mathrm{qz}[\textit{i}-1]` contains the value of the partial derivatives of the interpolant :math:`Q\left(x, y, z\right)` at :math:`\left(u_{\textit{i}}, v_{\textit{i}}, w_{\textit{i}}\right)`, for :math:`\textit{i} = 1,2,\ldots,n`. If any of these evaluation points lie outside the region of definition of the interpolant, the corresponding entries in :math:`\mathrm{qx}`, :math:`\mathrm{qy}` and :math:`\mathrm{qz}` are set to extrapolated approximations to the partial derivatives, and ``dim3_scat_shep_eval`` returns with :math:`\mathrm{errno}` = 3. .. _e01th-py2-py-errors: **Raises** **NagValueError** (`errno` :math:`1`) On entry, :math:`\textit{lrq}` is too small: :math:`\textit{lrq} = \langle\mathit{\boldsymbol{value}}\rangle`. (`errno` :math:`1`) On entry, :math:`\textit{liq}` is too small: :math:`\textit{liq} = \langle\mathit{\boldsymbol{value}}\rangle`. (`errno` :math:`1`) On entry, :math:`n = \langle\mathit{\boldsymbol{value}}\rangle`. Constraint: :math:`n\geq 1`. (`errno` :math:`1`) On entry, :math:`m = \langle\mathit{\boldsymbol{value}}\rangle`. Constraint: :math:`m\geq 10`. (`errno` :math:`2`) On entry, values in :math:`\mathrm{rq}` appear to be invalid. Check that :math:`\mathrm{rq}` has not been corrupted between calls to :meth:`dim3_scat_shep` and ``dim3_scat_shep_eval``. (`errno` :math:`2`) On entry, values in :math:`\mathrm{iq}` appear to be invalid. Check that :math:`\mathrm{iq}` has not been corrupted between calls to :meth:`dim3_scat_shep` and ``dim3_scat_shep_eval``. **Warns** **NagAlgorithmicWarning** (`errno` :math:`3`) On entry, at least one evaluation point lies outside the region of definition of the interpolant. At such points the corresponding values in :math:`\mathrm{q}` and :math:`\mathrm{qx}` contain extrapolated approximations. Points should be evaluated one by one to identify extrapolated values. .. _e01th-py2-py-notes: **Notes** ``dim3_scat_shep_eval`` takes as input the interpolant :math:`Q\left(x, y, z\right)` of a set of scattered data points :math:`\left(x_r, y_r, z_r, f_r\right)`, for :math:`\textit{r} = 1,2,\ldots,m`, as computed by :meth:`dim3_scat_shep`, and evaluates the interpolant and its first partial derivatives at the set of points :math:`\left(u_i, v_i, w_i\right)`, for :math:`\textit{i} = 1,2,\ldots,n`. ``dim3_scat_shep_eval`` must only be called after a call to :meth:`dim3_scat_shep`. This function is derived from the function QS3GRD described by Renka (1988). .. _e01th-py2-py-references: **References** Renka, R J, 1988, `Algorithm 661: QSHEP3D: Quadratic Shepard method for trivariate interpolation of scattered data`, ACM Trans. Math. Software (14), 151--152 """ raise NotImplementedError
[docs]def dim4_scat_shep(x, f, nw, nq): r""" ``dim4_scat_shep`` generates a four-dimensional interpolant to a set of scattered data points, using a modified Shepard method. .. _e01tk-py2-py-doc: For full information please refer to the NAG Library document for e01tk https://www.nag.com/numeric/nl/nagdoc_28.5/flhtml/e01/e01tkf.html .. _e01tk-py2-py-parameters: **Parameters** **x** : float, array-like, shape :math:`\left(4, m\right)` :math:`\mathrm{x}[0,\textit{r}-1],\ldots,\mathrm{x}[3,\textit{r}-1]` must be set to the Cartesian coordinates of the data point :math:`\mathbf{x}_{\textit{r}}`, for :math:`\textit{r} = 1,2,\ldots,m`. **f** : float, array-like, shape :math:`\left(m\right)` :math:`\mathrm{f}[\textit{r}-1]` must be set to the data value :math:`f_{\textit{r}}`, for :math:`\textit{r} = 1,2,\ldots,m`. **nw** : int The number :math:`N_w` of data points that determines each radius of influence :math:`R_w`, appearing in the definition of each of the weights :math:`w_{\textit{r}}`, for :math:`\textit{r} = 1,2,\ldots,m` (see :ref:`Notes <e01tk-py2-py-notes>`). Note that :math:`R_w` is different for each weight. If :math:`\mathrm{nw}\leq 0` the default value :math:`\mathrm{nw} = \mathrm{min}\left(32, {m-1}\right)` is used instead. **nq** : int The number :math:`N_q` of data points to be used in the least squares fit for coefficients defining the quadratic functions :math:`q_r\left(\mathbf{x}\right)` (see :ref:`Notes <e01tk-py2-py-notes>`). If :math:`\mathrm{nq}\leq 0` the default value :math:`\mathrm{nq} = \mathrm{min}\left(38, {m-1}\right)` is used instead. **Returns** **iq** : int, ndarray, shape :math:`\left(2\times m+1\right)` Integer data defining the interpolant :math:`Q\left(\mathbf{x}\right)`. **rq** : float, ndarray, shape :math:`\left(15\times m+9\right)` Real data defining the interpolant :math:`Q\left(\mathbf{x}\right)`. .. _e01tk-py2-py-errors: **Raises** **NagValueError** (`errno` :math:`1`) On entry, :math:`\mathrm{nw} = \langle\mathit{\boldsymbol{value}}\rangle` and :math:`m = \langle\mathit{\boldsymbol{value}}\rangle`. Constraint: :math:`\mathrm{nw}\leq \mathrm{min}\left(50, {m-1}\right)`. (`errno` :math:`1`) On entry, :math:`\mathrm{nq} = \langle\mathit{\boldsymbol{value}}\rangle` and :math:`m = \langle\mathit{\boldsymbol{value}}\rangle`. Constraint: :math:`\mathrm{nq}\leq \mathrm{min}\left(50, {m-1}\right)`. (`errno` :math:`1`) On entry, :math:`\mathrm{nq} = \langle\mathit{\boldsymbol{value}}\rangle`. Constraint: :math:`\mathrm{nq} \leq 0` or :math:`\mathrm{nq} \geq 14`. (`errno` :math:`1`) On entry, :math:`m = \langle\mathit{\boldsymbol{value}}\rangle`. Constraint: :math:`m\geq 16`. (`errno` :math:`2`) There are duplicate nodes in the dataset. :math:`{\mathrm{x}[i-1,k-1]} = {\mathrm{x}[j-1,k-1]}`, for :math:`i = \langle\mathit{\boldsymbol{value}}\rangle`, :math:`j = \langle\mathit{\boldsymbol{value}}\rangle` and :math:`k = 1,2,\ldots,4`. The interpolant cannot be derived. (`errno` :math:`3`) On entry, all the data points lie on the same three-dimensional hypersurface. No unique solution exists. .. _e01tk-py2-py-notes: **Notes** ``dim4_scat_shep`` constructs a smooth function :math:`Q\left(\mathbf{x}\right)`, :math:`\mathbf{x} \in \mathbb{R}^4` which interpolates a set of :math:`m` scattered data points :math:`\left(\mathbf{x}_r, f_r\right)`, for :math:`r = 1,2,\ldots,m`, using a modification of Shepard's method. The surface is continuous and has continuous first partial derivatives. The basic Shepard method, which is a generalization of the two-dimensional method described in Shepard (1968), interpolates the input data with the weighted mean .. math:: Q\left(\mathbf{x}\right) = \frac{{\sum_{{r = 1}}^mw_r\left(\mathbf{x}\right)q_r}}{{\sum_{{r = 1}}^mw_r\left(\mathbf{x}\right)}}\text{,} where :math:`q_r = f_r`, :math:`w_r\left(\mathbf{x}\right) = \frac{1}{{d_r^2}}` and :math:`d_r^2 = \left\lVert \mathbf{x}-\mathbf{x}_r\right\rVert_2^2`. The basic method is global in that the interpolated value at any point depends on all the data, but ``dim4_scat_shep`` uses a modification (see Franke and Nielson (1980) and Renka (1988a)), whereby the method becomes local by adjusting each :math:`w_r\left(\mathbf{x}\right)` to be zero outside a hypersphere with centre :math:`\mathbf{x}_r` and some radius :math:`R_w`. Also, to improve the performance of the basic method, each :math:`q_r` above is replaced by a function :math:`q_r\left(\mathbf{x}\right)`, which is a quadratic fitted by weighted least squares to data local to :math:`\mathbf{x}_r` and forced to interpolate :math:`\left(\mathbf{x}_r, f_r\right)`. In this context, a point :math:`\mathbf{x}` is defined to be local to another point if it lies within some distance :math:`R_q` of it. The efficiency of ``dim4_scat_shep`` is enhanced by using a cell method for nearest neighbour searching due to Bentley and Friedman (1979) with a cell density of :math:`3`. The radii :math:`R_w` and :math:`R_q` are chosen to be just large enough to include :math:`N_w` and :math:`N_q` data points, respectively, for user-supplied constants :math:`N_w` and :math:`N_q`. Default values of these arguments are provided by the function, and advice on alternatives is given in `Further Comments <https://www.nag.com/numeric/nl/nagdoc_28.5/flhtml/e01/e01tkf.html#fcomments2>`__. ``dim4_scat_shep`` is derived from the new implementation of QSHEP3 described by Renka (1988b). It uses the modification for high-dimensional interpolation described by Berry and Minser (1999). Values of the interpolant :math:`Q\left(\mathbf{x}\right)` generated by ``dim4_scat_shep``, and its first partial derivatives, can subsequently be evaluated for points in the domain of the data by a call to :meth:`dim4_scat_shep_eval`. .. _e01tk-py2-py-references: **References** Bentley, J L and Friedman, J H, 1979, `Data structures for range searching`, ACM Comput. Surv. (11), 397--409 Berry, M W, Minser, K S, 1999, `Algorithm 798: high-dimensional interpolation using the modified Shepard method`, ACM Trans. Math. Software (25), 353--366 Franke, R and Nielson, G, 1980, `Smooth interpolation of large sets of scattered data`, Internat. J. Num. Methods Engrg. (15), 1691--1704 Renka, R J, 1988, `Multivariate interpolation of large sets of scattered data`, ACM Trans. Math. Software (14), 139--148 Renka, R J, 1988, `Algorithm 661: QSHEP3D: Quadratic Shepard method for trivariate interpolation of scattered data`, ACM Trans. Math. Software (14), 151--152 Shepard, D, 1968, `A two-dimensional interpolation function for irregularly spaced data`, Proc. 23rd Nat. Conf. ACM, 517--523, Brandon/Systems Press Inc., Princeton """ raise NotImplementedError
[docs]def dim4_scat_shep_eval(x, f, iq, rq, xe): r""" ``dim4_scat_shep_eval`` evaluates the four-dimensional interpolating function generated by :meth:`dim4_scat_shep` and its first partial derivatives. .. _e01tl-py2-py-doc: For full information please refer to the NAG Library document for e01tl https://www.nag.com/numeric/nl/nagdoc_28.5/flhtml/e01/e01tlf.html .. _e01tl-py2-py-parameters: **Parameters** **x** : float, array-like, shape :math:`\left(4, m\right)` Note: the coordinates of :math:`x_r` are stored in :math:`\mathrm{x}[0,r-1]\ldots \mathrm{x}[3,r-1]`. **Must** be the same array supplied as argument :math:`\textit{x}` in the preceding call to :meth:`dim4_scat_shep`. It **must** remain unchanged between calls. **f** : float, array-like, shape :math:`\left(m\right)` **Must** be the same array supplied as argument :math:`\textit{f}` in the preceding call to :meth:`dim4_scat_shep`. It **must** remain unchanged between calls. **iq** : int, array-like, shape :math:`\left(2\times m+1\right)` **Must** be the same array returned as argument :math:`\textit{iq}` in the preceding call to :meth:`dim4_scat_shep`. It **must** remain unchanged between calls. **rq** : float, array-like, shape :math:`\left(15\times m+9\right)` **Must** be the same array returned as argument :math:`\textit{rq}` in the preceding call to :meth:`dim4_scat_shep`. It **must** remain unchanged between calls. **xe** : float, array-like, shape :math:`\left(4, n\right)` :math:`\mathrm{xe}[0,\textit{r}-1],\ldots,\mathrm{xe}[3,\textit{r}-1]` must be set to the evaluation point :math:`\mathbf{x}_{\textit{i}}`, for :math:`\textit{i} = 1,2,\ldots,n`. **Returns** **q** : float, ndarray, shape :math:`\left(n\right)` :math:`\mathrm{q}[\textit{i}-1]` contains the value of the interpolant, at :math:`\mathbf{x}_{\textit{i}}`, for :math:`\textit{i} = 1,2,\ldots,n`. If any of these evaluation points lie outside the region of definition of the interpolant the corresponding entries in :math:`\mathrm{q}` are set to an extrapolated approximation, and ``dim4_scat_shep_eval`` returns with :math:`\mathrm{errno}` = 3. **qx** : float, ndarray, shape :math:`\left(4, n\right)` :math:`\mathrm{qx}[j-1,i-1]` contains the value of the partial derivatives with respect to :math:`\mathbf{x}_j` of the interpolant :math:`Q\left(\mathbf{x}\right)` at :math:`\mathbf{x}_{\textit{i}}`, for :math:`\textit{i} = 1,2,\ldots,n`, and for each of the four partial derivatives :math:`j = 1,2,3,4`. If any of these evaluation points lie outside the region of definition of the interpolant, the corresponding entries in :math:`\mathrm{qx}` are set to extrapolated approximations to the partial derivatives, and ``dim4_scat_shep_eval`` returns with :math:`\mathrm{errno}` = 3. .. _e01tl-py2-py-errors: **Raises** **NagValueError** (`errno` :math:`1`) On entry, :math:`n = \langle\mathit{\boldsymbol{value}}\rangle`. Constraint: :math:`n\geq 1`. (`errno` :math:`1`) On entry, :math:`m = \langle\mathit{\boldsymbol{value}}\rangle`. Constraint: :math:`m\geq 16`. (`errno` :math:`2`) On entry, values in :math:`\mathrm{rq}` appear to be invalid. Check that :math:`\mathrm{rq}` has not been corrupted between calls to :meth:`dim4_scat_shep` and ``dim4_scat_shep_eval``. (`errno` :math:`2`) On entry, values in :math:`\mathrm{iq}` appear to be invalid. Check that :math:`\mathrm{iq}` has not been corrupted between calls to :meth:`dim4_scat_shep` and ``dim4_scat_shep_eval``. **Warns** **NagAlgorithmicWarning** (`errno` :math:`3`) On entry, at least one evaluation point lies outside the region of definition of the interpolant. At such points the corresponding values in :math:`\mathrm{q}` and :math:`\mathrm{qx}` contain extrapolated approximations. Points should be evaluated one by one to identify extrapolated values. .. _e01tl-py2-py-notes: **Notes** ``dim4_scat_shep_eval`` takes as input the interpolant :math:`Q\left(\mathbf{x}\right)`, :math:`x \in \mathbb{R}^4` of a set of scattered data points :math:`\left(\mathbf{x}_{\textit{r}}, f_{\textit{r}}\right)`, for :math:`\textit{r} = 1,2,\ldots,m`, as computed by :meth:`dim4_scat_shep`, and evaluates the interpolant and its first partial derivatives at the set of points :math:`\mathbf{x}_i`, for :math:`\textit{i} = 1,2,\ldots,n`. ``dim4_scat_shep_eval`` must only be called after a call to :meth:`dim4_scat_shep`. ``dim4_scat_shep_eval`` is derived from the new implementation of QS3GRD described by Renka (1988). It uses the modification for high-dimensional interpolation described by Berry and Minser (1999). .. _e01tl-py2-py-references: **References** Berry, M W, Minser, K S, 1999, `Algorithm 798: high-dimensional interpolation using the modified Shepard method`, ACM Trans. Math. Software (25), 353--366 Renka, R J, 1988, `Algorithm 661: QSHEP3D: Quadratic Shepard method for trivariate interpolation of scattered data`, ACM Trans. Math. Software (14), 151--152 """ raise NotImplementedError
[docs]def dim5_scat_shep(x, f, nw, nq): r""" ``dim5_scat_shep`` generates a five-dimensional interpolant to a set of scattered data points, using a modified Shepard method. .. _e01tm-py2-py-doc: For full information please refer to the NAG Library document for e01tm https://www.nag.com/numeric/nl/nagdoc_28.5/flhtml/e01/e01tmf.html .. _e01tm-py2-py-parameters: **Parameters** **x** : float, array-like, shape :math:`\left(5, m\right)` :math:`\mathrm{x}[0,\textit{r}-1],\ldots,\mathrm{x}[4,\textit{r}-1]` must be set to the Cartesian coordinates of the data point :math:`\mathbf{x}_{\textit{r}}`, for :math:`\textit{r} = 1,2,\ldots,m`. **f** : float, array-like, shape :math:`\left(m\right)` :math:`\mathrm{f}[\textit{r}-1]` must be set to the data value :math:`f_{\textit{r}}`, for :math:`\textit{r} = 1,2,\ldots,m`. **nw** : int The number :math:`N_w` of data points that determines each radius of influence :math:`R_w`, appearing in the definition of each of the weights :math:`w_{\textit{r}}`, for :math:`\textit{r} = 1,2,\ldots,m` (see :ref:`Notes <e01tm-py2-py-notes>`). Note that :math:`R_w` is different for each weight. If :math:`\mathrm{nw}\leq 0` the default value :math:`\mathrm{nw} = \mathrm{min}\left(32, {m-1}\right)` is used instead. **nq** : int The number :math:`N_q` of data points to be used in the least squares fit for coefficients defining the quadratic functions :math:`q_r\left(\mathbf{x}\right)` (see :ref:`Notes <e01tm-py2-py-notes>`). If :math:`\mathrm{nq}\leq 0` the default value :math:`\mathrm{nq} = \mathrm{min}\left(50, {m-1}\right)` is used instead. **Returns** **iq** : int, ndarray, shape :math:`\left(2\times m+1\right)` Integer data defining the interpolant :math:`Q\left(\mathbf{x}\right)`. **rq** : float, ndarray, shape :math:`\left(21\times m+11\right)` Real data defining the interpolant :math:`Q\left(\mathbf{x}\right)`. .. _e01tm-py2-py-errors: **Raises** **NagValueError** (`errno` :math:`1`) On entry, :math:`\mathrm{nw} = \langle\mathit{\boldsymbol{value}}\rangle` and :math:`m = \langle\mathit{\boldsymbol{value}}\rangle`. Constraint: :math:`\mathrm{nw}\leq \mathrm{min}\left(50, {m-1}\right)`. (`errno` :math:`1`) On entry, :math:`\mathrm{nq} = \langle\mathit{\boldsymbol{value}}\rangle` and :math:`m = \langle\mathit{\boldsymbol{value}}\rangle`. Constraint: :math:`\mathrm{nq}\leq \mathrm{min}\left(70, {m-1}\right)`. (`errno` :math:`1`) On entry, :math:`\mathrm{nq} = \langle\mathit{\boldsymbol{value}}\rangle`. Constraint: :math:`\mathrm{nq} \leq 0` or :math:`\mathrm{nq} \geq 20`. (`errno` :math:`1`) On entry, :math:`m = \langle\mathit{\boldsymbol{value}}\rangle`. Constraint: :math:`m\geq 23`. (`errno` :math:`2`) There are duplicate nodes in the dataset. :math:`{\mathrm{x}[i-1,k-1]} = {\mathrm{x}[j-1,k-1]}`, for :math:`i = \langle\mathit{\boldsymbol{value}}\rangle`, :math:`j = \langle\mathit{\boldsymbol{value}}\rangle` and :math:`k = 1,2,\ldots,5`. The interpolant cannot be derived. (`errno` :math:`3`) On entry, all the data points lie on the same four-dimensional hypersurface. No unique solution exists. .. _e01tm-py2-py-notes: **Notes** ``dim5_scat_shep`` constructs a smooth function :math:`Q\left(\mathbf{x}\right)`, :math:`\mathbf{x} \in \mathbb{R}^5` which interpolates a set of :math:`m` scattered data points :math:`\left(\mathbf{x}_r, f_r\right)`, for :math:`r = 1,2,\ldots,m`, using a modification of Shepard's method. The surface is continuous and has continuous first partial derivatives. The basic Shepard method, which is a generalization of the two-dimensional method described in Shepard (1968), interpolates the input data with the weighted mean .. math:: Q\left(\mathbf{x}\right) = \frac{{\sum_{{r = 1}}^mw_r\left(\mathbf{x}\right)q_r}}{{\sum_{{r = 1}}^mw_r\left(\mathbf{x}\right)}}\text{,} where :math:`q_r = f_r`, :math:`w_r\left(\mathbf{x}\right) = \frac{1}{{d_r^2}}` and :math:`d_r^2 = \left\lVert \mathbf{x}-\mathbf{x}_r\right\rVert_2^2`. The basic method is global in that the interpolated value at any point depends on all the data, but ``dim5_scat_shep`` uses a modification (see Franke and Nielson (1980) and Renka (1988a)), whereby the method becomes local by adjusting each :math:`w_r\left(\mathbf{x}\right)` to be zero outside a hypersphere with centre :math:`\mathbf{x}_r` and some radius :math:`R_w`. Also, to improve the performance of the basic method, each :math:`q_r` above is replaced by a function :math:`q_r\left(\mathbf{x}\right)`, which is a quadratic fitted by weighted least squares to data local to :math:`\mathbf{x}_r` and forced to interpolate :math:`\left(\mathbf{x}_r, f_r\right)`. In this context, a point :math:`\mathbf{x}` is defined to be local to another point if it lies within some distance :math:`R_q` of it. The efficiency of ``dim5_scat_shep`` is enhanced by using a cell method for nearest neighbour searching due to Bentley and Friedman (1979) with a cell density of :math:`3`. The radii :math:`R_w` and :math:`R_q` are chosen to be just large enough to include :math:`N_w` and :math:`N_q` data points, respectively, for user-supplied constants :math:`N_w` and :math:`N_q`. Default values of these arguments are provided, and advice on alternatives is given in `Further Comments <https://www.nag.com/numeric/nl/nagdoc_28.5/flhtml/e01/e01tmf.html#fcomments2>`__. ``dim5_scat_shep`` is derived from the new implementation of QSHEP3 described by Renka (1988b). It uses the modification for five-dimensional interpolation described by Berry and Minser (1999). Values of the interpolant :math:`Q\left(\mathbf{x}\right)` generated by ``dim5_scat_shep``, and its first partial derivatives, can subsequently be evaluated for points in the domain of the data by a call to :meth:`dim5_scat_shep_eval`. .. _e01tm-py2-py-references: **References** Bentley, J L and Friedman, J H, 1979, `Data structures for range searching`, ACM Comput. Surv. (11), 397--409 Berry, M W, Minser, K S, 1999, `Algorithm 798: high-dimensional interpolation using the modified Shepard method`, ACM Trans. Math. Software (25), 353--366 Franke, R and Nielson, G, 1980, `Smooth interpolation of large sets of scattered data`, Internat. J. Num. Methods Engrg. (15), 1691--1704 Renka, R J, 1988, `Multivariate interpolation of large sets of scattered data`, ACM Trans. Math. Software (14), 139--148 Renka, R J, 1988, `Algorithm 661: QSHEP3D: Quadratic Shepard method for trivariate interpolation of scattered data`, ACM Trans. Math. Software (14), 151--152 Shepard, D, 1968, `A two-dimensional interpolation function for irregularly spaced data`, Proc. 23rd Nat. Conf. ACM, 517--523, Brandon/Systems Press Inc., Princeton """ raise NotImplementedError
[docs]def dim5_scat_shep_eval(x, f, iq, rq, xe): r""" ``dim5_scat_shep_eval`` evaluates the five-dimensional interpolating function generated by :meth:`dim5_scat_shep` and its first partial derivatives. .. _e01tn-py2-py-doc: For full information please refer to the NAG Library document for e01tn https://www.nag.com/numeric/nl/nagdoc_28.5/flhtml/e01/e01tnf.html .. _e01tn-py2-py-parameters: **Parameters** **x** : float, array-like, shape :math:`\left(5, m\right)` **Must** be the same array supplied as argument :math:`\textit{x}` in the preceding call to :meth:`dim5_scat_shep`. It **must** remain unchanged between calls. **f** : float, array-like, shape :math:`\left(m\right)` **Must** be the same array supplied as argument :math:`\textit{f}` in the preceding call to :meth:`dim5_scat_shep`. It **must** remain unchanged between calls. **iq** : int, array-like, shape :math:`\left(2\times m+1\right)` **Must** be the same array returned as argument :math:`\textit{iq}` in the preceding call to :meth:`dim5_scat_shep`. It **must** remain unchanged between calls. **rq** : float, array-like, shape :math:`\left(21\times m+11\right)` **Must** be the same array returned as argument :math:`\textit{rq}` in the preceding call to :meth:`dim5_scat_shep`. It **must** remain unchanged between calls. **xe** : float, array-like, shape :math:`\left(5, n\right)` :math:`\mathrm{xe}[0,\textit{i}-1],\ldots,\mathrm{xe}[4,\textit{i}-1]` must be set to the evaluation point :math:`\mathbf{x}_{\textit{i}}`, for :math:`\textit{i} = 1,2,\ldots,n`. **Returns** **q** : float, ndarray, shape :math:`\left(n\right)` :math:`\mathrm{q}[\textit{i}-1]` contains the value of the interpolant, at :math:`\mathbf{x}_{\textit{i}}`, for :math:`\textit{i} = 1,2,\ldots,n`. If any of these evaluation points lie outside the region of definition of the interpolant the corresponding entries in :math:`\mathrm{q}` are set to an extrapolated approximation, and ``dim5_scat_shep_eval`` returns with :math:`\mathrm{errno}` = 3. **qx** : float, ndarray, shape :math:`\left(5, n\right)` :math:`\mathrm{qx}[j-1,i-1]` contains the value of the partial derivatives with respect to :math:`\mathbf{x}_j` of the interpolant :math:`Q\left(\mathbf{x}\right)` at :math:`\mathbf{x}_{\textit{i}}`, for :math:`\textit{i} = 1,2,\ldots,n`, and for each of the five partial derivatives :math:`j = 1,2,3,4,5`. If any of these evaluation points lie outside the region of definition of the interpolant, the corresponding entries in :math:`\mathrm{qx}` are set to extrapolated approximations to the partial derivatives, and ``dim5_scat_shep_eval`` returns with :math:`\mathrm{errno}` = 3. .. _e01tn-py2-py-errors: **Raises** **NagValueError** (`errno` :math:`1`) On entry, :math:`n = \langle\mathit{\boldsymbol{value}}\rangle`. Constraint: :math:`n\geq 1`. (`errno` :math:`1`) On entry, :math:`m = \langle\mathit{\boldsymbol{value}}\rangle`. Constraint: :math:`m\geq 23`. (`errno` :math:`2`) On entry, values in :math:`\mathrm{rq}` appear to be invalid. Check that :math:`\mathrm{rq}` has not been corrupted between calls to :meth:`dim5_scat_shep` and ``dim5_scat_shep_eval``. (`errno` :math:`2`) On entry, values in :math:`\mathrm{iq}` appear to be invalid. Check that :math:`\mathrm{iq}` has not been corrupted between calls to :meth:`dim5_scat_shep` and ``dim5_scat_shep_eval``. **Warns** **NagAlgorithmicWarning** (`errno` :math:`3`) On entry, at least one evaluation point lies outside the region of definition of the interpolant. At such points the corresponding values in :math:`\mathrm{q}` and :math:`\mathrm{qx}` contain extrapolated approximations. Points should be evaluated one by one to identify extrapolated values. .. _e01tn-py2-py-notes: **Notes** ``dim5_scat_shep_eval`` takes as input the interpolant :math:`Q\left(\mathbf{x}\right)`, :math:`\mathbf{x} \in \mathbb{R}^5` of a set of scattered data points :math:`\left(\mathbf{x}_{\textit{r}}, f_{\textit{r}}\right)`, for :math:`\textit{r} = 1,2,\ldots,m`, as computed by :meth:`dim5_scat_shep`, and evaluates the interpolant and its first partial derivatives at the set of points :math:`\mathbf{x}_{\textit{i}}`, for :math:`\textit{i} = 1,2,\ldots,n`. ``dim5_scat_shep_eval`` must only be called after a call to :meth:`dim5_scat_shep`. ``dim5_scat_shep_eval`` is derived from the new implementation of QS3GRD described by Renka (1988). It uses the modification for five-dimensional interpolation described by Berry and Minser (1999). .. _e01tn-py2-py-references: **References** Berry, M W, Minser, K S, 1999, `Algorithm 798: high-dimensional interpolation using the modified Shepard method`, ACM Trans. Math. Software (25), 353--366 Renka, R J, 1988, `Algorithm 661: QSHEP3D: Quadratic Shepard method for trivariate interpolation of scattered data`, ACM Trans. Math. Software (14), 151--152 """ raise NotImplementedError
[docs]def dimn_grid(narr, uniform, axis, v, point, method, k=1, wf=0.0): r""" ``dimn_grid`` interpolates data at a point in :math:`n`-dimensional space, that is defined by a set of gridded data points. It offers three methods to interpolate the data: Linear Interpolation, Cubic Interpolation and Weighted Average. .. _e01za-py2-py-doc: For full information please refer to the NAG Library document for e01za https://www.nag.com/numeric/nl/nagdoc_28.5/flhtml/e01/e01zaf.html .. _e01za-py2-py-parameters: **Parameters** **narr** : int, array-like, shape :math:`\left(d\right)` The number of data ordinates in each dimension, with :math:`\mathrm{narr}[\textit{i}-1] = n_i`, for :math:`\textit{i} = 1,2,\ldots,d`. **uniform** : bool States whether the data points are uniformly spaced. :math:`\mathrm{uniform} = \mathbf{True}` The data points are uniformly spaced. :math:`\mathrm{uniform} = \mathbf{False}` The data points are not uniformly spaced. **axis** : float, array-like, shape :math:`\left(\textit{lx}\right)` Defines the axis. If the data points are uniformly spaced (see argument :math:`\mathrm{uniform}`) :math:`\mathrm{axis}` should contain the start and end of each dimension :math:`\left(x_{{1 1}}, x_{{1 n_1}}, \ldots, x_{{d 1}}, x_{{d n_d}}\right)`. If the data points are not uniformly spaced, :math:`\mathrm{axis}` should contain all the data ordinates for each dimension :math:`\left(x_{{1 1}}, x_{{1 2}}, \ldots, x_{{1 n_1}}, \ldots, x_{{d 1}}, x_{{d 2}}, \ldots, x_{{d n_d}}\right)`. **v** : float, array-like, shape :math:`\left(\mathrm{prod}\left(\mathrm{narr}\right)\right)` Holds the values of the data points in such an order that the index of a data value with coordinates :math:`\left(z_1, z_2, \ldots, z_d\right)` is .. math:: \sum_{1}^{d}{z_{\textit{i}}\prod_{{n \in \mathbf{S}_i}}n}\text{,} where :math:`\mathbf{S}_i = \left\{\mathrm{narr}[l-1]:l = \left. 1, \ldots, {i-1}\right. \right\}` e.g., :math:`\left(\left(x_{{11}}, x_{{21}}, \ldots, x_{{d1}}\right), \left(x_{{12}}, x_{{21}}, \ldots, x_{{d1}}\right), \ldots, \left(x_{{1n_d}}, x_{{21}}, \ldots, x_{{d1}}\right), \left(x_{{11}}, x_{{22}}, \ldots, x_{{d1}}\right), \left(x_{{12}}, x_{{22}}, \ldots, x_{{d1}}\right), \ldots, \left(x_{{1n_d}}, x_{{2n_d}}, \ldots, x_{{dn_d}}\right)\right)`. **point** : float, array-like, shape :math:`\left(d\right)` :math:`\mathbf{x}`, the point at which the data value is to be interpolated. **method** : int The method to be used. :math:`\mathrm{method} = 1` Weighted Average. :math:`\mathrm{method} = 2` Linear Interpolation. :math:`\mathrm{method} = 3` Cubic Interpolation. **k** : int, optional If :math:`\mathrm{method} = 1`, :math:`\mathrm{k}` controls the number of data points used in the Weighted Average method, with :math:`\mathrm{k}` points used in each dimension, either side of the interpolation point. The total number of data points used for the interpolation will, therefore, be :math:`\left(2\mathrm{k}\right)^d`. If :math:`\mathrm{method} \neq 1`, then :math:`\mathrm{k}` is not referenced and need not be set. **wf** : float, optional The power used for the weighted average such that a high power will cause closer points to be more heavily weighted. **Returns** **ans** : float Holds the result of the interpolation. .. _e01za-py2-py-errors: **Raises** **NagValueError** (`errno` :math:`1`) On entry, :math:`d = \langle\mathit{\boldsymbol{value}}\rangle`. Constraint: :math:`d\geq 2`. (`errno` :math:`2`) On entry, :math:`\mathrm{narr}[\langle\mathit{\boldsymbol{value}}\rangle] = \langle\mathit{\boldsymbol{value}}\rangle`. Constraint: :math:`\mathrm{narr}[\textit{i}-1]\geq 2`. (`errno` :math:`4`) On entry, :math:`\mathrm{axis}` decreases in dimension :math:`\langle\mathit{\boldsymbol{value}}\rangle`. Constraint: :math:`\mathrm{axis}` definition must be strictly increasing. (`errno` :math:`5`) On entry, :math:`\textit{lx} = \langle\mathit{\boldsymbol{value}}\rangle`, sum of :math:`\mathrm{narr}`:math:`= \langle\mathit{\boldsymbol{value}}\rangle`. Constraint: if :math:`\mathrm{uniform} = \mathbf{False}`, :math:`\textit{lx} =` sum of :math:`\mathrm{narr}`. (`errno` :math:`5`) On entry, :math:`\textit{lx} = \langle\mathit{\boldsymbol{value}}\rangle`, :math:`d = \langle\mathit{\boldsymbol{value}}\rangle`. Constraint: if :math:`\mathrm{uniform} = \mathbf{True}`, :math:`\textit{lx} = 2d`. (`errno` :math:`7`) On entry, :math:`\mathrm{k} = \langle\mathit{\boldsymbol{value}}\rangle`. Constraint: if :math:`\mathrm{method} = 1`, :math:`\mathrm{k}\geq 1`. (`errno` :math:`8`) On entry, :math:`\mathrm{point}[\langle\mathit{\boldsymbol{value}}\rangle] = \langle\mathit{\boldsymbol{value}}\rangle` and data range :math:`= \left[\langle\mathit{\boldsymbol{value}}\rangle, \langle\mathit{\boldsymbol{value}}\rangle\right]`. Constraint: :math:`\mathrm{point}` must be within the data range. (`errno` :math:`9`) On entry, :math:`\mathrm{method} = 3` and :math:`\mathrm{uniform} = \mathbf{False}`. Constraint: if :math:`\mathrm{method} = 3`, :math:`\mathrm{uniform}` must be :math:`\mathbf{True}`. (`errno` :math:`9`) On entry, :math:`\mathrm{method} = \langle\mathit{\boldsymbol{value}}\rangle`. Constraint: :math:`\mathrm{method} = 1`, :math:`2` or :math:`3`. (`errno` :math:`10`) On entry, :math:`\mathrm{wf} = \langle\mathit{\boldsymbol{value}}\rangle`. Constraint: if :math:`\mathrm{method} = 1`, :math:`1.0\leq \mathrm{wf}\leq 15.0`. (`errno` :math:`101`) Cubic Interpolation method does not have enough data surrounding :math:`\mathrm{point}`; interpolation not possible. **Warns** **NagAlgorithmicWarning** (`errno` :math:`201`) Warning: the size of :math:`\mathrm{k}` has been reduced, due to too few data points around :math:`\mathrm{point}`. .. _e01za-py2-py-notes: **Notes** ``dimn_grid`` interpolates an :math:`n`-dimensional point within a set of gridded data points, :math:`\mathbf{Z} = \left\{z_{{1j_1}}, z_{{2j_2}}, \ldots, z_{{dj_d}}\right\}`, with corresponding data values :math:`\mathbf{F} = \left\{f_{{1j_1}}, f_{{2j_2}}, \ldots, f_{{dj_d}}\right\}` where, for the :math:`i`\ th dimension, :math:`j_i = 1,\ldots,n_i` and :math:`n_i` is the number of ordinates in the :math:`i`\ th dimension. A hypercube of :math:`\left(2k\right)^d` data points :math:`\left[\mathbf{h}_1, \mathbf{h}_2, \ldots, \mathbf{h}_{{\left(2k\right)^d}}\right]⊂\mathbf{Z}`, where :math:`\mathbf{h}_i = \left(h_{{i1}}, h_{{i2}}, \ldots, h_{{id}}\right)` and the corresponding data values are :math:`f\left(\mathbf{h}_i\right) \in \mathbf{F}`, around the given point, :math:`\mathbf{x} = \left(x_1, x_2, \ldots, x_d\right)`, is found and then used to interpolate using one of the following three methods. (i) Weighted Average, that is a modification of Shepard's method (Shepard (1968)) as used for scattered data in :meth:`dimn_scat_shep`. This method interpolates the data with the weighted mean .. math:: Q\left(\mathbf{x}\right) = \frac{{\sum_{{r = 1}}^{\left(2k\right)^d}w_r\left(\mathbf{x}\right)f_r}}{{\sum_{{r = 1}}^{\left(2k\right)^d}w_r\left(\mathbf{x}\right)}}\text{,} where :math:`f_r = f\left(\mathbf{h}_r\right)`, :math:`w_r\left(\mathbf{x}\right) = \frac{1}{{D\left(\left\lvert \mathbf{x}-\mathbf{h}_r\right\rvert \right)}}` and :math:`D\left(\mathbf{y}\right) = \left. y_1^{\rho }+y_2^{\rho } + \cdots +y_d^{\rho }\right.`, for a given value of :math:`\rho`. (#) Linear Interpolation, which takes :math:`2^d` surrounding data points (:math:`k = 1`) and performs two one-dimensional linear interpolations in each dimension on data points :math:`\mathbf{h}_a` and :math:`\mathbf{h}_b`, reducing the dimension every iteration until it has reached an answer. The formula for linear interpolation in dimension :math:`i` is simply .. math:: f = f_a+\left(x_i-h_{{ai}}\right)\frac{{f_b-f_a}}{{h_{{bi}}-h_{{ai}}}}\text{,} where :math:`f_r = f\left(\mathbf{h}_r\right)` and :math:`h_{{ai}} < x_i < h_{{bi}}`. (#) Cubic Interpolation, based on cubic convolution (Keys (1981)). In a similar way to the Linear Interpolation method, it performs the interpolations in one dimension reducing it each time, however it requires four surrounding data points in each dimension (:math:`k = 2`), two in each direction :math:`\left(\mathbf{h}_{-1}, \mathbf{h}_0, \mathbf{h}_1, \mathbf{h}_2\right)`. The following is used to calculate the one-dimensional interpolant in dimension :math:`i` .. math:: f = \frac{1}{2}\begin{pmatrix}1&t&t^2&t^3\end{pmatrix}\begin{pmatrix}0&2&0&0\\-1&0&1&0\\2&-5&4&-1\\-1&3&-3&1\end{pmatrix}\begin{pmatrix}f_{-1}\\f_0\\f_1\\f_2\end{pmatrix} where :math:`t = x_i-h_{{0i}}` and :math:`f_r = f\left(\mathbf{h}_r\right)`. .. _e01za-py2-py-references: **References** Keys, R, 1981, `Cubic Convolution Interpolation for Digital Image Processing`, IEEE Transactions on Acoutstics, Speech, and Signal Processing (Vol ASSP-29 No. 6), 1153--1160, http://hmi.stanford.edu/doc/Tech_Notes/HMI-TN-2004-004-Interpolation/Keys_cubic_interp.pdf Shepard, D, 1968, `A two-dimensional interpolation function for irregularly spaced data`, Proc. 23rd Nat. Conf. ACM, 517--523, Brandon/Systems Press Inc., Princeton """ raise NotImplementedError
[docs]def dimn_scat_shep(x, f, nw=-1, nq=-1): r""" ``dimn_scat_shep`` generates a multidimensional interpolant to a set of scattered data points, using a modified Shepard method. When the number of dimensions is no more than five, there are corresponding functions in submodule ``interp`` which are specific to the given dimensionality. :meth:`dim2_scat_shep` generates the two-dimensional interpolant, while :meth:`dim3_scat_shep`, :meth:`dim4_scat_shep` and :meth:`dim5_scat_shep` generate the three-, four - and five-dimensional interpolants respectively. .. _e01zm-py2-py-doc: For full information please refer to the NAG Library document for e01zm https://www.nag.com/numeric/nl/nagdoc_28.5/flhtml/e01/e01zmf.html .. _e01zm-py2-py-parameters: **Parameters** **x** : float, array-like, shape :math:`\left(d, m\right)` The :math:`\textit{d}` components of the first data point must be stored in elements :math:`0,1,\ldots,d-1` of :math:`\mathrm{x}`. The second data point must be stored in elements :math:`d,d+1,\ldots,2\times d-1` of :math:`\mathrm{x}`, and so on. In general, the :math:`\textit{m}` data points must be stored in :math:`\mathrm{x}[\textit{j},\textit{i}]`, for :math:`\textit{j} = 0,1,\ldots,d-1`, for :math:`\textit{i} = 0,1,\ldots,m-1`. **f** : float, array-like, shape :math:`\left(m\right)` :math:`\mathrm{f}[r-1]` must be set to the data value :math:`f_{\textit{r}}`, for :math:`\textit{r} = 1,2,\ldots,m`. **nw** : int, optional The number :math:`N_w` of data points that determines each radius of influence :math:`R_w`, appearing in the definition of each of the weights :math:`w_{\textit{r}}`, for :math:`\textit{r} = 1,2,\ldots,m` (see :ref:`Notes <e01zm-py2-py-notes>`). Note that :math:`R_w` is different for each weight. If :math:`\mathrm{nw}\leq 0` the default value :math:`\mathrm{nw} = \mathrm{min}\left({2\times \left(d+1\right)\times \left(d+2\right)}, {m-1}\right)` is used instead. **nq** : int, optional The number :math:`N_q` of data points to be used in the least squares fit for coefficients defining the quadratic functions :math:`q_r\left(\mathbf{x}\right)` (see :ref:`Notes <e01zm-py2-py-notes>`). If :math:`\mathrm{nq}\leq 0` the default value :math:`\mathrm{nq} = \mathrm{min}\left({\left(d+1\right)\times \left(d+2\right)\times 6/5}, {m-1}\right)` is used instead. **Returns** **iq** : int, ndarray, shape :math:`\left(2\times m+1\right)` Integer data defining the interpolant :math:`Q\left(\mathbf{x}\right)`. **rq** : float, ndarray, shape :math:`\left(\textit{lrq}\right)` Real data defining the interpolant :math:`Q\left(\mathbf{x}\right)`. .. _e01zm-py2-py-errors: **Raises** **NagValueError** (`errno` :math:`1`) On entry, :math:`\mathrm{nw} = \langle\mathit{\boldsymbol{value}}\rangle` and :math:`m = \langle\mathit{\boldsymbol{value}}\rangle`. Constraint: :math:`\mathrm{nw}\leq m-1`. (`errno` :math:`1`) On entry, :math:`\mathrm{nq} = \langle\mathit{\boldsymbol{value}}\rangle` and :math:`m = \langle\mathit{\boldsymbol{value}}\rangle`. Constraint: :math:`\mathrm{nq}\leq m-1`. (`errno` :math:`1`) On entry, :math:`\mathrm{nq} = \langle\mathit{\boldsymbol{value}}\rangle` and :math:`d = \langle\mathit{\boldsymbol{value}}\rangle`. Constraint: :math:`\mathrm{nq}\leq 0` or :math:`\mathrm{nq}\geq \left(d+1\right)\times \left(d+2\right)/2-1`. (`errno` :math:`1`) On entry, :math:`m = \langle\mathit{\boldsymbol{value}}\rangle` and :math:`d = \langle\mathit{\boldsymbol{value}}\rangle`. Constraint: :math:`m\geq \left(d+1\right)\times \left(d+2\right)/2+2`. (`errno` :math:`1`) On entry, :math:`\left(\left(d+1\right)\times \left(d+2\right)/2\right)\times m+2\times d+1` exceeds the largest machine integer. :math:`d = \langle\mathit{\boldsymbol{value}}\rangle` and :math:`m = \langle\mathit{\boldsymbol{value}}\rangle`. (`errno` :math:`1`) On entry, :math:`d = \langle\mathit{\boldsymbol{value}}\rangle`. Constraint: :math:`d\geq 2`. (`errno` :math:`2`) There are duplicate nodes in the dataset. :math:`{\mathrm{x}[k-1,i-1]} = {\mathrm{x}[k-1,j-1]}`, for :math:`i = \langle\mathit{\boldsymbol{value}}\rangle`, :math:`j = \langle\mathit{\boldsymbol{value}}\rangle` and :math:`k = 1,2,\ldots,d`. The interpolant cannot be derived. (`errno` :math:`3`) On entry, all the data points lie on the same hypersurface. No unique solution exists. .. _e01zm-py2-py-notes: **Notes** ``dimn_scat_shep`` constructs a smooth function :math:`Q\left(\mathbf{x}\right)`, :math:`\mathbf{x} \in \mathbb{R}^d` which interpolates a set of :math:`m` scattered data points :math:`\left(\mathbf{x}_r, f_r\right)`, for :math:`r = 1,2,\ldots,m`, using a modification of Shepard's method. The surface is continuous and has continuous first partial derivatives. The basic Shepard method, which is a generalization of the two-dimensional method described in Shepard (1968), interpolates the input data with the weighted mean .. math:: Q\left(\mathbf{x}\right) = \frac{{\sum_{{r = 1}}^mw_r\left(\mathbf{x}\right)q_r}}{{\sum_{{r = 1}}^mw_r\left(\mathbf{x}\right)}}\text{,} where :math:`q_r = f_r`, :math:`w_r\left(\mathbf{x}\right) = \frac{1}{\left\lVert \mathbf{x}-\mathbf{x}_r\right\rVert_2^2}`. The basic method is global in that the interpolated value at any point depends on all the data, but ``dimn_scat_shep`` uses a modification (see Franke and Nielson (1980) and Renka (1988a)), whereby the method becomes local by adjusting each :math:`w_r\left(\mathbf{x}\right)` to be zero outside a hypersphere with centre :math:`\mathbf{x}_r` and some radius :math:`R_w`. Also, to improve the performance of the basic method, each :math:`q_r` above is replaced by a function :math:`q_r\left(\mathbf{x}\right)`, which is a quadratic fitted by weighted least squares to data local to :math:`\mathbf{x}_r` and forced to interpolate :math:`\left(\mathbf{x}_r, f_r\right)`. In this context, a point :math:`\mathbf{x}` is defined to be local to another point if it lies within some distance :math:`R_q` of it. The efficiency of ``dimn_scat_shep`` is enhanced by using a cell method for nearest neighbour searching due to Bentley and Friedman (1979) with a cell density of :math:`3`. The radii :math:`R_w` and :math:`R_q` are chosen to be just large enough to include :math:`N_w` and :math:`N_q` data points, respectively, for user-supplied constants :math:`N_w` and :math:`N_q`. Default values of these parameters are provided, and advice on alternatives is given in `Further Comments <https://www.nag.com/numeric/nl/nagdoc_28.5/flhtml/e01/e01zmf.html#fcomments2>`__. ``dimn_scat_shep`` is derived from the new implementation of QSHEP3 described by Renka (1988b). It uses the modification for high-dimensional interpolation described by Berry and Minser (1999). Values of the interpolant :math:`Q\left(\mathbf{x}\right)` generated by ``dimn_scat_shep``, and its first partial derivatives, can subsequently be evaluated for points in the domain of the data by a call to :meth:`dimn_scat_shep_eval`. .. _e01zm-py2-py-references: **References** Bentley, J L and Friedman, J H, 1979, `Data structures for range searching`, ACM Comput. Surv. (11), 397--409 Berry, M W, Minser, K S, 1999, `Algorithm 798: high-dimensional interpolation using the modified Shepard method`, ACM Trans. Math. Software (25), 353--366 Franke, R and Nielson, G, 1980, `Smooth interpolation of large sets of scattered data`, Internat. J. Num. Methods Engrg. (15), 1691--1704 Renka, R J, 1988, `Multivariate interpolation of large sets of scattered data`, ACM Trans. Math. Software (14), 139--148 Renka, R J, 1988, `Algorithm 661: QSHEP3D: Quadratic Shepard method for trivariate interpolation of scattered data`, ACM Trans. Math. Software (14), 151--152 Shepard, D, 1968, `A two-dimensional interpolation function for irregularly spaced data`, Proc. 23rd Nat. Conf. ACM, 517--523, Brandon/Systems Press Inc., Princeton """ raise NotImplementedError
[docs]def dimn_scat_shep_eval(x, f, iq, rq, xe): r""" ``dimn_scat_shep_eval`` evaluates the multidimensional interpolating function generated by :meth:`dimn_scat_shep` and its first partial derivatives. .. _e01zn-py2-py-doc: For full information please refer to the NAG Library document for e01zn https://www.nag.com/numeric/nl/nagdoc_28.5/flhtml/e01/e01znf.html .. _e01zn-py2-py-parameters: **Parameters** **x** : float, array-like, shape :math:`\left(d, m\right)` Note: the :math:`i`\ th ordinate of the point :math:`x_j` is stored in :math:`\mathrm{x}[i-1,j-1]`. **Must** be the same array supplied as argument :math:`\textit{x}` in the preceding call to :meth:`dimn_scat_shep`. It **must** remain unchanged between calls. **f** : float, array-like, shape :math:`\left(m\right)` **Must** be the same array supplied as argument :math:`\textit{f}` in the preceding call to :meth:`dimn_scat_shep`. It **must** remain unchanged between calls. **iq** : int, array-like, shape :math:`\left(2\times m+1\right)` **Must** be the same array returned as argument :math:`\textit{iq}` in the preceding call to :meth:`dimn_scat_shep`. It **must** remain unchanged between calls. **rq** : float, array-like, shape :math:`\left(\textit{lrq}\right)` **Must** be the same array returned as argument :math:`\textit{rq}` in the preceding call to :meth:`dimn_scat_shep`. It **must** remain unchanged between calls. **xe** : float, array-like, shape :math:`\left(d, n\right)` Note: the :math:`i`\ th ordinate of the point :math:`x_j` is stored in :math:`\mathrm{xe}[i-1,j-1]`. :math:`\mathrm{xe}[0,\textit{j}-1],\ldots,\mathrm{xe}[d-1,\textit{j}-1]` must be set to the evaluation point :math:`\mathbf{x}_{\textit{j}}`, for :math:`\textit{j} = 1,2,\ldots,n`. **Returns** **q** : float, ndarray, shape :math:`\left(n\right)` :math:`\mathrm{q}[\textit{i}-1]` contains the value of the interpolant, at :math:`\mathbf{x}_{\textit{i}}`, for :math:`\textit{i} = 1,2,\ldots,n`. If any of these evaluation points lie outside the region of definition of the interpolant the corresponding entries in :math:`\mathrm{q}` are set to an extrapolated approximation, and ``dimn_scat_shep_eval`` returns with :math:`\mathrm{errno}` = 3. **qx** : float, ndarray, shape :math:`\left(d, n\right)` :math:`\mathrm{qx}[i-1,j-1]` contains the value of the partial derivatives with respect to the :math:`i`\ th independent variable (dimension) of the interpolant :math:`Q\left(\mathbf{x}\right)` at :math:`\mathbf{x}_{\textit{j}}`, for :math:`\textit{j} = 1,2,\ldots,n`, and for each of the partial derivatives :math:`i = 1,2,\ldots,d`. If any of these evaluation points lie outside the region of definition of the interpolant, the corresponding entries in :math:`\mathrm{qx}` are set to extrapolated approximations to the partial derivatives, and ``dimn_scat_shep_eval`` returns with :math:`\mathrm{errno}` = 3. .. _e01zn-py2-py-errors: **Raises** **NagValueError** (`errno` :math:`1`) On entry, :math:`n = \langle\mathit{\boldsymbol{value}}\rangle`. Constraint: :math:`n\geq 1`. (`errno` :math:`1`) On entry, :math:`m = \langle\mathit{\boldsymbol{value}}\rangle` and :math:`d = \langle\mathit{\boldsymbol{value}}\rangle`. Constraint: :math:`m\geq \left(d+1\right)\times \left(d+2\right)/2+2`. (`errno` :math:`1`) On entry, :math:`\left(\left(d+1\right)\times \left(d+2\right)/2\right)\times m+2\times d+1` exceeds the largest machine integer. :math:`d = \langle\mathit{\boldsymbol{value}}\rangle` and :math:`m = \langle\mathit{\boldsymbol{value}}\rangle`. (`errno` :math:`1`) On entry, :math:`d = \langle\mathit{\boldsymbol{value}}\rangle`. Constraint: :math:`d\geq 2`. (`errno` :math:`2`) On entry, values in :math:`\mathrm{rq}` appear to be invalid. Check that :math:`\mathrm{rq}` has not been corrupted between calls to :meth:`dimn_scat_shep` and ``dimn_scat_shep_eval``. (`errno` :math:`2`) On entry, values in :math:`\mathrm{iq}` appear to be invalid. Check that :math:`\mathrm{iq}` has not been corrupted between calls to :meth:`dimn_scat_shep` and ``dimn_scat_shep_eval``. **Warns** **NagAlgorithmicWarning** (`errno` :math:`3`) On entry, at least one evaluation point lies outside the region of definition of the interpolant. At such points the corresponding values in :math:`\mathrm{q}` and :math:`\mathrm{qx}` contain extrapolated approximations. Points should be evaluated one by one to identify extrapolated values. .. _e01zn-py2-py-notes: **Notes** ``dimn_scat_shep_eval`` takes as input the interpolant :math:`Q\left(\mathbf{x}\right)`, :math:`\mathbf{x} \in \mathbb{R}^d` of a set of scattered data points :math:`\left(\mathbf{x}_{\textit{r}}, f_{\textit{r}}\right)`, for :math:`\textit{r} = 1,2,\ldots,m`, as computed by :meth:`dimn_scat_shep`, and evaluates the interpolant and its first partial derivatives at the set of points :math:`\mathbf{x}_{\textit{i}}`, for :math:`\textit{i} = 1,2,\ldots,n`. ``dimn_scat_shep_eval`` must only be called after a call to :meth:`dimn_scat_shep`. ``dimn_scat_shep_eval`` is derived from the new implementation of QS3GRD described by Renka (1988). It uses the modification for high-dimensional interpolation described by Berry and Minser (1999). .. _e01zn-py2-py-references: **References** Berry, M W, Minser, K S, 1999, `Algorithm 798: high-dimensional interpolation using the modified Shepard method`, ACM Trans. Math. Software (25), 353--366 Renka, R J, 1988, `Algorithm 661: QSHEP3D: Quadratic Shepard method for trivariate interpolation of scattered data`, ACM Trans. Math. Software (14), 151--152 """ raise NotImplementedError