張中衛(wèi), 程 曉， 李 潮

(長(zhǎng)江大學(xué)城市建設(shè)學(xué)院，湖北荊州 434000)

?

PFC2D自定義模型使用研究

張中衛(wèi), 程 曉， 李 潮

(長(zhǎng)江大學(xué)城市建設(shè)學(xué)院，湖北荊州 434000)

論文研究了基于visual stdio 2013自定義接觸模型的實(shí)現(xiàn)方法，并嘗試將自定義接觸模型用在模擬中，分析內(nèi)置接觸模型與自定義接觸模型之間替代關(guān)系。得出結(jié)論：自定義接觸模型無(wú)法與常規(guī)的顆粒生成方法一起運(yùn)算，常規(guī)的顆粒生成方法指半徑擴(kuò)大法、explosive生成法；自定義接觸模型只能與固定坐標(biāo)生成的顆粒一同使用，巖土模擬中受限；自定義接觸模型在顆粒接觸剛度很小的情況下可作為內(nèi)置接觸模型的一種替代。

PFC2D； 接觸模型； 細(xì)觀參數(shù)； 數(shù)值模擬； c++

1 PFC2D自定義模型介紹

PFC2D(particle flow code for 2dimention)是一種對(duì)顆粒散體進(jìn)行模擬的有效數(shù)值手段。核心算法是離散元法，其中離散單元指的是每一顆小顆粒，顆粒和顆粒之間是按照粘結(jié)接觸聯(lián)系在一起的。一個(gè)完整的接觸模型包括接觸強(qiáng)度模塊、滑動(dòng)模塊和粘結(jié)模塊三部份。其中，接觸強(qiáng)度模塊用于給定顆粒和顆粒接觸的法向強(qiáng)度和切向強(qiáng)度，滑動(dòng)模塊是用來(lái)限定滑動(dòng)強(qiáng)度，粘結(jié)模塊用來(lái)規(guī)定顆粒和顆粒之間的粘結(jié)是點(diǎn)粘結(jié)還是平行粘結(jié)。現(xiàn)在大多數(shù)研究使用的都是PFC內(nèi)置的接觸模型。

PFC2D同時(shí)也支持用戶(hù)自己定義接觸模型，具體的做法是運(yùn)用c++將接觸模型相關(guān)的函數(shù)和參量編譯為dll文件，然后在PFC2D內(nèi)置命令流中調(diào)用dll文件。自定義接觸模型相對(duì)于PFC命令流來(lái)說(shuō)是相對(duì)獨(dú)立的，它有自己的名稱(chēng)、數(shù)據(jù)結(jié)構(gòu)、定義類(lèi)型。自定義模型與內(nèi)置接觸模型功用體現(xiàn)在被調(diào)用后可將用戶(hù)顆粒間的特征數(shù)據(jù)傳遞到顆粒粘結(jié)中并自己定義顆粒粘結(jié)破壞強(qiáng)度，以及粘結(jié)破壞后顆粒的軌跡。

2 自定義接觸模型國(guó)內(nèi)外研究現(xiàn)狀

PFC由于能從細(xì)觀方面研究宏觀的表現(xiàn)，從而被國(guó)內(nèi)外學(xué)者廣泛采用。但是幾乎沒(méi)有研究者用到自定義接觸模型。王剛[1]等利用PFC2D程序從宏觀和細(xì)觀角度探討節(jié)理在不同邊界條件下剪切過(guò)程中的力學(xué)演化規(guī)律和破壞機(jī)制；Axel Ruiken[2]等在指出前人關(guān)于由實(shí)驗(yàn)孔隙率推導(dǎo)2D情況下的孔隙率方法的不通用性的基礎(chǔ)上提出了自己由3D孔隙率得到2D孔隙率的方法；徐文杰[3]等提出了一種基于數(shù)字圖像處理的非均質(zhì)巖土材料細(xì)觀結(jié)構(gòu)PFC2D數(shù)值計(jì)算模型自動(dòng)生成方法。

莊德林[4]等首次采用了顆粒流程序PFC2D模擬垮落區(qū)充填中煤層采出、頂板垮落、充填的整個(gè)動(dòng)態(tài)過(guò)程；孔祥臣[5]等針對(duì)土石混合料振動(dòng)擊實(shí)過(guò)程中顆粒運(yùn)動(dòng)規(guī)律和密度形成機(jī)理問(wèn)題引入PFC2D顆粒離散元方法,模擬了一定級(jí)配條件下不同土石比的土石混合料的振動(dòng)擊實(shí)過(guò)程。

朱煥春[6]在PFC 及其在礦山崩落開(kāi)采研究中的應(yīng)用一文中對(duì)PFC2D進(jìn)行了功能介紹。其中提到PFC2D是幾乎用命令流操作的，所以上手難度很大，研究一般在高校中。另外也有提到“構(gòu)建 PFC 模型和進(jìn)行相應(yīng)的運(yùn)算準(zhǔn)備工作，必須使用 PFC 的二次開(kāi)發(fā)功能，可通過(guò)編程實(shí)現(xiàn)上述目的?！比欢诮鼉赡甑难芯恐?，PFC模擬都是使用內(nèi)置的接觸模型。論文對(duì)當(dāng)下學(xué)者對(duì)PFC研究中沒(méi)有使用自定義接觸模型的原因進(jìn)行研究。

3 構(gòu)建自定義接觸模型

自定義接觸是靠PFC2D主程序調(diào)用用戶(hù)編譯的dll文件來(lái)實(shí)現(xiàn)的。構(gòu)建自定義接觸的步驟有兩步：(1)取得dll文件關(guān)聯(lián)文件；(2)dll文件編譯。

3.1 取得dll文件關(guān)聯(lián)文件

編譯dll文件需要.cpp文件以及.h頭文件、外部關(guān)聯(lián)項(xiàng)。其中.cpp文件和.h頭文件可以從user manual中得到，直接嵌套運(yùn)行不會(huì)通過(guò)，需要進(jìn)行調(diào)試；外部關(guān)聯(lián)項(xiàng)是PFC2D程序自帶的cmodlewrv.dll和同名的object file library。

3.2 dll文件編譯

編譯dll文件的平臺(tái)需要為c++6.0以上版本，以vs2013為例，創(chuàng)建win32 dll工作項(xiàng)目，關(guān)聯(lián)以上提到的關(guān)聯(lián)文件，調(diào)試報(bào)錯(cuò)為0后，下拉菜單生成解決方案生成dll文件。

4 不同顆粒生成法與自定義接觸模型的嵌套

PFC2D生成模型的方法有：直接用gen命令生成少量顆粒、半徑擴(kuò)大法、explosive顆粒生成法、坐標(biāo)生成法。但是不同生成方法采用不同的參數(shù)值，因此可能導(dǎo)致無(wú)法在試樣中使用自定義模型。接下來(lái)對(duì)上述三種顆粒生成法與自定義接觸模型的耦合性進(jìn)行比較研究。

顆粒流建模首先要確定微觀參數(shù)，通過(guò)試錯(cuò)法，論文中的微觀參數(shù)的選取如表1所示[7]。

表1 微觀參數(shù)的選取

4.1 半徑擴(kuò)大法

半徑擴(kuò)大法的思路是先生成規(guī)定數(shù)目的顆粒填充限定的范圍，此時(shí)顆粒的半徑大小是在用戶(hù)設(shè)定半徑范圍內(nèi)的一個(gè)隨機(jī)值，顆粒與顆粒之間有重疊，然后定義半徑擴(kuò)大系數(shù)，也叫乘子。隨著cycle的進(jìn)行，顆粒的半徑逐漸增加，最終形成相互點(diǎn)接觸的顆粒試樣。

按照自定義模型的插入步驟調(diào)用編寫(xiě)好的dll文件，檢查命令流可以看出需要賦予顆粒的性狀已經(jīng)完成，但是在加載階段cycle計(jì)算開(kāi)始時(shí)，系統(tǒng)報(bào)錯(cuò)：property 0 not found in usercm1 model。證明半徑擴(kuò)大法生成的試樣缺失自定義模型需要的參數(shù)。

4.2 explosive 顆粒生成法

explosive顆粒生成法的思路是先生成規(guī)定數(shù)目的顆粒，此時(shí)生成的顆粒相對(duì)于半徑擴(kuò)大法來(lái)說(shuō)更加不均勻，顆粒間的重合更加明顯。explosive生成法的特點(diǎn)就是在cycle運(yùn)算的開(kāi)始，互相重疊的顆粒能迅速形成點(diǎn)接觸，并且這種生成方法可以控制生成試樣的孔隙率。重復(fù)半徑擴(kuò)大法中的自定義模型嵌套，同樣報(bào)錯(cuò)參數(shù)缺失。

4.3 坐標(biāo)生成法

坐標(biāo)生成方法是通過(guò)指定顆粒坐標(biāo)然后使用顆粒生成命令直接生成的。填充一份10×20試樣一般需要2 000左右顆粒，如果每個(gè)顆粒都要用戶(hù)自己定義坐標(biāo)會(huì)很繁瑣，所以這種顆粒生成方法適用范圍很小。

4.3.1 模型方案的建立

模型通過(guò)定義顆粒坐標(biāo)的方式生成3個(gè)顆粒，并且用4個(gè)墻體將顆粒圍起來(lái)，然后在顆粒接觸部分分別導(dǎo)入內(nèi)置接觸模型和自定義接觸模型。微觀參數(shù)的選取按照表2選取[8]，生成模型見(jiàn)圖1。

表2 微觀參數(shù)的選取

圖1 坐標(biāo)法建模

4.3.2 加壓與記錄

加壓采用賦予上部墻體固定速度0.001 m/s的方式，記錄顆粒間的平均非均衡力和平均接觸應(yīng)力，并繪出圖像(圖2)。

圖2 顆粒間應(yīng)力(kn=ks=1)

5 自定義接觸模型中參數(shù)影響

從圖2得出，內(nèi)置接觸模型和自定義接觸模型的均不平衡力-時(shí)步線(xiàn)和均接觸應(yīng)力-時(shí)步線(xiàn)近乎一致。當(dāng)kn=ks=1時(shí)，自定義接觸模型可以取代內(nèi)置接觸模型。

接下來(lái)參數(shù)進(jìn)行修改，來(lái)確定顆粒連接參數(shù)的選取是否會(huì)影響自定義接觸模型對(duì)內(nèi)置接觸模型的取代。

5.1 顆粒法向剛度和切向剛度的影響

保持摩擦系數(shù)為0.1的情況下，保持法向剛度和切向剛度相等，并分別取為0.3 Pa、0.5 Pa、1 Pa、3 Pa,繪出均不平衡力-時(shí)步曲線(xiàn)和均接觸應(yīng)力-時(shí)步曲線(xiàn)(圖3～圖5)。

圖3 顆粒間應(yīng)力(kn=ks=0.3 Pa)

圖4 顆粒間應(yīng)力(kn=ks=0.5 Pa)

圖5 顆粒間應(yīng)力(kn=ks=3 Pa)

從圖2～圖5可以得到以下結(jié)論：顆粒間的連接剛度kn與ks相等，并且維持在一個(gè)較小的數(shù)值范圍內(nèi)(0 Pa，1 Pa)，此時(shí)可以認(rèn)為自定義接觸模型可以完全代替內(nèi)置接觸模型。當(dāng)ks=ks≥1 Pa時(shí)，內(nèi)置接觸模型和自定義接觸模型的顆粒間應(yīng)力差距較大，此時(shí)不能用自定義接觸模型代替內(nèi)置接觸模型。

另外，數(shù)值模擬表明當(dāng)kn=ks≥10 Pa時(shí)，顆粒間的應(yīng)力大小會(huì)差5個(gè)數(shù)量級(jí)，所以未將圖像未放在一起比較。

5.2 顆粒間的摩擦系數(shù)的影響

接下來(lái)通過(guò)設(shè)定摩擦系數(shù)為0.1、1.2，控制kn=ks=1 Pa來(lái)確定顆粒間的摩擦系數(shù)是否會(huì)影響自定義接觸模型對(duì)內(nèi)置接觸模型的代替(圖6、圖7)。

圖6 顆粒間應(yīng)力(fric=0.1)

圖7 顆粒間應(yīng)力(fric=1.2)

圖6和圖7進(jìn)行對(duì)比，摩擦系數(shù)從0.1上升成1.2后，接觸應(yīng)力-時(shí)步曲線(xiàn)的峰值和線(xiàn)型都發(fā)生了明顯變化，證明摩擦系數(shù)的變動(dòng)也影響自定義接觸模型對(duì)內(nèi)置接觸模型的替代。

6 結(jié) 論

在PFC2D用戶(hù)手冊(cè)中，自定義接觸模型作為二次開(kāi)發(fā)手段被提出來(lái)，它的建立需要用戶(hù)掌握一定的c++語(yǔ)言，給PFC2D的建模增加了難度，并且經(jīng)過(guò)數(shù)據(jù)分析，自定義接觸取代內(nèi)置接觸模型的條件也較苛刻，所以在近幾年的PFC2D的研究中沒(méi)有關(guān)于自定義接觸模型的使用和論述。另外通過(guò)參數(shù)控制得出以下結(jié)論：

(1)使用自定義接觸模型替代內(nèi)置接觸模型的條件是：顆粒位置是通過(guò)坐標(biāo)定義生成；顆粒間的連接剛度kn與ks相等，并且維持在一個(gè)較小的數(shù)值范圍內(nèi)(0 Pa，1 Pa)。

(2)自定義接觸模型無(wú)法與半徑擴(kuò)大法生成的顆粒一起運(yùn)算也無(wú)法與explosive法生成的顆粒一同運(yùn)算。

(3)顆粒間摩擦系數(shù)的變化影響自定義接觸模型對(duì)內(nèi)置接觸模型的替代。

[1] 王剛，張學(xué)朋，賈宇靜.基于顆粒離散元法的巖石節(jié)理面剪切破壞細(xì)觀機(jī)理[J].中南大學(xué)學(xué)報(bào)， 2015， 46(4)：1-3.

[2] Zhijie Wang,Axel Ruiken ,Felix Jacobs and Martin Ziegler. A new suggestion for determining 2D porosities in DEM studies[J]. Geomechanics and Engineering, 2014,7(6): 665-678.

[3] 徐文杰,胡瑞林,王艷萍. 基于數(shù)字圖像的非均質(zhì)巖土材料細(xì)觀結(jié)構(gòu)PFC2D模型[J].煤炭學(xué)報(bào)， 2007， 32(4)：1-3.

[4] 莊德林,李興尚,許家林. 垮落區(qū)注漿充填壓實(shí)特性的PFC2D模擬試驗(yàn)[J].采礦與安全工程學(xué)報(bào)， 2008， 25(1)：1-3.

[5] 孔祥臣, 陳謙應(yīng),賈學(xué)明. 土石混合料振動(dòng)擊實(shí)試驗(yàn)的PFC2D模擬研究 [J].重慶交通學(xué)院學(xué)報(bào), 2005, 24(1)：1-4.

[6] 朱煥春. PFC及其在礦山崩落開(kāi)采研究中的應(yīng)用[J].巖石力學(xué)與工程學(xué)報(bào)， 2006， 25(9)：1-3.

[7] 彭博，陳玉明，袁利偉. 基于PFC2D馬鞍山鐵礦2號(hào)露天采場(chǎng)邊坡穩(wěn)定性分析[J].采礦工程， 2014， 35(12)：2-3.

[8] Cundall P A. PFC2D User's Manual(version3.0)[Z].Minnesota:Itasca Consulting Group Inc.

TU441+.3

A

[定稿日期]2016-04-17

[作者信息]張中衛(wèi)(1992～)， 男， 碩士研究生，研究方向?yàn)閹r土數(shù)值模擬。