Developer Reference for Intel® oneAPI Math Kernel Library for C

ID 766684
Date 11/07/2023
Public

A newer version of this document is available. Customers should click here to go to the newest version.

Document Table of Contents

Data Fitting Usage Examples

The examples below illustrate several operations that you can perform with Data Fitting routines.

You can get source code for similar examples in the .\examples\datafittingcsubdirectory of the Intel® oneAPI Math Kernel Library (oneMKL) installation directory.

The following example demonstrates the construction of a linear spline using Data Fitting routines. The spline approximates a scalar function defined on non-uniform partition. The coefficients of the spline are returned as a one-dimensional array:

Example of Linear Spline Construction

#include "mkl.h"
#define N 500                     /* Size of partition, number of breakpoints */
#define SPLINE_ORDER DF_PP_LINEAR /* Linear spline to construct */

int main()
{
    int status;          /* Status of a Data Fitting operation */
    DFTaskPtr task;      /* Data Fitting operations are task based */

    /* Parameters describing the partition */
    MKL_INT nx;          /* The size of partition x */
    double x[N];         /* Partition x */
    MKL_INT xhint;       /* Additional information about the structure of breakpoints */

    /* Parameters describing the function */
    MKL_INT ny;          /* Function dimension */
    double y[N];         /* Function values at the breakpoints */
    MKL_INT yhint;       /* Additional information about the function */
    
    /* Parameters describing the spline */
    MKL_INT  s_order;    /* Spline order */
    MKL_INT  s_type;     /* Spline type */
    MKL_INT  ic_type;    /* Type of internal conditions */
    double* ic;         /* Array of internal conditions */
    MKL_INT  bc_type;    /* Type of boundary conditions */
    double* bc;         /* Array of boundary conditions */

    double scoeff[(N-1)* SPLINE_ORDER];   /* Array of spline coefficients */
    MKL_INT scoeffhint;            /* Additional information about the coefficients */

    /* Initialize the partition */
    nx = N;
    /* Set values of partition x */
    ...
    xhint = DF_NO_HINT;    /* No additional information about the function is provided.
                              By default, the partition is non-uniform. */
    /* Initialize the function */
     ny = 1;               /* The function is scalar. */
   
    /* Set function values */
    ...
    yhint = DF_NO_HINT;    /* No additional information about the function is provided. */

    /* Create a Data Fitting task */
    status = dfdNewTask1D( &task, nx, x, xhint, ny, y, yhint );

    /* Check the Data Fitting operation status */
    ...

    /* Initialize spline parameters */
    s_order = DF_PP_LINEAR;    /* Spline is of the second order. */ 
    s_type = DF_PP_DEFAULT;    /* Spline is of the default type. */ 

    /* Define internal conditions for linear spline construction (none in this example) */
    ic_type = DF_NO_IC; 
    ic = NULL;

    /* Define boundary conditions for linear spline construction (none in this example) */
    bc_type = DF_NO_BC; 
    bc = NULL;
    scoeffhint = DF_NO_HINT;    /* No additional information about the spline. */ 

    /* Set spline parameters  in the Data Fitting task */
    status = dfdEditPPSpline1D( task, s_order, s_type, bc_type, bc, ic_type,
                                ic, scoeff, scoeffhint );
    
    /* Check the Data Fitting operation status */
    ...

    /* Use a standard computation method to construct a linear spline: */
    /* Pi(x) = ci,0+ci,1(x-xi), i=0,..., N-2 */
    /* The library packs spline coefficients to array scoeff. */ 
    /* scoeff[2*i+0]=ci,0 and scoeff[2*i+1]=ci,1, i=0,..., N-2 */
    status = dfdConstruct1D( task, DF_PP_SPLINE, DF_METHOD_STD );
    
    /* Check the Data Fitting operation status */
    ...

    /* Process spline coefficients */
    ...

    /* Deallocate Data Fitting task resources */
    status = dfDeleteTask( &task ) ;

    /* Check the Data Fitting operation status */
    ...
    return 0 ;
}

The following example demonstrates cubic spline-based interpolation using Data Fitting routines. In this example, a scalar function defined on non-uniform partition is approximated by Bessel cubic spline using not-a-knot boundary conditions. Once the spline is constructed, you can use the spline to compute spline values at the given sites. Computation results are packed by the Data Fitting routine in row-major format.

Example of Cubic Spline-Based Interpolation

#include "mkl.h"

#define NX 100                     /* Size of partition, number of breakpoints */
#define NSITE 1000                 /* Number of interpolation sites */
#define SPLINE_ORDER DF_PP_CUBIC   /* A cubic spline to construct */

int main()
{
    int status;          /* Status of a Data Fitting operation */
    DFTaskPtr task;      /* Data Fitting operations are task based */

    /* Parameters describing the partition */
    MKL_INT nx;          /* The size of partition x */
    double x[NX];         /* Partition x */
    MKL_INT xhint;       /* Additional information about the structure of breakpoints */

    /* Parameters describing the function */
    MKL_INT ny;          /* Function dimension */
    double y[NX];         /* Function values at the breakpoints */
    MKL_INT yhint;       /* Additional information about the function */
    
    /* Parameters describing the spline */
    MKL_INT  s_order;    /* Spline order */
    MKL_INT  s_type;     /* Spline type */
    MKL_INT  ic_type;    /* Type of internal conditions */
    double* ic;         /* Array of internal conditions */
    MKL_INT  bc_type;    /* Type of boundary conditions */
    double* bc;         /* Array of boundary conditions */

    double scoeff[(NX-1)* SPLINE_ORDER];   /* Array of spline coefficients */
    MKL_INT scoeffhint;            /* Additional information about the coefficients */

    /* Parameters describing interpolation computations */
    MKL_INT nsite;        /* Number of interpolation sites */
    double site[NSITE];   /* Array of interpolation sites */
    MKL_INT sitehint;     /* Additional information about the structure of 
                             interpolation sites */

    MKL_INT ndorder, dorder;    /* Parameters defining the type of interpolation */

    double* datahint;   /* Additional information on partition and interpolation sites */

    double r[NSITE];    /* Array of interpolation results */
    MKL_INT rhint;     /* Additional information on the structure of the results */
    MKL_INT* cell;      /* Array of cell indices */

    /* Initialize the partition */
    nx = NX;

    /* Set values of partition x */
    ...
    xhint = DF_NON_UNIFORM_PARTITION;  /* The partition is non-uniform. */

    /* Initialize the function */
     ny = 1;               /* The function is scalar. */
     
    /* Set function values */
    ...
    yhint = DF_NO_HINT;    /* No additional information about the function is provided. */

    /* Create a Data Fitting task */
    status = dfdNewTask1D( &task, nx, x, xhint, ny, y, yhint );

    /* Check the Data Fitting operation status */
    ...

    /* Initialize spline parameters */
    s_order = DF_PP_CUBIC;     /* Spline is of the fourth order (cubic spline). */ 
    s_type = DF_PP_BESSEL;     /* Spline is of the Bessel cubic type. */ 

    /* Define internal conditions for cubic spline construction (none in this example) */
    ic_type = DF_NO_IC; 
    ic = NULL;

    /* Use not-a-knot boundary conditions. In this case, the is first and the last 
     interior breakpoints are inactive, no additional values are provided. */
    bc_type = DF_BC_NOT_A_KNOT; 
    bc = NULL;
    scoeffhint = DF_NO_HINT;    /* No additional information about the spline. */ 

    /* Set spline parameters  in the Data Fitting task */
    status = dfdEditPPSpline1D( task, s_order, s_type, bc_type, bc, ic_type,
                                ic, scoeff, scoeffhint );
    
    /* Check the Data Fitting operation status */
    ...

    /* Use a standard method to construct a cubic Bessel spline: */
    /* Pi(x) = ci,0 + ci,1(x - xi) + ci,2(x - xi)2 + ci,3(x - xi)3, */
    /* The library packs spline coefficients to array scoeff: */
    /* scoeff[4*i+0] = ci,0, scoef[4*i+1] = ci,1,         */   
    /* scoeff[4*i+2] = ci,2, scoef[4*i+1] = ci,3,         */   
    /* i=0,...,N-2  */
    status = dfdConstruct1D( task, DF_PP_SPLINE, DF_METHOD_STD );

    /* Check the Data Fitting operation status */
    ...

    /* Initialize interpolation parameters */
    nsite = NSITE;

    /* Set site values */
    ...

    sitehint = DF_NON_UNIFORM_PARTITION; /* Partition of sites is non-uniform */

    /* Request to compute spline values */
    ndorder = 1;
    dorder = 1;
    datahint = DF_NO_APRIORI_INFO;  /* No additional information about breakpoints or
                                       sites is provided. */
    rhint = DF_MATRIX_STORAGE_ROWS; /* The library packs interpolation results 
                                       in row-major format. */
    cell = NULL;                    /* Cell indices are not required. */

    /* Solve interpolation problem using the default method: compute the spline values
       at the points site(i), i=0,..., nsite-1 and place the results to array r */ 
    status = dfdInterpolate1D( task, DF_INTERP, DF_METHOD_PP, nsite, site,
    sitehint, ndorder, &dorder, datahint, r, rhint, cell );
 
    /* Check Data Fitting operation status */ 
     ...

    /* De-allocate Data Fitting task resources */
     status = dfDeleteTask( &task );
    /* Check Data Fitting operation status */ 
     ...
    return 0;
}

The following example demonstrates how to compute indices of cells containing given sites. This example uses uniform partition presented with two boundary points. The sites are in the ascending order.

Example of Cell Search

#include "mkl.h"

#define NX 100                     /* Size of partition, number of breakpoints */
#define NSITE 1000                 /* Number of interpolation sites */

int main()
{    
    int status;          /* Status of a Data Fitting operation */
    DFTaskPtr task;      /* Data Fitting operations are task based */

    /* Parameters describing the partition */
    MKL_INT nx;          /* The size of partition x */
    float x[2];         /* Partition x is uniform and holds endpoints 
                            of interpolation interval [a, b] */
    MKL_INT xhint;       /* Additional information about the structure of breakpoints */

    /* Parameters describing the function */
    MKL_INT ny;          /* Function dimension */
    float   *y;          /* Function values at the breakpoints */
    MKL_INT yhint;       /* Additional information about the function */
    
    /* Parameters describing cell search */
    MKL_INT nsite;       /* Number of interpolation sites */
    float  site[NSITE]; /* Array of interpolation sites  */
    MKL_INT sitehint;    /* Additional information about the structure of sites */

    float* datahint;     /* Additional information on partition and interpolation sites */

    MKL_INT cell[NSITE]; /* Array for cell indices */ 

    /* Initialize a uniform partition */    
    nx = NX;
    /* Set values of partition x: for uniform partition,           */
    /* provide end-points of the interpolation interval [-1.0,1.0] */
    x[0] = -1.0f; x[1] = 1.0f;
    xhint = DF_UNIFORM_PARTITION; /* Partition is uniform */

    /* Initialize function parameters */ 
    /* In cell search, function values are not necessary and are set to zero/NULL values */
    ny    = 0;
    y     = NULL; 
    yhint = DF_NO_HINT;

    /* Create a Data Fitting task */
    status = dfsNewTask1D( &task, nx, x, xhint, ny, y, yhint );

    /* Check Data Fitting operation status */ 
    ...

    /* Initialize interpolation (cell search) parameters */
    nsite = NSITE;

   /* Set sites in the ascending order */
    ...
    sitehint = DF_SORTED_DATA;      /* Sites are provided in the ascending order. */
    datahint = DF_NO_APRIORI_INFO;  /* No additional information 
                                       about breakpoints/sites is provided.*/

    /* Use a standard method to compute indices of the cells that contain 
       interpolation sites. The library places the index of the cell containing     
       site(i) to the cell(i), i=0,...,nsite-1 */
    status = dfsSearchCells1D( task, DF_METHOD_STD, nsite, site, sitehint,
                            datahint, cell );
    /* Check Data Fitting operation status */ 
     ...

    /* Process cell indices */   
     ...  

    /* Deallocate Data Fitting task resources */
    status = dfDeleteTask( &task );

    /* Check Data Fitting operation status */ 
     ...
    return 0;
}