Is there something wrong in DataFitting interpolate with user defined extrapolate call back routine?

[Tianxiong_Lu]

Hi all,

I meet a strange problem in DataFitting Module while call df?InterpolateEx1D function. If I present the call back function for extrapolation , the result was wrong. But if set them to NULL, result is right.

My problem is:

x[] = {1, 1.5, 2.5, 3.5, 4.5, 5}, and y[] = {100, 120, 150, 170, 200, 200}, Interpolate method is stepwise constant interpolant. And I need found y value at x=4.

The following is the sample code. I hope that was something wrong in my code. Can you help me to explain it? thanks.

#include <stdio.h>
#include "mkl_df.h"


// fixed value extrapolate function parameters
struct fixvalue_extrap_params
{
    float value;
};

// fixed value extrapolate function
int fixvalue_extrap( MKL_INT64* n, MKL_INT64 cell[], float site[], float r[],
                     void *params )
{
    fixvalue_extrap_params *p = (fixvalue_extrap_params*)params;
    r[0] = p->value;
    return 0;
}


int main(int argc, char* argv[])
{
    float x[] = { 1., 1.5, 2.5, 3.5, 4.5, 5 };
    float y[] = { 100., 120., 150., 170., 200., 200. };
    int nx = 6;
    int ny = 1;
    const int nsite = 1;
    float site[nsite] = { 4 };
    float r[nsite];
    int cell[nsite];

    DFTaskPtr task;                     // Data Fitting task descriptor

    MKL_INT ndorder = 1;                    // size of array describing derivative
    MKL_INT dorder[] = { 1 };           // only value to calculate

    int errcode = 0;

    /***** Create Data Fitting task *****/
    errcode = dfsNewTask1D(&task, nx, x, DF_NON_UNIFORM_PARTITION, ny, y, DF_NO_HINT);

    /***** Edit task parameters for look up interpolant *****/
    errcode = dfsEditPPSpline1D(task, DF_PP_STD, DF_CL_STEPWISE_CONST_INTERPOLANT, 0, 0, 0, 0, 0, 0);

    /***** Interpolate using lookup method without extrapolate *****/
    errcode = dfsInterpolateEx1D(task, (DF_CELL | DF_INTERP), DF_METHOD_PP,
                                 nsite, site, 0, ndorder,
                                 dorder, 0, r, 0, cell,
                                 0, 0, 0, 0,
                                 0, 0, 0, 0);
    // this answer is right
    printf("interpolate reuslt for [%f] = [%f] in section %d\n", site[0], r[0], cell[0]);

    dfsInterpCallBack le_cb, re_cb;       // interpolation call backs
    void * le_params, *re_params;    // interpolation call backs parameters
    fixvalue_extrap_params le_fixparams, re_fixparams;

    le_cb = fixvalue_extrap;
    re_cb = fixvalue_extrap;
    le_fixparams.value = y[0];
    re_fixparams.value = y[nx-1];
    le_params = &le_fixparams;
    re_params = &re_fixparams;

    /***** Interpolate using lookup method *****/
    errcode = dfsInterpolateEx1D(task, (DF_CELL | DF_INTERP), DF_METHOD_PP,
                                 nsite, site, 0, ndorder,
                                 dorder, 0, r, 0, cell,
                                 le_cb, le_params, re_cb, re_params,
                                 0, 0, 0, 0);

    // this answer is wrong
    printf("interpolate reuslt with extrapolate for [%f] = [%f] in section %d\n", site[0], r[0], cell[0]);

    /***** Delete Data Fitting task *****/
    errcode = dfDeleteTask(&task);

    return 0;
}

 My output is :

interpolate reuslt for [4.000000] = [170.000000] in section 4
interpolate reuslt with extrapolate for [4.000000] = [150.000000] in section 4

[Tianxiong_Lu]

By the way, My MKL version is 11.3.

[Zhang_Z_Intel]

Thanks for providing a test case!

Yes, this looks weird. Even if you make the user-defined extrapolation function a "no-op" (delete everything and only return 0), you still don't get the same result as you would get without passing an extrapolation function at all. MK dev team is taking a look and will let us know soon.

[Tianxiong_Lu]

Thank you for your information. So I need to try another way to reach my objective and hope to hear good news from dev team.

[Andrey_N_Intel]

We reproduced this behavior, and it looks like a bug. We plan to analyze/fix it in one of the future releases of the library. Will keep you updated. Thanks, Andrey

[VictoriyaS_F_Intel]

Hello Tianxiong and Guix,

The issue is fixed in Intel MKL 2017 Beta Update 1. Can you please check if the fix works for you with this version of the library?