邢志中 郭衛(wèi) 張武剛 王淵 路正雄 李實軍

摘要：為分析前后立式粉碎機粉碎內(nèi)腔所受沖擊力的情況及其影響因素，對改進前后立式粉碎機、中黃39大豆進行實體測量，測定其物理特性與力學(xué)特性，建立其EDEM模型并進行模擬仿真，得知：改進前后粉碎機粉碎下蓋內(nèi)腔所受沖擊力均大于各自上蓋內(nèi)腔，改進后粉碎機下蓋內(nèi)腔所受沖擊力大于改進前，改進前后粉碎機上蓋內(nèi)腔所受沖擊力大小無明顯區(qū)別。以1. 01 -1. 20 s改進前后粉碎機粉碎下蓋內(nèi)腔所受沖擊力為評價指標，研究不同大豆數(shù)目以及不同粉碎刀具的轉(zhuǎn)速對下蓋內(nèi)腔所受沖擊力的影響。結(jié)果表明：隨著大豆顆粒從100粒增至800粒，改進前后粉碎機下蓋內(nèi)腔所受沖擊力增大值出現(xiàn)在大豆數(shù)目為600左右。隨著粉碎刀具轉(zhuǎn)速從6200增至7 000 rpm改進前后粉碎機下蓋內(nèi)腔所受沖擊力逐漸增大。且在相同轉(zhuǎn)速以及相同大豆顆數(shù)的情況下，改進后粉碎機下蓋內(nèi)腔所受沖擊力均大于改進前。改進前后粉碎機下蓋內(nèi)腔所受沖擊力最大值出現(xiàn)大豆數(shù)目為600，刀具轉(zhuǎn)速為7 000 rpm，其值分別為64. 97與69. 56 N，改進前后粉碎機下蓋內(nèi)腔所受沖擊力最小值出現(xiàn)大豆數(shù)目為400，刀具轉(zhuǎn)速為6 200 rpm，其值分別為39. 35與41. 27 N.試驗與軟件仿真結(jié)果所得結(jié)論一致，該結(jié)果為研究立式粉碎機的粉碎內(nèi)腔沖擊力提供理論參考與依據(jù)。

關(guān)鍵詞：離散元;中黃39大豆;立式粉碎機;沖擊力;粉碎內(nèi)腔

中圖分類號：TD 402

文獻標志碼：A

文章編號：1672 -9315（2019）06 -1041 -08

DOI：10. 13800/j.cnki.xakjdxxb. 2019. 0617

開放科學(xué)（資源服務(wù)）標識碼（OSID）：

收稿日期：2019 -04 -18

責任編輯：高佳

基金項目：國家重點研發(fā)計劃子課題（ 2017 YFC0804310）

第一作者：邢志中（1991 -]，男，山西運城人，博士研究生，E-mail：xustxingzz@ 163.com

通信作者：郭衛(wèi)（1955 -），男，陜西黃陵人，教授，博士生導(dǎo)師，E-mail：guow@ xust.edu.cn

0 引言

粉碎機的更新推動著各行各業(yè)的不斷進步，礦業(yè)、食品等行業(yè)的發(fā)展無一不需要粉碎機[1]。而粉碎機的磨損已經(jīng)是粉碎設(shè)備行業(yè)所常見的現(xiàn)象，同時它也是各制造商及學(xué)者非常重視的問題。粉碎內(nèi)腔所受沖擊力為粉碎機磨損的一個重要因素，因此，針對粉碎內(nèi)腔所受沖擊力展開研究，分析其影響因素，對延長粉碎機的使用壽命，避免不必要的損害，具有十分重要的意義[2-3]。

離散元法是求解與分析復(fù)雜離散系統(tǒng)的運動規(guī)律與力學(xué)特性的一種新型數(shù)值方法，而EDEM作為全球首個現(xiàn)代化離散元模型模擬的軟件，已在制造業(yè)、礦山、食品行業(yè)等領(lǐng)域發(fā)揮了重要作用[4-5]。中國礦業(yè)大學(xué)的梁艷坤通過離散元法對神東礦區(qū)哈拉溝煤礦垮落帶中破碎巖體潰砂的速度進行了深度研究，最終得知松散砂在潰砂通道內(nèi)流動的速度分布呈現(xiàn)出紡錘形[6]。清華大學(xué)的李睿采用離散元法對球床的堆積結(jié)構(gòu)開產(chǎn)研究，從仿真與試驗驗證結(jié)果得知：在球床的大空間區(qū)域顆粒堆積具有隨機性，邊緣部位則呈現(xiàn)有序堆積[7]。中國農(nóng)業(yè)大學(xué)的彭飛應(yīng)用EDEM軟件研究喂料器的喂料過程，得到了主軸直徑、螺距、物料轉(zhuǎn)速與喂料穩(wěn)定性的關(guān)系，并確定了這3個因素的最佳參數(shù)組合[8]。本研究應(yīng)用離散元軟件EDEM研究立式粉碎機粉碎內(nèi)腔所受沖擊力的大小與分布情況。同時由筆者已發(fā)表學(xué)術(shù)論文以及試驗研究，得知改進后的圓錐上蓋粉碎機粉碎效果優(yōu)于圓柱上蓋粉碎機（改進前粉碎機）[9]。在此基礎(chǔ)上，分析改進前后粉碎機粉碎內(nèi)腔所受沖擊力的差別，大豆顆數(shù)、粉碎刀具轉(zhuǎn)速對改進前后粉碎機內(nèi)腔所受沖擊力的影響，所得結(jié)論為今后的研究與實際生產(chǎn)、延長其使用壽命提供一定的參考依據(jù)。

1 研究對象及粉碎機模型的建立

1.1 立式粉碎機及其模型

研究對象為小型立式粉碎機與中黃39大豆，利用EDEM軟件模擬分析大豆對小型立式粉碎機內(nèi)腔的沖擊力。研究所用的立式粉碎機型號為HC-700（圖1），改進前后粉碎機上蓋內(nèi)腔如圖2所示，對立式粉碎機幾何尺寸進行測量，依此采用UG軟件建立其三維模型。在EDEM軟件仿真過程中，只需建立顆粒所接觸的幾何體，因此對其模型進行簡化[1O -12]，改進前后粉碎機內(nèi)腔三維模型如圖3（a）、（b）所示。

1.2 顆粒模型

本研究以中黃39大豆為樣本（圖4）。隨機選取800粒，測量其長、寬、高[13]。計算出每顆大豆的等效直徑D與球形率φ，從而求出大豆長、寬、

將前文UG軟件建立的粉碎機三維模型導(dǎo)入到EDEM中。根據(jù)粉碎機的工作轉(zhuǎn)速，設(shè)置粉碎刀具轉(zhuǎn)速為6 200 rpm.添加顆粒工廠于1×10 s內(nèi)生成100顆中黃39大豆，其大小采用正態(tài)分布生成，參照表1設(shè)置中黃39大豆直徑平均值為3. 452 mm.標準差為0.087.

在求解器中設(shè)置仿真時間總長4s，1-3s為粉碎機運行時間，除此之外，設(shè)置數(shù)據(jù)寫出頻率為0.01s，網(wǎng)格邊長為最小大豆顆粒半徑的2倍[22 -25]。

2.2 粉碎機內(nèi)腔所受沖擊力分析

在EDEM后處理的coloring模塊中選擇部件為粉碎機內(nèi)腔，為便于觀察，將粉碎機內(nèi)所有顆粒隱藏，分析改進前后粉碎機內(nèi)腔所受顆粒沖擊力的大小及分布情況。圖4（a）、（b）、（c）、（d）分別為1. 02及1.12 s時，改進前后粉碎機內(nèi)腔所受顆粒沖擊力的情況。

從圖4（a）、（b）、（c）、（d）上下腔體所受沖擊力云圖可以看出，顏色越趨近于深灰色沖擊力越小，顏色越趨近于紅色沖擊力越大，在粉碎機工作的過程中，改進前后粉碎機內(nèi)腔所受顆粒的沖擊力多集中于下腔體。通過EDEM后處理中的數(shù)據(jù)導(dǎo)出功能，將1.01 -1. 20 s時間段改進前后粉碎機上下蓋內(nèi)腔所受沖擊力導(dǎo)出進行比較，重復(fù)9次仿真實驗，求出各時間點上下蓋內(nèi)腔所受沖擊力均值。其中A代表改進前的粉碎機，B代表改進后的粉碎機（圖5）。

從圖5可以看出，1. 01 -1.20 s時間段內(nèi)，改進前后粉碎機下蓋內(nèi)腔所受沖擊力分別大于各自上蓋內(nèi)腔所受沖擊力。其中，改進前粉碎機上蓋內(nèi)腔的最大受力為10. 32 N，最小受力為1.17 N，受力波動較平緩;下蓋內(nèi)腔的最大受力為34. 06 N，最小受力為18. 69 N，受力波動較大。改進后粉碎機上蓋內(nèi)腔的最大受力為11. 74 N，最小受力為1.16N，受力波動較平緩，下蓋內(nèi)腔的最大受力為43. 98N，最小受力為25. 02 N，受力波動較大。因此在改進前后粉碎機上下蓋設(shè)計的過程中，應(yīng)著重考慮下蓋的強度。除此之外得知，改進前后粉碎機上蓋內(nèi)腔所受沖擊力無明顯大小區(qū)別，而改進后粉碎機下蓋內(nèi)腔所受沖擊力除在1. 12與1.17 s時小于改進前粉碎機，其余各時間點改進后粉碎機下蓋所受沖擊力均大于改進前。進而可以得出雖然改進后粉碎機粉碎效果有所提高，但其下蓋內(nèi)腔所受沖擊力也明顯增大。

3 大豆數(shù)量、刀具轉(zhuǎn)速對粉碎內(nèi)腔所受顆粒沖擊力的影響

3.1 大豆數(shù)量對粉碎機內(nèi)腔所受沖擊力的影響

按照上訴仿真步驟，分別在粉碎機內(nèi)生成100，200，300，400，500，600，700，800粒中黃39大豆進行仿真，其他仿真條件保持不變，每種大豆數(shù)目均進行9次重復(fù)試驗。同樣通過EDEM軟件將1.01-1. 20 s時間段內(nèi)，改進前后粉碎機下蓋內(nèi)腔所受沖擊力導(dǎo)出，并求出每次仿真不同大豆數(shù)目在各時間點所對應(yīng)沖擊力數(shù)值之和的平均值，再對不同大豆數(shù)目所對應(yīng)的9次重復(fù)仿真結(jié)果求均值，繪制相關(guān)折線圖（圖6）。其中A代表改進前的粉碎機，B代表改進后的粉碎機。

從圖6可以看出，1. 01 -1.20 s時間段內(nèi)，粉碎機內(nèi)生成不同大豆數(shù)量進行粉碎，下蓋內(nèi)腔所受到?jīng)_擊力均不相等。改進前后粉碎機下蓋內(nèi)腔所受沖擊力最小值，均出現(xiàn)在大豆數(shù)量為100顆時，其值分別為33. 83與37. 61 N;改進前后粉碎機下蓋內(nèi)腔所受沖擊力最大值，均出現(xiàn)在大豆數(shù)量為600顆時，其值分別為74. 59與83. 21 N.改進前后粉碎機下蓋內(nèi)腔所受沖擊力與大豆顆粒數(shù)目有關(guān)，并且當大豆數(shù)量小于600時，隨著大豆顆數(shù)的增加，下蓋內(nèi)腔所受沖擊力逐漸增大。

3.2 刀具轉(zhuǎn)速對粉碎機內(nèi)腔所受沖擊力的影響

同樣按照首次仿真步驟，分別設(shè)置粉碎刀具轉(zhuǎn)速為：6 200，6 400，6 600，6 800，7 000 rpm進行仿真，其他仿真條件不變，每種轉(zhuǎn)速下均進行9次重復(fù)試驗。將1. 01 -1. 20 s時間段內(nèi)，改進前后粉碎機下蓋內(nèi)腔所受沖擊力導(dǎo)出，并求出每次仿真不同轉(zhuǎn)速在各時間點所對應(yīng)沖擊力數(shù)值之和的平均值，再對不同刀具轉(zhuǎn)速所對應(yīng)的9次重復(fù)仿真結(jié)果求出均值，繪制相關(guān)折線圖（圖7）。其中A代表改進前的粉碎機，B代表改進后的粉碎機。

從圖7可以看出，在1. 01 -1.20 s時間段內(nèi)，當粉碎刀具轉(zhuǎn)速為6 200 rpm時，改進前后粉碎機的下蓋內(nèi)腔所受沖擊力最小，分別為38. 36和39. 01N.當粉碎刀具轉(zhuǎn)速為7 000 rpm時，改進前后粉碎機的下蓋內(nèi)腔所受沖擊力最大，分別為47. 03和47. 13 N.隨著粉碎刀具轉(zhuǎn)速的增大，改進前后粉碎機下蓋內(nèi)腔所受沖擊力均明顯增大。

3.3 大豆顆數(shù)、刀具轉(zhuǎn)速對粉碎機內(nèi)腔所受沖擊力的綜合影響

分別在中黃39大豆顆數(shù)為400，600，800時，設(shè)置改進前后粉碎機粉碎刀具轉(zhuǎn)速為6 200，6 400，6 600，6 800，7 000 rpm，其余條件不變，每種情況9次重復(fù)仿真，共進行270次仿真。分別求出1. 01 -1. 20 s時間段內(nèi)，各時間點對應(yīng)下蓋所受沖擊力數(shù)值之和的平均值，繪制表格如圖8所示。

從圖8可以看出，在相同大豆數(shù)目時，改進前后粉碎機隨著粉碎刀具轉(zhuǎn)速增加，下蓋內(nèi)腔所受沖擊力同樣明顯增大。其中，改進前粉碎機下蓋內(nèi)腔所受沖擊力最大值出現(xiàn)在粉碎刀具為7 000 rpm，大豆數(shù)目為600時，其值為：64. 97 N;最小值出現(xiàn)在粉碎刀具為6 200 rpm，大豆數(shù)目為400時，其值為：39. 35 N.改進后粉碎機下蓋內(nèi)腔所受沖擊力最大值出現(xiàn)在粉碎刀具為7 000 rpm，大豆數(shù)目為600時，其值為：69. 56 N;最小值出現(xiàn)在粉碎刀具為6 200 rpm，大豆數(shù)目為400時，其值為：41. 27 N.除此之外，在相同的轉(zhuǎn)速以及相同大豆顆數(shù)的情況下，改進后粉碎機下蓋內(nèi)腔所受沖擊力均大于改進前。

4 試驗驗證

試驗所用大豆為山東祥豐種業(yè)有限公司所生產(chǎn)的中黃39大豆，變頻器為上海本川自動化科技有限公司所生產(chǎn)的XYS-AT2型變頻器（圖9）。采用BSF120-3AA-T常溫應(yīng)變片（南京聚航科技有限公司生產(chǎn)）布置在改進前后粉碎機上下蓋內(nèi)腔，并在粉碎機上下蓋進行打孔處理，應(yīng)變片的連接導(dǎo)線通過孔洞與外部ASMCl-9型電阻應(yīng)變儀連接（每通道最高采樣率為2 000 S/s，南京聚航科技有限公司生產(chǎn)），最終由計算機對所收集的數(shù)據(jù)進行處理（圖10）。在中黃39大豆顆數(shù)為400.600，800時，通過變頻器調(diào)至改進前后粉碎機粉碎刀具轉(zhuǎn)速為6200，6 400，6 600，6 800，7 000 rpm，其余條件不變，每種情況重復(fù)9次試驗。結(jié)果如圖11，圖12所示。

通過對比圖5與圖11，圖8與圖12.可以看出，在相同的轉(zhuǎn)速以及相同大豆顆數(shù)的情況下，改進前后粉碎機上下蓋內(nèi)腔所受沖擊力的仿真與試驗值有略微偏差。但粉碎機上下蓋內(nèi)腔受力試驗結(jié)果與仿真結(jié)果的總體改變趨勢一致。

5 結(jié)論

1）改進前后粉碎機下蓋內(nèi)腔所受沖擊力均大于各自上蓋內(nèi)腔所受沖擊力，即下蓋內(nèi)腔比上蓋內(nèi)腔的磨損更為嚴重，在粉碎機設(shè)計過程中，應(yīng)著重考慮下蓋的強度。改進后粉碎機下蓋內(nèi)腔所受沖擊力大于改進前粉碎機，改進后粉碎機上蓋內(nèi)腔所受沖擊力與改進前上蓋內(nèi)腔所受沖擊力無明顯大小差別。

2）在相同粉碎機內(nèi)，當中黃39大豆數(shù)目小于600時，隨著中黃39大豆數(shù)目的增加，下蓋內(nèi)腔所受沖擊力明顯增大。即中黃39大豆數(shù)目為600左右時，下蓋內(nèi)腔的磨損隨之越嚴重。

3）在相同粉碎機內(nèi)，隨著粉碎刀具轉(zhuǎn)速的增加，下蓋內(nèi)腔所受沖擊力明顯增大。即粉碎刀具轉(zhuǎn)速越高，下蓋內(nèi)腔磨損隨之越嚴重。

4）在相同的轉(zhuǎn)速以及相同大豆顆數(shù)的情況下，改進后粉碎機下蓋內(nèi)腔所受沖擊力均大于改進前。即雖然改進后粉碎機提高了粉碎效果，但其下蓋內(nèi)腔磨損程度比改進前粉碎機更為嚴重。

5） EDEM軟件仿真所得改進前后粉碎機內(nèi)腔受力結(jié)果與試驗驗證結(jié)果基本一致，表明了采用EDEM軟件分析粉碎機內(nèi)腔受力的可行性。

參考文獻（ References）：

[1] 陳樹召.大型露天煤礦他移式破碎站半連續(xù)工藝系統(tǒng)優(yōu)化與應(yīng)用研究[J].煤炭學(xué)報，2012，37（3）：533 - 534.

CHEN Shu-zhao. Study on optimization and application ofsemi-continuous process system of other moving crushingstation in large open pit coal mine[J].Journal of CoalMine，2012 ，37 （3）：533 - 534.

[2]李潤軍，單仁亮，李潤圣，等，盾構(gòu)機主驅(qū)動密封維修改造關(guān)鍵技術(shù)[J].西安科技大學(xué)學(xué)報，2014，34（5）：

579 -584.

LI Run-jun， SHAN Ren-Iiang， LI Run-sheng， et al.Keytechnologies of shield machine main drive seal repair andtransformation[J].Journal of Xi' an University of Sci-ence and Technology ，2014 ，34（5）：579 - 584.

[3]劉瑜，杜長龍，付林，等，煤塊沖擊破碎速度研究[J].振動與沖擊，2011，30（3）：18 -21.

LIU Yu， DU Chang-long， FU Lin，et al.Study on the im-pact crushing speed of coal blocks[J].Study on the im-pact crushing speed of coal blocks[J].Vibration and im-pact，2011 ，30（3）：18 - 21.

[4] 李守巨，李德，于申.壓頭作用下巖石破碎過程的細觀模擬[J].西安科技大學(xué)學(xué)報，2016，36（6）：769 - 774.

LI Shou-ju， LI De， YU Shen. Microscopic simulation ofrock fragmentation process under pressure head[ J].Journal of Xi' an University of Science and Technology，2016，36（6）：769 -774.

[5]孫學(xué)陽，安孝會，苗霖田，等.煤礦井工開采對上覆反向滑坡擾動的模擬研究[J].西安科技大學(xué)學(xué)報，2017 ，37 （1）：71 - 77.

SUN Xue-yang， AN Xiao-hui， MIAO Li-tian， et al.A sim-ulation study on the disturbance of overlying reverselandslide caused by coal mining[J].Journal of Xi' an U—niversitv of Science and Technology， 2017， 37（1）：71 - 77.

[6] 梁艷坤，隋旺華，朱濤，等.哈拉溝煤礦垮落帶破碎巖體潰砂的離散元數(shù)值模擬研究[J].煤炭學(xué)報，2017 ，42（2）：470 - 476.

LIANC Yan-kun， SUI Wang-hua， ZHU Tao， et al.Dis-crete element numerical simulation of sand burst in frac-tured rock mass of collapse zone in Haragou Coal Mine[J]. Journal of Coal Mine ，2017 ，42（2）：470 - 476.

[7]李睿，任成，楊星團，等.基于單元體理論的球床堆積結(jié)構(gòu)與傳熱算法研究[J].核動力工程，2016，37（3）：39 -42.

LI Rui， REN Cheng ，YANC Xing-tuan，et al.An analysisof structures and heat transfer for packed beds[J].Nu-clear Power Engineering， 2016 ，37（3）：39 - 42.

[8]彭飛，李騰飛，康宏彬，等.小型制粒機喂料器參數(shù)優(yōu)化與試驗[J].農(nóng)業(yè)機械學(xué)報，2016，47（2）：51 -58.

PENG Fei， LI Teng-fei， KANC Hong-bin， et al.Parameteroptimization and test of feeder for small granulator[J].Journal of Agricultural Machinery， 2016， 47（ 2）： 51 - 58.

[9]邢志中，李貴榮，張海東，等.基于EDEM的小型立式粉碎機數(shù)值模擬與優(yōu)化研究[J].食品工業(yè)，2017 （1）：206 - 209.

XING Zhi-zhong， LI Gui-rong， ZHANG Hai-dong， et al.Study on numerical simulation and optimization of smallvertical mill based on EDEM[J].Food IndustW， 2017（1）：206 - 209.

[10]毛君，劉歆妍，陳洪月，等.基于EDEM的采煤機滾筒工作性能的仿真研究[J].煤炭學(xué)報，2017，42（4）：1069 - 1077.

MAO Jun， LIU Xin-yan， CHEN Hong-yue，et al.Simula-tion study on working performance of shearer drum basedon EDEM[J].Journal of Coal Science， 2017， 42（ 4）：1069 - 1077.

[11]

Ouyang Y， Feng G，Read J J，et al.Estimating the ratio ofpond size to irrigated soybean land in Mississippi：a casestudy[J].Water Science and Technology： Water Supply，2016 ，16（6）：1639 - 1647.

[12] Thakur S C，Morrissey J P，Sun J， et al.Micromechanicalanalysis of cohesive granular materials using the discreteelement method with an adhesive elasto-plastic contactmodel[J].Cranular Matter ，2014， 16（3）：383 - 400.

[13] Woo S M， Uyeh D D， Sagong M S，et al.Development ofseeder for mixed planting of corn and soybeans[J].In-ternational Journal of Agricultural and Biological Engi-neering，2017 ，10（3）：95 - 101.

[14]陳海濤，李桐輝，王洪飛，等.氣吸滾筒式壟上三行大豆密植排種器設(shè)計與參數(shù)優(yōu)化[J].農(nóng)業(yè)工程學(xué)報，2018 ，34（ 17）：16 - 24.

CHEN Hai-tao， LI Tong-hui， WANG Hong-fei，et al.De-sign and parameter optimization of three rows of soybeanseed metering device on air-suction drum ridge[J].Jour-nal of Agricultwal Engineering ，2018 ，34（ 17）： 16 - 24.

[15]常杰，王菊，劉勇，等，低共熔離子液體催化大豆油環(huán)氧化制備增塑劑[J].華南理工大學(xué)學(xué)報（自然科學(xué)版），2018 ，46（2）：93 - 101.

CHANG Jie， WANG Ju， LIU Yong， et al.Preparation ofplasticizers by epoxidation of soybean oil catalyzed bylow eutectic ionic liquids[J].Joumal of South China U一niversity of Technology（ Natural Science Edition）， 2018，46（2）：93 - 101.

[16] Constantine D A，Wang Y， Terrell E J.Effect of recipro-cation frequency on friction and wear of vibrating con-tacts lubricated with soybean-based BlOO biodiesel[J].Tribology Letters ，2013 ，50（2）：279 - 285.

[17]

Lisjak A，Liu Q，Zhao Q，et al.Numerical simulation of a_coustic emission in brittle rocks by two-dimensional fi—nite-discrete element analysis[J].Ceophysical JournalIntemational ，2013， 195（1）：423 - 443.

[18] Mukhopadhyay G，Palit P，Bhattacharvya S. Developmentof AISI A2 tool steel beater head for an impact crusher ina sinter plant[J].Metallography， Microstructure， and A-nalysis，2015 ，4（2）：114 - 121.

[19]王德福，王沫，李利橋.錘片式粉碎機粉碎玉米秸稈機理分析與參數(shù)優(yōu)化[J].農(nóng)業(yè)機械學(xué)報，2017 ，48（ 11）：165 - 171.

WANG De-fu， WANG Mo， LI Li-qiao. Mechanism analy-sis and parameter optimization of com straw crushing byhammer mill[ J]. Jomnal of Agricultural Machinery，2017 ，48（ 11）：165 - 171.

[20]金磊，曾亞武，葉陽，等.不規(guī)則顆粒及其集合體三維離散元建模方法的改進[J].巖土I程學(xué)報，2017 ，39（7）：1273 - 1281.

JIN Lei， ZENG Ya-wu， YE Yang， et al.Improvement of3D discrete element modeling method for irregular parti-cles and their aggregates[J].Journal of Geotechnical En-gineering ，2017 ，39（7）：1273 - 1281.

[21]王揚，呂鳳妍，徐天月，等，大豆籽粒形狀和尺寸分析及其建模[J].吉林大學(xué)學(xué)報（工學(xué)版），2018，48（2）：507 -517.

WANG Yang， LV Feng-yan， XU Tian-yue， et al.Shapeand size analysis and modeling of soybean seeds[J].Journal of Jilin University（ Engineering Edition）， 2018，48（2）：507 - 517.

[22]孫逸飛，宋順翔，高玉峰.級配散粒體循環(huán)荷載下變形特性的離散元模擬[J].中國礦業(yè)大學(xué)學(xué)報，2018 ，47（4）：874 - 878.

SUN Yi-fei， SONC Shun-xiang， CAO Yu-feng. DEM sim-ulation of deformation characteristics of graded granularmaterials under cyclic loading[J].Journal of China Uni-versity of Mining and Technology， 2018， 47（ 4）： 874 -878.

[23]

Celigueta M A，Latorre S， Arrufat F， et al.Accurate mod-elling of the elastic behavior of a continuum with the Dis-crete Element Method[J].Computational Mechanics，2017 ，60（6）：1-14.

[24] Lee S H， Choi M， Kim T T，et al_Switching terahertz waves with gate-controlled active graphene metamaterials[J]. Nature Materials ，2012 ，11（ 11）：936 - 941.

[25]龍雪，劉社文，季順迎，水位變化對正倒錐體冰載荷影響的離散元分析[J].力學(xué)學(xué)報，2019 ，51（1）：74 - 84.

LONC Xue， LIU She-wen， JI Shun-ying. Discrete elementanalysis of the impact of water level change on ice load offorward and inverted cones[J]. Joumal of Mechanics，2019 ，51（1）：74 - 84.