郭 鴻，陳 茜，李 軍，陳棟梁

(1.陜西理工大學(xué)土木工程與建筑學(xué)院，陜西 漢中 723000；2. 西北農(nóng)林科技大學(xué)水利水電科學(xué)研究院，陜西 楊凌 7121001)

二維理想顆粒受剪行為的細(xì)觀機理*

郭 鴻1,2，陳 茜1，李 軍1，陳棟梁1

(1.陜西理工大學(xué)土木工程與建筑學(xué)院，陜西 漢中 723000；2. 西北農(nóng)林科技大學(xué)水利水電科學(xué)研究院，陜西 楊凌 7121001)

通過對二維理想顆粒(圓盤)的直接剪切細(xì)觀力學(xué)機理分析，研究了顆粒受剪時的運動規(guī)律，建立了能夠描述等徑顆粒情況下最密排列顆粒的最大剪應(yīng)力公式。并通過對另外四種不同級配的圓盤顆粒的直接剪切離散元模擬，探討了所建立的最大剪應(yīng)力公式的應(yīng)用。結(jié)果表明：在豎向荷載一定時，最密排列顆粒的最大剪應(yīng)力取決于邊界顆粒影響參數(shù)和非邊界顆粒影響參數(shù)，且這兩個參數(shù)在最密排列等徑顆粒的情況下數(shù)值保持一致；對于一般顆粒級配，邊界影響參數(shù)基本保持不變，非邊界影響參數(shù)隨著d10，d30和d60的增加而線性減小。

顆粒材料；離散元法；直接剪切；顆粒級配

顆粒材料廣泛的存在于自然界中，在巖土工程、道路工程以及鐵路工程中的有著極廣泛的應(yīng)用，滑坡、地基剪切失穩(wěn)等現(xiàn)象都與顆粒材料的剪切特性有直接關(guān)系[1-4]。目前的對于顆粒材料剪切特性的研究大致有兩種方法，即宏觀室內(nèi)試驗以及細(xì)觀數(shù)值模擬[5-9]。近年來興起的離散元方法對研究顆粒材料起到了很大的推動作用，其基本方法是針對材料定義不同的顆粒接觸模型，利用離散元模擬軟件來探索其細(xì)觀指標(biāo)諸如力鏈、配位數(shù)、孔隙率、局部剪應(yīng)力和位移場等對宏觀相應(yīng)的影響[10-13]。然而，由于目前還少有文獻在顆粒材料在受剪時，將其單個顆粒的受力及運動特征與宏觀層面剪切應(yīng)力建立定量關(guān)系，所以針對不同的顆粒材料，仍然需要試驗測定或離散元模擬的方法確定其剪切特性，諸如內(nèi)摩擦角等。理想圓形顆粒是顆粒材料最簡化的物理模型，分析其受力規(guī)律有望會對實體顆粒材料的研究提供有效的思路，目前已有研究針對理想顆粒的力學(xué)特性研究[14-19]，主要探討顆粒形狀、級配等因素對剪切特性的影響，但是存在的問題仍然是尚未建立顆粒細(xì)觀受力與宏觀相應(yīng)之間的定量關(guān)系，基本都是從力鏈演化等角度定性描述。

基于此，本文以二維等徑顆粒為研究對象，通過分析單顆粒在剪切過程中的受力平衡，推導(dǎo)了二維等徑顆粒規(guī)則排列時的剪力模型。并且運用基于PFC2D的數(shù)值剪切試驗，對所建立的模型，針對其他不同級配的二維理想顆粒，進行了應(yīng)用探討和驗證。

1 顆粒剪切運動受力分析

1.1 邊界顆粒

顆粒受剪時，同時存在法向力和切向剪力。分析單顆粒在受剪時的受力規(guī)律，如圖1(a)、(b)和(c)所示。顆粒在邊界時的初始狀態(tài)，上部顆粒處于平衡穩(wěn)定狀態(tài)，隨著剪力Fi的增大，顆粒發(fā)生緩慢逆時針旋轉(zhuǎn)，同時逐漸離開右下方所接觸的顆粒，當(dāng)剛好離開右下方的顆粒時，剪力達(dá)到最大值，隨后剪力減小，顆?！疤S”至前方另外兩個顆粒之間，此時顆粒材料已經(jīng)剪壞。豎向荷載作用在邊界單個顆粒上的分力Pi仍然是豎直分向，其值為

Pi=pD/m

(1)

式中p為豎向總應(yīng)力，D為試樣水平方向總寬度，m為水平方向的顆?？倲?shù)。

當(dāng)單個顆粒在同時受到Pi和Fi作用下，它的運動從初始穩(wěn)定點逐漸運動至另一個穩(wěn)定點，在該運動過程中，顆粒被剪壞。如圖1(d)受力分析，這種情況時其水平方向的動態(tài)平衡方程為

(2)

式中N為下部顆粒對其產(chǎn)生的法向支持力，μ為顆粒之間的摩擦系數(shù)，θ為該顆粒和下部顆粒幾何中心連線與水平方向的夾角。

豎直方向的平衡方程為

(3)

式中m為單個顆粒的質(zhì)量，g為重力加速度。

由式(2)和式(3)得出邊界顆粒情況下，單個顆粒移動所需要的水平剪力為

(4)

圖1 邊界顆粒受剪運動過程及受力分析Fig.1 Movement and force analysis of singe particle in the boundary under shear force

1.2 非邊界顆粒顆粒

另一種情況就是顆粒在非邊界上，其剪切運動過程如圖2(a)、(b)和(c)所示。豎向荷載作用在非邊界單個顆粒上的分力P′不在豎直方向，假設(shè)左上角相鄰顆粒與該單個顆粒同步剪切運動，則其值為

(5)

圖2 非邊界顆粒受剪運動過程及受力分析Fig.2 Movement and force analysis of singe particle within the boundary under shear force

對于此種情況來說，當(dāng)單個顆粒在同時受到P′和Fi作用時，它的運動從初始穩(wěn)定點逐漸運動至另一個穩(wěn)定點，受力分析如圖2(d)示，此時，水平方向的動態(tài)平衡方程為

(6)

豎直方向的平衡方程為

(7)

由式(6)和式(7)得出非邊界顆粒情況下，單個顆粒移動所需要的水平剪力為

(8)

1.3 一般情況

推廣到一般的等徑顆粒最密排列情景，如圖3所示，若每層m個顆粒，考慮到顆粒數(shù)目的影響效應(yīng)，則要使該層的所有顆粒都同步運動，則對于每層來說需要的的剪力為

F=fjFi

(9)

式中fj(j=1，2)是與邊界條件有關(guān)的系數(shù)，f1表示邊界層顆粒影響系數(shù)，f2表示非邊界層顆粒影響系數(shù)。

假設(shè)有一共有2n層顆粒，所有層總剪力為

(10)

得到的剪應(yīng)力τ為

(11)

理想情況下，規(guī)則排列的剛性顆粒在受剪時其運動是同步的，然而考慮到顆粒之間非剛性接觸，其運動并非同步進行，且每層不同顆粒運動相應(yīng)存在一定的相位差。由于非同步現(xiàn)象涉及到顆粒之間的接觸本構(gòu)，且影響因素諸多，故為了簡單起見，本文只考慮邊界層顆粒影響系數(shù)f1和非邊界層顆粒影響系數(shù)f2。

2 二維等徑顆粒DEM剪切試驗

建立二維等徑DEM(離散元)數(shù)值剪切試樣，如圖4所示，所有試樣都是等徑顆粒最密排列。保持理想試樣的尺寸為6cm×2cm，設(shè)置不同粒徑的二維顆粒，即5，2.5，1.667，1.25，1，0.833，0.714和0.625mm。設(shè)定顆粒間的摩擦系數(shù)為1，顆粒接觸為最簡線彈性模型，墻體剛度和顆粒剛度分別為500MN/m和50MN/m。

圖3 m×2n理想剪切試樣Fig.3 m×2nideal shear sample

圖4 不同粒徑的二維離散元剪切試樣Fig.4 2D DEM shear sample of different particle sizes

對不同粒徑的試樣做離散元數(shù)值剪切試驗，編寫自動加壓伺服程序段，使得在試樣剪切前顆粒體上部的豎向壓應(yīng)力為100kPa，然后固定數(shù)值剪切盒的下半部分，向左移動數(shù)值剪切盒的上半部分。需要說明的是，為了提高模擬計算速度，虛擬剪切速度采用2cm/s，在下文中所有的數(shù)值試驗均采用該剪切速度。

離散元數(shù)值剪切試驗得到的內(nèi)摩擦角和半徑的關(guān)系如圖5所示。由此圖可以看出，粒徑對摩擦角的影響十分顯著，且根據(jù)擬合公式，二者基本呈線性關(guān)系。當(dāng)豎直荷載為100kPa時，數(shù)值剪切試驗得到的最大剪應(yīng)力和顆粒半徑的關(guān)系如圖6所示，菱形點為DEM試驗值，方形點數(shù)據(jù)為本文式(11)的計算值。在本文所設(shè)定的計算參數(shù)和DEM參數(shù)的情況下，f1=2.4,f2=3.3。說明非邊界顆粒影響系數(shù)對最大剪切力的影響較大，這一方面是因為非邊界顆粒的平均配位數(shù)較大，另一方面是由于內(nèi)部顆粒在剪切的過程中顆粒相互咬合更多(而邊界顆粒的相對滑動更多)。換言之就是非邊界顆粒相對邊界顆粒的運動阻力更大，而各個顆粒的定向運動阻力之和和宏觀上的剪切力是平衡的，因此f2的數(shù)值大于f1。

圖5 顆粒半徑和內(nèi)摩擦角的關(guān)系Fig.5 Rrelation of particle radius vs internal friction angle

圖6 顆粒半徑和最大剪應(yīng)力的關(guān)系Fig.6 Relation of particle radius vs maximum shear stress

3 非規(guī)則顆粒數(shù)值試驗

實際工程中存在的材料大都屬于非等徑且具有復(fù)雜級配的顆粒介質(zhì)。據(jù)此，本節(jié)以四種不同級配的顆粒材料為對象，進行數(shù)值剪切試驗，進一步探索顆粒大小級配對最大剪應(yīng)力的影響規(guī)律，同時討論二維非等徑理想顆粒情況下式(11)的應(yīng)用。如圖7所示，為四種級配的二維顆粒(圓盤)，采用Jiang氏方法[20]分四層生成。顆粒級配曲線如圖8所示，Cc和Cu值如表1所示。

表1 四種級配砂土的Cc和Cu值

在處理非等徑顆粒情況時，先計算出每種顆粒級配情況下的加權(quán)平均粒徑

圖7 四種顆粒級配的離散元模型Fig.7 Discrete model of the 4 types sands

圖8 四種級配試樣Fig.8 Four DEM samples of different particle size distributions

(12)

式中Ri代表某種級配情況下第i級的顆粒半徑，Ni代表該i級顆粒的數(shù)目，n為不同粒徑的總級數(shù)。

然后將處理后的等徑顆粒按照矩形方式排列為2n行×m列，再用式(10)進行計算最大剪切力。在以上四種顆粒級配情況下，f1=1，f2與d10，d30和d60有關(guān)，如圖9所示?？梢钥闯?，f2基本都隨d10，d30和d60的增大而減小，且其線性關(guān)系比較明顯。需要進一步說明的是，f1是與邊界顆粒有關(guān)的系數(shù)，f2是與非邊界顆粒有關(guān)的系數(shù)，所以當(dāng)數(shù)值試樣內(nèi)的顆粒數(shù)量增大時，f1的影響相對會減弱，而f2的影響相對會加強。

圖9 d10，d30和d60和f2的關(guān)系曲線Fig.9 Relation of d10，d30 and d60 vs f2

這是因為在數(shù)值試樣內(nèi)部，隨著顆粒數(shù)量的增加，水平方向的力鏈數(shù)量更多，形狀更細(xì)，更彎曲，如圖10所示，而對于邊界上的顆粒來說，由于上部豎向均布壓應(yīng)力，顆粒的排列相對規(guī)則，形態(tài)變化也相對不明顯。f2和最大剪應(yīng)力的關(guān)系如圖10所示，可以看出最大剪應(yīng)力隨著f2的增大而線性增大，表明了力鏈越密集，剪應(yīng)力越大(見圖11)。

圖10 f2和τ的關(guān)系曲線Fig.10 Relation of f2 vs τ

圖11 四種級配離散元剪切過程力鏈圖Fig.11 Force chains of four samples with different particle size distributions in DEM shear test

4 結(jié) 論

本文基于二維等徑顆粒體，分析了顆粒材料在受剪時的力學(xué)機理，提出了顆粒規(guī)則最密規(guī)則排列時的最大剪應(yīng)力公式。通過對四種不同級配顆粒體的數(shù)值剪切試驗，發(fā)現(xiàn)本文提出的最大剪應(yīng)力公式仍然適用，所不同的只是f1和f2，其中，兩者在非等徑顆粒情況下的數(shù)值都較等徑顆粒情況下的數(shù)值小，而且后者隨著d10，d30和d60的增大而減小。f2和最大剪應(yīng)力的關(guān)系表明，力鏈形態(tài)越密集，剪應(yīng)力越大。

[1] 郭鴻,駱亞生,程大偉,等. 傾斜推移體激發(fā)顆粒流的連續(xù)介質(zhì)模型[J]. 中山大學(xué)學(xué)報(自然科學(xué)版),2015,54(6):94-98.

[2] 李洋溢. 條形基礎(chǔ)加筋砂土地基室內(nèi)模型試驗的分析研究[D].上海：同濟大學(xué),2006.

[3] 曹延波,彭芳樂,小竹望,等. 加筋砂土邊坡漸進性變形破壞的數(shù)值分析[J]. 巖石力學(xué)與工程學(xué)報,2010,29(S2):3905-3915.

[4] 吳劍. 滑帶剪切過程的離散元模擬研究[D].北京：中國科學(xué)院研究生院,2007.

[5] 嚴(yán)穎,季順迎. 顆粒形態(tài)對離散介質(zhì)剪切強度的影響[J]. 巖土力學(xué),2009,30(S1):225-230，249.

[6] 邊學(xué)成,李偉,李公羽,等. 基于顆粒真實幾何形狀的鐵路道砟剪切過程三維離散元分析[J]. 工程力學(xué),2015,32(5):64-75，83.

[7] 秦建敏,張洪武. 顆粒材料的微觀臨界狀態(tài)理論模型[J]. 巖土力學(xué),2010,31(12):3697-3703.

[8] 秦建敏,遲璐璐. 顆粒材料剪脹性的微觀力學(xué)分析[J]. 巖土力學(xué),2013,34(5):1508-1514.

[9] 孫其誠,程曉輝,季順迎,等. 巖土類顆粒物質(zhì)宏-細(xì)觀力學(xué)研究進展[J]. 力學(xué)進展,2011,41(3):351-371.

[10] 張宇,余飛,陳善雄,等.CAS-1模擬月壤動剪切模量與阻尼比的試驗研究[J]. 巖土力學(xué),2014,35(1):74-82.

[11] 蔣明鏡,彭鏑,申志福,等. 深海能源土剪切帶形成機理離散元分析[J]. 巖土工程學(xué)報,2014,36(9):1624-1630.

[12] 蔣明鏡,張望城,王劍鋒. 密實散粒體剪切破壞能量演化的離散元模擬[J]. 巖土力學(xué),2013,34(2):551-558.

[13] 劉動,陳曉平. 滑帶土環(huán)剪剪切面的微觀觀測與分析[J]. 巖石力學(xué)與工程學(xué)報,2013,32(9):1827-1834.

[14] 李罡,劉映晶,尹振宇,等. 粒狀材料臨界狀態(tài)的顆粒級配效應(yīng)[J]. 巖土工程學(xué)報,2014,36(3):452-457.

[15] 孔亮,季亮亮,曹杰峰. 應(yīng)力路徑和顆粒級配對砂土變形影響的細(xì)觀機制[J]. 巖石力學(xué)與工程學(xué)報,2013,32(11):2334-2341.

[16] 陸勇,周國慶,賴澤金. 砂土剪切特性的荷載與粒徑效應(yīng)研究[J]. 中國礦業(yè)大學(xué)學(xué)報,2014,43(2):195-202.

[17]DABEETA,WIJEWICKREMED,BYRNEP.Applicationofdiscreteelementmodelingforsimulationofcyclicdirectsimpleshearresponseofgranularmaterials[C]//NationalConferenceonEarthquakeEngineering， 2014: 21-25.

[18]CAVARRETTAI,O’SULLIVANC,IBRAIME,etal.CharacterizationofartificialsphericalparticlesforDEMvalidationstudies[J].Particuology, 2012, 10(2): 209-220.

[19]O’SULLIVANC.Particulatediscreteelementmodelling[M].Taylor&Francis, 2011.

[20]JIANGMingjing,KONRADJM,LEROUEILS.AnefficienttechniqueforgeneratinghomogeneousspecimensforDEMstudies[J].ComputersandGeotechnics, 2003, 30 (5): 579-597.

The microscopic mechanism of the shear properties of two dimensional ideal particle

GUOHong1,2,CHENXi1，LIJun1,CHENDongliang1

(1. School of Civil Engineering and Architecture, Shaanxi University of Technology,Hanzhong 723000, China；2. Institute of Water Resources and Hydropower Research, Northwest A&F University,Yangling 712100, China)

By analyzing the micromechanics of the two-dimensional ideal particles (disc) under direct shear test, the movement during the shearing is investigated. The formula described the maximum shear stress in the regular dense packing condition is established. In addition, following the direct shear simulations of the four samples with different particle size distributions, application of the formula which could predict the maximum shear stress was discussed. The results show that the maximum shear stress is depending on two parameters namely boundary effect parameter and the non-boundary effect parameters, and these two parameters keep the same under the regular dense packing when the vertical loading is constant. For the general packing, the boundary effect parameters remained unchanged, and the non-boundary effect parameter is linearly decreases with thed10,d30andd60.

particulate material; discrete element method; direct shear; particle size distribution

10.13471/j.cnki.acta.snus.2016.04.008

2015-10-26

國家自然科學(xué)基金資助項目(51178392) ；陜西理工大學(xué)科研資助項目(SLGKY15-24)

郭鴻(1984年生)，男;研究方向：離散元數(shù)值模擬；通訊作者： 駱亞生;E-mail:lyas1967@nwsuaf.edu.cn

TU442

A

0529-6579(2016)04-0047-06