[oneMKL] add Data Fitting domain to oneMKL

source/elements/oneMKL/source/architecture/api_design.inc.rst

@@ -16,27 +16,28 @@ oneMKL namespaces
 
 The oneMKL library uses C++ namespaces to organize routines by mathematical domain.  All oneMKL objects and routines shall be contained within the ``oneapi::mkl`` base namespace.  The individual oneMKL domains use a secondary namespace layer as follows:
 
-========================  =======================================================================================================
-namespace                 oneMKL domain or content
-========================  =======================================================================================================
-``oneapi::mkl``           oneMKL base namespace, contains general oneMKL data types, objects, exceptions and routines
-``oneapi::mkl::blas``     Dense linear algebra routines from BLAS and BLAS like extensions. The oneapi::mkl::blas namespace should contain two namespaces column_major and row_major to support both matrix layouts. See :ref:`onemkl_blas`
-``oneapi::mkl::lapack``   Dense linear algebra routines from LAPACK and LAPACK like extensions. See :ref:`onemkl_lapack`
-``oneapi::mkl::sparse``   Sparse linear algebra routines from Sparse BLAS and Sparse Solvers. See :ref:`onemkl_sparse_linear_algebra`
-``oneapi::mkl::dft``      Discrete and fast Fourier transformations. See :ref:`onemkl_dft`
-``oneapi::mkl::rng``      Random number generator routines. See :ref:`onemkl_rng`
-``oneapi::mkl::vm``       Vector mathematics routines, e.g. trigonometric, exponential functions acting on elements of a vector. See :ref:`onemkl_vm`
-========================  =======================================================================================================
+=========================================== =========================================================================
+namespace                                   oneMKL domain or content
+=========================================== =========================================================================
+``oneapi::mkl``                             oneMKL base namespace, contains general oneMKL data types, objects, exceptions and routines
+``oneapi::mkl::blas``                       Dense linear algebra routines from BLAS and BLAS like extensions. The oneapi::mkl::blas namespace should contain two namespaces column_major and row_major to support both matrix layouts. See :ref:`onemkl_blas`
+``oneapi::mkl::lapack``                     Dense linear algebra routines from LAPACK and LAPACK like extensions. See :ref:`onemkl_lapack`
+``oneapi::mkl::sparse``                     Sparse linear algebra routines from Sparse BLAS and Sparse Solvers. See :ref:`onemkl_sparse_linear_algebra`
+``oneapi::mkl::dft``                        Discrete and fast Fourier transformations. See :ref:`onemkl_dft`
+``oneapi::mkl::rng``                        Random number generator routines. See :ref:`onemkl_rng`
+``oneapi::mkl::vm``                         Vector mathematics routines, e.g. trigonometric, exponential functions acting on elements of a vector. See :ref:`onemkl_vm`
+``oneapi::mkl::experimental::data_fitting`` Data fitting routines, e.g. interpolation. See :ref:`data_fitting`
+=========================================== =========================================================================
 
 .. note::
    :name: Implementation Requirement
 
-   Inside each oneMKL domain, there are many routines, classes, enums and objects defined which constitute the breadth and scope of that oneMKL domain.  
+   Inside each oneMKL domain, there are many routines, classes, enums and objects defined which constitute the breadth and scope of that oneMKL domain.
    It is permitted for a library implementation of the oneMKL specification to implement either all, one or more than one of the domains in oneMKL. However, within an implementation of a specific domain, all relevant routines, classes, enums and objects (including those relevant enums and objects which live outside a particular domain in the general ``oneapi::mkl`` namespace must be both declared and defined in the library so that an application that uses that domain could build and link against that library implementation successfully.
 
-   It is however acceptable to throw the runtime exception :ref:`oneapi::mkl::unimplemented<onemkl_exception_unimplemented>` inside of the routines or class member functions in that domain that have not been fully implemented.  
-   For instance, a library may choose to implement the oneMKL BLAS functionality and in particular may choose to implement only the :ref:`onemkl_blas_gemm` api for their library, in which case they must also include all the other blas namespaced routines and throw the :ref:`oneapi::mkl::unimplemented<onemkl_exception_unimplemented>` exception inside all the others.  
-   
+   It is however acceptable to throw the runtime exception :ref:`oneapi::mkl::unimplemented<onemkl_exception_unimplemented>` inside of the routines or class member functions in that domain that have not been fully implemented.
+   For instance, a library may choose to implement the oneMKL BLAS functionality and in particular may choose to implement only the :ref:`onemkl_blas_gemm` api for their library, in which case they must also include all the other blas namespaced routines and throw the :ref:`oneapi::mkl::unimplemented<onemkl_exception_unimplemented>` exception inside all the others.
+
    In such a case, the implemented routines in such a library should be communicated clearly and easily understood by users of that library.
 
 
@@ -73,7 +74,7 @@ oneMKL uses the following DPC++ data types:
 oneMKL defined datatypes
 ++++++++++++++++++++++++
 
-oneMKL dense and sparse linear algebra routines use scoped enum types as type-safe replacements for the traditional character arguments used in C/Fortran implementations of BLAS and LAPACK. These types all belong to the ``oneapi::mkl`` namespace.  
+oneMKL dense and sparse linear algebra routines use scoped enum types as type-safe replacements for the traditional character arguments used in C/Fortran implementations of BLAS and LAPACK. These types all belong to the ``oneapi::mkl`` namespace.
 
 Each enumeration value comes with two names: A single-character name (the traditional BLAS/LAPACK character) and a longer, more descriptive name. The two names are exactly equivalent and may be used interchangeably.
 
@@ -238,7 +239,7 @@ Each enumeration value comes with two names: A single-character name (the tradit
          :name: layout
          :class: sectiontitle
 
-      The ``layout`` type specifies how a dense matrix ``A`` with leading dimension ``lda`` is stored as one dimensional array in memory. 
+      The ``layout`` type specifies how a dense matrix ``A`` with leading dimension ``lda`` is stored as one dimensional array in memory.
       The layouts are traditionally provided in one of two forms: C/C++-style using ``row_major`` layout,
       or Fortran-style using ``column_major`` layout. The ``layout`` type can take the following values:
 

source/elements/oneMKL/source/domains/data_fitting/cubic.rst

@@ -0,0 +1,88 @@
+.. SPDX-FileCopyrightText: 2022 Intel Corporation
+..
+.. SPDX-License-Identifier: CC-BY-4.0
+
+.. _cubic:
+
+Cubic Splines
+=============
+
+Cubic splines are splines whose degree is equal to 3.
+
+Cubic splines are described by the following polynomial
+
+.. math::
+  P_i\left( x \right) = c_{1,i}+ c_{2,i}\left( x - x_i \right) + c_{3,i}{\left( x - x_i \right)}^2+ c_{4,i}{\left( x - x_i \right)}^3,
+
+where
+
+.. math::
+  x \in \left[ x_i, x_{i+1} \right),
+
+.. math::
+  i = 1,\dots , n-1.
+
+There are many different kinds of cubic splines: Hermite, natural, Akima, Bessel.
+However, the current version of DPC++ API supports only one type: Hermite.
+
+Header File
+-----------
+
+.. code:: cpp
+
+  #include<oneapi/mkl/experimental/data_fitting.hpp>
+
+Namespace
+---------
+
+.. code:: cpp
+
+  oneapi::mkl::experimental::data_fitiing
+
+Hermite Spline
+--------------
+
+The coefficients of Hermite spline are calculated using the following formulas:
+
+.. math::
+  c_{1,i} = f\left( x_i \right),
+
+.. math::
+  c_{2,i} = s_i,
+
+.. math::
+  c_{3,i} = \left( \left[ x_i, x_{i+1} \right]f - s_i \right)  / \left( \Delta x_i \right) - c_{4,i}\left( \Delta x_i \right),
+
+.. math::
+  c_{4,i} = \left( s_i + s_{i+1} - 2\left[ x_i, x_{i+1} \right]f \right) / {\left( \Delta x_i \right)}^2,
+
+.. math::
+  s_i = f^{\left( 1 \right)}\left( x_i \right).
+
+The following boundary conditions are supported for Hermite spline:
+
+ - Free end (:math:`f^{(2)}(x_1) = f^{(2)}(x_n) = 0`).
+ - Periodic.
+ - First derivative.
+ - Second Derivative.
+
+Syntax
+^^^^^^
+
+.. code:: cpp
+
+  namespace cubic_spline {
+    struct hermite {};
+  }
+
+Example
+^^^^^^^
+
+To create a cubic Hermite spline object use the following:
+
+.. code:: cpp
+
+  spline<float, cubic_spline::hermite> val(
+    /*SYCL queue object*/q,
+    /*number of spline functions*/ny
+  );

source/elements/oneMKL/source/domains/data_fitting/data-fitting.inc.rst

@@ -0,0 +1,33 @@
+.. SPDX-FileCopyrightText: 2022 Intel Corporation
+..
+.. SPDX-License-Identifier: CC-BY-4.0
+
+.. _data_fitting:
+
+Data Fitting
+============
+
+oneMKL provides spline-based interpolation capabilities that can be used for
+spline construction (Linear, Cubic, Quadratic etc.),
+to perform cell-search operations, and to approximate functions,
+function derivatives, or integrals.
+
+APIs are experimental. It means that no API or ABI backward compatibility are guaranteed.
+
+APIs are based on SYCL USM (Unfied Shared Memory) input data.
+
+Routines
+--------
+
+Splines:
+    :ref:`splines`
+Interpolate function:
+    :ref:`interpolate`
+
+.. toctree::
+   :maxdepth: 1
+   :hidden:
+
+   data_fitting/terms
+   data_fitting/splines
+   data_fitting/interpolate

source/elements/oneMKL/source/domains/data_fitting/interpolate.rst

@@ -0,0 +1,107 @@
+.. SPDX-FileCopyrightText: 2022 Intel Corporation
+..
+.. SPDX-License-Identifier: CC-BY-4.0
+
+.. _interpolate:
+
+Interpolate Function
+====================
+
+Interpolate function performs computations of function and derivatives values at interpolation sites.
+
+If the sites do not belong to interpolation interval ``[a, b)``, the library uses:
+
+  - interpolant :math:`I_0` coefficients computed for interval :math:`[x_0, x_1)` for the
+    computations at the sites to the left of ``a``.
+  - interpolant :math:`I_{n-2}` coefficients computed for interval
+    :math:`[x_{n-2}, x_{n-1})` for the computations at the sites to the right of ``b``.
+
+Interpolation algorithm depends on interpolant's type (e.g., for a cubic spline,
+the evaluation of a third-order polynomial is performed to obtain the resulting interpolant value).
+
+Header File
+-----------
+
+.. code:: cpp
+
+  #include<oneapi/mkl/experimental/data_fitting.hpp>
+
+Namespace
+---------
+
+.. code:: cpp
+
+  oneapi::mkl::experimental::data_fitting
+
+Syntax
+------
+
+.. code:: cpp
+
+  template <typename Interpolant>
+  sycl::event interpolate(
+      Interpolant& interpolant,
+      typename Interpolant::value_type* sites,
+      std::int64_t n_sites,
+      typename Interpolant::value_type* results,
+      const std::vector<sycl::event>& dependencies,
+      interpolate_hint ResultHint = interpolate_hint::funcs_sites_ders,
+      site_hint SiteHint = site_hint::non_uniform); // (1)
+
+  template <typename Interpolant>
+  sycl::event interpolate(
+      Interpolant& interpolant,
+      typename Interpolant::value_type* sites,
+      std::int64_t n_sites,
+      typename Interpolant::value_type* results,
+      std::bitset<32> derivatives_indicator,
+      const std::vector<sycl::event>& dependencies = {},
+      interpolate_hint ResultHint = interpolate_hint::funcs_sites_ders,
+      site_hint SiteHint = site_hint::non_uniform); // (2)
+
+  template <typename Interpolant>
+  sycl::event interpolate(
+      sycl::queue& queue,
+      const Interpolant& interpolant,
+      typename Interpolant::value_type* sites,
+      std::int64_t n_sites,
+      typename Interpolant::value_type* results,
+      const std::vector<sycl::event>& dependencies,
+      interpolate_hint ResultHint = interpolate_hint::funcs_sites_ders,
+      site_hint SiteHint = site_hint::non_uniform); // (3)
+
+  template <typename Interpolant>
+  sycl::event interpolate(
+      sycl::queue& queue,
+      const Interpolant& interpolant,
+      typename Interpolant::value_type* sites,
+      std::int64_t n_sites,
+      typename Interpolant::value_type* results,
+      std::bitset<32> derivatives_indicator,
+      const std::vector<sycl::event>& dependencies = {},
+      interpolate_hint ResultHint = interpolate_hint::funcs_sites_ders,
+      site_hint SiteHint = site_hint::non_uniform); // (4)
+
+For all functions users can provide ``SiteHint`` and ``ResultHint`` to specify
+the layout of ``sites`` and ``results`` respectively.
+If ``results`` layout doesn't satisfy ``ResultHint`` and/or
+``sites`` layout doesn't satisfy ``SiteHint``, behavior is undefined.
+Returns the SYCL event of the submitted task.
+
+Interfaces (3) and (4) can be useful when users try to interpolate over different parts of ``sites``
+using different ``queue`` objects.
+
+#. Performs computations of function values only using the SYCL queue
+   associated with ``interpolant``.
+#. Performs computations of certain derivatives
+   (function values is considered as a zero derivative) which are indicated in
+   ``derivatives_indicator`` (each bit corresponds to certain derivative starting from lower bit)
+   using the SYCL queue associated with ``interpolant``.
+#. Performs computations of function values only using ``queue`` as an input argument
+   that should be created from the same context and device as the SYCL queue
+   associated with ``interpolant``.
+#. Performs computations of certain derivatives
+   (function values is considered as a zero derivative) which are indicated in
+   ``derivatives_indicator`` (each bit corresponds to certain derivative starting from lower bit)
+   using ``queue`` as an input argument that should be created from
+   the same context and device as the SYCL queue associated with ``interpolant``.

source/elements/oneMKL/source/domains/data_fitting/linear.rst

@@ -0,0 +1,64 @@
+.. SPDX-FileCopyrightText: 2022 Intel Corporation
+..
+.. SPDX-License-Identifier: CC-BY-4.0
+
+.. _linear:
+
+Linear Spline
+=============
+
+Linear spline is a spline whose degree is equal to 1.
+
+It's described by the following polynomial
+
+.. math::
+  P_i\left( x \right) = c_{1,i} + c_{2,i}\left( x - x_i \right),
+
+where
+
+.. math::
+  x \in \left[ x_i, x_{i+1} \right),
+
+.. math::
+  c_{1,i} = f\left( x_i \right),
+
+.. math::
+  c_{2,i} = \left[ x_i, x_{i+1} \right]f,
+
+.. math::
+  i = 1, \dots, n-1.
+
+Header File
+-----------
+
+.. code:: cpp
+
+  #include<oneapi/mkl/experimental/data_fitting.hpp>
+
+Namespace
+---------
+
+.. code:: cpp
+
+  oneapi::mkl::experimental::data_fitiing
+
+Syntax
+------
+
+.. code:: cpp
+
+  namespace linear_spline {
+    struct default_type {};
+  }
+
+Example
+-------
+
+To create a linear spline object use the following:
+
+.. code:: cpp
+
+  spline<float, linear_spline::default_type> val(
+    /*SYCL queue object*/q,
+    /*number of spline functions*/ny
+  );