李 爽，劉 洋，吳可嘉

(北京科技大學(xué) 土木與資源工程學(xué)院，北京 100083)

砂土直剪試驗(yàn)離散元數(shù)值模擬與細(xì)觀變形機(jī)理研究

李 爽，劉 洋，吳可嘉

(北京科技大學(xué) 土木與資源工程學(xué)院，北京 100083)

基于離散元法進(jìn)行砂土直剪試驗(yàn)?zāi)M與剪切宏細(xì)觀響應(yīng)研究。采用試樣分條帶法研究砂土剪切過(guò)程中的變形發(fā)育情況；記錄應(yīng)力-位移關(guān)系，分析砂土剪切應(yīng)力、應(yīng)變特性與剪脹特性。通過(guò)分析接觸力鏈、接觸組構(gòu)分布、配位數(shù)及滑動(dòng)接觸百分比等細(xì)觀行為演化研究顆粒細(xì)觀接觸狀態(tài)演變；基于局部孔隙率變化分析顆粒細(xì)觀分布狀態(tài)演變；采用顆粒累積轉(zhuǎn)動(dòng)量和速度場(chǎng)分析顆粒細(xì)觀運(yùn)動(dòng)狀態(tài)及剪切帶發(fā)育的細(xì)觀機(jī)制；對(duì)離散元軟件PFC2D進(jìn)行二次開(kāi)發(fā)，實(shí)現(xiàn)主應(yīng)力的可視化。研究表明：應(yīng)變局部化集中于顆粒轉(zhuǎn)動(dòng)、接觸滑動(dòng)、體脹明顯的剪切帶內(nèi)，剪切帶內(nèi)外顆粒不同的運(yùn)動(dòng)狀態(tài)導(dǎo)致細(xì)觀參量的變化，剪切帶厚度為(10～15)d50；力鏈網(wǎng)絡(luò)、接觸組構(gòu)和主應(yīng)力場(chǎng)本質(zhì)上是一致的，三者主方向在剪切過(guò)程中發(fā)生了近似一致的偏轉(zhuǎn)。

砂土；離散元法；直剪試驗(yàn)；細(xì)觀力學(xué)；剪切帶

1 研究背景

土作為離散顆粒聚集體，具有離散性、非各向同性、細(xì)觀隨機(jī)性等特殊性質(zhì)，造成其剪切機(jī)制的復(fù)雜性。直剪試驗(yàn)是一種簡(jiǎn)單有效的土強(qiáng)度室內(nèi)試驗(yàn)測(cè)定方法[1]，廣泛應(yīng)用于測(cè)定散粒介質(zhì)的抗剪強(qiáng)度和流動(dòng)性質(zhì)。土體試樣在直剪時(shí)受剪切作用會(huì)發(fā)生變形，伴隨顆粒的轉(zhuǎn)動(dòng)、滑移、破碎，當(dāng)變形逐漸發(fā)育并集中于某一區(qū)域內(nèi)時(shí)，試樣發(fā)生失穩(wěn)并破壞，通常將此局部變形區(qū)域稱為剪切帶。剪切帶的形成和演變導(dǎo)致土體的失穩(wěn)和破壞，是決定土體強(qiáng)度的薄弱區(qū)域[2-3]。

目前剪切帶理論的研究成果主要有：局部分叉理論、梯度塑性理論、非局部理論、Cosserat理論及復(fù)合體理論等。常規(guī)直剪儀雖能得到抗剪強(qiáng)度指標(biāo)，但不能分析試樣在剪切過(guò)程中的應(yīng)力-應(yīng)變關(guān)系和微細(xì)觀力學(xué)響應(yīng)，對(duì)剪切帶的研究也受到限制?；诖耍瑒⑽陌椎萚4]開(kāi)發(fā)了半膜式細(xì)觀結(jié)構(gòu)觀測(cè)直剪儀，研究了砂土試樣在剪切過(guò)程中的應(yīng)力應(yīng)變場(chǎng)演變，但研究中未從顆粒運(yùn)動(dòng)角度對(duì)相關(guān)微細(xì)觀力學(xué)行為進(jìn)行分析。由于物理試驗(yàn)對(duì)于微細(xì)觀力學(xué)研究來(lái)說(shuō)限制較多，且耗時(shí)費(fèi)力，數(shù)值試驗(yàn)手段可作為其有效補(bǔ)充。其中，Potts等[5]、Tejchman等[6]分別采用基于彈塑性模型和亞塑性模型的有限元法研究了砂土直剪力學(xué)性質(zhì)及剪切帶的應(yīng)變局部化演變，并分析了邊界效應(yīng)和顆粒直徑對(duì)剪切帶厚度的影響。但由于有限單元法無(wú)法針對(duì)微細(xì)觀(顆粒)層面剪切應(yīng)力與變形演變的機(jī)理進(jìn)行探究，而離散元法在土體細(xì)觀機(jī)理研究方面體現(xiàn)出獨(dú)特的優(yōu)勢(shì)。Thornton等[7]基于離散元法對(duì)散粒體直剪試驗(yàn)進(jìn)行了模擬，Liu等[8]基于離散元法研究了粒狀體的直剪力學(xué)特性，Masson等[9]基于離散元法對(duì)松密試樣的剪切力學(xué)響應(yīng)進(jìn)行了宏細(xì)觀分析與對(duì)比研究，Wang等[10]采用離散元法對(duì)散粒體直剪試驗(yàn)的尺寸效應(yīng)進(jìn)行了探究，蔣明鏡等[11]通過(guò)離散元法研究了單粒組密砂剪切帶的形成過(guò)程。因此，離散元法作為土微細(xì)觀力學(xué)研究的有效手段，是對(duì)物理試驗(yàn)研究的有益補(bǔ)充，離散元數(shù)值試驗(yàn)有重要的研究與指導(dǎo)意義。

本文利用二維顆粒流程序?qū)ι皹又奔粼囼?yàn)進(jìn)行數(shù)值模擬，并進(jìn)行細(xì)觀組構(gòu)統(tǒng)計(jì)開(kāi)發(fā)和主應(yīng)力可視化開(kāi)發(fā)，選定特征狀態(tài)對(duì)砂樣直剪過(guò)程中的受力及變形宏細(xì)特性進(jìn)行全面分析，揭示砂土剪切的本質(zhì)機(jī)理。著重研究以下問(wèn)題：①砂樣直剪中的變形發(fā)育情況及剪切帶形成過(guò)程；②從顆粒接觸狀態(tài)、顆粒分布狀態(tài)和顆粒運(yùn)動(dòng)狀態(tài)3方面分析直剪過(guò)程中砂樣的細(xì)觀力學(xué)特性演變；③主應(yīng)力可視化開(kāi)發(fā)與直剪過(guò)程中砂樣的主應(yīng)力矢量分布演變研究。

2 直剪試驗(yàn)離散元分析過(guò)程

砂土離散元試樣的尺寸為70 mm×30 mm，試樣的粒徑范圍為0.5～0.7 mm均勻分布，平均粒徑d50為0.6 mm。試樣初始孔隙率0.175，顆?？倲?shù)6 127個(gè)。離散元模型采用6道剛性墻體邊界模擬直剪試驗(yàn)的剪切盒，并生成2道水平剛性墻體防止剪切過(guò)程中的顆粒溢出。

粒間接觸模型采用接觸剛度模型，符合無(wú)黏性土的特性。試樣主要的細(xì)觀參數(shù)包括顆粒間及顆粒與墻體間法向接觸剛度、切向接觸剛度、摩擦系數(shù)，取值列于表1中。其中，顆粒間的接觸剛度是決定散粒試樣彈性模量的重要參數(shù)，根據(jù)散粒體彈性理論，當(dāng)顆粒接觸剛度的取值在108～109N/m之間，法向和切向剛度比值在1～3之間時(shí)，模型收斂性良好[12]。本次數(shù)值試驗(yàn)設(shè)定顆粒間及顆粒與墻體間的法向接觸剛度、切向接觸剛度分別為1.5×108，1.0×108N/m,二者比值為1.5，顆粒間摩擦系數(shù)設(shè)為0.5，顆粒與墻體間的摩擦系數(shù)設(shè)為0。

表1 離散元模型參數(shù)Table 1 Parameters of discrete element model

圖1 直剪試驗(yàn)數(shù)值模型及測(cè)量圈Fig.1 Model of direct shear simulation test and measurement circles

試樣生成后，通過(guò)伺服程序施加100 kPa 法向力，使試樣固結(jié)穩(wěn)定。隨后固定下盒，保持法向應(yīng)力為100 kPa的常應(yīng)力條件，并給上盒施加一個(gè)向左的速度實(shí)現(xiàn)上下盒的相對(duì)運(yùn)動(dòng)，以模擬砂土直剪試驗(yàn)。數(shù)值試驗(yàn)的示意圖如圖1(a)所示。

為了監(jiān)測(cè)試樣內(nèi)局部孔隙率、配位數(shù)、滑動(dòng)接觸百分比等細(xì)觀參量的變化，并實(shí)現(xiàn)主應(yīng)力偏轉(zhuǎn)可視化二次開(kāi)發(fā)，在試驗(yàn)前建立了如圖1(b)所示的監(jiān)測(cè)測(cè)量圈。

3 宏觀力學(xué)特性分析

3.1 宏觀應(yīng)力-應(yīng)變特性

土工直剪儀上下盒各自的側(cè)面邊界并不能發(fā)生相對(duì)運(yùn)動(dòng)，故直剪試驗(yàn)過(guò)程中試樣的體積只能沿豎直方向發(fā)生膨脹或收縮，試驗(yàn)時(shí)可記錄試樣的高度變化反映其剪脹性。同時(shí)設(shè)置監(jiān)測(cè)記錄剪切過(guò)程中的剪切應(yīng)力、剪切位移，得到砂土試樣剪切5 mm過(guò)程中的剪切應(yīng)力、剪切位移、體積變化三者之間的特征曲線如圖2所示。

圖2 剪切應(yīng)力-剪切位移及體積變化-剪切位移曲線Fig.2 Curves of shear stress vs. shear displacement and volume change vs. shear displacement

由圖2中的剪切應(yīng)力-剪切位移曲線可知，試樣符合應(yīng)變軟化型特征，達(dá)到峰值應(yīng)力狀態(tài)后會(huì)發(fā)生應(yīng)變軟化現(xiàn)象，最終發(fā)育到殘余狀態(tài)。試樣的峰值強(qiáng)度約為49 kPa，峰值應(yīng)力比為0.49；殘余強(qiáng)度約為31 kPa，殘余應(yīng)力比為0.31。根據(jù)剪切應(yīng)力-剪切位移、體積變化確定圖中A，B，C，D，E 5個(gè)狀態(tài)作為特征狀態(tài)以研究試樣在剪切過(guò)程中的宏細(xì)觀力學(xué)行為演變。這5個(gè)狀態(tài)對(duì)應(yīng)的剪切位移、體變情況、剪切應(yīng)力及描述列于表2中。

表2 特征狀態(tài)Table 2 Characteristic states

根據(jù)體變曲線可知，砂土試樣在剪切過(guò)程中先發(fā)生體積收縮后發(fā)生體積膨脹，剪切完成后整體表現(xiàn)出剪脹性，剪切位移為5 mm時(shí)的試樣沿豎向膨脹0.4 mm。進(jìn)入殘余狀態(tài)時(shí)，試樣體積基本不發(fā)生變化。結(jié)合剪切應(yīng)力-剪切位移特征曲線可知，該砂土試樣表現(xiàn)出密砂的性狀。

3.2 剪切中的變形演變

土體力學(xué)特性研究中，非連續(xù)、局部化發(fā)育的變形是重點(diǎn)研究問(wèn)題。為了研究剪切過(guò)程中試樣內(nèi)各個(gè)位置的變形演變情況，在本次直剪試驗(yàn)?zāi)M前，于試樣中設(shè)置了9個(gè)黑色縱列，監(jiān)測(cè)直剪過(guò)程中縱列的變形演變過(guò)程如圖3。

圖3 試樣變形演變Fig.3 Evolution of specimen’s deformation

由圖3(c)可以看出，剪切開(kāi)始時(shí)，試樣內(nèi)統(tǒng)一純剪切是變形的主導(dǎo)模式， 縱列基本不發(fā)生變形，顆粒的位移連續(xù)且均勻，處于純剪切狀態(tài)。達(dá)到峰值應(yīng)力狀態(tài)時(shí)，僅靠近側(cè)邊界的縱列1和縱列9發(fā)生了輕微的變形。峰值應(yīng)力狀態(tài)后，非線性的應(yīng)力應(yīng)變狀態(tài)成為主導(dǎo)，隨著剪切位移的增大，縱列變形由外向內(nèi)逐漸發(fā)育，應(yīng)變局部化現(xiàn)象越來(lái)越明顯。剪切完成后，變形集中于一個(gè)條帶狀區(qū)域內(nèi)，即為應(yīng)變局部化的剪切帶區(qū)域。剪切帶邊界區(qū)域變形最大，向內(nèi)逐漸減小。因此，試樣在峰值前以純剪切變形模式為主導(dǎo)，總體變形很小；峰值后以非線性應(yīng)力應(yīng)變?yōu)橹鲗?dǎo)，由剪切面兩側(cè)向內(nèi)逐漸發(fā)育變形，最終形成應(yīng)變局部化條帶區(qū)。

圖4 試樣剪切變形情況Fig.4 Shear deformation of specimen

圖4所示的是剪切后縱列的變形情況，可以看出，圖示陰影條帶區(qū)域內(nèi)的縱列發(fā)生了相對(duì)更為明顯的彎曲變形。可認(rèn)為該變形局部化條帶區(qū)域即為剪切帶，其厚度近似為(13～15)d50。

4 細(xì)觀力學(xué)特性分析

4.1 顆粒接觸狀態(tài)分析

4.1.1 力鏈網(wǎng)絡(luò)

力鏈網(wǎng)絡(luò)是土骨架體系上粒間接觸力的宏觀形式，可體現(xiàn)顆粒體系受外荷載作用時(shí)受力響應(yīng)機(jī)制的敏感性。圖5所示的是試樣在不同剪切階段的力鏈網(wǎng)絡(luò)演變。試樣在法向力作用下固結(jié)穩(wěn)定后，強(qiáng)力鏈基本沿豎直方向均勻分布。剪切開(kāi)始后，峰值應(yīng)力狀態(tài)前，在右上側(cè)墻體施加的剪切力作用下，力鏈網(wǎng)絡(luò)的分布開(kāi)始發(fā)生變化，變化的區(qū)域集中于右上角和左下角區(qū)域，試樣中部的力鏈基本還是沿著豎直方向。達(dá)到峰值應(yīng)力狀態(tài)時(shí)，力鏈網(wǎng)絡(luò)的均勻性和分布特征完全改變，力鏈集中形成幾條強(qiáng)力鏈，強(qiáng)力鏈在左下角和右上角沿著水平方向，在試樣中部沿著傾斜30°～45°的方向。到達(dá)峰值應(yīng)力狀態(tài)后，隨著剪切位移的增大，力鏈網(wǎng)絡(luò)進(jìn)一步演化，強(qiáng)力鏈更為集中，剪切5 mm后的傾角也近似為30°～45°。

圖5 接觸力鏈網(wǎng)絡(luò)演變Fig.5 Contact force chains in specimen

4.1.2 組構(gòu)分布

土的接觸組構(gòu)指標(biāo)[13]是反映其微細(xì)觀各向異性發(fā)展的重要指標(biāo)，其可有效量化土體微細(xì)觀特性對(duì)宏觀力學(xué)響應(yīng)的影響，以建立宏觀力學(xué)與細(xì)觀響應(yīng)的內(nèi)在聯(lián)系。接觸組構(gòu)中接觸法向組構(gòu)(E(θ))、法向接觸力組構(gòu)(fn(θ))和切向接觸力組構(gòu)(ft(θ))可采用Rothenburg等[14]建議的傅里葉函數(shù)分別表示為：

(1)

式中：θ是顆粒間接觸法向的角度；f0是所有顆粒間法向接觸力平均值；θa，θn，θt分別是顆粒間接觸法向、法向接觸力和切向接觸力各向異性分布主方向的角度；a，an，at是傅里葉擬合系數(shù)，分別反映顆粒間接觸法向、法向接觸力和切向接觸力分布的各向異性程度。

圖6 接觸法向分布Fig.6 Normal direction distribution of contact

圖7 法向接觸力分布Fig.7 Distribution of normal contact force

本次數(shù)值試驗(yàn)得到剪切過(guò)程中試樣的接觸法向分布、法向接觸力和切向接觸力分布見(jiàn)圖6至圖8所示。其中：圖6中的數(shù)字為接觸數(shù)目；圖7、圖8中的

圖8 切向接觸力分布Fig.8 Distribution of tangential contact force

數(shù)字均為接觸力(N);接觸法向指標(biāo)為接觸切面法線在對(duì)應(yīng)方向上的接觸數(shù)目；法向接觸力和切向力為對(duì)應(yīng)接觸法向上法向接觸力和切向力的平均值；粗實(shí)線即為式(1)中傅里葉分布函數(shù)的擬合曲線?？砂l(fā)現(xiàn)剪切前試樣內(nèi)粒間接觸法向分布較為均勻，而剪切破壞了試樣接觸法向分布的均勻性。隨著剪切位移增加，0°～60°內(nèi)的接觸數(shù)逐漸增多，而120°～180°內(nèi)的后粒間的接觸數(shù)逐漸減少。達(dá)到峰值應(yīng)力狀態(tài)時(shí)，接觸法向分布的形狀由圓形變?yōu)闄E圓形，主軸偏轉(zhuǎn)至45°左右。峰值后的軟化階段與殘余階段，接觸法向的主軸方向基本不發(fā)生變化，形狀只有微小的改變。

接觸力可分解為法向和切向接觸力。剪切過(guò)程中切向接觸力比法向接觸力小一個(gè)量級(jí)，二者矢量疊加所得接觸力的大小和方向?qū)嶋H與法向接觸力非常相近。因此，理論上可認(rèn)為法向接觸力分布的演變與力鏈網(wǎng)絡(luò)的演變本質(zhì)上是非常相近的。由模擬結(jié)果可知，剪切過(guò)程中法向接觸力分布近似呈“花生狀”或橢圓形，主軸方向的演變近似為：90°→60°→25°→35°→40°，與力鏈網(wǎng)絡(luò)的演變情況基本一致，也可發(fā)現(xiàn)法向接觸力在達(dá)到峰值前逐漸增大的現(xiàn)象。切向接觸力分布呈對(duì)稱的“花瓣?duì)睢?，在剪切過(guò)程中也發(fā)生了大小變化與類似的主軸偏轉(zhuǎn)現(xiàn)象。

圖9 配位數(shù)分布演變Fig.9 Evolution of coordinate number distribution

4.1.3 配位數(shù)

配位數(shù)指顆粒的平均接觸數(shù)，配位數(shù)越大，顆粒間的受力傳遞越充分[15]。砂土剪脹性本質(zhì)上是由砂土顆粒在荷載作用下發(fā)生重排列而導(dǎo)致的，顆粒的滑移、轉(zhuǎn)動(dòng)和重組必然會(huì)導(dǎo)致試樣內(nèi)配位數(shù)的變化。剪切前后的配位數(shù)分布如圖9所示。

若將試樣分成如圖10所示的6個(gè)寬度相等、距盒底高度不同的條帶，由下到上分別編號(hào)為Z1，Z2， Z3，Z4，Z5，Z6，可得到剪切過(guò)程中各個(gè)條帶區(qū)域內(nèi)的平均配位數(shù)變化情況。

圖10 6個(gè)條帶的配位數(shù)變化曲線Fig.10 Curves of coordination number in six zones

由圖9可以看出，剪切作用下試樣內(nèi)配位數(shù)的大小及分布發(fā)生變化，剪切后試樣內(nèi)多數(shù)區(qū)域內(nèi)配位數(shù)都有所降低。由圖10可以看出，各條帶內(nèi)的配位數(shù)在剪切開(kāi)始后隨剪切位移的增大而逐漸減小，當(dāng)剪切位移達(dá)到1.5 mm后基本保持不變，穩(wěn)定在3.0～3.3之間。同時(shí)也可看出，剪切帶內(nèi)配位數(shù)相對(duì)于帶外降低更顯著。

4.1.4 滑動(dòng)接觸百分比

滑動(dòng)接觸百分比是顆粒間發(fā)生滑動(dòng)的接觸數(shù)占總接觸數(shù)的比例。剪切前后的滑動(dòng)接觸百分比分布如圖11。

圖11 滑動(dòng)接觸百分比分布Fig.11 Evolution of sliding contact distribution

圖11中出現(xiàn)了滑動(dòng)接觸百分比為-1%的情況，這是由于云圖網(wǎng)格計(jì)算方法導(dǎo)致的，實(shí)際并不會(huì)出現(xiàn)此種情況。對(duì)比剪切前后的滑動(dòng)接觸百分比分布圖可知，剪切完成后剪切面附近區(qū)域接觸滑動(dòng)比例明顯增大，意味著剪切作用下剪切面附近區(qū)域內(nèi)顆粒間的接觸發(fā)生了相對(duì)更為劇烈的滑動(dòng)，可認(rèn)為此條帶即為剪切帶。

4.2 局部孔隙率演化分析

剪切過(guò)程中試樣孔隙率分布的變化情況如圖12所示。從圖12中可以看出，剪切開(kāi)始后，試樣由剪切面兩端附近的局部區(qū)域開(kāi)始發(fā)生孔隙率增大現(xiàn)象，并隨著剪切位移的增大逐漸向內(nèi)發(fā)育直至貫通。最終，試樣在剪切面附近形成了明顯的孔隙率增大的條帶區(qū)域，可認(rèn)為其即為剪切帶，孔隙率增大的現(xiàn)象本質(zhì)上是由密砂的剪脹性導(dǎo)致的。由圖12(c)得到剪切帶厚度近似為(10～15)d50。

圖12 孔隙率分布Fig.12 Distribution of porosity

按與圖10一致的形式進(jìn)行分區(qū)，得到剪切過(guò)程中各個(gè)條帶內(nèi)的孔隙率演化情況如圖13所示?？梢钥闯?，剪切面附近的條帶Z3和Z4發(fā)生了明顯的孔隙率增大現(xiàn)象，其他條帶內(nèi)的孔隙率變化量相對(duì)較小，驗(yàn)證了上述的孔隙演變分析。因此，密砂在剪切過(guò)程中，剪切面附近區(qū)域由外向內(nèi)逐漸發(fā)育孔隙率增大區(qū)域，貫通形成剪切帶，剪切帶內(nèi)發(fā)生剪脹，導(dǎo)致試樣整體的剪脹；帶剪切帶外發(fā)生輕微剪脹或剪縮，對(duì)試樣的剪脹性質(zhì)影響不大。

圖13 6個(gè)條帶內(nèi)的孔隙率變化Fig.13 Porosity changes in six zones

4.3 顆粒運(yùn)動(dòng)狀態(tài)分析

4.3.1 顆粒速度場(chǎng)

砂樣在B和C這2個(gè)特征狀態(tài)的瞬時(shí)速度場(chǎng)如圖14所示。由峰值前B狀態(tài)的速度場(chǎng)可以看出，剛開(kāi)始剪切時(shí)，試樣除了4個(gè)角的區(qū)域外，整體速度較均勻，符合峰值前統(tǒng)一純剪變形模式。達(dá)到圖14(b)所示的峰值應(yīng)力狀態(tài)時(shí)，試樣內(nèi)速度出現(xiàn)上盒大下盒小的現(xiàn)象，剪切面附近尤其是兩端的顆粒大多都有豎向速度分量，此時(shí)顆粒發(fā)生上下錯(cuò)動(dòng)，會(huì)促進(jìn)剪脹現(xiàn)象的發(fā)生。認(rèn)為剪切面附近的剪切帶為速度過(guò)渡區(qū)域，不僅發(fā)生速度大小的過(guò)渡，方向也發(fā)生變化。試樣在剪切面左側(cè)邊界和右側(cè)邊界局部區(qū)域內(nèi)出現(xiàn)了速度渦旋現(xiàn)象。

圖14 速度場(chǎng)Fig.14 Velocity fields

4.3.2 顆粒累計(jì)轉(zhuǎn)動(dòng)量

基于顆粒流程序的累計(jì)轉(zhuǎn)動(dòng)量函數(shù)統(tǒng)計(jì)顆粒的累計(jì)轉(zhuǎn)動(dòng)量，以特定轉(zhuǎn)動(dòng)量(本數(shù)值試驗(yàn)設(shè)為15°)為界限，將轉(zhuǎn)動(dòng)量大于特定值和不大于特定值的顆粒分別設(shè)定為黃色和紅色。得到剪切過(guò)程中的顆粒累計(jì)轉(zhuǎn)動(dòng)情況如圖15。可發(fā)現(xiàn)，峰值前試樣內(nèi)顆粒累計(jì)轉(zhuǎn)動(dòng)量很小，達(dá)到峰值應(yīng)力狀態(tài)時(shí)，在剪切面兩側(cè)的邊界形成了劇烈轉(zhuǎn)動(dòng)區(qū)域。峰值應(yīng)力狀態(tài)后，隨著剪切位移的增大，劇烈變形區(qū)域由兩側(cè)逐漸向內(nèi)部發(fā)育，貫通后形成初始的帶狀區(qū)域，進(jìn)一步發(fā)育，條帶的厚度逐漸增大。

圖15 顆粒累計(jì)轉(zhuǎn)動(dòng)量演變Fig.15 Evolution of particle’s cumulative rotation

結(jié)合變形和孔隙演變可知，剛開(kāi)始剪切時(shí)，剪切帶并未形成，隨著剪切位移的增大，變形不斷增大，顆粒轉(zhuǎn)動(dòng)不斷加劇，變形、轉(zhuǎn)動(dòng)和膨脹由剪切面兩側(cè)的局部區(qū)域逐漸向內(nèi)發(fā)育貫通，最終形成了變形、顆粒轉(zhuǎn)動(dòng)和體脹集中的剪切帶區(qū)域。

4.4 主應(yīng)力場(chǎng)分析

基于二向應(yīng)力狀態(tài)解析法計(jì)算各測(cè)量圈內(nèi)的主應(yīng)力和主方向，并以測(cè)量圈的圓心為基點(diǎn)計(jì)算主應(yīng)力矢量的端點(diǎn)，建立表示2個(gè)主應(yīng)力相對(duì)大小和分布方向的可視化墻單元，以實(shí)現(xiàn)主應(yīng)力的可視化。得到砂樣在直剪過(guò)程中主應(yīng)力分布如圖16所示。

圖16 主應(yīng)力場(chǎng)分布Fig.16 Distribution of principal stress field

剪切前試樣的初始應(yīng)力場(chǎng)分布較為均勻，主應(yīng)力分布反映的試樣內(nèi)部受力狀態(tài)與力鏈網(wǎng)絡(luò)一致。隨著剪切位移的增大，主應(yīng)力場(chǎng)開(kāi)始出現(xiàn)明顯的偏轉(zhuǎn)并且變得不均勻，左上角和右下角的主應(yīng)力相對(duì)減小。峰值應(yīng)力狀態(tài)時(shí)，主應(yīng)力的偏轉(zhuǎn)最明顯且統(tǒng)一。結(jié)合力鏈網(wǎng)絡(luò)及接觸組構(gòu)分布變化可知，峰值應(yīng)力狀態(tài)時(shí)，接觸組構(gòu)的主軸方向及強(qiáng)力鏈的方向與主應(yīng)力偏轉(zhuǎn)后的方向基本一致，此時(shí)砂樣的剪切抗力能得到最大發(fā)揮。這表明，主應(yīng)力場(chǎng)、接觸組構(gòu)分布和力鏈網(wǎng)絡(luò)反映的細(xì)觀受力機(jī)制在本質(zhì)上是一致的，剪切中三者發(fā)生了近似一致的主方向偏轉(zhuǎn)。

5 結(jié) 論

(1) 直剪過(guò)程中砂樣的應(yīng)力、應(yīng)變分布都存在不均勻性。剪切過(guò)程中，試樣變形集中于剪切帶內(nèi)；主應(yīng)力和力鏈分布區(qū)域化且發(fā)生了主方向的偏轉(zhuǎn)。強(qiáng)力鏈、接觸組構(gòu)主軸及主應(yīng)力方向發(fā)生了一致的偏轉(zhuǎn)，三者本質(zhì)上具有一致性，反映了顆粒體系在外荷載作用下的內(nèi)部響應(yīng)。

(2) 受剪切力和法向力綜合作用，試樣由剪切面兩端開(kāi)始發(fā)育局部變形，顆粒發(fā)生劇烈轉(zhuǎn)動(dòng)及上下錯(cuò)動(dòng)，導(dǎo)致孔隙率的增大和接觸的滑動(dòng)，隨著剪切位移增大逐漸向內(nèi)貫通形成剪切帶。帶內(nèi)外顆粒處于不同的運(yùn)動(dòng)狀態(tài)，導(dǎo)致了變形局部化區(qū)域的發(fā)育，剪切帶的厚度近似為(10～15)d50。

(3) 直剪試驗(yàn)條件下剪切帶內(nèi)外顆粒處于不同的運(yùn)動(dòng)狀態(tài)，帶內(nèi)顆粒發(fā)生轉(zhuǎn)動(dòng)和上下錯(cuò)移，帶外顆粒發(fā)生平動(dòng)，導(dǎo)致試樣組構(gòu)、孔隙率、配位數(shù)、滑動(dòng)接觸百分比等細(xì)觀參量分布的變化，這些細(xì)觀參量的變化也反映了剪切帶宏觀應(yīng)變局部化現(xiàn)象與剪脹現(xiàn)象。

[1] 史旦達(dá), 周 健, 劉文白, 等. 砂土直剪力學(xué)性狀的非圓粒模擬與宏細(xì)觀機(jī)理研究[J]. 巖土工程學(xué)報(bào), 2010, 32(10): 1557-1565.

[2] 蔡正銀, 李相崧. 無(wú)黏性土中剪切帶的形成過(guò)程[J]. 巖土工程學(xué)報(bào), 2003, 25(2) : 129-134.

[3] 張 兆, 劉元雪. 應(yīng)變局部化理論及在巖土工程中的應(yīng)用[J]. 后勤工程學(xué)院學(xué)報(bào), 2004, (2):12-15.

[4] 劉文白, 張 輝, 鄧一兵. 基于DPDM技術(shù)的砂土直剪試驗(yàn)剪切過(guò)程中的應(yīng)力場(chǎng)分析[J]. 中國(guó)水運(yùn), 2008, 8(7): 235-239.

[5] POTTS D M, DOUNIAS G T, VAUGHAN P R. Finite Element Analysis of the Direct Shear Box Test[J]. Geotechnique, 1987, 37(1): 11-23.

[6]TEJCHMAN J, BAUER E. FE-simulation of a Direct Shear and a True Simple Shear Test within a Polar Hypoplasticity[J]. Computers and Geotechnics, 2005, 32(1):1-16.

[7] THORNTON C, ZHANG L. Numerical Simulation of the Direct Shear Test[J]. Chemical Engineering and Technology, 2003, 26(2): 153-156.

[8] LIU S H, SUN D A, MATSUOKA H. On the Interface Friction in Direct Shear Test[J]. Computers and Geotechnics, 2005, 32(5):317-325.

[9] MASSON S, MARTINEZ J. Micromechanical Analysis of the Shear Behavior of a Granular Material[J]. Journal of Engineering Mechanics, 2001, 127(10): 1007-1016.

[10]WANG J, GUTIERREZ M. Discrete Element Simulation of Direct Shear Specimen Scale Effects[J]. Geotechnique, 2010, 60(5): 359-409.

[11]蔣明鏡, 王富周, 朱合華. 單粒組密砂剪切帶的直剪試驗(yàn)離散元數(shù)值分析[J]. 巖土力學(xué), 2010, 31(1): 253-257.

[12]MISHRA B K, MEHROTRA S P. Modeling of Particle Stratification in Jigs by the Discrete Element Method[J]. Minerals Engineering, 1998, 11(6): 511-522.

[13]ODA M. Initial Fabrics and Their Relations to Mechanical Properties of Granular Material[J]. Soils and Foundations, 1972, 12(1):17-36.

[14]ROTHENBRUG L, BATHURST J. Analytical Study of Induced Anisotropy in Idealized Granular Materials[J]. Geotechnique, 1989, 39(4):601-614.

[15]劉 洋, 吳順川, 周 健. 單調(diào)荷載下砂土變形過(guò)程數(shù)值模擬及細(xì)觀機(jī)制研究[J]. 巖土力學(xué), 2008, 29(12): 3199-3204.

(編輯：姜小蘭)

Exploring Mesoscopic Deformation Mechanism of Sand in Direct ShearTest by Numerical Simulation Using Discrete Element Method

LI Shuang, LIU Yang, WU Ke-jia

(School of Civil and Resource Engineering, University of Science and Technology Beijing,Beijing 100083, China)

Direct shear test of sand was simulated using discrete element method and the macro-and-mesoscopic responses of sand during shearing were analyzed in this research. The development of sand deformation was observed through strip partition method; and the relationship between shear stress and shear displacement was recorded to analyze the macroscopic behaviors of sand. In mesoscopic scale, the development of contact status was researched through analyzing contact force chain, contact fabric distribution, coordination number and sliding contact percentage; the evolution of granular distribution under shearing was observed through the analysis of porosity variation;the granular motion state and the mesoscopic formation mechanism of shear band were studied through observing the cumulative rotation and velocity of particles. Furthermore, principal stresses and their inclinations in specimen were visualized through secondary development of PFC2D. The study shows that strain is concentrated in the shear band in which particle rotation, contact sliding and obvious dilatancy are observed. Different motion state of particles inside and outside the shear band leads to the changes of mesoscopic parameters, and the thickness of shear band is 10-15 times ofd50. Contact force chain network, contact fabric distribution and principal stress distribution are essentially consistent, and the principal direction changes consistently during shearing.

sand; discrete element method; direct shear test; meso-mechanics; shear band

2015-11-27；

2015-12-21

國(guó)家自然科學(xué)基金項(xiàng)目(51178044)；北京高校青年英才計(jì)劃資助項(xiàng)目(YETP0340)

李 爽(1993-)，男，安徽滁州人，碩士研究生，主要從事散粒體細(xì)觀力學(xué)研究與數(shù)值模擬，(電話)15201451398(電子信箱)ls19931014@126.com。

10.11988/ckyyb.20151006

2017,34(4):104-110,116

TU43

A

1001-5485(2017)04-0104-07