專利名稱：一種基于cpu和gpu協(xié)作的分子動力學(xué)加速方法
技術(shù)領(lǐng)域：
本發(fā)明涉及一種分子動力學(xué)加速方法，尤指基于CPU和GPU協(xié)作的分子動力學(xué)加速方法。
背景技術(shù)：
分子動力學(xué)方法是一種依靠牛頓力學(xué)定律來模擬分子體系的運動，從而計算分子體系宏觀性質(zhì)的方法。它是研究納米尺度物理現(xiàn)象的重要手段，廣泛應(yīng)用于材料科學(xué)、生物物理和藥物設(shè)計等領(lǐng)域。隨著模擬精度的提高和模擬規(guī)模的增大，分子動力學(xué)模擬時間也隨之增長，從而對其實際應(yīng)用產(chǎn)生了很大影響。目前通常的做法有以下幾種1、通過組建集群系統(tǒng)來并行的進行模擬。但這種方法硬件成本較高，同時也需要開發(fā)人員不斷改進并行算法來挖掘集群系統(tǒng)性能。因此在考慮硬件、人力成本時，不適宜采用該方法。2、通過定制專用計算機來加速模擬。針對分子動力學(xué)模擬的不同階段定制相應(yīng)的專用硬件進行計算，從而達到很高的模擬速度。但是這種機器硬件成本很高，目前該方法很少采用。3、利用通用加速單元GPU (Graphics Processing Unit，圖形處理單元)來協(xié)作普通計算機進行模擬。目前GPU廠商為用戶提供了通用編程模型如CUDA(統(tǒng)一計算架構(gòu)模型)等對GPU進行編程，利用GPU浮點計算能力比CPU強的特點，使用GPU完成分子動力學(xué)中計算任務(wù)，從而提高模擬速度。但是該方法目前存在一些局限性，主要表現(xiàn)在使用GPU計算時，計算機中CPU處于閑置狀態(tài)，計算機的總體利用率較低。第三種方法由于成本較低，是目前前景較大的分子動力學(xué)加速方法。但是目前該方法仍存在一些問題一方面CPU和GPU的計算模型不同，GPU適合一次性處理大量數(shù)據(jù)， 因此它要求分配的計算任務(wù)不能太小，而大塊的計算任務(wù)又不利于CPU進行調(diào)度，為了使 CPU和GPU都能高效計算，必須設(shè)計大小適中的計算任務(wù)；另一方面CPU和GPU完成計算任務(wù)的時間難以估量，靜態(tài)地分配計算任務(wù)會造成CPU或GPU出現(xiàn)閑置等待的情況，應(yīng)根據(jù)兩者的忙閑動態(tài)調(diào)整計算任務(wù)的分配，以保證兩者皆不出現(xiàn)閑置等待的情況。但目前還沒有公開的技術(shù)方案能較好解決這兩方面的問題。

發(fā)明內(nèi)容
本發(fā)明要解決的技術(shù)問題在于針對當(dāng)前基于GPU的分子動力學(xué)加速方法效率較低的問題，提出一種通過協(xié)同使用CPU和GPU的加速方法，使系統(tǒng)中CPU和GPU都能高效計算，從而在保持較低成本的情況下實現(xiàn)分子動力學(xué)的加速。為了解決上述技術(shù)問題，本發(fā)明提出的技術(shù)方案為使用配置有多核CPU和GPU的計算機，對待模擬的分子系統(tǒng)建立“網(wǎng)格-單元”結(jié)構(gòu)存儲分子信息，使用CPU和GPU動態(tài)的申請和處理網(wǎng)格，在處理網(wǎng)格時，CPU以單元為單位處理網(wǎng)格，GPU則直接處理整個網(wǎng)格， 從而實現(xiàn)CPU和GPU高效協(xié)作地完成分子動力學(xué)模擬。
具體技術(shù)方案為第一步、構(gòu)建配置有多核CPU和GPU的計算機系統(tǒng)。選擇包含Q個O! > = 4)核的CPU作為計算機系統(tǒng)的處理器，選擇支持CUDA的GPU通過PCI-E接口與CPU相連。CPU 的Q個核順序編號為0，...，Q-I。第二步、從存儲于硬盤的文件中讀入分子系統(tǒng)相關(guān)參數(shù)，包括截斷半徑I 。和所有分子的數(shù)據(jù)量。其中截斷半徑I 。表示每個分子與其他分子發(fā)生作用的最大距離，如果兩個分子之間的距離超過R。，則它們之間不發(fā)生作用。分子數(shù)據(jù)量包括分子總數(shù)目MN以及每個分子的編號、位置向量和速度向量。在CPU的內(nèi)存中申請二維雙精度浮點型數(shù)組Pos存儲分子的位置向量，申請二維雙精度浮點型數(shù)組Vel存儲分子的速度向量。第三步、設(shè)置模擬參數(shù)，包括目標(biāo)網(wǎng)格權(quán)重TW，目標(biāo)模擬步數(shù)T和實際模擬步數(shù)t。 TW和T的值由模擬需求確定，Tff的值為期望網(wǎng)格包含的分子數(shù)目，T的值為需要模擬的步數(shù)，實際模擬步數(shù)t初始化為0。第四步、使用“網(wǎng)格-單元”結(jié)構(gòu)表示分子系統(tǒng)，“網(wǎng)格-單元”結(jié)構(gòu)的建立過程如下4. 1使用“單元”結(jié)構(gòu)表示分子系統(tǒng)，步驟為4. 1. 1將分子系統(tǒng)劃分為均勻組合的邊長為Rc的η個正方體區(qū)域，稱這些正方體區(qū)域為單元，將這些單元順序編號為0，. . .，η-1。4. 1. 2在CPU的內(nèi)存中初始化結(jié)構(gòu)體數(shù)組SCell保存單元中的分子編號以及單元的鄰接單元編號，稱該結(jié)構(gòu)體數(shù)組為“單元”結(jié)構(gòu)，使用符號SCell[k] (0 ^ k < η)表示序號為k的單元。每個“單元”擁有一個變量成員和兩個數(shù)組成員，分別為單元內(nèi)分子數(shù)目AN、 單元內(nèi)分子編號數(shù)組AIdx和鄰接單元編號數(shù)組NB。AIdx順序存儲單元中分子編號。NB順序存儲單元在三維空間中鄰接沈個單元的編號。用符號SCell [k]. AN表示序號為k的單元包含的分子數(shù)目；用符號SCell [k]. AIdx表示序號為k的單元中分子編號數(shù)組；用符號 SCell [k]. NB表示序號為k的單元鄰接單元編號數(shù)組。4.2根據(jù)單元內(nèi)分子數(shù)目SCelUk].AN(0<k<n)將多個單元組成區(qū)域，使每個區(qū)域中的分子數(shù)目相等，稱這些區(qū)域為網(wǎng)格。在CPU的內(nèi)存中用鏈表保存網(wǎng)格中分子信息以及網(wǎng)格處理標(biāo)識，鏈表的每一個節(jié)點稱為一個“網(wǎng)格”結(jié)構(gòu)。每個“網(wǎng)格”結(jié)構(gòu)擁有一個指針、一個整型變量、一個布爾變量、一個整型數(shù)組和四個雙精度浮點型數(shù)組，分別為后續(xù)“網(wǎng)格”指針Next，網(wǎng)格內(nèi)單元數(shù)目CN，網(wǎng)格處理標(biāo)識R，網(wǎng)格內(nèi)單元編號數(shù)組CIdx，網(wǎng)格內(nèi)分子位置信息數(shù)組Ppos，網(wǎng)格內(nèi)分子速度信息數(shù)組Pvel，網(wǎng)格內(nèi)分子受力信息數(shù)組Pforce， 網(wǎng)格外分子位置信息數(shù)組Gpos。為了便于描述，用SPatch表示指向鏈表首節(jié)點的指針，稱這個節(jié)點為SPatch指向的網(wǎng)格；用符號“一”引用節(jié)點中的元素。例如，用SPatch —Next 表示指向下一個網(wǎng)格的指針；用SPatch —CN表示SPatch指向的網(wǎng)格中單元數(shù)目；用 SPatch — R表示SPatch指向的網(wǎng)格處理標(biāo)識，若SPatch — R的值為iTrue則表示SPatch 指向的網(wǎng)格已經(jīng)被處理，若SPatch — R的值為False則表示SPatch指向的網(wǎng)格未被處理； 用SPatch — CIdx表示SPatch指向的網(wǎng)格中單元編號數(shù)組，用SPatch — CIdx[χ] (O^x < SPatch — CN)存儲網(wǎng)格中序號為χ的單元編號；用SPatch — Ppos表示SPatch指向的網(wǎng)格內(nèi)分子位置信息數(shù)組；用SPatch — Pvel表示SPatch指向的網(wǎng)格內(nèi)分子速度信息數(shù)組；用SPatch — Pforce表示SPatch指向的網(wǎng)格內(nèi)分子受力信息數(shù)組；用SPatch — Gpos表示SPatch指向的網(wǎng)格外分子位置信息數(shù)組。生成網(wǎng)格結(jié)構(gòu)的步驟為4. 2. 1計算每個單元的權(quán)重W。，W1, ...，Wlri，其中Wk (0彡k < η)的值為單元k中分子數(shù)目SCell [k]. AN。初始化當(dāng)前網(wǎng)格負(fù)載PW為0，起始單元編號St為0，結(jié)束單元編號 Ed為0。初始化指針CP,使指針CP = SPatch04. 2. 2將指向下一個網(wǎng)格的指針CP — Next初始化為NULL，CP指向的網(wǎng)格中單元數(shù)目CP — CN初始化為0，CP指向的網(wǎng)格處理標(biāo)識CP — R初始化為i^alse，CP指向的數(shù)組 CP — CIdx, CP — Ppos, CP — Pvel.CP — Pforce, CP — Gpos 全部初始化為 O。4. 2. 3判斷PW是否小于目標(biāo)網(wǎng)格權(quán)重TW，若是，執(zhí)行4. 2. 4，否則轉(zhuǎn)4. 2. 6。4. 2. 4判斷&是否等于n-1，若是，轉(zhuǎn)4. 2. 6，否則執(zhí)行4. 2. 5。4. 2.bPW = PW + WEd+l，& = Ed+1，轉(zhuǎn) 4· 2· 3。4. 2. 6 將編號 St，St+1，...，4順序填入 CP — CIdx 數(shù)組中。根據(jù) SCell [a]. AIdx (St ^ a < Ed)中的分子編號將Pos和Vel數(shù)組對應(yīng)分子的位置向量和速度向量順序填入CP — Ppos和CP — Pvel中。順序檢查SCell [a]. AIdx (St彡a < Ed)數(shù)組中的單元編號 b，若b < St或b > Ed，則按SCell [b]. AIdx中的分子編號將Pos和Vel數(shù)組對應(yīng)分子的位置向量順序填入CP —Gpos中。將CP —CN賦值為4-St+l。判斷4是否等于n_l，若是， 轉(zhuǎn)第五步，否則執(zhí)行4. 2.7。4. 2. 7 使 PW = 0，St = Ed+1, Ed = Ed+1, CP = CP — Next。轉(zhuǎn) 4. 2. 2。第五步、啟動兩個線程控制GPU和CPU并行的處理網(wǎng)格，稱控制GPU的線程為O號線程，控制CPU的線程為1號線程，分別運行在CPU編號為O和1的核上。并行處理網(wǎng)格的方法為0號線程初始化指針Ntl = SPatchdf N。一 R置為True，執(zhí)行第六步；1號線程初始化指針N1 = SPatch — Next，將N1 — R置為True，轉(zhuǎn)第九步。第六步、O號線程控制GPU對Ntl指向的網(wǎng)格進行處理，具體做法為6. 1將Ntl指向的網(wǎng)格的數(shù)據(jù)拷貝到GPU的內(nèi)存中。6. 2啟動GPU線程采用牛頓力學(xué)方程計算N。中每個分子所受的力，根據(jù)力和初始速度計算每個分子當(dāng)前速度和位置，將計算結(jié)果存放在GPU的內(nèi)存中。6. 3將計算結(jié)果從GPU的內(nèi)存拷貝到位于CPU的內(nèi)存中的N。一 Ppos和N。一 Pvel 數(shù)組中。第七步判斷Ntl — Next是否為NULL，若是，轉(zhuǎn)第十二步，否則執(zhí)行第八步。第八步、使N。= N。一 Next，判斷N。一 R是否為True，若是轉(zhuǎn)第七步，否則將N。一 R 置為True并轉(zhuǎn)第六步。第九步、1號線程根據(jù)CPU的核數(shù)Q啟動Q-2個子線程Ttl，T1，...，TQ_3運行于CPU的編號為2，3，...，Q-I的核上，Tc (O ^ C < Q-3)采用牛頓力學(xué)方程計算單元
權(quán)利要求
1. 一種基于CPU和GPU協(xié)作的分子動力學(xué)加速方法，其特征在于包括以下步驟 第一步、構(gòu)建配置有Q核CPU和GPU的計算機系統(tǒng)，Q > = 4，采用的GPU支持CUDA，通過PCI-E接口與CPU相連，CPU的Q個核順序編號為0，...，GL1 ；第二步、從存儲于硬盤的文件中讀入分子系統(tǒng)相關(guān)參數(shù)，包括截斷半徑I 。和所有分子的數(shù)據(jù)量，其中截斷半徑I 。表示每個分子與其他分子發(fā)生作用的最大距離，分子數(shù)據(jù)量包括分子總數(shù)目MN以及每個分子的編號、位置向量和速度向量；在CPU的內(nèi)存中申請二維雙精度浮點型數(shù)組Pos存儲分子的位置向量，申請二維雙精度浮點型數(shù)組Vel存儲分子的速度向量；第三步、設(shè)置模擬參數(shù)，包括目標(biāo)網(wǎng)格權(quán)重TW，目標(biāo)模擬步數(shù)T和實際模擬步數(shù)t ； 第四步、使用“網(wǎng)格-單元”結(jié)構(gòu)表示分子系統(tǒng)，“網(wǎng)格-單元”結(jié)構(gòu)的建立過程如下 4. 1使用“單元”結(jié)構(gòu)表示分子系統(tǒng)，步驟為-4. 1. 1將分子系統(tǒng)劃分為均勻組合的邊長為Rc的η個正方體區(qū)域，稱這些正方體區(qū)域為單元，將這些單元順序編號為0，. . .，η-1 ；-4. 1.2在CPU的內(nèi)存中初始化結(jié)構(gòu)體數(shù)組SCell保存單元中的分子編號以及單元的鄰接單元編號，稱該結(jié)構(gòu)體數(shù)組為“單元”結(jié)構(gòu)，使用符號SCell [k]表示序號為k的單元， O彡k < η ；每個“單元”擁有一個變量成員和兩個數(shù)組成員，分別為單元內(nèi)分子數(shù)目AN、 單元內(nèi)分子編號數(shù)組Aldx和鄰接單元編號數(shù)組NB ；Aldx順序存儲單元中分子編號，NB順序存儲單元在三維空間中鄰接沈個單元的編號；用符號SCell [k]. AN表示序號為k的單元包含的分子數(shù)目，用符號SCell [k]. AIdx表示序號為k的單元中分子編號數(shù)組，用符號 SCelUk]. NB(表示序號為k的單元鄰接單元編號數(shù)組；-4. 2根據(jù)單元內(nèi)分子數(shù)目SCell [k]. AN將多個單元組成區(qū)域，稱這些區(qū)域為網(wǎng)格，生成網(wǎng)格結(jié)構(gòu)的步驟為-4. 2. 1計算每個單元的權(quán)重Wtl，W1,...，Wn+1，其中Wk的值為單元k中分子數(shù)目SCelUk]. AN,初始化當(dāng)前網(wǎng)格負(fù)載PW為0，起始單元編號St為0，結(jié)束單元編號&為0，初始化指針 CP，使指針CP = SPatch ；在CPU的內(nèi)存中用鏈表保存網(wǎng)格中分子信息以及網(wǎng)格處理標(biāo)識， 鏈表的每一個節(jié)點稱為一個“網(wǎng)格”結(jié)構(gòu)；每個“網(wǎng)格”結(jié)構(gòu)擁有一個指針、一個整型變量、一個布爾變量、一個整型數(shù)組和四個雙精度浮點型數(shù)組，分別為后續(xù)“網(wǎng)格”指針Next，網(wǎng)格內(nèi)單元數(shù)目CN，網(wǎng)格處理標(biāo)識R，網(wǎng)格內(nèi)單元編號數(shù)組CIdx，網(wǎng)格內(nèi)分子位置信息數(shù)組Ppos， 網(wǎng)格內(nèi)分子速度信息數(shù)組Pvel，網(wǎng)格內(nèi)分子受力信息數(shù)組Pforce，網(wǎng)格外分子位置信息數(shù)組Gpos ；SPatch表示指向鏈表首節(jié)點的指針，稱這個節(jié)點為SPatch指向的網(wǎng)格；用符號 “一”引用節(jié)點中的元素；-4. 2. 2將指向下一個網(wǎng)格的指針CP — Next初始化為NULL，CP指向的網(wǎng)格中單元數(shù)目CP — CN初始化為0，CP指向的網(wǎng)格處理標(biāo)識CP — R初始化為i^alse，CP指向的數(shù)組 CP — CIdx、CP — Ppos、CP — Pvel、CP — Pforce、CP — Gpos 全部初始化為 O ； 4. 2. 3判斷Pff是否小于目標(biāo)網(wǎng)格權(quán)重Tff,若是，執(zhí)行4. 2.4,否則轉(zhuǎn)4. 2. 6 ； ---4. 2. 4判斷&是否等于η-1，若是，轉(zhuǎn)4. 2. 6，否則執(zhí)行4. 2. 5 ； 4. 2.^PW = PW + WEd+1，& = Ed+1，轉(zhuǎn) 4· 2· 3 ；-4. 2. 6將編號St，St+1，. . · ,Ed順序填入CP — CIdx數(shù)組中，根據(jù)SCell [a]. AIdx中的分子編號將Pos和Vel數(shù)組對應(yīng)分子的位置向量和速度向量順序填入CP — Ppos和CP — Pvel中，St彡a < Ed ；順序檢查SCell[a] · AIdx數(shù)組中的單元編號b，若b < St或b > Ed，則按 SCelltb]. AIdx中的分子編號將Pos和Vel數(shù)組對應(yīng)分子的位置向量順序填入CP — Gpos 中，將CP — CN賦值為；判斷&是否等于n-1，若是，轉(zhuǎn)第五步，否則執(zhí)行4. 2. 7 ； 4. 2. 7 使 PW = 0，St = Ed+1, Ed = Ed+1, CP = CP — Next，轉(zhuǎn) 4. 2. 2 ； 第五步、啟動兩個線程控制GPU和CPU并行的處理網(wǎng)格，稱控制GPU的線程為O號線程，控制CPU的線程為1號線程，分別運行在CPU編號為O和1的核上；并行處理網(wǎng)格的方法為0號線程初始化指針Ntl = SPatchJf N。一 R置為True，執(zhí)行第六步；1號線程初始化指針N1 = SPatch — Next，將N1 — R置為True，轉(zhuǎn)第九步；第六步、O號線程控制GPU對Ntl指向的網(wǎng)格進行處理，具體做法為 6. 1將Ntl指向的網(wǎng)格的數(shù)據(jù)拷貝到GPU的內(nèi)存中；·6. 2啟動GPU線程采用牛頓力學(xué)方程計算Ntl中每個分子所受的力，根據(jù)力和初始速度計算每個分子當(dāng)前速度和位置，將計算結(jié)果存放在GPU的內(nèi)存中；·6. 3將計算結(jié)果從GPU的內(nèi)存拷貝到位于CPU的內(nèi)存中的N。一 Ppos和N。一 Pvel數(shù)組中；第七步判斷Ntl — Next是否為NULL，若是，轉(zhuǎn)第十二步，否則執(zhí)行第八步； 第八步、使N。= N。一 Next，判斷N。一 R是否為True，若是轉(zhuǎn)第七步，否則將N。一 R置為True并轉(zhuǎn)第六步；第九步、1號線程根據(jù)CPU的核數(shù)Q啟動Q-2個子線程Ttl，T1, ...，TQ_3運行于CPU 的編號為2，3，...，Q-I的核上，Τ。采用牛頓力學(xué)方程計算單元 ，
2.如權(quán)利要求1所述的一種基于CPU和GPU協(xié)作的分子動力學(xué)加速方法，其特征在于所述目標(biāo)網(wǎng)格權(quán)重TW和目標(biāo)模擬步數(shù)T的值由模擬需求確定，Tff的值為期望網(wǎng)格包含的分子數(shù)目，T的值為需要模擬的步數(shù)，實際模擬步數(shù)t初始化為0。
全文摘要
本發(fā)明公開了一種基于CPU和GPU協(xié)作的分子動力學(xué)加速方法，目的是提出一種通過協(xié)同使用CPU和GPU的加速方法，在較低成本的情況下實現(xiàn)分子動力學(xué)的加速。技術(shù)方案是使用配置有多核CPU和GPU的計算機，對待模擬的分子系統(tǒng)建立“網(wǎng)格-單元”結(jié)構(gòu)存儲分子信息，使用CPU和GPU動態(tài)申請和處理網(wǎng)格，在處理網(wǎng)格時，CPU以單元為單位處理網(wǎng)格，GPU則直接處理整個網(wǎng)格，實現(xiàn)CPU和GPU高效協(xié)作地完成分子動力學(xué)模擬。采用本發(fā)明可使得CPU核之間負(fù)載均衡，避免靜態(tài)分配計算任務(wù)時CPU和GPU可能出現(xiàn)閑置等待的情況，CPU和GPU都能充分發(fā)揮自己的計算性能，提高了整個計算機系統(tǒng)的利用率，以低成本實現(xiàn)了分子動力學(xué)的加速。
文檔編號G06F17/50GK102411658SQ20111038200
公開日2012年4月11日 申請日期2011年11月25日 優(yōu)先權(quán)日2011年11月25日
發(fā)明者吳強, 左克, 廖湘科, 彭林, 李春江, 杜云飛, 楊燦群, 石志才, 陳娟 申請人:中國人民解放軍國防科學(xué)技術(shù)大學(xué)