離散元方法

1971年Cundall提出此方法時採用distinct element method是為了和連續介質力學中的finite element method相區別。後來用discrete element method取代了distinct element method，以反映系統是離散的這一本質特徵。

1971年Cundall提出適於岩石力學的離散元法， 1979年Cundall和Strack又提出適於土力學的離散元法，並推出二維圓盤(disc)程序BALL和三維圓球程序TRUBAL(后發展成商業軟體PFC-2D/3D)，形成較系統的模型與方法，被稱為軟顆粒模型；

離散元與分子動力學的比較

從本質上來講，離散元和分子動力學方法類似（molecular dynamics)，以至於有些作者在文獻中不加區別的使用MD和DEM兩個名字。然而離散元和分子動力學相似性只體現在形式上的相似（顆粒和牛頓定理）。二者還是有很大差別，在於分子動力學計算原子如何在給定相互作用勢下如何運動，而離散元計算的顆粒通常為微米及毫米量級。此外，離散元方法中需要考慮顆粒體在外力作用下的旋轉運動，顆粒的形狀，顆粒尺寸分佈，以及顆粒之間填充氣體，液體對顆粒材料宏觀性能都有很大的影響。總之，即使計算模擬一個最簡單的顆粒系統，單一尺寸的球形顆粒考慮摩擦作用下的運動問題都涉及到許多需要仔細考慮的細節，然而正如其他模擬方法一樣，這些細節往往不會被作者在文章中出版，大多靠自己在實踐中去不斷領悟。

第一步：建立所需要的幾何模型併產生顆粒

幾何模型可以根據實際計算模型需要建立，顆粒產生通常為隨機產生，及在給定的幾何空間內隨機產生所需要的顆粒。產生顆粒時需要實時監測新產生的顆粒和已有顆粒之間的位置關係，任意兩顆粒之間不能有重疊，否則顆粒之間相互作用力可能很大而導致系統崩潰。所以如果幾何模型尺寸，顆粒尺寸以及顆粒數目之間關係不合適，有可能導致顆粒產生失敗。顆粒的初始速度需要根據模擬需要而給定。

第二步：接觸探測

計算顆粒之間的相互距離，如果顆粒之間存在相互接觸（顆粒之間的距離小於兩個顆粒之間的半徑和），則要採用接觸模型計算相互作用力。

第三步：確定接觸模型

接觸模型是離散元計算的核心。所謂接觸模型就是確定顆粒接觸時的相互作用力。離散元計算中首先把相互作用力分解為法向力和切向力（法向指的是兩接觸顆粒中心之間的連線），所以接觸模型一般包含法向相互作用和切向相互作用。

在目前的離散元計算模擬中，一般來講，可以把所涉及到的接觸模型分為兩大類：

1）非結合性接觸力模型。此種接觸模型不考慮顆粒之間的相互吸引力，顆粒之間的相互作用以彈簧-粘壺模型來近似。彈簧代表顆粒之間的彈性相互作用，粘壺代表顆粒之間由於碰撞而引起的能量耗散。切向相互作用還要考利庫侖最大摩擦力約束。

2）結合性接觸力模型。微觀上來講，任何顆粒材料都由原子分子構成，然後原子分子之間存在van der Waals相互吸引力。當顆粒尺寸比較大時，van der Waals相互吸引力和顆粒自身重力相比較而言，對顆粒本身運動的影響可以忽略不計；然後，隨著顆粒尺寸減小，顆粒之間的相互吸引力變的和重力相當；當顆粒尺寸進一步減小時，顆粒之間相互吸引力要遠遠大於顆粒自身的重力。此時，就必須要考慮顆粒之間的相互吸引力。

目前來講，為了描述顆粒之間的相互引力作用力，有以下幾個模型可以選擇：

A. Johnson-Kendall-Roberts(JKR) 模型

JKR模型首先是基於能量平衡原理，從赫茲的非結合性接觸力模型而建立的。在JKR 模型中，顆粒之間的相互吸引力用表面能代表。並且顆粒之間的相互吸引作用只存在於接觸面積之類，不考慮顆粒之間其他部分的吸引相互作用。因此，JKR模型對大的，比較軟的顆粒材料描述比較好。

B, Derjaguin-Muller-Toporov(DMT)模型

DMT模型建設兩個顆粒之間的接觸變形和赫茲模型相同，但在接觸面積之外，存在一個額外的相互吸引作用力。DMT模型對小的比較硬的顆粒材料描述比較好。

早期，當JKR和DMT模型被提出時，人們爭論哪一個模型是正確的，後來研究發現，兩個模型都正確，他們考慮的顆粒之間接觸的兩個極端情況，基於材料性質，應該有區別的加以選取。1977年，人們提出Tabor 係數作為選取依據，Tabor係數是表面能，等效模量等材料性質的函數。JKR模型對大的Tabor係數適用，而DMT模型則反之。

C. van der Walls模型

此模型是基於顆粒之間的van der Walls相互作用直接推導得出。van der Waals相互作用能可以知道，然後對距離求導就直接得到力的表達式。

第四步：考慮其他相互作用力

有時，在顆粒系統中，取決於外部條件，需要考慮其他類型的相互作用力。例如，如果環境比較潮濕，則需要考慮顆粒之間的液橋相互作用力；如果顆粒表面本身帶有電荷，則要引入顆粒之間的靜電相互作用；如果有外加磁場並且顆粒材料本身有磁矩，則要考慮外磁場引起的磁性力。所有這些類型的相互作用力，都可以嵌入在離散元模型之中。

第五步：考慮顆粒和邊界之間的相互作用

計算模型的邊界可以是全部周期性邊界條件（x,y,z方向，全部採用周期性邊界），也可以是部分周期性邊界。在非周期性邊界條件時，則要給出顆粒本身和邊界之間如何相互作用。接觸模型和顆粒之間相互作用的接觸模型類似。

第六步：計算總的受力、加速度

綜合顆粒之間的相互作用力以及其他需要考慮的特殊相互作用力、以及顆粒與邊界之間的相互作用力，可以得到顆粒本身受的合外力，以此可以求得加速度

第七步：更新顆粒速度，坐標

根據加速度，更新顆粒速度，角速度，坐標等變數。坐標更新后，再進行接觸探測，整個計算流程進入下一個循環。

第八步：保存數據

第九步：分析處理

離散元方法的編碼思想十分簡單。集合中每一個單元都是獨立的，每個單元都具有相應的尺寸、質量、轉動慣量和接觸參數等屬性。它以牛頓第二定律和力一位移定理為基礎，對每一個單元首先確定與之接觸的單元，根據單元之間的重疊量，運用力一位移定理計算單元之間的接觸力，從而得到單元的合力和合力矩，之後用牛頓第二定律確定單元的運動規律，如此循環計算，直到系統中所有顆粒都計算完畢。