How do you properly call the ARPACK library wrapped by MATLAB from MEX?

15 views (last 30 days)
I have a code base written in C++ that I'm trying to turn into a MEX file under the 64-bit Linux version of Matlab R2013a. The code base uses the ARPACK (ARnoldi PACKage), LAPACK, and BLAS libraries. All three of these libraries have been modified and wrapped by Matlab (libmwarpack, libmwlapack, libmwblas). It is my understanding that these libraries are the same as their native counterparts except that all integers passed to the underlying Fortran functions are 64-bit (signed) integers (for 64-bit installations) instead of 32-bit in the native libraries.
I attempted to MEX the code and link to the native libraries (libarpack, liblapack, libblas). The MEX file compiles, but at runtime, Matlab crashes. Turns out that since the symbols in these libraries were already loaded from Matlab's versions, the program does not dynamically link to the native libraries, but instead to the Matlab ones. Since my MEX passes 32-bit integers, it crashes.
I have changed my code to use these 64-bit integers (mwSignedIndex, i.e. ptrdiff_t) when calling LAPACK and BLAS explicitly. However, my code uses the object-oriented implementation of ARPACK, ARPACK++, which Matlab has not wrapped (to my knowledge). Since the native ARPACK++ uses 32-bit integers for its calls to the ARPACK library (which internally calls LAPACK and BLAS), I took the liberty of downloading the source code from Debian and modifying all integers to be ptrdiff_t (I tried int64_t, too, FWIW) and rebuilding the library in order for it to be compatible with libmwarpack.
MEXing against this new ARPACK++ library and the 3 Matlab-wrapped libraries, the code compiles and links. Running the MEX, everything runs as it should, but at the very end when the destructor for the ARPACK++ object is called, the program crashes. Following the trail with the debugger, it fails when it tries to free dynamic memory at an illegal address (0xc). The class defines a pointer which, when needed, is used to dynamically allocate memory for a work array to pass to the ARPACK Fortran subroutines. The pointer that was passed to free() was supposed to be NULL since it was never allocated (and therefore nothing would happen upon freeing the null pointer), but the pointer was modified to be 0xc during one of the calls to one of the Fortran subroutines in libmwarpack. The GDB debugger shows:
Old value = 0
New value = 1
0x00007fffcf3cbe3b in dsaup2_ () from /usr/local/bin/Matlab/R2013a/bin/glnxa64/libmwarpack.so
The watchpoint put on the pointer's value gets hit 25 times before changing in value, which happens to be the number of eigenvalues I am requesting. So, I know the watchpoint is being hit with every iteration of some loop, which it shouldn't be being accessed at all because the pointer in question is NOT being used by this function at all. My guess is that since some of these 64-bit integers were declared near this pointer in code, and therefore they may be neighboring in memory. And for some reason, they have been cast as the wrong size when passed as pointers to the Fortran routines. However, I can't figure out what is going on in the library that would cause this because I don't know the details of the implementation by Matlab.
Does anyone have any knowledge of Matlab's implementation of the ARPACK library (libmwarpack) which could help me diagnose this problem? Or could someone direct me to a resource which details it? Did Matlab change anything in the underlying Fortran code (seems to me that you would HAVE to in order to support 64-bit integers and pointers).
I would essentially have to upload my entire project here to provide code samples, so I haven't done that. And for the record, if I set my environment variables BLAS_VERSION and LAPACK_VERSION to the paths to the native libraries, Matlab loads these libraries in place of its own versions, which causes dynamic linking of these (correct) libraries and the MEX works perfectly without ANY of the 64-bit modifications I've made. The only reason I don't run with this particular solution is because I don't know how using these libraries will affect the performance of Matlab outside of my MEX routine, and there's nothing that I can see where I can switch the loaded libraries at run time.
Any help is appreciated.
EDIT: Found the fix, but not the error. IPARAM and IPTR should be size 11 and 15 integer arrays, respectively. However, I suspect that since Matlab is a 1-based index language and C++ is a 0-based index, the Matlab-wrapped classes compensate for this and probably index into the arrays less one, or the pointer is passed less one index position. If C++ calls those same functions, then the pointer position is incorrect (or the indexing) because the compensation is not needed, so the Fortran subroutine was writing to one of these arrays in the wrong part of memory and inadvertently overwriting RWORK. I found extending each array by one in size alleviates the issue. I would need to know the implementation of these functions in the class to make a less kludgy fix.
  1 Comment
Stephen
Stephen on 9 Apr 2014
I have a similar problem. With mwblas and mwlapack, there are two nice things: first, blas.h and lapack.h are provided, so I know how to call them (i.e., that they do indeed want ptrdiff_t instead of int), and secondly, they also have 32-bit compatible libraries like libmwblascompat32 (and as you mention, I can also change Matlab's blas and lapack versions to other libraries). With arpack, I have none of these options... too bad!

Sign in to comment.

Answers (0)

Categories

Find more on Write C Functions Callable from MATLAB (MEX Files) in Help Center and File Exchange

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!