Are the trigonometric transforms thread-safe?

[Rodney_H_]

I'm using the trig transforms (MKL_STAGGERED_COSINE_TRANSFORM) within an OpenMP parallel region on Windows.

Each OpenMP thread has it's own array of data to be transformed.

When I use one OpenMP thread my program works perfectly. When I use more than one, results are unpredictable.

If I enclose the trig transform code within a #pragma omp critical{} it's OK.

It seems to me that the trigonometric transforms are not thread-safe. Results are similar with sequential and parallel MKL libraries.

My code is a bit too complex to post it all here.

Thanks

Rodney

[Evgueni_P_Intel]

Dear Rodney H.,

Trigonometric transforms are thread-safe. Please let us know which version of MKL you are using and the sequence of MKL calls related to trigonometric transforms.

Thank you.

Evgueni.

[Rodney_H_]

Must be a programming error on my part then.

My sequence of MKL calls is made by calling the following function.

void dct1(double * data, MKL_INT num_cols, int direction) {

DFTI_DESCRIPTOR_HANDLE handle;
MKL_INT stat = 0;
MKL_INT ipar[128];

double * par = new double[(3*num_cols/2) + 2]();

// use function pointers to handle whether float is a float or double

// and direction of the transform

void (*d_trig_transform)(double *, DFTI_DESCRIPTOR_HANDLE *, MKL_INT *, double *, MKL_INT *);

if(direction < 0)
d_trig_transform = d_forward_trig_transform;
else
d_trig_transform = d_backward_trig_transform;

MKL_INT tt_type = MKL_STAGGERED_COSINE_TRANSFORM;

d_init_trig_transform(&num_cols, &tt_type, ipar, par, &stat);

if(stat != 0){
cerr << "Error initializing row transform: " << stat << endl;
return;
}

d_commit_trig_transform(data, &handle, ipar, par, &stat);
if(stat != 0){
cerr << "Error committing row transform: " << stat << endl;
return;
}

d_trig_transform(data, &handle, ipar, par, &stat);
if(stat != 0){
cerr << "Error in row transform: " << stat << endl;
return;
}

free_trig_transform(&handle, ipar, &stat);
if(stat != 0){
cerr << "Error freeing row transform: " << stat << endl;
return;

}
　delete[] par;
}

Each OpenMP thread calls it like this, Q is the size of the data

double * temp = new double[Q+1]();

for(int n1 = 0; n1 < Q; ++n1){
temp[n1]=transdata[n1]*2.0/((double)Q);
}

dct1(temp, Q, 1);

I'm using MKL 11.1  Update 1 on Windows 7 64-bit installed with C++ 2013.1.1.139, compiling under 64 bit mode on an Intel Core i7 laptop.

If I use more than one OpenMP thread, then the program behaves unpredictably. Sometimes it runs, often fails around calling the trig transform, as far as I can tell a heap corruption error occurs.

If I surround the call as follows, it runs OK.
   #pragma omp critical (DCT)
{
                dct1(temp, Q, 1);
}

If I substitute a different DCT function, for example if I call the DCT from the code at http://www.kurims.kyoto-u.ac.jp/~ooura/fft.html, and call it in the same manner, then my program runs fine with multiple openMP threads.

Note that I'm trying to get each thread to run it's own single-threaded DCT on it's own data, not run a DCT in parallel on one set of data.

[Evgueni_P_Intel]

Dear Rodney H.,

Please check MKL documentation -- par should be (5*num_cols/2) + 2 elements long (you allocate only (3*num_cols/2) + 2 elements.)

Thank you,

Evgueni.

[Rodney_H_]

Dear Evgueni,

I don't think it makes much difference. From MKL documentation: "If tt_type=MKL_STAGGERED_COSINE_TRANSFORM, the transform uses only the first 3n/2 elements, which contain tabulated sine and cosine values. "

I tried your suggestion and still get heap corruption problems when calling from multiple threads.

Thanks,

Rodney

[Kai1]

Hi, Rodney,

I want to know if you have solved the problem, because I have met exactly the same issue, and suspect that the TT in Intel MKL is not thread-safe. 

Can any maester from Intel help look at the problem? At least, what are the appropriate options should we set to make TT thread-safe when used with OpenMP?

Thanks

[Kai1]

and below is my TT routines:

        call s_init_trig_transform(m, transform_type, ipar, spar, dct_status)
        call s_commit_trig_transform(w, handle, ipar, spar, dct_status)
        call s_backward_trig_transform(w, handle, ipar, spar, dct_status)
        call free_trig_transform(handle, ipar, dct_status)

where m = size(w), w is a 1D real array, spar has length nint(3*m/2.0). trasnform_type = mkl_staggered_cosine_transform. 

Thanks

[Gennady_F_Intel]

that's not enough to reproduce the issue. You may submit the issue to https://supporttickets.intel.com/?lang=en-US  and  provide us the reproducer of the problem.