【並行計算】用MPI進行分布式記憶體編程（一）

標籤：

    通過上一篇關於並行計算準備部分的介紹，我們知道MPI（Message-Passing-Interface 訊息傳遞介面）實現並行是進程層級的，通過通訊在進程之間進行訊息傳遞。MPI並不是一種新的開發語言，它是一個定義了可以被C、C++和Fortran程式調用的函數庫。這些函數庫裡面主要涉及的是兩個進程之間通訊的函數。MPI可以在Windows和linux環境中都有相應的庫，本篇以Windows10作為示範開發環境。

1、Windows10+VS 2015上搭建MPI開發環境

• 下載 mpi for windows

Windows為了相容MPI，自己做了一套基於一般個人電腦的MPI實現。如果要安裝正真意義上的MPI的話，請直接去www.mpich.org下載，裡面根據對應的系統下載相應的版本。

 

我由於是需要在個人筆記本上進行實驗，就用的是微軟的 HPC Pack 2008 R2 MS-MPI Redistributable Package with Service Pack 4 - 中文(簡體)，： http://www.microsoft.com/zh-cn/download/details.aspx?id=14737 。

• 安裝 mpi

   我的電腦是64位的，所以安裝的是 mpi_x64.msi ，預設安在C:\Program Files\Microsoft HPC Pack 2008 R2，在此，為了之後調試代碼方便，最好設定一下環境變數：在使用者變數PATH中，加入：C:\Program Files\Microsoft HPC Pack 2008 R2\Bin\。

• 配置2015

配置目錄，即載入Include和Lib庫

 

載入依賴項

• 編譯

幾乎所有人的第一個程式是從“hello，world”程式開始學習的，我也寫了一個這樣測試例子：

#include "mpi.h"  #include <stdio.h>  int main(int argc, char* argv[]){    int rank, numproces;    int namelen;    char processor_name[MPI_MAX_PROCESSOR_NAME];    MPI_Init(&argc, &argv);    MPI_Comm_rank(MPI_COMM_WORLD, &rank);//獲得進程號    MPI_Comm_size(MPI_COMM_WORLD, &numproces);//返回通訊子的進程數    MPI_Get_processor_name(processor_name, &namelen);    fprintf(stderr, "hello world! process %d of %d on %s\n", rank, numproces, processor_name);    MPI_Finalize();    return 0;}

上述代碼中，第1行中的#include "mpi.h" 標頭檔必須包含，在VS2015下編譯產生exe檔案（產生在debug檔案中），通過cmd命令，進入debug檔案夾目錄中，敲入：mpiexec –n 4 TestForMPI.exe。其中命令中-n 4 表示使用4個進程進行並行計算，具體結果：

 

2、理論知識

    通過上述的例子，我們對MPI編寫的並行計算有了一個初步的認識，但是我們還不知道如何去真正編寫一個MPI的並行程式，這需要我們學習一定的理論知識。在上面的例子中，有幾個函數對於初學MPI的人來說並不明白是什麼意思，下面就從這些函數入手。

MPI_Init：告知MPI系統進行所有必要的初始化設定。它是寫在啟動MPI並行計算的最前面的。具體的文法結構為：

MPI_Init(     int* argc_p,     char*** argv_p    );

    參數argc_p和argv_p分別指向main函數中的指標參數，為了弄明白這部分，還得從main函數的參數說起：Ｃ語言規定main函數的參數只能有兩個，習慣上這兩個參數寫為argc和argv。因此，main函數的函數頭可寫為： main (argc,argv)。Ｃ語言還規定argc(第一個形參)必須是整型變數,argv( 第二個形參)必須是指向字串的指標數組。其中argc參數表示了命令列中參數的個數(注意：檔案名稱本身也 算一個參數)，argc的值是在輸入命令列時由系統按實際參數的個數自動賦予的。例如有命令列為： C:">E6 24 BASIC dbase FORTRAN由於檔案名稱E6 24本身也算一個參數，所以共有4個參數，因此argc取得的值為4。argv參數是字串指標數組，其各元素值為命令列中各字串(參數均按字串處 理)的首地址。

然而在MPI_Init函數中，並不一定都需要設定argc_p和argv_p這兩個參數的，不需要的時候，將它們設定為NULL即可。

通訊子（communicator）：MPI_COMM_WORLD表示一組可以互相發送訊息的進程集合。

MPI_Comm_rank:用來擷取正在調用進程的通訊子中的進程號的函數。

MPI_Comm_size:用來得到通訊子的進程數的函數。

這兩個函數的具體結構如下：

int MPIAPI MPI_Comm_rank(    __in MPI_Comm comm,    __out int* rank    );int MPIAPI MPI_Comm_size(    __in MPI_Comm comm,    __out int* size    );

它們的第一個參數都傳入通訊子作為參數，第二參數都用傳出參數分別把正在調用通訊子的進程號和通訊的個數。

MPI_Finalize：告知MPI系統MPI已經使用完畢。它總是放到做並行計算的功能塊的最後面，在此函數之後就不能再出現任何有關MPI相關的東西了。

以上只是表達了作為一個MPI並行計算的基本結構，並沒有真正涉及進程之間的通訊，為了更好的進行並行，必然需要進程間的通訊，下面介紹兩個處理序間通訊的函數，它們就是MPI_Send和MPI_Recv，分別用於訊息的發送和接收。

MPI_Send：阻塞型訊息發送。其結構為：

int MPI_Send (void *buf, int count, MPI_Datatype datatype,int dest, int tag,MPI_Comm comm)

參數buf為發送緩衝區；count為發送的資料個數；datatype為發送的資料類型；dest為訊息的目的地址(進程號)，其取值範圍為0到np－1間的整數(np代表通訊器comm中的進程數) 或MPI_PROC_NULL；tag為訊息標籤，其取值範圍為0到MPI_TAG_UB間的整數；comm為通訊器。

MPI_Recv：阻塞型訊息接收。

int MPI_Recv (void *buf, int count, MPI_Datatype datatype,int source, int tag, MPI_Comm comm,MPI_Status *status)

參數buf為接收緩衝區；count為資料個數，它是接收資料長度的上限，具體接收到的資料長度可通過調用MPI_Get_count函數得到；datatype為接收的資料類型；source為訊息源地址(進程號)，其取值範圍為0到np－1 間的整數(np代表通訊器comm 中的進程數)，或MPI_ANY_SOURCE，或MPI_PROC_NULL；tag為訊息標籤，其取值範圍為0到MPI_TAG_UB間的整數或MPI_ANY_TAG；comm為通訊器；status返回接收狀態。

MPI_Status：返回訊息傳遞的完成情況。資料結構的相關變數的意義就比較多了，具體可以參數使用手冊。

typedef struct {... ...int MPI_SOURCE;             /*訊息源地址*/int MPI_TAG;                /*訊息標籤*/int MPI_ERROR;              /*錯誤碼*/... ...} MPI_Status;

3、舉例

介紹了最基本的處理序間通訊的函數，我們就能編寫一個更為複雜也更為有意義的程式，通過編寫程式來實現資料積分中的梯形積分法。

梯形積分法的基本思想是：將x軸上的區間劃分為n個等長的子區間。然後計運算元區間的和。

假設子區間的端點為xi和xi+1，那麼子區間的長度為：h=xi+1-xi。那麼梯形的面積就為：

由於n個子區間是等分的，邊界分別為xi=a和x=b，則：

這片地區的所有梯形的面積和為：

變換為：

因此，串列的程式碼就可以這樣寫：

h = (b - a) / h;approx = (f(a) + f(b)) / 2;for (int i = 1; i < n-1; i++){    x_i = a + i*h;    approx += f(x_i);}approx = h*approx;

通過對串列程式的分析，對於這個例子，我們可以識別出兩種任務：第一種擷取單個矩形地區的面積，另一種是計算這些地區的面積和。

假設求f(x)=x3將梯形劃分為1024個子領域計算[0,3]地區內的積分。

#include "mpi.h"  #include <stdio.h>  #include <cmath>double Trap(double left_endpt, double right_endpt, double trap_count, double base_len);double f(double x);int main(int argc, char* argv[]){    int my_rank = 0, comm_sz = 0, n = 1024, local_n = 0;    double a = 0.0, b = 3.0, h = 0, local_a = 0, local_b = 0;    double local_int = 0, total_int = 0;    int source;    MPI_Init(&argc, &argv);    MPI_Comm_rank(MPI_COMM_WORLD, &my_rank);    MPI_Comm_size(MPI_COMM_WORLD, &comm_sz);    h = (b - a) / n;       /*  h is the same for all processes  */    local_n = n / comm_sz; /*  So is the number of trapezoids */    local_a = a + my_rank*local_n*h;    local_b = local_a + local_n*h;    local_int = Trap(local_a, local_b, local_n, h);    if (my_rank != 0)    {        MPI_Send(&local_int, 1, MPI_DOUBLE, 0, 0, MPI_COMM_WORLD);    }    else    {        total_int = local_int;        for (source = 1; source < comm_sz; source++)        {            MPI_Recv(&local_int, 1, MPI_DOUBLE, source, 0,                MPI_COMM_WORLD, MPI_STATUS_IGNORE);            total_int += local_int;        }    }    if (my_rank == 0)    {        printf("With n = %d trapezoids, our estimate\n", n);        printf("of the integral from %f to %f = %.15e\n", a, b, total_int);            }    MPI_Finalize();    return 0;}//子領域的積分函數double Trap(double left_endpt, double right_endpt, double trap_count, double base_len){    double estimate = 0, x = 0;    int i;    estimate = (f(left_endpt) + f(right_endpt)) / 2.0;    for (i = 1; i <= trap_count - 1; i++)    {        x = left_endpt + base_len;        estimate += f(x);    }    estimate = estimate*base_len;    return estimate;}//數學函數double f(double x){    return pow(x, 3);}

上述代碼中，運行結果

 

這段程式碼中的意思是，通過輸入的進程數，將1024個劃分的子區間等分的分配到控制台輸入的進程（100個）進行子任務求解，求解完成之後，1-99進程計算的結果值通過MPI_Send函數發送出去，而0號進程使用MPI_Recv函數接收匯總，將每個進程的結果求和，得到區間的積分值。

本次並行計算訊息之間的通訊如：

至此，我們已經能使用MPI_Send訊息發送函數和MPI_Recv訊息接收函數進行簡單的並行程式計算了，但我們想一下，最後的求和都是用0號進程去做，為了更加的提高效能，還需更進一步使用集合通訊，下一章節將深的進行講解。

【並行計算】用MPI進行分布式記憶體編程（一）