BLAS in mex file -- what all needs to be updated?

Hi, I wrote some software a while ago, and it uses BLAS in a mex file. Everything compiles well in 2008b, but in the newer versions the files with BLAS crash when running (newer version could only be tested on 64-bit unix). I already replaced all integers by mxSignedIndex - it did not help. Any suggestions what else changed and needs to be updated? Below is an excerpt from the code. Thanks!

#ifdef _WIN32
  double ddot(mwSignedIndex*, double*, mwSignedIndex*, double*, mwSignedIndex*);
#else
  double ddot_(mwSignedIndex*, double*, mwSignedIndex*, double*, mwSignedIndex*);
#endif

void mexFunction(int nlhs, mxArray *plhs[], 
     int nrhs, const mxArray *prhs[])
{

double *k, *dK, *t, *w, *tdiff;
double *inds;
double sigma, dprod, oneD=-1.0;
mwSignedIndex n, m, nw, i, j, colind, oneI = 1, nm;

k = mxGetPr(prhs[0]);
inds = (double*) mxGetPr(prhs[1]);
w = mxGetPr(prhs[2]);
t = mxGetPr(prhs[3]);
sigma = *mxGetPr(prhs[4]);
sigma *= -sigma;
sigma = 1/sigma;
n = mxGetM(prhs[0]);  /* number of samples */
m = mxGetN(prhs[0]);  /* number of indices (called d elsewhere) */
nw = mxGetN(prhs[2]); /* length of w: number of sources */
nm = n*m;

plhs[0] = mxCreateDoubleMatrix(nm, nw, mxREAL);
dK = mxGetPr(plhs[0]);

tdiff = mxCalloc(nw, sizeof(double));

for (i=0; i<m; ++i){ /* column index */
  colind = (mwSignedIndex)(inds[i]) - 1;
  for (j=0; j<n; ++j){

      /* tdiff = tj - ti */
      dcopy_(&nw, &t[j*nw], &oneI, tdiff, &oneI);
      daxpy_(&nw, &oneD, &t[colind*nw], &oneI, tdiff, &oneI);
      /* dprod = w * (ti - tj) */
      dprod = ddot_(&nw, w, &oneI, tdiff, &oneI);
      /* dK( (i,j),:) = -K(i,j)/sigma * w * (tj - ti)(tj - ti)' */
      dcopy_(&nw, tdiff, &oneI, &dK[i*n + j], &nm);
      dprod *= k[i*n + j] * sigma;
      dscal_(&nw, &dprod, &dK[i*n + j], &nm);
    }
  }
mxFree(tdiff);
return;
}

[Jan]

Either use a debugger or set some mexPrintf's in the code to find out, which line causes the crash. It will be nmuch more efficient, if you locate the problem exactly.

[Kaustubha Govind]

Here are some tips for debugging: http://www.mathworks.com/help/techdoc/matlab_external/f32489.html

[James Tursa]

Can you provide a small example input & output set?

[someone]

It crashes in dcopy. I am pretty sure it is something with the BLAS, since the exact same code compiles and works with Matlab 2008b.

Example input/output:

out=func(A, b, c, D, e)

A = [ 1 4;

2 5;

3 6];

b = [1, 2];

c = [1, 2];

D = [ 1 3 5;

2 4 6];

e = 1;

out will be

out =

0 0

-24 -24

-144 -144

-48 -48

0 0

-72 -72

[Jan]

There are 2 dcopy calls. Which causes the crash?

[someone]

sorry, the first one ends in a segmentation fault.

[James Tursa]

Which BLAS library are you linking with? The one with MATLAB or another one? Also, are you using prototypes for all the BLAS functions you use? I don't see all of them in your listing above.

Answer 1

[Jan]

A strange idea, but try ptrdiff_t instead of mxSignedIndex.

[someone]

thanks a lot, weird, but that seems to help :-)