This package contains the command "lapack", which provides a simple interface to call any LAPACK or BLAS routine from inside Matlab, as well as the command "lapackhelp", which brings up detailed information about any LAPACK or BLAS routine.
A summary of the features include:
· Uses Matlab's internal LAPACK and BLAS libraries.
· Automatic argument type conversion (ie. converts arguments from Matlab types to FORTRAN types).
· Easy access to detailed descriptions about the behavior and arguments for every LAPACK and BLAS function
· Allows use of LAPACK and BLAS without any knowledge of MEX, C, or FORTRAN.
Qianyong Chen, I think the segfault was because lwork was too small; doing as follows doesn't seem to crash matlab. The fortranError as well as fortranError_SVD remain both numerically zero.
A = rand(118,106);
% 1st: results from direct matlab
[UM, SM, VM] = svd(A);
matlabError = max(max(abs(A-UM*SM*VM')))
% 2nd: results when using 'lapack' by Tim Toolan
[m,n] = size(A);
Dude, your code rocks! You're awesome! :) I have been able to use the dsyev routine on a very large matrix (size 200k^2). It actually worked just as well, if not better, than Matlab 2011b's eig.
I really like this idea. But couldn't make it work on my machine: 64bit Win7 with SDK 7.1 with Matlab 2011b.
The first run is always fine. But when I run the same code a few more times, it either stops due to segmentation error,
or gives unreliable results. Also, the difference for U and V (between Matlab results and using lapack.c provided here)
is large when the matrix is as small as 118 by 106. The error for the matrix itself (A - U * S * V) is very small though.
What I did is: 1) compile with mex lapack.c
2) run the following code:
% test a couple subroutines
clear all;
%--------------------
%1st: dgesvd.f
%--------------------
A = rand(118,106);
% 1st: results from direct matlab
[UM, SM, VM] = svd(A);
matlabError = max(max(abs(A-UM*SM*VM')))
% 2nd: results when using 'lapack' by Tim Toolan
[m,n] = size(A);
Tim, I'm using a slightly modified version of your lapack package in my trilin package (#31125). I'd be glad if you'd consider incorporating the modifications in my submission into your original submission. Thanks.
I was suuuuper excited about this package. I am considering ditching all my LAPACK-calling C mex code in favor of it. The compilation is a bliss, no hassle in finding out where either or both LAPACK and BLAS are installed, depending on MATLAB version, computer, and compiler!
Before that, though, I was wondering if could you please comment on the following issues.
- Compiling without the -largeArrayDims flag in a 64-bit MATLAB seems to give unreliable results. Any idea why? It should be optional.
- The instructions on top of lapack.c say:
To build on Windows, run: mex lapack.c
Since lapack.m already beautifully takes care of that using largeArrayDims, you might want to change lapack.c to instruct users to simply run lapack(.m).
- Using lapack.mexw64 compiled on MATLAB 2009 64-bit gives unreliable results when used on MATLAB 2007 64-bit; it seems related to casting integer types, e.g., ipiv in a call to dgetrf has some very large values:
17179869189
21474836483
0
0
- I see that all arguments are copied to new memory. For me that's a concern, because I need to handle very large matrices. I understand doing that for arguments that will be overwritten by LAPACK, but that's clearly unnecessary for, e.g., the matrix factor in solving a linear system. I was wondering if you'd agree with the inclusion of a new argument, something like this:
else if (c == 'U'){
/* USER DEFINED **********************************************
/ User guarantees type and also whether or not LAPACK will overwrite this argument.
/ lowercase: input-only: need not copy.
/ UPPERCASE: input/output: need to copy.
/ TODO: complex-valued variables need at least to have their real/imag parts reorganized. */
if (mxIsComplex(args[i])){
sprintf(msg, "Can't convert argument %d to FORTRAN type "
"DOUBLE PRECISION, use double.", i+1);
mexErrMsgTxt(msg);
}
if (isupper(argType[i])){
value = mxDuplicateArray(args[i]);
a[i] = mxGetPr(value);
if (plhs[0] && mxIsCell(plhs[0]))
mxSetCell(plhs[0], indArgOut[i], value);
else
plhs[indArgOut[i]] = value;
}
else{
a[i] = mxGetPr(args[i]);
}
}
- e.g., lapackhelp('dtrcon') leads to an unreachable website; the most stable URL seems to be:
http://www.netlib.org/lapack/double/dtrcon.f
(where the double/single subdir would have to be chosen automatically)
I'm trying to edit my previous review. I am now able to get this to run! To compile on 2006b under Gentoo linux:
mex -D_GNU_SOURCE -ldl -lmwlapack lapack.c
I have heard that -lmwlapack is not needed for >2006 versions of matlab, tho.
I look forward to using this! (the initial timings look good, too.)
03 Jan 2008
worker b
Unable to get this to run on my machine, but this is probably due to a general matlab/lapack issue. However, I had to use the -D_GNU_SOURCE flag in compiling (RTLD_DEFAULT is a GNU extension evidently. whatever.)
Updates
22 Oct 2009
Improved ease of use by adding internal prototypes and the ability to return the results in a cell array.
Added the command lapackhelp, which will bring up detailed information about any BLAS/LAPACK routine and its arguments from inside Matlab.