2014-12-24 41 views
0

儘管我已經閱讀了Mathworks文檔,但是我發現構造一個調用dgesv LAPACK例程的mex文件非常困難。有人能幫我一下嗎?我在matlab中是好的,在C中是初學者,所以我不能解決這個問題。調用dgesv的mex文件LAPACK例程

我有一段可以玩的代碼,但是我對如何繼續進行了一點線索。

#include "mex.h" 
// Include headers for your library 

void mexFunction(int nlhs,mxArray *plhs[],int nrhs,const mxArray *prhs[]) 
{ 
    void* x = mxGetData(prhs[0]); // Assume one input. Check nrhs 
    plhs[0] = mxCreateDoubleMatrix(10,10,mxREAL); // Create 10x10 double matrix for output 
    void* y = mxGetData(plhs[0]); 
    yourLibraryFunction(x, y); // Read from x and write to y. Pass sizes in if needed 
} 

回答

0

我不是專家,但我會盡力解釋你。從mathworks中看到這個例子。

參數prhs有輸入參數(表示它們有n個),plhs有輸出參數(表示有nlhs)。

首先,您必須「轉換」您可以從C中使用它們的類型中的參數.C中的參數指向矩陣,因此我們使用mxGetPr作爲數據。在文檔中說,返回一個指向實際數據的第一個元素的指針,在C中,數組的名稱是第一個元素的指針。

然後用mxGetM和mxGetN得到矩陣A和矢量b的維數,在bnrhs中保存矩陣b的列數,因爲您需要將它作爲dgesv的第二個參數傳遞。但是,因爲你在C中工作,而且你有指向真實數據的指針(通常在函數中我們不想使用指向真實數據的指針),所以你必須複製矩陣,如果你想複製一個數組,一元一元地做,因爲如果你做AT = A他們指向相同的位置。在mex文件中,我們使用mxCalloc而不是malloc。

然後,你必須調用dgesv函數,最後你必須通過使指針plhs [0]指向dgesv作爲解決方案返回的向量來返回結果。

/*========================================================= 
* dgesv.c 
* example for illustrating how to use LAPACK within a C 
* MEX-file on HP700, IBM_RS and PC. 
* 
* Solves the symmetric indefinite system of linear 
* equations A*X=B for X. 
* X = dgesv(A,B) computes LU decomposition of A and returns 
* the result X such that A*X is B. 
* B must have as many rows as A. A is n-by-n matrix and b is n-by-bnrsh matrix 
* http://www.netlib.org/lapack/double/dgesv.f 
* This is a MEX-file for MATLAB. 
* 
* Compiling command in Matlab 
* mex -v -g dgesv.c <matlab>\extern\examples\refbook\fort.c -I<matlab>\extern\examples\refbook <matlab>\extern\lib\win32\microsoft\msvc60\libmwlapack.lib 
*=======================================================*/ 
#include "mex.h" 
#include "c:\\matlab6p5\\extern\\examples\\refbook\\fort.h" 

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

    /* mex interface to LAPACK functions dsysvx and zhesvx */ 

    int m, n, bnrhs, lda, *ipiv, ldb, info; 
    double *A,*b, *AT, *bT, *ret; 
    int i; 

    if (nrhs != 2) { 
    mexErrMsgTxt("Expect 2 input arguments and return 1 output argument"); 
    } 

    A = mxGetPr(prhs[0]); 
    b = mxGetPr(prhs[1]); 
    m = mxGetM(prhs[0]); 
    n = mxGetN(prhs[0]); 
    bnrhs = mxGetN(prhs[1]); 

    AT = (double *)mxCalloc(m*n,sizeof(double)); 
    bT = (double *)mxCalloc(m*bnrhs,sizeof(double)); 
    ipiv = (int *)mxCalloc(m,sizeof(int)); 

    for(i=0;i<m*bnrhs;i++) { 
     bT[i] = b[i]; 
    } 
    mexPrintf("\n"); 
    for(i=0;i<m*n;i++) AT[i] = A[i]; 
    lda = m; 
    ldb = m; 
    dgesv(&m, &bnrhs, AT, &lda, ipiv, bT, &ldb, &info);  

    for(i=0;i<m;i++) { 
     mexPrintf("%f\t",bT[i]); 
    } 

    plhs[0] = mxCreateDoubleMatrix(m,bnrhs,mxREAL); 
    ret = mxGetPr(plhs[0]); 
    for(i=0;i<m*bnrhs;i++) ret[i] = bT[i]; 

    mxFree(AT); 
    mxFree(bT); 
    mxFree(ipiv); 
} 

源:http://www.mathworks.com/matlabcentral/fx_files/8006/1/dgesv.chttp://www.mathworks.com/help/matlab/mex-library.html