1. 離散單元法誕生的背景

巖土，經受長期的地質構造作用，在一定的地質環境中形成一定的結構，顯現出寬廣和多變的材料響應范圍。與一般的工程材料相比，它顯現出結構上的不連續性、不均勻性和各向異性，且在物理力學性質上存在非線性。巖土材料的這些特性促使了許多數值模擬方法的發展以研究它的力學行為，如有限差分法、有限單元法和離散單元法。能夠模擬連續和非連續材料各力學行為的數值模擬工具已成為了研究者們追求的目標。

在巖土工程的早期研究階段，太沙基、比奧等先賢們讓碎散的土擁有了和其他連續介質一樣的“方程”，使得連續介質的理論也能夠為其所用。隨著計算機技術的發展以及科學研究的深入，基于連續介質理論的計算方法難以滿足研究者們對于計算精度的要求。受到分子動力學原理的啟發，在20世紀70年代，Cundall P. A. 教授[1]首次提出了一種顆粒離散體材料的分析方法，即離散單元法（Discrete Element Method），并將其應用于巖石塊體力學問題的分析。為了研究顆粒尺度上顆粒集合體的力學特性，1979年Cundall和Strack[2]又提出了適用于土力學的離散單元法。與常規有限單元法不同的是，離散單元法允許單元間的相對運動，不一定滿足位移連續和變形協調條件，計算速度快且所需的存儲空間較小，特別適用于巖土體材料的大變形/位移的分析。在隨后的幾十年中，離散單元法的應用領域不斷拓寬，逐漸被應用于散狀物料、粉體工程等領域。

2. 離散單元法的基本原理

世界上所有的物體都是由原子組成的，原子之間相互作用，進而構成分子、實體，并在外界作用下發生運動。理論上，如果知道了每一個原子的運動狀態，那么由這些原子所構成的實體的運動狀態便是確定的。離散單元法的原理與之類似，其最核心的思想便是通過大量的顆粒單元來模擬實際的研究對象，通過求解每一個顆粒的運動狀態來反映實體結構或者微觀結構的力學行為。

離散單元法的一般求解過程為：

將求解空間離散為離散元單元陣,并根據實際問題采用連接元件（即接觸模型）將相鄰單元連接起來；單元間相對位移是基本變量，由力與相對位移的關系可得到單元間法向和切向的作用力；對單元在各個方向上與其它單元間的作用力以及其它物理場對單元作用所引起的外力求合力和合力矩，根據牛頓運動第二定律求得單元的加速度；對其進行時間積分，得到單元的速度和位移。從而得到所有單元在任意時刻的速度、加速度、角速度、線位移和轉角等物理量。

在離散單元法中，接觸模型用來計算接觸力，進而計算顆粒的運動信息，是離散元法的理論核心。Cundall等最先提出的是簡單的彈簧-阻尼器接觸模型，如圖1 (a) 和 (b) 所示。圖中，kn、ks分別為法向和切向剛度，dn和ds分別為法向和切向阻尼。

圖1 接觸模型

上述接觸模型未考慮接觸上的顆粒滾動效應，顆粒容易發生轉動，導致數值模擬結果與實際情況有較大出入。Iwashita和Oda[3,4]引入接觸力矩，提出了抗滾動接觸模型，見圖1 (c)。圖中，kr和dr分別為滾轉剛度和滾轉阻尼。

3. 離散單元法的應用領域

隨著離散單元法理論的完善，該方法逐漸被人們所熟悉并應用于各個科學領域。以下是一些運用離散元方法求解科學問題的成功案例：

隧道滲漏

地質滑坡（http://matdem.com/content/?721.html）

地震波傳播（http://matdem.com/content/?698.html）

單軸壓縮

4. 離散單元法的局限性

雖然離散單元法能夠有效模擬巖土體的非連續性、不均勻性以及大變形破壞，在地質、巖土工程和能源開采等領域有非常廣泛的應用價值，但是我們也要意識到，離散元法在應用于實際工程問題時也面臨著許多困難：

（1）離散元的計算量巨大。以常用的離散元軟件PFC為例，其模擬的對象通常是細微觀的物質，所以它所建立的模型尺寸非常有限，而且在數值模擬的時候通常需要迭代計算至少一百萬次，迭代計算量非常大。

（2）定量建模困難。離散元模擬過程當中，通常都是賦予顆粒接觸參數和接觸模型來進行計算，使得其表現出的宏觀力學性質和實際材料的宏觀力學性質相匹配，而這種堆積模型的宏觀力學性質和單元力學參數間的關系是不明確的，需要進行大量的調試實現兩者的匹配。

（3）多場耦合理論尚未完善。現代工程中常常面臨多場和多相耦合的問題，如流固耦合、熱固耦合，而這些問題在離散單元法中實現的理論還不夠完善，一定程度上限制了離散單元法應用的范圍。

當然，針對上述問題，廣大科研工作者們也在離散單元方法的基礎上，開發出了一系列高效的數值計算軟件，如采用GPU進行并行計算求解的MatDEM大大提高了離散元問題的計算效率；又如開源的離散元程序LIGGGHTS，能夠引入OpenFoam來實現流固耦合問題的精確求解。

5. 結語

總體來說，作為一種非連續介質的數值計算方法，離散單元法從被創立至今一直承載著國內外科研工作者們的關注和重視，并在原有基礎上不斷完善和發展。隨著計算機技術的發展，基于離散單元法的數值模擬研究正不斷地向精細化、規模化和多尺度的方向發展，并將逐漸成為各工程領域不可缺少的數值計算方法之一。

參考文獻：

[1] Cundall P A. A computer model for simulating progressive, large-scale movements in blocky rock systems[C]. In: Proceedings of the International Symposium on Rock Mechanics. Nancy, France: [s.n.], 1971: 129–136.

[2] Cundall P A, Strack O D L. A discrete element model for granular assemblies[J]. Géotechnique, 1979, 29(1): 47–65.

[3] Iwashita K, Oda M. Rolling resistance at contacts in simulation of shear band development by DEM[J]. Journal of engineering mechanics, 1998, 124(3): 285-292.

[4] Iwashita K, Oda M. Micro-deformation mechanism of shear banding process based on modified distinct element[J]. Powder Technol, 2000, 109(1-3):192–205.

作者簡介：

skyhunter（筆名），同濟大學碩士在讀。

專注于巖土材料的數值模擬研究。

-End