1. 引言
1.1 起因
這周去見應(yīng)用數(shù)學(xué)的畢設(shè)老師Z当悔，閑聊時(shí)Z說另一個(gè)同學(xué)代碼里某個(gè)for循環(huán)要跑幾小時(shí)，但這個(gè)for循環(huán)里操作并不復(fù)雜，主要就是轉(zhuǎn)換和賦值。然后我就想悴务，可不可以使用C/C++加速呢？
作為堅(jiān)定的C/C++黨譬猫，正好也想學(xué)學(xué)python中使用C/C++，C/C++中使用asm匯編這些技術(shù)羡疗，索性在github上建了一個(gè)新倉庫Extend-Embed染服，專門學(xué)習(xí)怎樣在一種語言編寫的代碼中使用另一種語言。
這篇文章記錄怎樣通過MEX方式叨恨，在MATLAB中使用C/C++柳刮，最主要的參考資料是：Writing MATLAB C/MEX Code
1.2 環(huán)境
操作系統(tǒng)：Windows 7
MATLAB：R2013a
Visual Studio：Community 2013
2. 準(zhǔn)備
2.1 基本概念
簡(jiǎn)單粗暴地講，在MATLAB中使用C/C++就是三步走的流程：
- 編寫C/C++文件
- 用mex命令把C/C++文件編譯成MEX(Matlab Executable)文件
- 像普通MATLAB函數(shù)一樣直接使用
舉個(gè)例子痒钝，我們寫了一個(gè)hello.c秉颗，功能是輸出一句"Hello, world!"，那么使用：
mex hello.c
后根據(jù)平臺(tái)不同會(huì)得到hello.mex或者h(yuǎn)ello.mexw64送矩，然后使用：
hello
就會(huì)輸出一句”Hello蚕甥， world！“
2.2 配置MEX
不得不說配置過程還是挺坑的栋荸。自己試驗(yàn)加上網(wǎng)上查資料試了無數(shù)遍菇怀，終于發(fā)現(xiàn)配置過程要經(jīng)常反其道而行之，也就是說按正常思路配置基本都要GG晌块。最終選了一堆no反而配好了環(huán)境時(shí)：
輸入mex -setup爱沟，會(huì)詢問是否讓mex自己找編譯器，注意選擇y是不行的：
這時(shí)要選擇n匆背，得到各種支持的編譯器名稱：
我安裝了Visual Studio Community 2013呼伸，所以選擇21。因?yàn)樗o的默認(rèn)路徑不對(duì)钝尸，所以選擇n后輸入正確路徑：
輸入路徑后選擇y確認(rèn)括享，輸出信息里有Done就是配置成功了：
3. 混合編程
在MATLAB安裝目錄下有一個(gè)extern\examples\mex文件夾搂根，里面有很多例子；也可以參考Writing MATLAB C/MEX Code奶浦。不過必須知道的新概念也就幾條兄墅，其它都可以邊用邊學(xué)，所以實(shí)際上手非常容易澳叉。
3.1 新概念
頭文件方面隙咸，使用MEX文件時(shí)都需要include “mex.h”。
main函數(shù)方面成洗，入口不再是main函數(shù)五督，而是mexFunction，聲明格式如下：
void mexFunction(int nlhs, mxArray *plhs[],
int nrhs, const mxArray *prhs[])
其中mxArray是mex.h里自己定義的數(shù)據(jù)結(jié)構(gòu)瓶殃，直接看四個(gè)奇怪變量名的含義：
lhs = left hand side充包，rhs = right hand side
nlhs = 輸出變量的數(shù)量
plhs = 指向輸出變量的指針數(shù)組
nrhs = 輸入變量的數(shù)量
prhs = 指向輸入變量的指針數(shù)組（注意const）
用MEX方式實(shí)現(xiàn)了函數(shù)test，在MATLAB里面使用：
F = test(X, Y)
那么輸入方面遥椿，nrhs為2基矮，prhs[0]是X，prhs[1]是Y冠场；輸出方面家浇，nlhs為1，plhs[0]是F碴裙。
需要注意的是兩件事：第一钢悲，要自己分配輸出的存儲(chǔ)空間；第二舔株，要自己檢查輸入是否合法莺琳。這兩件事直接講不好講清楚，看后面的例子就很容易理解了载慈。
此外惭等，mex.h里提供了不少輔助函數(shù)，輸出信息用mexPrintf娃肿、mexErrMsgTxt咕缎，得到矩陣尺寸用mxGetM和mxGetN，得到指向矩陣數(shù)據(jù)的指針用mxGetPr料扰，etc
3.2 例子
寫了一個(gè)discrete_cpp.cpp凭豪，功能是傳入一個(gè)矩陣A返回一個(gè)矩陣B，B和A的形狀相同晒杈，B中元素B(i, j)和A中元素A(i, j)的對(duì)應(yīng)關(guān)系是：
- 若 A(i, j) 小于1/3嫂伞，則 B(i, j) = 0
- 若 A(i, j) 大于等于1/3 且小于2/3，則 B(i, j) = 1
- 若 A(i, j) 大于2/3，則 B(i, j) = 2
類似于做了離散化帖努。
#include "mex.h"
// 使用MEX必須包含的頭文件
int discrete(double d);
void mexFunction(int nlhs, mxArray *plhs[],
int nrhs, const mxArray *prhs[]){
if (nrhs != 1)
mexErrMsgTxt("Wrong number of input arguments.\n");
// 檢查輸入變量數(shù)量是否正確撰豺，否則報(bào)錯(cuò)
if (nlhs > 1)
mexErrMsgTxt("Too many output argumnents.\n");
// 檢查輸出變量數(shù)量是否正確，否則報(bào)錯(cuò)
#define A_IN prhs[0]
#define B_OUT plhs[0]
int M = mxGetM(A_IN);
int N = mxGetN(A_IN);
// 得到輸入矩陣A的行數(shù)和列數(shù)
B_OUT = mxCreateDoubleMatrix(M, N, mxREAL);
// 為輸出矩陣B分配存儲(chǔ)空間
double *A = mxGetPr(A_IN);
double *B = mxGetPr(B_OUT);
// 取得輸入矩陣A和輸出矩陣B的數(shù)據(jù)指針
for (int i = 0; i < M * N; ++i)
B[i] = discrete(A[i]);
// 調(diào)用discrete拼余，根據(jù)A(i, j)計(jì)算B(i, j)
}
int discrete(double d) {
if (d < 1.0 / 3.0)
return 0;
else if (d < 2.0 / 3.0)
return 1;
return 2;
}
4. 測(cè)試
編寫discrete的MATLAB版本discrete_matlab.m：
function [ B ] = discrete_matlab( A )
[M, N] = size(A);
B = zeros(M, N);
for i = 1:M
for j = 1:N
if A(i, j) < 1/3
B(i, j) = 0;
elseif A(i, j) < 2/3
B(i, j) = 1;
else
B(i, j) = 2;
end
end
end
end
編寫測(cè)試腳本test_discrete：
MATRIX_M = 100;
MATRIX_N = 100;
MATRIX_NUMBER = 1000;
data = rand(MATRIX_M, MATRIX_N, MATRIX_NUMBER);
tic
for i = 1:MATRIX_NUMBER
discrete_cpp(data(:, :, i));
end
toc
tic
for i = 1:MATRIX_NUMBER
discrete_matlab(data(:, :, i));
end
toc
也就是生成MATRIX_NUMBER個(gè)MATRIX_M * MATRIX_N的隨機(jī)矩陣污桦，然后分別用mex版本和普通matlab函數(shù)做discrete，看各自消耗的時(shí)間匙监。
測(cè)試結(jié)果……第一行是MEX版本的時(shí)間凡橱，第二行是MATLAB版本的時(shí)間，可以看出MEX方式在數(shù)據(jù)量大時(shí)還是挺有優(yōu)勢(shì)的亭姥。
所以這圖還真是各種生動(dòng)形象：
5. 參考資料
- Writing MATLAB C-MEX Code：
http://www.docin.com/p-684614383.html - 解決"mex setup或mbuild setup找不到編譯器"問題
http://blog.csdn.net/sanqima/article/details/49514889
