最近在看多核編程。簡(jiǎn)單來(lái)說(shuō)，由于現(xiàn)在電腦CPU一般都有兩個(gè)核，4核與8核的CPU也逐漸走入了尋常百姓家，傳統(tǒng)的單線程編程方式難以發(fā)揮多核CPU的強(qiáng)大功能，于是多核編程應(yīng)運(yùn)而生。按照我的理解，多核編程可以認(rèn)為是對(duì)多線程編程做了一定程度的抽象，提供一些簡(jiǎn)單的API，使得用戶不必花費(fèi)太多精力來(lái)了解多線程的底層知識(shí)，從而提高編程效率。這兩天關(guān)注的多核編程的工具包括openMP和TBB。按照目前網(wǎng)上的討論，TBB風(fēng)頭要蓋過(guò)openMP，比如openCV過(guò)去是使用openMP的，但從2.3版本開(kāi)始拋棄openMP，轉(zhuǎn)向TBB。但我試下來(lái)，TBB還是比較復(fù)雜的，相比之下，openMP則非常容易上手。因?yàn)榫蜁r(shí)間有限，沒(méi)辦法花費(fèi)太多時(shí)間去學(xué)習(xí)TBB，就在這里分享下這兩天學(xué)到的openMP的一點(diǎn)知識(shí)，和大家共同討論。
openMP支持的編程語(yǔ)言包括C語(yǔ)言、C++和Fortran，支持OpenMP的編譯器包括Sun Studio，Intel Compiler，Microsoft Visual Studio，GCC。我使用的是Microsoft Visual Studio 2008，CPU為Intel i5 四核，首先講一下在Microsoft Visual Studio 2008上openMP的配置。非常簡(jiǎn)單，總共分2步：
(1) 新建一個(gè)工程。這個(gè)不再多講。
(2) 建立工程后，點(diǎn)擊 菜單欄->Project->Properties，彈出菜單里，點(diǎn)擊 Configuration Properties->C/C++->Language->OpenMP Support，在下拉菜單里選擇Yes。
至此配置結(jié)束。下面我們通過(guò)一個(gè)小例子來(lái)說(shuō)明openMP的易用性。這個(gè)例子是 有一個(gè)簡(jiǎn)單的test()函數(shù)，然后在main()里，用一個(gè)for循環(huán)把這個(gè)test()函數(shù)跑8遍。
1
#include <iostream>
2
#include <time.h>
3
void
test()
4
{
5
int
a =
0
;
6
for
(
int
i=
0
;i<
100000000
;i++)
7
a++;
8
}
9
int
main()
10
{
11
clock_t t1 = clock();
12
for
(
int
i=
0
;i<
8
;i++)
13
test();
14
clock_t t2 = clock();
15
std::cout<<
"
time:
"
<<t2-t1<<std::endl;
16
}
編譯運(yùn)行后，打印出來(lái)的耗時(shí)為：1.971秒。下面我們用一句話把上面代碼變成多核運(yùn)行。
1
#include <iostream>
2
#include <time.h>
3
void
test()
4
{
5
int
a =
0
;
6
for
(
int
i=
0
;i<
100000000
;i++)
7
a++;
8
}
9
int
main()
10
{
11
clock_t t1 = clock();
12
#pragma
omp parallel for
13
for
(
int
i=
0
;i<
8
;i++)
14
test();
15
clock_t t2 = clock();
16
std::cout<<
"
time:
"
<<t2-t1<<std::endl;
17
}
編譯運(yùn)行后，打印出來(lái)的耗時(shí)為：0.546秒，幾乎為上面時(shí)間的1/4。
由此我們可以看到openMP的簡(jiǎn)單易用。在上面的代碼里，我們一沒(méi)有額外include頭文件，二沒(méi)有額外link庫(kù)文件，只是在for循環(huán)前加了一句#pragma omp parallel for。而且這段代碼在單核機(jī)器上，或者編譯器沒(méi)有將openMP設(shè)為Yes的機(jī)器上編譯也不會(huì)報(bào)錯(cuò)，將自動(dòng)忽略#pragma這行代碼，然后按照傳統(tǒng)單核串行的方式編譯運(yùn)行！我們唯一要多做的一步，是從C:\Program Files\Microsoft Visual Studio 9.0\VC\redist\x86\Microsoft.VC90.OPENMP和C:\Program Files\Microsoft Visual Studio 9.0\VC\redist\Debug_NonRedist\x86\Microsoft.VC90.DebugOpenMP目錄下分別拷貝vcomp90d.dll和vcomp90.dll文件到工程文件當(dāng)前目錄下。
對(duì)上面代碼按照我的理解做個(gè)簡(jiǎn)單的剖析。
當(dāng)編譯器發(fā)現(xiàn)#pragma omp parallel for后，自動(dòng)將下面的for循環(huán)分成N份，(N為電腦CPU核數(shù))，然后把每份指派給一個(gè)核去執(zhí)行，而且多核之間為并行執(zhí)行。下面的代碼驗(yàn)證了這種分析。
1
#include <iostream>
2
int
main()
3
{
4
#pragma
omp parallel for
5
for
(
int
i=
0
;i<
10
;i++)
6
std::cout<<i<<std::endl;
7
return
0
;
8
}
會(huì)發(fā)現(xiàn)控制臺(tái)打印出了0 3 4 5 8 9 6 7 1 2。注意：因?yàn)槊總€(gè)核之間是并行執(zhí)行，所以每次執(zhí)行時(shí)打印出的順序可能都是不一樣的。
下面我們來(lái)了談?wù)劯?jìng)態(tài)條件(race condition)的問(wèn)題，這是所有多線程編程最棘手的問(wèn)題。該問(wèn)題可表述為，當(dāng)多個(gè)線程并行執(zhí)行時(shí)，有可能多個(gè)線程同時(shí)對(duì)某變量進(jìn)行了讀寫操作，從而導(dǎo)致不可預(yù)知的結(jié)果。比如下面的例子，對(duì)于包含10個(gè)整型元素的數(shù)組a，我們用for循環(huán)求它各元素之和，并將結(jié)果保存在變量sum里。
1
#include <iostream>
2
int
main()
3
{
4
int
sum =
0
;
5
int
a[
10
] = {
1
,
2
,
3
,
4
,
5
,
6
,
7
,
8
,
9
,
10
};
6
#pragma
omp parallel for
7
for
(
int
i=
0
;i<
10
;i++)
8
sum = sum + a[i];
9
std::cout<<
"
sum:
"
<<sum<<std::endl;
10
return
0
;
11
}
如果我們注釋掉#pragma omp parallel for，讓程序先按照傳統(tǒng)串行的方式執(zhí)行，很明顯，sum = 55。但按照并行方式執(zhí)行后，sum則會(huì)變成其他值，比如在某次運(yùn)行過(guò)程中，sum = 49。其原因是，當(dāng)某線程A執(zhí)行sum = sum + a[i]的同時(shí)，另一線程B正好在更新sum，而此時(shí)A還在用舊的sum做累加，于是出現(xiàn)了錯(cuò)誤。
那么用openMP怎么實(shí)現(xiàn)并行數(shù)組求和呢？下面我們先給出一個(gè)基本的解決方案。該方案的思想是，首先生成一個(gè)數(shù)組sumArray，其長(zhǎng)度為并行執(zhí)行的線程的個(gè)數(shù)(默認(rèn)情況下，該個(gè)數(shù)等于CPU的核數(shù))，在for循環(huán)里，讓各個(gè)線程更新自己線程對(duì)應(yīng)的sumArray里的元素，最后再將sumArray里的元素累加到sum里，代碼如下
1
#include <iostream>
2
#include <omp.h>
3
int
main(){
4
int
sum =
0
;
5
int
a[
10
] = {
1
,
2
,
3
,
4
,
5
,
6
,
7
,
8
,
9
,
10
};
6
int
coreNum = omp_get_num_procs();
//
獲得處理器個(gè)數(shù)
7
int
* sumArray =
new
int
[coreNum];
//
對(duì)應(yīng)處理器個(gè)數(shù)，先生成一個(gè)數(shù)組
8
for
(
int
i=
0
;i<coreNum;i++)
//
將數(shù)組各元素初始化為0
9
sumArray[i] =
0
;
10
#pragma
omp parallel for
11
for
(
int
i=
0
;i<
10
;i++)
12
{
13
int
k = omp_get_thread_num();
//
獲得每個(gè)線程的ID
14
sumArray[k] = sumArray[k]+a[i];
15
}
16
for
(
int
i =
0
;i<coreNum;i++)
17
sum = sum + sumArray[i];
18
std::cout<<
"
sum:
"
<<sum<<std::endl;
19
return
0
;
20
}
需要注意的是，在上面代碼里，我們用omp_get_num_procs()函數(shù)來(lái)獲取處理器個(gè)數(shù)，用omp_get_thread_num()函數(shù)來(lái)獲得每個(gè)線程的ID，為了使用這兩個(gè)函數(shù)，我們需要include <omp.h>。
上面的代碼雖然達(dá)到了目的，但它產(chǎn)生了較多的額外操作，比如要先生成數(shù)組sumArray，最后還要用一個(gè)for循環(huán)將它的各元素累加起來(lái)，有沒(méi)有更簡(jiǎn)便的方式呢？答案是有，openMP為我們提供了另一個(gè)工具，歸約(reduction)，見(jiàn)下面代碼：
1
#include <iostream>
2
int
main(){
3
int
sum =
0
;
4
int
a[
10
] = {
1
,
2
,
3
,
4
,
5
,
6
,
7
,
8
,
9
,
10
};
5
#pragma
omp parallel for reduction(+:sum)
6
for
(
int
i=
0
;i<
10
;i++)
7
sum = sum + a[i];
8
std::cout<<
"
sum:
"
<<sum<<std::endl;
9
return
0
;
10
}
上面代碼里，我們?cè)?pragma omp parallel for 后面加上了 reduction(+:sum)，它的意思是告訴編譯器：下面的for循環(huán)你要分成多個(gè)線程跑，但每個(gè)線程都要保存變量sum的拷貝，循環(huán)結(jié)束后，所有線程把自己的sum累加起來(lái)作為最后的輸出。
reduction雖然很方便，但它只支持一些基本操作，比如+,-,*,&,|,&&,||等。有些情況下，我們既要避免race condition，但涉及到的操作又超出了reduction的能力范圍，應(yīng)該怎么辦呢？這就要用到openMP的另一個(gè)工具，critical。來(lái)看下面的例子，該例中我們求數(shù)組a的最大值，將結(jié)果保存在max里。
1
#include <iostream>
2
int
main(){
3
int
max =
0
;
4
int
a[
10
] = {
11
,
2
,
33
,
49
,
113
,
20
,
321
,
250
,
689
,
16
};
5
#pragma
omp parallel for
6
for
(
int
i=
0
;i<
10
;i++)
7
{
8
int
temp = a[i];
9
#pragma
omp critical
10
{
11
if
(temp > max)
12
max = temp;
13
}
14
}
15
std::cout<<
"
max:
"
<<max<<std::endl;
16
return
0
;
17
}
上例中，for循環(huán)還是被自動(dòng)分成N份來(lái)并行執(zhí)行，但我們用#pragma omp critical將 if (temp > max) max = temp 括了起來(lái)，它的意思是：各個(gè)線程還是并行執(zhí)行for里面的語(yǔ)句，但當(dāng)你們執(zhí)行到critical里面時(shí)，要注意有沒(méi)有其他線程正在里面執(zhí)行，如果有的話，要等其他線程執(zhí)行完再進(jìn)去執(zhí)行。這樣就避免了race condition問(wèn)題，但顯而易見(jiàn)，它的執(zhí)行速度會(huì)變低，因?yàn)榭赡艽嬖诰€程等待的情況。
有了以上基本知識(shí)，對(duì)我來(lái)說(shuō)做很多事情都足夠了。下面我們來(lái)看一個(gè)具體的應(yīng)用例，從硬盤讀入兩幅圖像，對(duì)這兩幅圖像分別提取特征點(diǎn)，特征點(diǎn)匹配，最后將圖像與匹配特征點(diǎn)畫出來(lái)。理解該例子需要一些圖像處理的基本知識(shí)，我不在此詳細(xì)介紹。另外，編譯該例需要opencv，我用的版本是2.3.1，關(guān)于opencv的安裝與配置也不在此介紹。我們首先來(lái)看傳統(tǒng)串行編程的方式。
1
#include
"
opencv2/highgui/highgui.hpp
"
2
#include
"
opencv2/features2d/features2d.hpp
"
3
#include <iostream>
4
#include <omp.h>
5
int
main( ){
6
cv::SurfFeatureDetector detector(
400
);
7
cv::SurfDescriptorExtractor extractor;
8
cv::BruteForceMatcher<cv::L2<
float
> > matcher;
9
std::vector< cv::DMatch > matches;
10
cv::Mat im0,im1;
11
std::vector<cv::KeyPoint> keypoints0,keypoints1;
12
cv::Mat descriptors0, descriptors1;
13
double
t1 = omp_get_wtime( );
14
//
先處理第一幅圖像
15
im0 = cv::imread(
"
rgb0.jpg
"
, CV_LOAD_IMAGE_GRAYSCALE );
16
detector.detect( im0, keypoints0);
17
extractor.compute( im0,keypoints0,descriptors0);
18
std::cout<<
"
find
"
<<keypoints0.size()<<
"
keypoints in im0
"
<<std::endl;
19
//
再處理第二幅圖像
20
im1 = cv::imread(
"
rgb1.jpg
"
, CV_LOAD_IMAGE_GRAYSCALE );
21
detector.detect( im1, keypoints1);
22
extractor.compute( im1,keypoints1,descriptors1);
23
std::cout<<
"
find
"
<<keypoints1.size()<<
"
keypoints in im1
"
<<std::endl;
24
double
t2 = omp_get_wtime( );
25
std::cout<<
"
time:
"
<<t2-t1<<std::endl;
26
matcher.match( descriptors0, descriptors1, matches );
27
cv::Mat img_matches;
28
cv::drawMatches( im0, keypoints0, im1, keypoints1, matches, img_matches );
29
cv::namedWindow(
"
Matches
"
,CV_WINDOW_AUTOSIZE);
30
cv::imshow(
"
Matches
"
, img_matches );
31
cv::waitKey(
0
);
32
return
1
;
33
}
很明顯，讀入圖像，提取特征點(diǎn)與特征描述子這部分可以改為并行執(zhí)行，修改如下：
1
#include
"
opencv2/highgui/highgui.hpp
"
2
#include
"
opencv2/features2d/features2d.hpp
"
3
#include <iostream>
4
#include <vector>
5
#include <omp.h>
6
int
main( ){
7
int
imNum =
2
;
8
std::vector<cv::Mat> imVec(imNum);
9
std::vector<std::vector<cv::KeyPoint>>keypointVec(imNum);
10
std::vector<cv::Mat> descriptorsVec(imNum);
11
cv::SurfFeatureDetector detector(
400
); cv::SurfDescriptorExtractor extractor;
12
cv::BruteForceMatcher<cv::L2<
float
> > matcher;
13
std::vector< cv::DMatch > matches;
14
char
filename[
100
];
15
double
t1 = omp_get_wtime( );
16
#pragma
omp parallel for
17
for
(
int
i=
0
;i<imNum;i++){
18
sprintf(filename,
"
rgb%d.jpg
"
,i);
19
imVec[i] = cv::imread( filename, CV_LOAD_IMAGE_GRAYSCALE );
20
detector.detect( imVec[i], keypointVec[i] );
21
extractor.compute( imVec[i],keypointVec[i],descriptorsVec[i]);
22
std::cout<<
"
find
"
<<keypointVec[i].size()<<
"
keypoints in im
"
<<i<<std::endl;
23
}
24
double
t2 = omp_get_wtime( );
25
std::cout<<
"
time:
"
<<t2-t1<<std::endl;
26
matcher.match( descriptorsVec[
0
], descriptorsVec[
1
], matches );
27
cv::Mat img_matches;
28
cv::drawMatches( imVec[
0
], keypointVec[
0
], imVec[
1
], keypointVec[
1
], matches, img_matches );
29
cv::namedWindow(
"
Matches
"
,CV_WINDOW_AUTOSIZE);
30
cv::imshow(
"
Matches
"
, img_matches );
31
cv::waitKey(
0
);
32
return
1
;
33
}
兩種執(zhí)行方式做比較，時(shí)間為：2.343秒v.s. 1.2441秒
在上面代碼中，為了改成適合#pragma omp parallel for執(zhí)行的方式，我們用了STL的vector來(lái)分別存放兩幅圖像、特征點(diǎn)與特征描述子，但在某些情況下，變量可能不適合放在vector里，此時(shí)應(yīng)該怎么辦呢？這就要用到openMP的另一個(gè)工具，section，代碼如下：
1
#include
"
opencv2/highgui/highgui.hpp
"
2
#include
"
opencv2/features2d/features2d.hpp
"
3
#include <iostream>
4
#include <omp.h>
5
int
main( ){
6
cv::SurfFeatureDetector detector(
400
); cv::SurfDescriptorExtractor extractor;
7
cv::BruteForceMatcher<cv::L2<
float
> > matcher;
8
std::vector< cv::DMatch > matches;
9
cv::Mat im0,im1;
10
std::vector<cv::KeyPoint> keypoints0,keypoints1;
11
cv::Mat descriptors0, descriptors1;
12
double
t1 = omp_get_wtime( );
13
#pragma
omp parallel sections
14
{
15
#pragma
omp section
16
{
17
std::cout<<
"
processing im0
"
<<std::endl;
18
im0 = cv::imread(
"
rgb0.jpg
"
, CV_LOAD_IMAGE_GRAYSCALE );
19
detector.detect( im0, keypoints0);
20
extractor.compute( im0,keypoints0,descriptors0);
21
std::cout<<
"
find
"
<<keypoints0.size()<<
"
keypoints in im0
"
<<std::endl;
22
}
23
#pragma
omp section
24
{
25
std::cout<<
"
processing im1
"
<<std::endl;
26
im1 = cv::imread(
"
rgb1.jpg
"
, CV_LOAD_IMAGE_GRAYSCALE );
27
detector.detect( im1, keypoints1);
28
extractor.compute( im1,keypoints1,descriptors1);
29
std::cout<<
"
find
"
<<keypoints1.size()<<
"
keypoints in im1
"
<<std::endl;
30
}
31
}
32
double
t2 = omp_get_wtime( );
33
std::cout<<
"
time:
"
<<t2-t1<<std::endl;
34
matcher.match( descriptors0, descriptors1, matches );
35
cv::Mat img_matches;
36
cv::drawMatches( im0, keypoints0, im1, keypoints1, matches, img_matches );
37
cv::namedWindow(
"
Matches
"
,CV_WINDOW_AUTOSIZE);
38
cv::imshow(
"
Matches
"
, img_matches );
39
cv::waitKey(
0
);
40
return
1
;
41
}
上面代碼中，我們首先用#pragma omp parallel sections將要并行執(zhí)行的內(nèi)容括起來(lái)，在它里面，用了兩個(gè)#pragma omp section，每個(gè)里面執(zhí)行了圖像讀取、特征點(diǎn)與特征描述子提取。將其簡(jiǎn)化為偽代碼形式即為：
1
#pragma
omp parallel sections
2
{
3
#pragma
omp section
4
{
5
function1();
6
}
7
#pragma
omp section
8
{
9
function2();
10
}
11
}
意思是：parallel sections里面的內(nèi)容要并行執(zhí)行，具體分工上，每個(gè)線程執(zhí)行其中的一個(gè)section，如果section數(shù)大于線程數(shù)，那么就等某線程執(zhí)行完它的section后，再繼續(xù)執(zhí)行剩下的section。在時(shí)間上，這種方式與人為用vector構(gòu)造for循環(huán)的方式差不多，但無(wú)疑該種方式更方便，而且在單核機(jī)器上或沒(méi)有開(kāi)啟openMP的編譯器上，該種方式不需任何改動(dòng)即可正確編譯，并按照單核串行方式執(zhí)行。
以上分享了這兩天關(guān)于openMP的一點(diǎn)學(xué)習(xí)體會(huì)，其中難免有錯(cuò)誤，歡迎指正。另外的一點(diǎn)疑問(wèn)是，看到各種openMP教程里經(jīng)常用到private,shared等來(lái)修飾變量，這些修飾符的意義和作用我大致明白，但在我上面所有例子中，不加這些修飾符似乎并不影響運(yùn)行結(jié)果，不知道這里面有哪些講究。
在寫上文的過(guò)程中，參考了包括以下兩個(gè)網(wǎng)址在內(nèi)的多個(gè)地方的資源，不再一 一列出，在此一并表示感謝。
http://blog.csdn.net/drzhouweiming/article/details/4093624
http://software.intel.com/zh-cn/articles/more-work-sharing-with-openmp
更多文章、技術(shù)交流、商務(wù)合作、聯(lián)系博主
微信掃碼或搜索：z360901061
微信掃一掃加我為好友
QQ號(hào)聯(lián)系： 360901061
您的支持是博主寫作最大的動(dòng)力，如果您喜歡我的文章，感覺(jué)我的文章對(duì)您有幫助，請(qǐng)用微信掃描下面二維碼支持博主2元、5元、10元、20元等您想捐的金額吧，狠狠點(diǎn)擊下面給點(diǎn)支持吧，站長(zhǎng)非常感激您！手機(jī)微信長(zhǎng)按不能支付解決辦法：請(qǐng)將微信支付二維碼保存到相冊(cè)，切換到微信，然后點(diǎn)擊微信右上角掃一掃功能，選擇支付二維碼完成支付。
【本文對(duì)您有幫助就好】元
