砂土直剪試驗離散元數(shù)值模擬與細觀變形機理研究

李 爽，劉 洋，吳可嘉

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

砂土直剪試驗離散元數(shù)值模擬與細觀變形機理研究

李 爽，劉 洋，吳可嘉

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

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

砂土；離散元法；直剪試驗；細觀力學；剪切帶

1 研究背景

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

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

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

2 直剪試驗離散元分析過程

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

粒間接觸模型采用接觸剛度模型，符合無黏性土的特性。試樣主要的細觀參數(shù)包括顆粒間及顆粒與墻體間法向接觸剛度、切向接觸剛度、摩擦系數(shù)，取值列于表1中。其中，顆粒間的接觸剛度是決定散粒試樣彈性模量的重要參數(shù)，根據(jù)散粒體彈性理論，當顆粒接觸剛度的取值在108～109N/m之間，法向和切向剛度比值在1～3之間時，模型收斂性良好[12]。本次數(shù)值試驗設(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 直剪試驗數(shù)值模型及測量圈Fig.1 Model of direct shear simulation test and measurement circles

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

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

3 宏觀力學特性分析

3.1 宏觀應力-應變特性

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

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

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

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

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

3.2 剪切中的變形演變

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

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

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

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

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

4 細觀力學特性分析

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

4.1.1 力鏈網(wǎng)絡

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

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

4.1.2 組構(gòu)分布

土的接觸組構(gòu)指標[13]是反映其微細觀各向異性發(fā)展的重要指標，其可有效量化土體微細觀特性對宏觀力學響應的影響，以建立宏觀力學與細觀響應的內(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ù)值試驗得到剪切過程中試樣的接觸法向分布、法向接觸力和切向接觸力分布見圖6至圖8所示。其中：圖6中的數(shù)字為接觸數(shù)目；圖7、圖8中的

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

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

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

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

4.1.3 配位數(shù)

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

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

圖10 6個條帶的配位數(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ù)在剪切開始后隨剪切位移的增大而逐漸減小，當剪切位移達到1.5 mm后基本保持不變，穩(wěn)定在3.0～3.3之間。同時也可看出，剪切帶內(nèi)配位數(shù)相對于帶外降低更顯著。

4.1.4 滑動接觸百分比

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

圖11 滑動接觸百分比分布Fig.11 Evolution of sliding contact distribution

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

4.2 局部孔隙率演化分析

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

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

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

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

4.3 顆粒運動狀態(tài)分析

4.3.1 顆粒速度場

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

圖14 速度場Fig.14 Velocity fields

4.3.2 顆粒累計轉(zhuǎn)動量

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

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

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

4.4 主應力場分析

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

圖16 主應力場分布Fig.16 Distribution of principal stress field

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

5 結(jié) 論

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

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

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

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

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

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

[4] 劉文白, 張 輝, 鄧一兵. 基于DPDM技術(shù)的砂土直剪試驗剪切過程中的應力場分析[J]. 中國水運, 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]蔣明鏡, 王富周, 朱合華. 單粒組密砂剪切帶的直剪試驗離散元數(shù)值分析[J]. 巖土力學, 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)荷載下砂土變形過程數(shù)值模擬及細觀機制研究[J]. 巖土力學, 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

國家自然科學基金項目(51178044)；北京高校青年英才計劃資助項目(YETP0340)

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

10.11988/ckyyb.20151006

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

TU43

A

1001-5485(2017)04-0104-07