許衛(wèi)明, 王建華
(1.馬鞍山師范高等專(zhuān)科學(xué)校,安徽 馬鞍山 243041; 2.菏澤市城市管理局,山東 菏澤 274000)
?
OpenCL離散元顆粒系統(tǒng)的優(yōu)化設(shè)計(jì)*
許衛(wèi)明1, 王建華2
(1.馬鞍山師范高等專(zhuān)科學(xué)校,安徽 馬鞍山 243041; 2.菏澤市城市管理局,山東 菏澤 274000)
針對(duì)傳統(tǒng)離散元顆粒系統(tǒng)進(jìn)行大量顆粒運(yùn)算時(shí),出現(xiàn)的運(yùn)算時(shí)間長(zhǎng)、運(yùn)算效果不理想的問(wèn)題,通過(guò)OpenCL算法對(duì)其進(jìn)行優(yōu)化改進(jìn).將大量顆粒運(yùn)算分為顆粒間碰撞運(yùn)算和顆粒與邊界的碰撞運(yùn)算.顆粒間碰撞的運(yùn)算利用OpenCL編程實(shí)現(xiàn),與邊界碰撞的運(yùn)算利用CPU計(jì)算實(shí)現(xiàn),同時(shí)在OpenCL算法運(yùn)行時(shí),根據(jù)不同的數(shù)據(jù)類(lèi)型,劃分不同的存儲(chǔ)器,在內(nèi)存訪問(wèn)上提高算法的運(yùn)行速度.通過(guò)測(cè)試,優(yōu)化后的算法運(yùn)算效率明顯提升.
顆粒系統(tǒng);OpenCL;GPU;內(nèi)存訪問(wèn)
離散元方法的基本思想是把不連續(xù)體分離成剛性元素集合,利用運(yùn)動(dòng)方程來(lái)描述剛性元素運(yùn)動(dòng)情況,通過(guò)迭代的方法對(duì)運(yùn)動(dòng)方程進(jìn)行求解,從而得出不連續(xù)體的運(yùn)動(dòng)形態(tài).最初該方法是針對(duì)非連續(xù)介質(zhì)的力學(xué)行為的[1].在離散元的顆粒系統(tǒng)中,由于顆粒數(shù)量多,存在著大量的顆粒間的碰撞和顆粒與邊界的碰撞,對(duì)系統(tǒng)的力學(xué)進(jìn)行計(jì)算,需處理、跟蹤大量的信息,僅憑CPU的運(yùn)算,效果太差.近年來(lái)隨著GPU科學(xué)計(jì)算能力的逐步提高,尤其是GPU的并行計(jì)算能力的發(fā)展,將計(jì)算機(jī)的計(jì)算能力大幅提升,使得在對(duì)顆粒系統(tǒng)中數(shù)億的顆粒模擬計(jì)算成為可能,但整體運(yùn)算時(shí)間還是較長(zhǎng),本文對(duì)基于OpenCL離散元顆粒系統(tǒng)進(jìn)行研究,通過(guò)在算法上適應(yīng)GPU構(gòu)架,同時(shí)在內(nèi)存訪問(wèn)上進(jìn)行優(yōu)化,提升離散元的仿真效率.
OpenCL是由蘋(píng)果公司和AMD、NVIDIA、Intel等公司制定的一個(gè)工業(yè)標(biāo)準(zhǔn)和規(guī)范,是面向CPU、GPU與其他處理器組合構(gòu)成的計(jì)算機(jī)編程的行業(yè)標(biāo)準(zhǔn).該標(biāo)準(zhǔn)可通過(guò)公布硬件來(lái)提高編程的可移植性.平臺(tái)模型如圖1所示,平臺(tái)包括各種OpenCL設(shè)備及一個(gè)宿主機(jī),宿主機(jī)與各個(gè)OpenCL設(shè)備之間的交互通信,通過(guò)宿主機(jī)定義內(nèi)核,建立邏輯關(guān)系,對(duì)各個(gè)設(shè)備下達(dá)命令的方式來(lái)完成[2].
圖1 OpenCL平臺(tái)模型
OpenCL定義了五種內(nèi)存區(qū)域,首先為宿主機(jī)內(nèi)存,該內(nèi)存只有宿主機(jī)可見(jiàn);全局內(nèi)存,所有設(shè)備可見(jiàn),并允許讀寫(xiě)操作;常量?jī)?nèi)存,用于存儲(chǔ)常量,只允許讀;局部?jī)?nèi)存,存儲(chǔ)不同設(shè)備的公用常量;私有內(nèi)存,計(jì)算中一個(gè)工作組的專(zhuān)用內(nèi)存.在利用OpenCL進(jìn)行編程時(shí)一般遵循以下步驟:1)查支持OpenCL的硬件設(shè)備的參數(shù)、特性;2)創(chuàng)建OpenCL平臺(tái)的內(nèi)核函數(shù);3)在設(shè)備上創(chuàng)建上下邏輯關(guān)系,管理程序執(zhí)行時(shí)的內(nèi)核對(duì)象;4)傳遞所需的參數(shù);5)執(zhí)行內(nèi)核函數(shù).
離散元顆粒系統(tǒng)仿真是將顆粒系統(tǒng)中的顆粒假象為獨(dú)立的球形或橢球型顆粒,計(jì)算這些顆粒之間碰撞的受力情況,計(jì)算位移與速度,得出下一時(shí)步的位置和速度,循環(huán)此過(guò)程直到計(jì)算所有顆粒停止運(yùn)動(dòng)為止,或到達(dá)規(guī)定運(yùn)算時(shí)步,仿真結(jié)束.本次的優(yōu)化運(yùn)算首先將顆粒的碰撞分為兩類(lèi):顆粒與顆粒之間的碰撞、顆粒與邊界之間的碰撞.
2.1顆粒之間的碰撞計(jì)算
顆粒間的碰撞計(jì)算流程如圖2所示.在OpenCL平臺(tái)上創(chuàng)建上下邏輯關(guān)系,選擇硬件設(shè)備,創(chuàng)建命令隊(duì)列,并對(duì)命令在隊(duì)列中排隊(duì),從.cl文件中加載OpenCL內(nèi)核代碼.在內(nèi)存中為內(nèi)核函數(shù)的參數(shù)分配空間[3].對(duì)于顆粒系統(tǒng),函數(shù)的主要參數(shù)為:力學(xué)常量,顆粒的數(shù)量、坐標(biāo)、大小、速度和角速度等信息.由于內(nèi)核函數(shù)不能命名指針,只能將顆粒的信息轉(zhuǎn)換成一維數(shù)組,核函數(shù)最終得出的結(jié)果為該時(shí)步下的顆粒碰撞受力情況,然后通過(guò)分析顆粒的受力情況來(lái)計(jì)算顆粒的運(yùn)動(dòng)狀態(tài),得出顆粒的速度,位置信息等狀態(tài)的變化.
圖2 顆粒間碰撞計(jì)算流程
2.1.1接觸判斷與位移增量計(jì)算
首先建立局部坐標(biāo)系,以遍歷球顆粒i的球心為坐標(biāo)原點(diǎn),與顆粒j的球心連線作為x軸,選取x軸的垂直平面與X-Y平面平行的直線為y軸,再由右手螺旋法則確定z軸,建立空間坐標(biāo)系,如圖3所示.
圖3 顆粒間局部坐標(biāo)系
全局坐標(biāo)與局部坐標(biāo)的轉(zhuǎn)換矩陣為
在轉(zhuǎn)換矩陣中Lij表示兩球i和j的球心距,(Xi,Yi,Zi)、(Xj,Yj,Zj)表示兩球在全局坐標(biāo)系中的坐標(biāo),(l,m,n)是x軸在全局坐標(biāo)內(nèi)的方向向量.進(jìn)行接觸判斷時(shí),要將該顆粒與比該顆粒編號(hào)大的顆粒進(jìn)行判斷,計(jì)算過(guò)程通過(guò)GPU并行進(jìn)行計(jì)算,判定球體是否接觸的方法為:計(jì)算兩顆球的球心距與兩個(gè)球的半徑和的大小,當(dāng)球心距大于半徑和時(shí),兩球沒(méi)有接觸,當(dāng)球心距小于半徑和時(shí),兩球接觸.
2.1.2力學(xué)計(jì)算
兩球發(fā)生碰撞時(shí),根據(jù)受力分析,顆粒的受力包括法向接觸力、切向接觸力和回轉(zhuǎn)力矩.力學(xué)的計(jì)算過(guò)程在計(jì)算接觸判斷的線程中完成.
在t時(shí)步下,發(fā)生碰撞的的顆粒i在x軸方向的法向接觸力為
(1)
在公式中kx為剛性系數(shù),Cij為兩球疊加量,hx為法向阻尼系數(shù).
在t時(shí)步下沿局部坐標(biāo)系的y、z的作用力為切向彈性力與切向阻尼力的和,則沿y軸的切向接觸力為
(2)
根據(jù)相同的原理可得沿z軸的切向接觸力
(3)
在t時(shí)步下,回轉(zhuǎn)力矩為彈性分量和阻尼分量的和,則回轉(zhuǎn)力矩的計(jì)算公式為
M(t)=Mt-Δt-krΔφ-hrΔφ/Δt
(4)
kr,hr分別為回轉(zhuǎn)度系數(shù)和阻尼系數(shù).
以上的接觸判斷及力學(xué)計(jì)算,使用OpenCL算法,在GPU中完成的,選取平臺(tái)并創(chuàng)建命令隊(duì)列,主程序通過(guò)命令隊(duì)列對(duì)核函數(shù)進(jìn)行任務(wù)分配和操作運(yùn)算.GPU的多線程并行運(yùn)算與CPU的不同,每個(gè)線程訪問(wèn)內(nèi)存都是連續(xù)的.
2.2顆粒邊界的碰撞計(jì)算
顆粒系統(tǒng)的邊界分為好多種,如遠(yuǎn)行平面、扇形平面等,在進(jìn)行顆粒與邊界碰撞時(shí)要區(qū)分對(duì)待.
2.2.1接觸判斷
判斷過(guò)程中使用虛函數(shù)P2BContact()返回接觸標(biāo)志,以此來(lái)判定與邊界是否接觸.判定的方法為比較球心到邊界的距離和球顆粒的半徑,當(dāng)半徑大于到邊界的距離時(shí),則判定為接觸,否則不接觸.
設(shè)平面方程ax+by+cz+e=0,邊界的半徑為R,則球心到邊界的距離為
(5)
2.2.2顆粒與邊界的力學(xué)計(jì)算
與邊界的力學(xué)計(jì)算同顆粒間的力學(xué)計(jì)算大致相同,主要的區(qū)別在于顆粒相對(duì)于邊界自身的質(zhì)量可以忽略不計(jì),因此邊界不會(huì)因?yàn)轭w粒的碰撞而產(chǎn)生運(yùn)動(dòng)[1].計(jì)算之前要根據(jù)不同的材質(zhì)建立數(shù)據(jù)庫(kù),方便在程序錄入的情況下調(diào)用相對(duì)應(yīng)的參數(shù).
(6)
(7)
計(jì)算結(jié)束后對(duì)每個(gè)顆粒的作用力和力矩進(jìn)行相加之后再與核函數(shù)返回計(jì)算結(jié)果相加,根據(jù)受力情況計(jì)算顆粒的運(yùn)動(dòng)狀態(tài).
在OpenCL算法運(yùn)行時(shí),好的內(nèi)存訪問(wèn)模式可以加快讀、寫(xiě)緩沖區(qū)的速度,從而提升計(jì)算運(yùn)行速度和效率.傳統(tǒng)的CPU優(yōu)化內(nèi)存的方式不適應(yīng)OpenCL設(shè)備內(nèi)存,必須探尋新的內(nèi)存優(yōu)化技術(shù),來(lái)提高OpenCL設(shè)備的運(yùn)行速度.
3.1全局內(nèi)存訪問(wèn)
在OpenCL的核函數(shù)中應(yīng)用比較廣泛的為全局內(nèi)存,決定全局內(nèi)存性能的因素有訪存次數(shù)和訪存模式.訪存次數(shù)就是核函數(shù)中全局內(nèi)存被讀、寫(xiě)的次數(shù),次數(shù)過(guò)度就會(huì)降低程序整體的運(yùn)行效率,而在這些內(nèi)存中,被訪問(wèn)的次數(shù)最多的數(shù)據(jù)一般為工作組的共享數(shù)據(jù),解決的方法為利用局部?jī)?nèi)存訪問(wèn)速度比全局內(nèi)存快的特點(diǎn),將共享數(shù)據(jù)寫(xiě)入該工作組的局部?jī)?nèi)存中去[5].訪存模式,GPU中一般采用多線程分段處理,CPU的分段遍歷方式為每個(gè)線程都順序訪問(wèn)數(shù)組的每一個(gè)元素,如果采用這種方式,在訪問(wèn)第一個(gè)元素時(shí),破壞了GPU數(shù)據(jù)訪問(wèn)的局部要求,會(huì)對(duì)整個(gè)運(yùn)算效率造成影響,會(huì)降低運(yùn)算性能,解決的方法是采用合并讀取的方式,這樣可以將數(shù)組按線程分為不同的運(yùn)算區(qū)域,讀取的都是數(shù)組相鄰的元素,這樣既符合了OpenCL的內(nèi)存局部性,又有效的提高了緩沖區(qū)數(shù)據(jù)的讀、寫(xiě)速度.
3.2局部?jī)?nèi)存訪問(wèn)
局部?jī)?nèi)存是GPU每個(gè)線程組上都能快速訪問(wèn)的內(nèi)存區(qū)域,訪問(wèn)速度要比全局內(nèi)存速度快,利用局部?jī)?nèi)存可以減少對(duì)全局內(nèi)存的讀、寫(xiě)操作.程序運(yùn)行中如果無(wú)法將數(shù)據(jù)分段合并時(shí),將這些數(shù)據(jù)放入局部?jī)?nèi)存中,這樣可以提高緩沖區(qū)的讀寫(xiě)速度,從而提高運(yùn)行的效率.存儲(chǔ)體沖突時(shí)降低局部?jī)?nèi)存訪存性能下降的主要因素,本文通過(guò)采用均勻訪存的方法盡量避免沖突的發(fā)生.
均勻訪存的方法基本思想為,任意兩個(gè)相鄰的線程同時(shí)訪問(wèn)局部?jī)?nèi)存中位置間隔相同的區(qū)域,核心就是訪問(wèn)間隔,間隔不同,訪問(wèn)的區(qū)域就不同,這樣的訪問(wèn)也會(huì)造成不同的線程訪問(wèn)同一個(gè)內(nèi)存區(qū)域[6].OpenCL對(duì)局部?jī)?nèi)存的讀寫(xiě)映射到GPU的LDS上,OpenCL可以將32組線程為一組并行訪問(wèn).間隔為2的線程訪問(wèn)如圖4所示.
圖4 間隔為2的均勻訪問(wèn)
從圖4可看出存儲(chǔ)0和2處發(fā)生沖突,如果將間隔設(shè)為3,訪問(wèn)情況如圖5所示.
圖5 間隔為3的均勻訪問(wèn)
從圖中可以得出每一個(gè)線程都有不同的存儲(chǔ)體,整個(gè)周期內(nèi)訪問(wèn)均沒(méi)出現(xiàn)存儲(chǔ)體沖突.因此只要選擇的間隔合適,采用均勻訪存的方法可以避免存儲(chǔ)體沖突的發(fā)生,提高整個(gè)程序運(yùn)行的效率.
通過(guò)OpenCL算法對(duì)傳統(tǒng)離散元顆粒系統(tǒng)進(jìn)行優(yōu)化改進(jìn),將大量顆粒運(yùn)算分為顆粒間碰撞運(yùn)算和顆粒與邊界的碰撞運(yùn)算,并利用OpenCL編程實(shí)現(xiàn)顆粒間碰撞的運(yùn)算,利用CPU計(jì)算實(shí)現(xiàn)邊界碰撞的運(yùn)算,同時(shí)在OpenCL算法運(yùn)行時(shí),根據(jù)不同的數(shù)據(jù)類(lèi)型,劃分不同的存儲(chǔ)器,在內(nèi)存訪問(wèn)上提高算法的運(yùn)行速度.極大的發(fā)揮了GPU的并行計(jì)算能力,提升離散元顆粒系統(tǒng)的仿真效率.
[1]趙啦啦,趙躍民,劉初升,等.濕顆粒堆力學(xué)特性的離散元法模擬研究[J]. 物理學(xué)報(bào), 2014(03):1-9.
[2]季順迎,趙金鳳,狄少丞,等.面向環(huán)境力學(xué)的離散元分析軟件研發(fā)和工程應(yīng)用[J]. 計(jì)算機(jī)輔助工程, 2014(1):69-75.
[3] 詹云,趙新?tīng)N,譚同德.基于OpenCL的異構(gòu)系統(tǒng)并行編程[J]. 計(jì)算機(jī)工程與設(shè)計(jì), 2012(11):4191-4195.
[4] 陳鋼,吳百鋒.面向OpenCL模型的GPU性能優(yōu)化[J]. 計(jì)算機(jī)輔助設(shè)計(jì)與圖形學(xué)學(xué)報(bào), 2011(04):571-581.
[5] 韓博,周秉鋒.GPGPU性能模型及應(yīng)用實(shí)例分析[J]. 計(jì)算機(jī)輔助設(shè)計(jì)與圖形學(xué)學(xué)報(bào), 2009(09).
[6] 徐泳,孫其誠(chéng),張凌,等.顆粒離散元法研究進(jìn)展[J]. 力學(xué)進(jìn)展, 2003(02):251-259.
On Optimization Design of OPENCL Discrete Element Particle System
XU Wei-ming1, WANG Jian-hua2
(1.Maanshan Teacher’s College, Maanshan Anhui 243041, China 2. Heze City Administration Bureau, Heze Shandong 274000, China)
Due to the long operation time and unsatisfactory effect of traditional discrete element particle system carrying out a large number of granular computing, OPENCL algorithm is used for the optimization. The large number of particles are divided into inter particle collision and the collision between the particles and the boundary. The operation of inter particle collision is realized by OPENCL programming and the operation of the boundary collision is realized by CPU. At the same time, different memory is partitioned according to different data types in OPENCL algorithm. By testing, the optimized algorithm is obviously improved.
particle system; OPENCL; GPU; memory access
1673-2103(2016)02-0023-05
2015-12-07
安徽省一般教學(xué)研究項(xiàng)目(2013jyxm300 )
許衛(wèi)明,男(1981-),安徽懷寧人,講師,碩士,研究方向:計(jì)算機(jī)網(wǎng)絡(luò)
F416.67
A