我想学习如何在MATLAB中调用内置的LAPACK/BLAS例程。我在MATLAB和mex文件方面有经验，但实际上我不知道如何调用LAPACK或BLAS库。我在文件交换中发现了网关例程，它简化了调用，因为我不必为任何函数（如本函数）编写mex文件。我需要任何玩具示例来学习MATLAB和这些内置库之间的基本消息传递。任何玩具示例，如矩阵乘法或LU分解都是受欢迎的。
如果您查看上述FEX提交的lapack.m文件内部，您将看到几个关于如何使用该函数的示例：
示例：使用DGESVD进行SVD分解：
X = rand(4,3);
[m,n] = size(X);
C = lapack('dgesvd', ...
'A', 'A', ... % compute ALL left/right singular vectors
m, n, X, m, ... % input MxN matrix
zeros(n,1), ... % output S array
zeros(m), m, ... % output U matrix
zeros(n), n, .... % output VT matrix
zeros(5*m,1), 5*m, ... % workspace array
0 ... % return value
);
[s,U,VT] = C{[7,8,10]}; % extract outputs
V = VT';
（注意：我们使用这些伪变量作为输出变量的原因是Fortran函数希望所有参数都通过引用传递，但MATLAB中的MEX函数不允许修改它们的输入，因此它被编写为返回单元格数组中所有输入的副本，并进行任何修改）
我们得到：
U =
-0.44459 -0.6264 -0.54243 0.3402
-0.61505 0.035348 0.69537 0.37004
-0.41561 -0.26532 0.10543 -0.86357
-0.50132 0.73211 -0.45948 -0.039753
s =
2.1354
0.88509
0.27922
V =
-0.58777 0.20822 -0.78178
-0.6026 -0.75743 0.25133
-0.53981 0.61882 0.57067
这相当于MATLAB自己的SVD函数：
[U,S,V] = svd(X);
s = diag(S);
给出：
U =
-0.44459 -0.6264 -0.54243 0.3402
-0.61505 0.035348 0.69537 0.37004
-0.41561 -0.26532 0.10543 -0.86357
-0.50132 0.73211 -0.45948 -0.039753
s =
2.1354
0.88509
0.27922
V =
-0.58777 0.20822 -0.78178
-0.6026 -0.75743 0.25133
-0.53981 0.61882 0.57067
编辑：
为了完整起见，我在下面展示了一个MEX函数的例子，它直接调用DGESVD例程的Fortran接口。
好消息是，MATLAB提供了libmwlapack和libmwblas库以及我们可以使用的两个相应的头文件blas.h和lapack.h。事实上，文档中有一页解释了从MEX文件调用BLAS/LAPACK函数的过程。
在我们的案例中，lapack.h定义了以下原型：
extern void dgesvd(char *jobu, char *jobvt,
ptrdiff_t *m, ptrdiff_t *n, double *a, ptrdiff_t *lda,
double *s, double *u, ptrdiff_t *ldu, double *vt, ptrdiff_t *ldvt,
double *work, ptrdiff_t *lwork, ptrdiff_t *info);
svd_lapack.c
#include "mex.h"
#include "lapack.h"
void mexFunction(int nlhs, mxArray *plhs[], int nrhs, const mxArray *prhs[])
{
mwSignedIndex m, n, lwork, info=0;
double *A, *U, *S, *VT, *work;
double workopt = 0;
mxArray *in;
/* verify input/output arguments */
if (nrhs != 1) {
mexErrMsgTxt("One input argument required.");
}
if (nlhs > 3) {
mexErrMsgTxt("Too many output arguments.");
}
if (!mxIsDouble(prhs[0]) || mxIsComplex(prhs[0])) {
mexErrMsgTxt("Input matrix must be real double matrix.");
}
/* duplicate input matrix (since its contents will be overwritten) */
in = mxDuplicateArray(prhs[0]);
/* dimensions of input matrix */
m = mxGetM(in);
n = mxGetN(in);
/* create output matrices */
plhs[0] = mxCreateDoubleMatrix(m, m, mxREAL);
plhs[1] = mxCreateDoubleMatrix((m<n)?m:n, 1, mxREAL);
plhs[2] = mxCreateDoubleMatrix(n, n, mxREAL);
/* get pointers to data */
A = mxGetPr(in);
U = mxGetPr(plhs[0]);
S = mxGetPr(plhs[1]);
VT = mxGetPr(plhs[2]);
/* query and allocate the optimal workspace size */
lwork = -1;
dgesvd("A", "A", &m, &n, A, &m, S, U, &m, VT, &n, &workopt, &lwork, &info);
lwork = (mwSignedIndex) workopt;
work = (double *) mxMalloc(lwork * sizeof(double));
/* perform SVD decomposition */
dgesvd("A", "A", &m, &n, A, &m, S, U, &m, VT, &n, work, &lwork, &info);
/* cleanup */
mxFree(work);
mxDestroyArray(in);
/* check if call was successful */
if (info < 0) {
mexErrMsgTxt("Illegal values in arguments.");
} else if (info > 0) {
mexErrMsgTxt("Failed to converge.");
}
}
在我的64位Windows上，我将MEX文件编译为：mex -largeArrayDims svd_lapack.c "C:Program FilesMATLABR2013aexternlibwin64microsoftlibmwlapack.lib"
这里有一个测试：
>> X = rand(4,3);
>> [U,S,VT] = svd_lapack(X)
U =
-0.5964 0.4049 0.6870 -0.0916
-0.3635 0.3157 -0.3975 0.7811
-0.3514 0.3645 -0.6022 -0.6173
-0.6234 -0.7769 -0.0861 -0.0199
S =
1.0337
0.5136
0.0811
VT =
-0.6065 -0.5151 -0.6057
0.0192 0.7521 -0.6588
-0.7949 0.4112 0.4462
与。
>> [U,S,V] = svd(X);
>> U, diag(S), V'
U =
-0.5964 0.4049 0.6870 0.0916
-0.3635 0.3157 -0.3975 -0.7811
-0.3514 0.3645 -0.6022 0.6173
-0.6234 -0.7769 -0.0861 0.0199
ans =
1.0337
0.5136
0.0811
ans =
-0.6065 -0.5151 -0.6057
0.0192 0.7521 -0.6588
-0.7949 0.4112 0.4462
（请记住，U和V中特征向量的符号是任意的，因此比较两者可能会得到翻转符号）
这里解释了如何在matlab中使用SVD的一个好的/实用的例子：将捕获的坐标转换为屏幕坐标
关于如何使用lapack在objective-c中计算svd的更多信息，请参阅此处。使用xcode 中的lapack从svd计算objective-c中的V=USVt