振動形成的單一尺寸硬球晶的結(jié)構(gòu)表征及其動力學(xué)的離散元模擬

安希忠，李長興，楊潤宇，余艾冰，曲德毅

(1.東北大學(xué) 材料與冶金學(xué)院，沈陽 110004;2.新南威爾士大學(xué) 材料科學(xué)與工程學(xué)院顆粒系統(tǒng)模擬中心，澳大利亞 悉尼新南威爾士 2052;3.沈陽鼓風(fēng)機集團股份有限公司，沈陽 110869)

振動形成的單一尺寸硬球晶的結(jié)構(gòu)表征及其動力學(xué)的離散元模擬

安希忠1，李長興1，楊潤宇2，余艾冰2，曲德毅3

(1.東北大學(xué) 材料與冶金學(xué)院，沈陽 110004;2.新南威爾士大學(xué) 材料科學(xué)與工程學(xué)院顆粒系統(tǒng)模擬中心，澳大利亞 悉尼新南威爾士 2052;3.沈陽鼓風(fēng)機集團股份有限公司，沈陽 110869)

采用離散元(DEM)數(shù)值方法結(jié)合物理實驗，分別實現(xiàn)了單一尺寸硬球在三維(3D)間歇振動及批量加料條件下堆積的有序化(結(jié)晶)，通過對晶體的結(jié)構(gòu)表征、數(shù)值仿真中的微觀性能如配位數(shù)(CN)、徑向分布函數(shù)(RDF)、角分布函數(shù)(ADF)、Voronoi/Delaunay孔的尺寸分布、力的網(wǎng)絡(luò)分析以及晶體形成過程中的動力學(xué)分析，得出如下結(jié)論:(1)數(shù)值仿真和物理實驗中所獲得的均為{111}取向FCC晶體，晶體中粒子的最大堆積密度可達0.74;(2)晶體中粒子的配位數(shù)分布在CN=12處出現(xiàn)峰值，表明每個粒子具有12個近鄰;RDF和ADF分布表明堆積粒子之間分別在徑向具有長程相關(guān)性及角度上具有相關(guān)性，這是晶態(tài)結(jié)構(gòu)的典型特征;Voronoi/Delaunay孔的尺寸分布與非晶態(tài)相比表現(xiàn)為高且窄的曲線，且向左移，表明所獲得的硬球晶中孔隙很小且分布均勻;粒子間力的網(wǎng)絡(luò)更直觀地表明了晶體的結(jié)構(gòu)，同時還指出在FCC硬球晶中分別存在著面內(nèi)的力和面間的力，其中后者小于前者;(3)晶體形成過程中的速度場和力場等動力學(xué)信息的變化表明粒子在堆積過程中一旦形成有序，它將作為一個整體運動，并為后續(xù)粒子的有序化提供模板，起到了結(jié)晶的晶核的作用.

顆粒堆積;致密化;振動;離散元模擬;硬球晶

無論是在工業(yè)生產(chǎn)還是在科學(xué)研究中，對硬球堆積致密化的研究一直是一個重要的課題［1，2］.因為對于許多物理學(xué)家和材料學(xué)家而言，硬球的堆積可用于研究多種系統(tǒng)結(jié)構(gòu)(直接研究這些系統(tǒng)的結(jié)構(gòu)是相當(dāng)困難的)的有效出發(fā)點.在過 去 的 幾 十 年 里，從 Bernal［3］，Scott［4］和Finney［5］等人首次利用硬球的堆積模擬簡單液體結(jié)構(gòu)取得成功以來，許多科學(xué)工作者又成功實現(xiàn)了硬球堆積非晶態(tài)(隨機堆積)的轉(zhuǎn)變［6～13］，即從初始的松裝密度約0.60轉(zhuǎn)變?yōu)榉蔷B(tài)結(jié)構(gòu)中硬球的最大堆積密度約0.64，在這個過程中，外部施加的振動起到了至關(guān)重要的作用.采用振動及合理地控制工藝參數(shù)也可以實現(xiàn)粒子的晶態(tài)堆積結(jié)構(gòu).這一點無論是從物理實驗上還是從數(shù)值計算中都得到證實［14～23］.但由于數(shù)值及物理實驗中的特殊性，如現(xiàn)有的物理實驗中對硬球結(jié)晶的動力學(xué)分析的難度，大多數(shù)數(shù)值仿真都是基于幾何的模擬而忽略了動力學(xué)的范疇［21］，從而導(dǎo)致所獲得的結(jié)果很難再現(xiàn)實際的過程，故此對于完美硬球晶的形成過程及其動力學(xué)還有待于進一步的研究.

使用基于分子動力學(xué)的離散元(DEM)數(shù)值仿真，可以克服其他數(shù)值仿真方法的不足，從而從動力學(xué)的角度真實再現(xiàn)粒子堆積的致密化，DEM方法的有效性已經(jīng)在許多科學(xué)研究中得到了證實［12，21，23～25］.在以前的工作中［21］，使用 DEM 模擬，我們已經(jīng)成功實現(xiàn)了振動條件下完美硬球晶的獲得.本文將在此基礎(chǔ)上，對該硬球晶的宏觀和微觀性能及其形成過程的動力學(xué)進行分析，以期找出硬球振動堆積過程中有序化的規(guī)律，為實際完美硬球晶的獲得提供參照.

1 模擬方法及條件

1.1 數(shù)值仿真方法

在DEM模擬中，每個粒子有兩種運動形式，即平動和轉(zhuǎn)動.根據(jù)牛頓第二定律，平動和轉(zhuǎn)動方程可分別表示為

其中 mi，vi，ωi和 Ii分別表示粒子 i的質(zhì)量、平動速度、角速度及轉(zhuǎn)動慣量.Ri是由粒子i的中心指向接觸點的矢量，其大小等于粒子i的半徑Ri.法向接觸力和切向接觸力分別表示為

1.2 模擬條件

模擬由在一個長方體容器中隨機產(chǎn)生一批無重疊的等同球形粒子開始，容器底面為邊長為12 d(d表示粒子直徑，模擬中所用的顆粒子直徑d=1 cm)的正方形，為了避免容器壁的影響，在水平方向上采用了周期性的邊界條件.初始時，隨機產(chǎn)生的粒子在重力作用下自由下落，向下運動的過程中粒子與其近鄰之間互相碰撞并來回彈跳躍.該動力學(xué)過程一直持續(xù)到所有粒子都達到穩(wěn)定狀態(tài)，此時由于能量完全耗散的結(jié)果，所有粒子的運動速度為零，這一過程在一秒鐘內(nèi)完成.一秒鐘后，容器在相互垂直的3個方向上進行三維(3D)周期性正弦振動，每個方向上振動的控制方程為:Z(t)=Asin［ω(t－t0)］，其中 Z 表示容器的位移，A和ω分別表示振動的振幅和頻率.t0表示振動開始時間(本工作中t0=1.0 s)，如果是1D振動，則只在豎直方向上施加周期性正弦振動.從位移方程很容易推導(dǎo)出振動的加速度方程:a(t)= －Aω2sin［ω(t－t0)］.本研究中采用的是間歇振動、批量加料，即從1 s開始，粒子先振動1 s，此過程中以一定的批量加料，然后停止1 s使粒子達到穩(wěn)定狀態(tài)，接著再振動1 s，此時再加入下一批料，之后再停止1 s，……，這樣的過程一直持續(xù)到振動結(jié)束.實驗中所用的粒子為玻璃球，為了簡化起見，取3D振動中3個方向上的振幅和頻率分別相同.表1中列出了模擬中所用的參數(shù)及其所對應(yīng)的取值，其中假定粒子與容器具有相同的性能.

表1 模擬中所用的參數(shù)Table 1 The parameters used in the simulation

2 研究結(jié)果及分析

2.1 硬球晶的實現(xiàn)

圖1給出了硬球在3D間歇振動、批量加料的條件下，結(jié)晶過程中所對應(yīng)的宏觀性能如堆積密度ρ隨振動時間t的演變及最終所獲得的局部完美硬球晶體結(jié)構(gòu).圖中的內(nèi)插圖是對密度變化曲線的局部放大，其中振幅A=0.1 d，振動頻率ω=200 rad/s，每次加料的批量為Nb=98個粒子/批［21］.可以看出在ρ－t曲線上存在著兩個區(qū)域，第一個區(qū)域一直持續(xù)到28.0 s，在這期間堆積密度隨時間增加.其中曲線上的起伏表示振動過程中對粒子堆積周期性壓縮與松馳的結(jié)果.在第二個區(qū)域，振動使堆積密度沿著一個相對穩(wěn)定的值變化.當(dāng)振動停止時，球逐漸停下來形成穩(wěn)定堆積，最終的堆積密度可達約0.74，結(jié)構(gòu)表征發(fā)現(xiàn)其對應(yīng)的是{111}取向的FCC晶體.圖1中的DEM數(shù)值仿真結(jié)果由圖2(其中的內(nèi)插圖為去除表面堆積粒子后所得到的內(nèi)部有序結(jié)構(gòu))中的物理實驗結(jié)果所驗證［22］.圖2中的結(jié)果是在5個不同尺寸容器中展開的單一尺寸玻璃球在3D間歇振動(振幅A=0.4 d，振動頻率ω=55 rad/s)及批量加料(Nb=1層/批)條件下得到的，為了消除容器壁的影響，我們對每種情況下不同容器中的實驗結(jié)果進行了外推，結(jié)果得到在無限大容器中粒子的最大堆積密度可達約0.74，且對最終的堆積結(jié)構(gòu)表征發(fā)現(xiàn)為{111}取向的FCC晶體，即在豎直的{111}方向上六方密排面的堆垛順序為ABCABC……，體現(xiàn)了數(shù)值仿真結(jié)果和物理實驗結(jié)果很好的一致性.

圖1 3D間歇振動、批量加料時，晶態(tài)結(jié)構(gòu)(內(nèi)插圖為原圖的局部放大)的實現(xiàn)(左)及最終所對應(yīng)的部分晶體(右)，振動條件為:A=0.1 d，ω=200 rad/s，加料批量Nb=98粒子/批，圖中A1為每一個振動周期內(nèi)的振動開始點(在此處開始加料)，A2為振動結(jié)束點［21］Fig.1 The realization of crystalline structure(left，the inset figure is the partial enlargement)and part of the final crystal(right)obtained under 3D interval vibration and batch－wised feeding，where the vibration condition is:A=0.1 d，ω =200 rad/s，batch of feeding Nb=98 particles/batch;A1represents the vibration starting point in each vibration period where the particles are fed at this moment，and A2is the end of vibration［21］

2.2 硬球晶的微觀結(jié)構(gòu)表征

圖2 物理實驗中，3D間歇振動、批量加料且A=0.4 d，ω=55 rad/s，Nb=1層/批時容器尺寸D(容器直徑)對堆積密度的影響，內(nèi)插圖為去除表面粒子而得到的有序結(jié)構(gòu)［22］Fig.2 The effects of container size on the packing density in physical experiments with 3D interval vibration and batch－wised feeding when A=0.4 d，ω =55 rad/s，and Nb=1 layer/batch，where the inset figures are the ordered structure obtained when the particles on the packing surface are removed［22］

圖3 晶態(tài)結(jié)構(gòu)的配位數(shù)分布Fig.3 Coordination number distribution of the crystalline structure

除了硬球晶宏觀性能的演變外，我們也對所獲得的硬球晶的微觀性能進行了表征，這些微觀性能包括:配位數(shù)(CN)、徑向分布函數(shù)(RDF)、角分布函數(shù)(ADF)、Voronoi/Delaunay孔的尺寸分布以及晶體內(nèi)力的網(wǎng)絡(luò)結(jié)構(gòu)，它們都是形成穩(wěn)態(tài)的硬球晶后得到的.圖3給出了硬球晶的配位數(shù)的分布，CN被定義為與給定粒子“接觸”的粒子數(shù)，本工作中取臨界距離為1.005 d，即如果兩個粒子之間的距離小于該值，就被認(rèn)為是接觸.與非晶態(tài)結(jié)構(gòu)CN分布的峰值為7不同地是［5，12，23，24］，圖3 中的 CN 分布的峰值對應(yīng)于12，這是FCC晶體對應(yīng)的配位數(shù)，表明我們所獲得的結(jié)構(gòu)是有序的.

與配位數(shù)CN相比，徑向分布函數(shù)RDF和角分布函數(shù)ADF是表征晶體結(jié)構(gòu)的更直接的方式.RDF定義為堆積中與特定粒子中心給定距離范圍內(nèi)發(fā)現(xiàn)一個粒子中心的概率，它包含了堆積粒子間徑向相關(guān)性的有用信息;相對地，ADF定義為堆積中與特定粒子中心給定角度范圍內(nèi)發(fā)現(xiàn)一個粒子中心的概率，它包含了堆積粒子間角度相關(guān)性的有用信息.圖4給出了我們所獲得的最終硬球晶的RDF和ADF分布，與非晶態(tài)(隨機密排)RDF 分布上的分叉第二峰不同［5，12，23，24］，晶態(tài)結(jié)構(gòu)的RDF分布還在r＞2d的許多相關(guān)位置上出現(xiàn)了各種孤立的峰，這些在較遠(yuǎn)位置上的峰的存在表明了堆積體內(nèi)粒子之間位置的長程相關(guān)性，這是晶態(tài)結(jié)構(gòu)的典型特征［21］.同樣，我們所獲得的晶態(tài)結(jié)構(gòu)的ADF也與非晶態(tài)結(jié)構(gòu)的ADF分布完全不同，后者在圖上表現(xiàn)在角度在60(°)～180(°)之間的一條連續(xù)的曲線，而前者的ADF只在60(°)、90(°)、120(°)和180(°)的角度上存在明顯的峰，這表明堆積體內(nèi)粒子之間存在著角度的相關(guān)性，這進一步表明了我們所獲得的硬球堆積結(jié)構(gòu)為晶態(tài)而不是非晶態(tài)［23］.

對所獲得硬球晶局部結(jié)構(gòu)更具體的分析來自于 Voronoi和 Delaunay 對空間的劃分［29～33］，它比圖4中的一維RDF和ADF能夠提供更豐富的信息［34］.Voronoi空間劃分是通過 Voronoi多面體來實現(xiàn)的，它定義為堆積體內(nèi)垂直平分給定粒子中心與其他粒子中心之間線段的面所圍成的最小封閉凸多面體［30］，如果連接堆積體內(nèi)每個粒子與其近鄰粒子中心，則會構(gòu)成一個空間網(wǎng)絡(luò)，該網(wǎng)絡(luò)則稱為Delaunay的空間劃分.網(wǎng)絡(luò)中的四面體單體形成了對空間的完全填充.圖5中給出了Voronoi(內(nèi)插圖為Delaunay)的孔的尺寸分布，我們將Voronoi或Delaunay單體內(nèi)孔的體積V除以粒子的體積Vp而使其無量綱化.由圖可見，無論是Voronoi還是Delaunay孔的尺寸分布，都呈現(xiàn)出高且窄的趨勢，并且明顯小于粒子的體積.這表明晶態(tài)結(jié)構(gòu)中孔的尺寸很小且分布很均勻;這一點與非晶態(tài)結(jié)構(gòu)中呈現(xiàn)右移的寬且矮的Voronoi/Delaunay孔的尺寸分布完全不同［12］.

為了更直觀的表征硬球晶的內(nèi)部結(jié)構(gòu)，我們也給出了最終部分硬球晶的力的網(wǎng)絡(luò)結(jié)構(gòu)，如圖6所示.圖中的“球”表示晶體中堆積的粒子，“棍”表示粒子之間的力，“棍”的粗細(xì)與力的大小成正比.從圖中除了可以很清晰地看出所形成的是{111}取向FCC晶體以外，更重要地是我們可以得到力在晶體中傳播的信息，即在晶體中存在著面內(nèi)的力和面間的力.在堆積中的每個球都可以看作是一個“外部載荷”，在這些“外部載荷”的作用下，力可以沿著三角架的方向傳遞，從而形成面間的網(wǎng)絡(luò)，這與Mueggenburg等人［35］的實驗結(jié)果相一致，所不同的是，在我們的結(jié)果有存在著較大的面內(nèi)的力.

2.3 硬球晶形成的動力學(xué)分析

除了前面對最終硬球晶的宏觀性能和微觀性能進行表征外，我們還從動力學(xué)的角度研究了硬球晶的形成過程.硬球結(jié)晶過程中的速度場和力場的變化分別如圖7和8所示.這兩個圖中顯示地是三維矢量場中的二維視場，其中視角的方向為沿著X軸的方向.從這兩個圖中我們可以發(fā)現(xiàn)堆積過程中出現(xiàn)的幾個特征.首先，是振動過程所帶來的對粒子的激發(fā)及松馳過程，其中對粒子的激發(fā)發(fā)生在t=53.318 6s，在這個過程中粒子的運動速度和力都很大.而粒子的松馳發(fā)生在t=47.044 2 s，此時粒子的速度和力逐漸減弱，漸漸形成穩(wěn)定結(jié)構(gòu).不管是激發(fā)還是松馳，粒子運動最劇烈地是在表面，這可以從速度和力的變化反映出來，這是因為表面粒子在振動堆積致密化過程中由于缺少約束而處于流化狀態(tài)，隨著外部輸入能量的停止，它們將找到各自的穩(wěn)定位置從而形成有序的結(jié)構(gòu).另一方面，在有序堆積過程中，先堆積的粒子的運動狀態(tài)基本相同，也就是說在底層形成有序結(jié)構(gòu)后，其中粒子的速度和力的變化不大，因此，在振動過程中其是作為一個整體運動的，它作為后續(xù)粒子有序化堆積的模板，起到了結(jié)晶的晶核的作用，這一點已經(jīng)在我們的數(shù)值仿真［23］和物理實驗［22］中得到了證實.

圖4 晶態(tài)結(jié)構(gòu)的徑向分布函數(shù)(左)和角分布函數(shù)(右)Fig.4 Radial distribution function(left)and angular distribution function(Right)for the crystalline structure

3 結(jié)論

采用離散元(DEM)數(shù)值方法模擬了單一尺寸球形粒子在三維(3D)間歇振動及批量加料條件下堆積的有序化(結(jié)晶)，并對所獲得的硬球晶的宏觀性能如堆積密度及微觀性能如配位數(shù)CN、徑向分布函數(shù) RDF、角分布函數(shù) ADF、Voronoi/Delaunay孔的尺寸分布及力的網(wǎng)絡(luò)結(jié)構(gòu)進行了表征，同時從動力學(xué)角度分析了硬球堆積有序化過程中的速度場和力場的變化，得出結(jié)論如下:

(1)對單一尺寸粒子在3D間歇振動及批量加料條件下堆積的DEM模擬可以獲得最大堆積密度達0.74的硬球晶，同時，該過程已得到了物理實驗的驗證.

(2)硬球晶的配位數(shù)(CN)分布在CN=12處出現(xiàn)峰值，結(jié)構(gòu)表征發(fā)現(xiàn)對應(yīng)的為FCC晶體;徑向分布函數(shù)(RDF)和角分布函數(shù)(ADF)的表征發(fā)現(xiàn)晶體中的粒子分別存在徑向的長程相關(guān)性及角度的相關(guān)性，這是晶態(tài)結(jié)構(gòu)區(qū)別于非晶態(tài)結(jié)構(gòu)的典型特征;Voronoi/Delaunay孔的尺寸分布表現(xiàn)為左移的高且窄的峰，表明堆積結(jié)構(gòu)中的孔的尺寸很小且分布均勻;力的網(wǎng)絡(luò)直觀的表現(xiàn)了晶體的結(jié)構(gòu)，同時表明所形成的{111}取向FCC硬球晶中，在面內(nèi)和面間均有力的存在，且面內(nèi)的力大于面間的力.

(3)通過對粒子振動堆積有序化過程中的動力學(xué)信息如速度場和力場變化的分析，發(fā)現(xiàn)硬球在結(jié)晶的過程中，一旦有序結(jié)構(gòu)形成，該有序結(jié)構(gòu)中的粒子在后面的振動過程中保持相同的運動狀態(tài)，即作為一個整體運動，它為后面有序化的繼續(xù)提供了模板，起到了晶核的作用.

［1］German R M. Particle packing characteristics［M］.Princeton，New Jersey:Metal Powder Industries，1989:1－18.

［2］Bideau D，Hansen A.Disorder and granular media，random materials and process［M］.Edited by Stanley H E，Guyon E，North － Holland.Amsterdam:Elsevier Science Publishers，1993:1－25.

［3］Bernal J D.A geometrical approach to structure of liquids［J］.Nature，1959，183(4655):141－147.

［4］Scott G D.Packing of spheres［J］.Nature，1960，188(4754):908－909.

［5］Finney J L.Random packings and the structure of simple liquids:I the geometry of random close packing ［J］.Proc Roy Soc Lond A，1970，319:479－493.

［6］Woodcock L V.Glass transition in the hard－sphere mode［J］.J Chem Soc Faraday II，1976，72:1667－1672.

［7］Finney J L.Modelling the structures of amorphous metals and alloys［J］.Nature，1977，266:309 －314.

［8］Zallen R.The physics of amorphous solids［M］.New York:Wiley，1983:1－15.

［9］Knight J B，F(xiàn)andrich C G，Lau C N，et al.Density relaxation in a vibrated granular material［J］.Phys Rev E，1995，51(5):3957－3963.

［10］Torquato S.Glass transition:hard knock for thermodynamics［J］.Nature，2000，405:521 －523.

［11］Stachurski Z H.Definition and properties of ideal amorphous solids［J］.Phys Rev Lett，2003，90(15):155502.

［12］An X Z，Yang R Y，Dong K J，et al.Micromechanical simulation and analysis of one－dimensional vibratory sphere packing［J］.Phys Rev Lett，2005，95:205502 －1 －4.

［13］An X Z，Li C X，Yang R Y，et al.Experimental study of the packing of mono－sized spheres subjected to one－dimensional vibration［J］.Powder Technol，2009，196:50 －55.

［14］Owe Berg T G，Mcdonald R L，Trainor R J.The packing of spheres［J］.Powder Technol，1969/1970，3:183 －188.

［15］Pouliquen O，Nicolas M，Weidman P D.Crystallization of non－Brownian spheres under horizontal shaking［J］.Phys Rev Lett，1997，79:3640 －3643.

［16］Van Blaaderen A，Ruel R，Wiltzius P.Template－directed colloidal crystallization［J］.Nature，1997，385:321 －324.

［17］Blair D L，Mueggenburg N W，Marshall A H，et al.Force distributions in three－dimensional granular assemblies:effects of packing order and interparticle friction［J］.Phys Rev E，2001，63:041304.

［18］Nahmad－Molinari Y，Ruiz－Suárez J C.Epitaxial growth of granular single crystals ［J］. Phys Rev Lett， 2002，89:264302.

［19］Kansal A R，Torquato S，Stillinger F H.Diversity of order and densities in jammed hard － particle packings［J］.Phys Rev E，2002，66:041109－1－8.

［20］Spannuth M J，Mueggenburg N W，Jaeger H M，et al.Stress transmission through three－dimensional granular crystals with stacking faults［J］.Granular Matter，2004，6:215 －219.

［21］Yu A B，An X Z，Zou R P，et al.Self－assembly of particles for densest packing by mechanical vibration［J］.Phys Rev Lett，2006，97:265501.

［22］Li C X，An X Z，Yang R Y，et al.Experimental study on the packing of uniform spheres under three－dimensional vibration［J］.Powder Technol，2011，208(3):617 － 622.

［23］An X Z，Yang R Y，Dong K J，et al.DEM study of crystallization of monosized spheres under mechanical vibrations［J］.Computer Physics Communications，2011，182:1989 －1994.

［24］Yang R Y，Zou R P，Yu A B.Computer simulation of the packing of fine particles［J］.Phys Rev E，2000，62(3):3900－3908.

［25］An X Z，Yang R Y，Zou R P，et al.Effect of vibration condition and inter－particle frictions on the packing of uniform spheres［J］.Powder Technol，2008，188(2):102－109.

［26］Martin C L， Bouvard D， Shima S. Studyofparticle rearrangement during powder compaction bythe discrete element method［J］.J Mech Phys Solids，2003，51(4):667－693.

［27］Schwager T，P?schel T.Coefficient of normal restitution of viscous particles and cooling rate of granular gases［J］.Phys Rev E，1998，57(1):650－654.

［28］Langston P A，Tüzün U，Heyes D M.Discrete element simulation of granular flow in 2D and 3D hoppers:dependence of discharge rate and wall stress on particle interactions［J］.Chem Eng Sci，1995，50(6):967 －987.

［29］Oger L，Gervois A，Troadec J P，et al.Voronoi tessellation of packing of spheres:topological correlation and statistics［J］.Phil Mag B，1996，74:177－197.

［30］Yang R Y，Zou R P，Yu A B.Voronoi tessellation of the packing of fine uniform spheres［J］.Phys Rev E，2002，65:041302.

［31］Aste T，Saadatfar M，Senden T J.Geometrical structure of disorderedspherepackings ［J］. Phys Rev E，2005，71:061302.

［32］Yang R Y，Zou R P，Yu A B，et al.Pore structure of the packing of fine particles［J］.J Colloid Interface Sci，2006，299:719－725.

［33］An X Z.Evolution of Voronoi/Delaunay characterized micro structure with transition from loose to dense sphere packing［J］.Chin Phys Lett，2007，24(8):2327 －2330.

［34］Montoro J C G，Abascal J L F.The Voronoi polyhedra as tools for structure determination in simple disordered systems［J］.J Phys Chem，1993，97:4211－4215.

［35］Mueggenburg N W，Jaeger H M，Nagel S R.Stress transmission through three－dimensional ordered granular arrays［J］.Phys Rev E，2002，66:031304－1.

Structure characterization on the vibration induced mono－sized hard sphere crystal and its dynamics by DEM modeling

AN Xi－zhong1，LI Chang－xing1，YANG Run－yu2，YU Ai－bing2，QU De－yi3

(1.School of Materials and Metallurgy，Northeastern University，Shenyang 110004，China;2.Laboratory for Simulation and Modeling of Particulate System，School of Materials Science and Engineering，University of New South Wales，Sydney 2052，Australia;3.Shenyang Blower Works Group，Co.，Ltd.，Shenyang 110869，China)

In this paper，the ordering(Crystallization)of mono－sized sphere packing under three－dimensional(3D)interval vibration and batch－wised feeding was studied with numerical discrete element method(DEM)combined with the physical experiments.Through the characterization on the crystalline structure and the analysis on the micro－properties obtained numerically such as coordination number(CN)，radial distribution function(RDF)，angular distribution function(ADF)，Voronoi/Delaunay void size distribution，force network and the dynamics formed during ordering，following conclusions can be drawn:(1)The crystalline obtained both numerically and physically all indicates{111} －oriented FCC structure，where the maximum packing density can reach 0.74;(2)There is a peak at CN=12 in the coordination number distribution curve，which implies that each particle in the hard sphere crystal has 12 neighbors in contact;the distributions for RDF and ADF show the long range correlation in radial distance and angle correlation among particles，which are the typical characteristics of the crystalline structure;the Voronoi/Delaunay void size distribution indicates the higher and narrower curve shifting to the left compared with non－crystal state，which means that the voids in crystalline structure are small and uniform;the force network among particles can identify the crystalline structure more intuitively，meanwhile，it’s also indicated that the in － plane forces and inter－ plane forces coexist in the crystal，and the latter is smaller than the former;(3)The variation of dynamic information such as the velocity field and the force field all indicates that once the ordered structure formed during crystallization，it will move as a whole part，which served as the template for the other particles，i.e.as the nucleus inthe crystallization.

particle packing;densification;vibration;DEM simulation;hard sphere crystal

TF 122;TB 44;O 756

A

1671－6620(2011)04－0318－07

2011－09－20.

國家自然基金 (50974040)，中央高?；究蒲袠I(yè)務(wù)費 (N100402009)，澳大利亞研究委員會基金 (DP0665276).

安希忠 (1973—)，男，遼寧蓋州人，東北大學(xué)副教授，E－mail:anxz@mail.neu.edu.cn