基于離散元法的雪崩效應(yīng)的仿真研究

馮靖禹 韓韌 張宇峰 辛昱鋆

摘? 要： 伴隨著全球氣候變暖，雪崩災(zāi)害的風(fēng)險(xiǎn)逐漸增加。因此，調(diào)查雪崩的發(fā)展規(guī)律和內(nèi)在機(jī)理十分迫切和必要。而使用實(shí)驗(yàn)的方法來(lái)研究往往效率低下，并且只能獲得宏觀數(shù)據(jù)。為了更深入地了解雪崩效應(yīng)，提出結(jié)合計(jì)算機(jī)仿真技術(shù)和轉(zhuǎn)鼓雪崩實(shí)驗(yàn)，通過(guò)建立基于離散元法（DEM）的仿真模型來(lái)模擬發(fā)生雪崩時(shí)的顆粒流態(tài)。通過(guò)運(yùn)用Hertz-Mindlin接觸理論和牛頓第二定律塑造的仿真模型以及在相同配置下的物理實(shí)驗(yàn)的對(duì)比，證明了離散元法研究雪崩效應(yīng)的有效性?；谠撃Ｐ?，我們跟蹤每個(gè)顆粒的運(yùn)動(dòng)狀態(tài)，并研究了自由表面上參與雪崩的顆粒數(shù)量的變化。參與雪崩的顆粒數(shù)量的增長(zhǎng)速度經(jīng)歷了先增加，然后恒定，最后減少的過(guò)程。我們還發(fā)現(xiàn)不同區(qū)域內(nèi)的顆粒平均速度隨時(shí)間都呈現(xiàn)先增大后減小的規(guī)律。離散元仿真方法為雪崩顆粒流的研究提供了參考。

關(guān)鍵詞： 水平轉(zhuǎn)鼓;雪崩;仿真;離散元法

中圖分類號(hào)： TP391. 9 ???文獻(xiàn)標(biāo)識(shí)碼： A??? DOI：10.3969/j.issn.1003-6970.2020.07.055

本文著錄格式：馮靖禹，韓韌，張宇峰，等. 基于離散元法的雪崩效應(yīng)的仿真研究[J]. 軟件，2020，41（07）：264-268

Simulation of Avalanche Effect Based on Discrete Element Method

FENG Jing-yu， HAN Ren^*， ZHANG Yu-feng， XIN Yu-yun

（School of Optical-Electrical and Computer Engineering， University of Shanghai for Science and Technology， Shanghai 200093， China）

【Abstract】： With global warming， the risk of avalanche disasters is increasing. Therefore， it is urgent and necessary ?to investigate the development rule and internal mechanism of avalanches. Studing by experimental methods is often inefficient and can only obtain macro-data. In order to learn more about the avalanche effect， a combination of computer simulation technology and drum avalanche experiments is proposed. The particle flow pattern during an avalanche is simulated by establishing a simulation model based on the discrete element method （DEM）. By comparing the simulation model using Hertz-Mindlin contact theory and Newton's second law with the physical experiments under the same configuration， the validity of studing avalanche effect by discrete element method is proved. Based on this model， we track the motion of each particle and investigate the quantitative change of particles involved in avalanches on the free surface. The growth rate of avalanching particles number experiences a process of increasing， then stabilizing， and finally reducing. We also find that the average velocity of particles in different regions increase at first and then decrease. The discrete element method provides a reference for the research of avalanche particle flow.

【Key words】： Horizontal rotating drum; Avalanche; Simulation; Discrete element method

0 ?引言

顆粒物質(zhì)廣泛存在于自然界，并與人類的日常生活以及生產(chǎn)活動(dòng)密切相關(guān)。一般認(rèn)為，單個(gè)顆粒的典型尺度在10^–6～10 m的范圍內(nèi)，其運(yùn)動(dòng)規(guī)律服從牛頓運(yùn)動(dòng)定律^[1]。顆粒系統(tǒng)在受到外力或內(nèi)部應(yīng)力狀況變化時(shí)發(fā)生運(yùn)動(dòng)，表現(xiàn)出流體的屬性，從而形成顆粒流。了解和研究顆粒流不僅對(duì)預(yù)防自然災(zāi)害，例如雪崩、山體滑坡和泥石流^[2-3]，而且對(duì)提高工業(yè)生產(chǎn)效率，包括制粒、研磨、干燥和混合過(guò)程^[4]都具有重要意義。

水平轉(zhuǎn)鼓由于其簡(jiǎn)單性和可控性，是研究顆粒流動(dòng)特性的理想工具^[5-6]。轉(zhuǎn)鼓在低速旋轉(zhuǎn)的情況下，轉(zhuǎn)鼓中的顆粒流會(huì)表現(xiàn)出一系列離散的周期性雪崩^[7]。不少學(xué)者針對(duì)轉(zhuǎn)鼓內(nèi)雪崩問(wèn)題做了大量實(shí)驗(yàn)研究，甚至借助正電子發(fā)射跟蹤（PEPT）和磁共振成像（MRI）。雖然這些研究大大提高了我們對(duì)顆粒流運(yùn)動(dòng)行為的理解，但由于條件有限，此類實(shí)驗(yàn)往往是基于宏觀觀察，研究?jī)?nèi)容不夠全面和深入。因此，在微觀層面上缺乏實(shí)驗(yàn)數(shù)據(jù)不但意味著對(duì)實(shí)驗(yàn)結(jié)果的分析往往要憑經(jīng)驗(yàn)確定，例如每個(gè)顆粒的位置和速度等，而且實(shí)驗(yàn)效率低，耗費(fèi)時(shí)間長(zhǎng)。

近年來(lái)，由于計(jì)算機(jī)的高速發(fā)展，把計(jì)算機(jī)仿真技術(shù)^[8-12]運(yùn)用在分析顆粒流動(dòng)行為成為一種可行且高效的研究方法。由于顆粒流的行為取決于粒子間相互作用的結(jié)果，單個(gè)顆粒的信息將有助于揭示整個(gè)顆粒系統(tǒng)的基本原理。因此，英國(guó)皇家工程院院士Peter Cundall博士提出的直接對(duì)固體顆粒進(jìn)行跟蹤的離散元法（DEM）得到了關(guān)注。這種方法特別適用于模擬離散顆粒組合體在準(zhǔn)靜態(tài)或動(dòng)態(tài)條件下的運(yùn)動(dòng)及變形過(guò)程。運(yùn)用上述這種建模思想，借助計(jì)算機(jī)強(qiáng)大的計(jì)算能力^[13，14]，能極大地提高研究速度。把離散元法仿真和轉(zhuǎn)鼓雪崩相結(jié)合來(lái)研究具體有如下幾點(diǎn)優(yōu)勢(shì)：

（1）使我們能夠準(zhǔn)確而快速地模擬出真實(shí)世界的顆粒系統(tǒng)。類似于雪崩、泥石流等大型自然災(zāi)害如果通過(guò)實(shí)驗(yàn)的方法復(fù)現(xiàn)，費(fèi)用極其高昂，且存在諸多危險(xiǎn);對(duì)于生產(chǎn)環(huán)節(jié)，制備不同的顆粒材料和生產(chǎn)機(jī)器需要耗費(fèi)大量的人力物力，且受到各種復(fù)雜環(huán)境因素的影響，實(shí)驗(yàn)難度很大。通過(guò)仿真，可以高效地搭建實(shí)驗(yàn)環(huán)境，迅速開展研究工作。

（2）離散元仿真能直接獲得離散物質(zhì)大量復(fù)雜行為信息以及不易測(cè)得的顆粒尺度行為的信息，為分析顆粒流的運(yùn)動(dòng)、受力、熱量和能量傳遞等信息提供高級(jí)的解決途徑。甚至對(duì)于連續(xù)介質(zhì)理論無(wú)法解釋的物質(zhì)力學(xué)行為，也可以進(jìn)行準(zhǔn)確的預(yù)測(cè)和分析。

（3）對(duì)于實(shí)際的工作生產(chǎn)，通過(guò)計(jì)算機(jī)仿真的直觀感受，能拓展工程師設(shè)計(jì)思路，加速在設(shè)計(jì)階段找出設(shè)計(jì)的潛在問(wèn)題，使用數(shù)據(jù)精確突出設(shè)計(jì)方案的優(yōu)勢(shì)，加速產(chǎn)品創(chuàng)新。

綜上所示，運(yùn)用基于離散元法的仿真模型來(lái)研究雪崩相較于單純地實(shí)驗(yàn)觀測(cè)有著巨大的優(yōu)勢(shì)，為研究轉(zhuǎn)鼓內(nèi)的雪崩現(xiàn)象開辟了一條全新的研究方向。

因此，本文基于離散元法進(jìn)行雪崩效應(yīng)的仿真研究，以進(jìn)一步剖析雪崩顆粒流的運(yùn)動(dòng)特性。我們基于離散元法基本原理搭建迭代式計(jì)算流程框架?？紤]到轉(zhuǎn)鼓以及內(nèi)部顆粒的物理特性，運(yùn)用Hertz-Mindlin接觸模型和軟球離散元模型進(jìn)行仿真建模。并根據(jù)理論基礎(chǔ)對(duì)模型的材料屬性進(jìn)行參數(shù)優(yōu)化。通過(guò)仿真結(jié)果和實(shí)驗(yàn)數(shù)據(jù)在表面角和顆粒速度方面的對(duì)比，驗(yàn)證仿真模型的可靠性和有效性。我們還對(duì)參與雪崩的顆粒數(shù)量以及不同區(qū)域的顆粒速度進(jìn)行統(tǒng)計(jì)。為高效預(yù)測(cè)雪崩發(fā)展過(guò)程以及科學(xué)設(shè)計(jì)轉(zhuǎn)鼓內(nèi)的顆粒系統(tǒng)提供了一定的參考價(jià)值。

1 ?離散元法基本原理

離散元法自提出以來(lái)，在巖土工程和計(jì)算機(jī)仿真應(yīng)用領(lǐng)域發(fā)揮了其他數(shù)值算法不可替代的作用。其基本原理是把不連續(xù)體分離為剛性元素的集合，使得各個(gè)剛性元素滿足運(yùn)動(dòng)方程，用時(shí)步迭代的方法求解各剛性元素的運(yùn)動(dòng)方程，繼而求得不連續(xù)體的整體運(yùn)動(dòng)形態(tài)。該方法允許單元間的相對(duì)運(yùn)動(dòng)，不一定要滿足位移連續(xù)和變形協(xié)調(diào)條件，計(jì)算速度塊，所需存儲(chǔ)空間小，尤其適合求解大位移和非線性的問(wèn)題。

在建模過(guò)程中，把物料中的每個(gè)顆粒單獨(dú)作為一個(gè)粒子單元建立數(shù)學(xué)模型，并給定粒子單元的尺寸和物理性質(zhì)，如質(zhì)量、剛度和阻尼等，各個(gè)粒子之間存在接觸和分離這兩種關(guān)系。

在進(jìn)行數(shù)值計(jì)算時(shí)，把總的仿真計(jì)算時(shí)間劃分為一個(gè)個(gè)時(shí)間步長(zhǎng)。選取的時(shí)間步長(zhǎng)足夠小，使得在一個(gè)單獨(dú)的時(shí)間步長(zhǎng)內(nèi)，除了與選定單元直接接觸的單元外，來(lái)自其他任何單元的擾動(dòng)都不能傳播開來(lái)，即在任意時(shí)刻單元所受到的作用力只取決于該單元本身及與之直接接觸的其他單元。

當(dāng)接觸發(fā)生時(shí)，接觸點(diǎn)處就會(huì)產(chǎn)生接觸力和力矩，其大小可以根據(jù)接觸力學(xué)模型求出。利用牛頓第二定律建立每個(gè)粒子的運(yùn)動(dòng)方程，并用中心差分法求解，整個(gè)介質(zhì)的變形和演化由各個(gè)單元的運(yùn)動(dòng)和相互位置來(lái)描述。通過(guò)相鄰顆粒之間的碰撞產(chǎn)生的接觸力、力矩和不平衡力，計(jì)算每個(gè)顆粒在特定時(shí)刻的運(yùn)動(dòng)特征（速度、加速度等）。每一個(gè)時(shí)間步長(zhǎng)進(jìn)行一次迭代，并根據(jù)前一次計(jì)算所得到的單元的位置，作為下一次迭代的出發(fā)點(diǎn)，用以求出單元的新位置。通過(guò)遍歷跟蹤計(jì)算每個(gè)單元的微觀運(yùn)動(dòng)，即能獲得整個(gè)顆粒系統(tǒng)的宏觀運(yùn)動(dòng)規(guī)律。其計(jì)算流程圖如圖1所示^[15]。

2 ?三維模型搭建

2.1? 接觸模型

接觸模型是離散單元法模擬散體物料運(yùn)動(dòng)的重要基礎(chǔ)，其本質(zhì)是在準(zhǔn)靜態(tài)下分析散體顆粒接觸時(shí)產(chǎn)生的接觸力和彈塑性變形。接觸模型中的分析計(jì)算對(duì)顆粒所受的力和力矩的大小起決定性的作用。在本研究中，為了最準(zhǔn)確地仿真出現(xiàn)實(shí)世界的轉(zhuǎn)鼓內(nèi)雪崩現(xiàn)象，使用基于Hertz-Mindlin接觸模型的軟球離散元模型來(lái)進(jìn)行研究。

2.2? 模型搭建

仿真過(guò)程分為三個(gè)階段：轉(zhuǎn)鼓內(nèi)填充顆粒、顆粒自然沉降和轉(zhuǎn)鼓運(yùn)動(dòng)。在第一階段，生成轉(zhuǎn)鼓和顆粒。設(shè)置顆粒和幾何模型的物理參數(shù)是利用離散單元法模擬顆粒運(yùn)動(dòng)的重要環(huán)節(jié)之一，參數(shù)設(shè)置的正確與否對(duì)仿真結(jié)果產(chǎn)生直接影響。本研究假定轉(zhuǎn)鼓和顆粒為同種玻璃材料。鼓沿x軸方向施加周期性邊界條件。如圖2所示。這種處理可以在生成足夠的數(shù)據(jù)進(jìn)行分析的同時(shí)，使用較少的顆粒以縮短仿真時(shí)間。之后，在生成的鼓內(nèi)填充52000個(gè)半徑為1 mm的顆粒。填充方法為在鼓中隨機(jī)位置判斷有無(wú)填充空間，如果有則生成顆粒。在第二階段，懸浮在轉(zhuǎn)鼓中的顆粒在重力g= 9.8 m/s²的作用下完成自由落體，填充轉(zhuǎn)鼓下半部分50%的空間。在第三階段，在顆粒沉降穩(wěn)定后，轉(zhuǎn)鼓開始以0.19 rpm的速度旋轉(zhuǎn)，帶動(dòng)顆粒床運(yùn)動(dòng)。模擬的時(shí)間步長(zhǎng)設(shè)置為3×10^–7秒?？傔\(yùn)行時(shí)間為300秒，模擬計(jì)算時(shí)間約為100小時(shí)。

3? 仿真結(jié)果

3.1? 仿真準(zhǔn)確性驗(yàn)證

為了驗(yàn)證離散元法仿真的準(zhǔn)確性，采用空間濾波測(cè)速法（SFV）獲得實(shí)驗(yàn)數(shù)據(jù)^[16]。實(shí)驗(yàn)數(shù)據(jù)和仿真結(jié)果在相同配置的裝置下進(jìn)行。對(duì)于此轉(zhuǎn)速，實(shí)驗(yàn)和仿真中都顯示了離散的雪崩。

首先，我們?cè)诜抡婧蛯?shí)驗(yàn)中研究表面角θ的變化。圖3顯示了θ隨時(shí)間t的變化（取50個(gè)連續(xù)雪崩的平均值），其中紅線表示仿真結(jié)果，藍(lán)點(diǎn)表示實(shí)驗(yàn)數(shù)據(jù)。兩條曲線都顯示出相同的變化趨勢(shì)，兩者之間的最大差異也在1°以內(nèi)。因此，仿真結(jié)果與實(shí)驗(yàn)數(shù)據(jù)吻合良好。

顆粒速度是雪崩的另一個(gè)重要特征。圖4顯示了在雪崩期間顆粒床表面的顆粒平均速度v沿y軸的分布（取50個(gè)連續(xù)雪崩的平均值）。仿真曲線與實(shí)驗(yàn)曲線一致，它們之間的平均差為0.78 mm/s，并且都在表面徑向中心出現(xiàn)最大顆粒速度。此外，高速區(qū)的標(biāo)準(zhǔn)差值大于低速區(qū)的標(biāo)準(zhǔn)差。Lin等人也在在較高轉(zhuǎn)速的滾動(dòng)狀態(tài)下觀察到了相似的結(jié)果^[17]。

基于上述兩種不同的全局屬性的對(duì)比，我們驗(yàn)證了所提出的仿真系統(tǒng)可以準(zhǔn)確地模擬旋轉(zhuǎn)鼓中的雪崩現(xiàn)象。因此，離散元法仿真模型可以提供一種有效的方法來(lái)研究雪崩過(guò)程中的顆粒流行為。

3.2? 參與雪崩的顆粒數(shù)量研究

參與雪崩的顆粒數(shù)量是表明雪崩規(guī)模的關(guān)鍵指標(biāo)。我們觀察到部分顆粒在滾動(dòng)一定距離后會(huì)停止運(yùn)動(dòng)。因此，我們使用通過(guò)跟蹤每個(gè)顆粒的運(yùn)動(dòng)軌跡來(lái)判斷顆粒是否被卷入雪崩。雪崩開始時(shí)顆粒所處位置被標(biāo)記為起點(diǎn)。在雪崩期間，當(dāng)其位置距離起點(diǎn)大于一個(gè)粒徑時(shí)被認(rèn)為參與了雪崩。圖5顯示了一個(gè)粒徑厚度的表面層在雪崩過(guò)程中參與雪崩的顆粒數(shù)量N隨時(shí)間t的變化（取50個(gè)連續(xù)雪崩前5秒的平均值）。顆粒的數(shù)量在持續(xù)增長(zhǎng)，這意味著雪崩的規(guī)模在持續(xù)擴(kuò)大?；诓煌瑫r(shí)段的不同增長(zhǎng)速率，我們將它分為3個(gè)時(shí)段，分別用T1-T3標(biāo)記。在T1（0-0.2 s）期間，顆粒數(shù)呈現(xiàn)出加速增長(zhǎng)的趨勢(shì)。在T2（0.2-0.3 s）期間，顆粒數(shù)以均速增加（在圖5中用實(shí)心黑線標(biāo)記）。在T3（0.3-0.5 s）期間，增長(zhǎng)速度減慢并最終變成0，顆粒數(shù)量增加到最大。

3.3? 參與雪崩的顆粒速度研究

床表面不同位置的參與雪崩的顆粒速度是另一項(xiàng)重要參數(shù)。我們把顆粒床一個(gè)粒徑厚度的表面層按y軸方向等分成5個(gè)區(qū)域，按所處y軸由小到大的位置標(biāo)號(hào)為1-5號(hào)區(qū)域。圖5顯示了這5個(gè)區(qū)域內(nèi)的顆粒平均速度v隨時(shí)間t的變化規(guī)律。5個(gè)區(qū)域的顆粒平均速度都呈現(xiàn)先增大后減小的規(guī)律。把同一時(shí)刻不同區(qū)域的速度進(jìn)行對(duì)比可以發(fā)現(xiàn)，相對(duì)處于中心位置的2-4號(hào)區(qū)域的顆粒平均速度更大，這意味著此區(qū)域內(nèi)的顆粒運(yùn)動(dòng)最劇烈。而把同一區(qū)域在不同時(shí)刻的速度進(jìn)行對(duì)比可以發(fā)現(xiàn)，相對(duì)處于轉(zhuǎn)鼓更高位置的4-5號(hào)區(qū)域的最大速度出現(xiàn)的時(shí)刻更晚。

4? 結(jié)語(yǔ)

雪崩嚴(yán)重影響著人們的生命安危，且傳統(tǒng)的實(shí)驗(yàn)的宏觀分析十分低效，因此通過(guò)計(jì)算機(jī)仿真的方法來(lái)研究雪崩勢(shì)在必行。本文提出使用基于離散元法的仿真來(lái)研究轉(zhuǎn)鼓中雪崩狀態(tài)下的顆粒流動(dòng)行為。詳細(xì)介紹了所使用的離散元計(jì)算流程，顆粒接觸時(shí)采用的Hertz-Mindlin模型以及模型的生成過(guò)程。通過(guò)將模擬結(jié)果與實(shí)驗(yàn)數(shù)據(jù)（如表面角和表面顆粒速度）進(jìn)行比較，證明了離散元法仿真能實(shí)現(xiàn)轉(zhuǎn)鼓內(nèi)雪崩的有效仿真。此外，運(yùn)用離散元法能跟蹤每一個(gè)顆粒，對(duì)每個(gè)顆粒運(yùn)動(dòng)狀態(tài)進(jìn)行分析。結(jié)果表明雪崩規(guī)模的增長(zhǎng)速度在每個(gè)時(shí)間段中是不同的：在T1期間加速，在T2期間以均速增長(zhǎng)，最后在T3期間減速直至變成0。最后，我們還發(fā)現(xiàn)不同區(qū)域內(nèi)的顆粒平均速度隨時(shí)間都呈現(xiàn)先增大后減小的規(guī)律。這充分展示了離散元仿真技術(shù)在分析雪崩效應(yīng)的可靠性和優(yōu)越性，為雪崩的預(yù)測(cè)喝分析提供了有效的途徑。

參考文獻(xiàn)

1. 胡國(guó)明. 顆粒系統(tǒng)的離散元素法分析仿真： 離散元素法的工業(yè)應(yīng)用與EDEM軟件簡(jiǎn)介[M]. 武漢： 武漢理工大學(xué)出版社， 2010： 3-7.

2. H. Yang， B. F. Zhang， R. Li， G. Zheng， V. Zivkovic. Particle dynamics in avalanche flow of irregular sand particles in the

slumping regime of a rotating drum[J]， Powder Technology， 2017， 311（1）： 439-448.

1. E. Marteau， J. E. Andrade. A model for decoding the life cycle of granular avalanches in a rotating drum[J]. Acta Geotechnica， 2017， 13： 549-555.

2. J. Rajchenbach. Development of grain avalanches[J]. Phys Rev Lett， 2002， 89： 074301.

3. D. A. Santos， C. R. Duarte， M. A. S. Barrozo. Segregation phenomenon in a rotary drum： Experimental study and CFD simulation[J]. Powder Technology， 2016， 294： 1-10.

4. N. Sepúlveda， G. Krstulovic， S. Rica. Scaling laws in granular continuous avalanches in a rotating drum[J]. Physica A： Statistical Mechanics and its Applications， 2005， 356： 178-183.

5. R. Li， H. Yang， G. Zheng， Q. C. Sun. Granular avalanches in slumping regime in a 2D rotating drum[J]. Powder Technology， 2018， 326： 322-326.

6. 劉波， 郭平， 沈岳， 等. 農(nóng)業(yè)文本特征重組向量矩陣算法優(yōu)化與仿真[J]. 軟件， 2015， 36（4）： 6-13

7. 曹妍， 陳偉， 徐森. 圖像去噪方法研究與仿真[J]. 軟件， 2015， 36（4）： 33-36

8. 聶敬云， 李春青， 李威威， 等. 關(guān)于遺傳算法優(yōu)化的最小二乘支持向量機(jī)在MBR 仿真預(yù)測(cè)中的研究[J]. 軟件， 2015， 36（5）： 40-44.

9. 劉果， 陳凡， 李劍鋒， 等. 構(gòu)建SDN 仿真實(shí)驗(yàn)平臺(tái)的探討與實(shí)踐[J]. 軟件， 2015， 36（6）： 103-108.

10. 王燕妮， 吳文輝. Cloudsim 3. 0 仿真流程分析[J]. 軟件， 2015， 36（6）： 109-113.

11. S. H. Lin， H. Yang， R. Li， G. Zheng， V. Zivkovic， Velocities of irregular particles in a continuously avalanching surface flow within a rotating drum[J]. Powder Technology， 2018， 338： 376-382.

12. 何牧君， 郭力， 嚴(yán)歷. 大規(guī)模并行粒子模擬系統(tǒng)代碼級(jí)優(yōu)化研究與實(shí)現(xiàn)[J]. 計(jì)算機(jī)與應(yīng)用化學(xué)， 2008， 09： 65-68.

13. 王國(guó)強(qiáng)， 郝萬(wàn)軍， 王繼新. 離散單元法及其在EDEM上的實(shí)踐[M]. 西北工業(yè)大學(xué)出版社， 2010： 1-2.

14. 付宏， 烏蘭， 黃萬(wàn)風(fēng). 基于圖元的三維離散元法邊界建模方法[J]. 計(jì)算機(jī)集成制造系統(tǒng)， 2008， 014（12）： 2328-2333.

15. 徐偉， 李玉忱， 王麗. 分子動(dòng)力學(xué)并行算法研究[J]. 計(jì)算機(jī)工程與應(yīng)用， 2002， 14： 218-221.