顆粒土中剪切帶臨界狀態數學描述及其完全解*

黃文雄, 崔 賢

(河海大學 力學與材料學院, 南京 210098)

0 引 言

臨界狀態在土力學中指土體剪切失效時所趨近的無剪脹塑性流動狀態．土體的臨界狀態只與土的物理特性及作用于土體的平均壓力有關,不依賴土體的初始密度,因此是建立土體本構模型的重要參考狀態[1-2]．

土體均勻趨于臨界狀態的失效模式只是一種理想的情形,實際土體破壞往往伴隨著剪切帶的發展．這種集中于帶狀區域的應變局部化的形成和演化是認識土體結構失效與破壞的關鍵,也一直是土力學與顆粒材料力學所關注的研究課題．相對于土體結構的宏觀尺度,剪切帶的厚度是很小的,實驗所觀察到的剪切帶厚度約為土體顆粒平均粒徑的10～15倍左右[3-4]．因此,剪切帶兩側有限的相對剪切位移就可在帶內引起較大的剪切變形,使得帶內材料很快達到臨界狀態．

經典連續介質本構模型不包含任何與材料細觀結構相關的特征長度,預測的剪切帶理論厚度為零．因此,常規的土體本構模型只能預測剪切帶的產生[5-6],一般不能規范剪切帶的厚度并正確預測土體涉及剪切帶演化的后失效行為．在采用有限元等方法分析土體的后失效行為時,數值解具有網格尺寸依賴性．

為克服上述困難,可采用高階連續介質力學理論,如微極場理論[7]、高階梯度理論[8]等．高階連續介質力學模型允許引入材料細觀特征長度,可以在一定程度上刻畫材料中應力和變形在細觀尺度上的變化．其中Cosserat介質屬于微極連續介質力學理論中較簡單的一種, 通過在連續體中引入質點的轉動自由度和偶應力,可以在一定程度上描述材料細觀尺度上的應力和變形的變化,從而能恰當模擬剪切帶的發展．Mühlhaus和Vardoulakis[9]采用Cosserat介質彈塑性模型討論了土體中剪切帶厚度與顆粒平均粒徑的關系;Oda[10]研究了顆粒材料細觀結構、偶應力與剪切帶的關系;Huang等[11]基于Cosserat介質理論,將Gudehus[12]和Bauer[13]所建立的模擬砂土等顆粒土的亞塑性模本構型推廣為微極亞塑性模型,并用于砂土中局部化應變發展模式的數值分析[14]．

基于文獻[11]建立的微極亞塑性模型,Huang和其合作者[15]分析了顆粒材料在平面Couette剪切中,條形區域內的應變局部化的形成與演化機制,并且針對充分發展的平直剪切帶,推導得到了臨界狀態下的控制方程,但針對該控制方程并未得到完全解．本文主要目的是在前期工作[15]的基礎上,探討剪切帶在臨界狀態下的關鍵變量,即Cosserat轉動角速度所滿足的非線性常微分方程的主要特性及求解途徑,并給出問題的完全解．

1 微極亞塑性模型簡述

(1)

Cosserat介質中的應力張量σ和偶應力張量μ滿足平衡方程:

divσ+b=0, divμ-∶σ=0,

(2)

其中,b為體力矢量;grad代表矢量和張量的梯度運算,div代表散度運算;為置換符號;點乘積運算代表兩個張量相鄰的一對下標的縮并,雙點積運算代表兩個張量前后兩對對應下標的縮并．式(1)和(2)也表明Cosserat介質中的應力、偶應力、應變率和微曲率率張量在一般情況下是非對稱的．

Gudehus[12]和Bauer[13]以Cauchy應力σ和孔隙比e為狀態變量,建立了一個有效實用的臨界狀態亞塑性模型(G-B模型),能較好地模擬砂土等無黏性顆粒土的力學響應．為了能客觀模擬顆粒土中剪切帶的形成與發展,Huang等[11]基于 Cosserat介質理論將G-B模型推廣為微極亞塑性模型,用下列張量方程描述Cauchy應力及偶應力張量的客觀時間導數與應變率、微曲率率張量之間的關系:

(3)

(4)

(5)

(6)

式(5)和(6)實際為建立模型時預設的臨界狀態流動法則和強度條件[11]．

2 剪切帶內臨界狀態的數學描述

圖1 剪切帶脫離體及受力示意圖及坐標系選取Fig. 1 Schematic diagram for the shear band with applied forces and coordinates

(7)

根據上述應力狀態,剪切帶內臨界狀態強度條件(5)簡化為

(8)

(9)

這里引入了無量綱長度因子ξ=x2/lc,并記θ′=dθ/dξ,其中

(10)

稱為材料的細觀特征長度,它與土體平均粒徑d50呈正比,并與材料細觀和宏觀強度參數之比相關．

將流動法則(9)代入平衡方程(7)的第三式,可得到關于θ的控制方程(參考附錄):

(11)

其中,ξd=d/lc為剪切帶厚度的無量綱參數或剪切帶的厚度因子;

(12)

3 剪切帶臨界狀態完全解

上節推導表明,微極亞塑性模型所模擬顆粒土中的剪切帶趨于臨界狀態(無剪脹流動狀態)時,帶內質點的Cosserat轉動角速度θ滿足非線性常微分方程(11)．本文的主要目標是通過求解方程(11)獲得剪切帶厚度因子及帶內各變量的分布規律．具體討論如下．

3.1 方程的特點

本文所要求解的控制方程(11)具有以下幾個主要特點:

2) 方程(11)關于基本變量θ具有齊次性, 若θ(ξ)是方程的解, 則對于任意常數A,Aθ(ξ)也是方程的解．

3) 剪切帶的基本特征是其內部材料處于流動狀態,外部材料作不變形的剛性運動,因此在剪切帶內部θ≠0,在邊界ξ=±ξd/2處應有θ=0．因此剪切帶的厚度因子ξd由θ的一個完整變化周期決定,并且與θ的絕對值無關．

3.2 問題的形式解

為求解方程(11),引入如下的變量替換:

(13)

容易驗證

考慮先在半帶厚(0≤ξ≤ξd/2)內求解θ和其余變量,再利用對稱性和反對稱性獲得完全解．根據圖1坐標系選擇及剪切帶中各變量正負號規定,應有θ≤0,θ′≥0,因此有y≤0,z≥0,從而

(14)

此外,在剪切帶中心和邊界處分別有

ξ=0:z=0,y=-1;ξ=ξd/2:y=0,z=1．

(15)

將替換變量(13)和(14)代入式(11),問題歸結為在半帶厚內求解z=z(ξ),滿足方程

(16)

中間變量z=z(ξ)應對式(16)直接積分獲得．另一中間變量y=y(ξ)可按式(14)求出．同時,利用ξ=ξd/2處邊界條件可得厚度因子表達式:

(17)

求得中間變量z(ξ)和y(ξ)后,利用流動法則(9)可求出正則化的應力和偶應力分量:

(18)

(19)

(20)

(21)

令ξ=ξd/2,可以得到剪切帶中心Cosserat轉動角速度與邊界處剪切速率的線性關系式:

(22)

3.3 參數的確定和問題的完全解

由于未能求出式(16)的解析積分,上述剪切帶中各變量的表達式僅是問題的形式解．本文中擬采用數值積分求出z～ξ關系及其余各變量的數值解．為此,需要確定方程(11)中參數β的具體數值．

(23)

根據平衡方程(7)的第三式知,(σ12-σ21)|ξ=0=(?μ32/?x2)|ξ=0>0．

(24)

上述分析雖然給出了參數β的范圍,但尚不能確定β的具體值．通過假定r0(或β)的不同值,按數值解可求得半剪切帶厚度因子ξd/2的對應值．結果表明(如圖3所示),剪切帶厚度因子ξd對參數r0(或β)具有強烈的依賴關系．

圖3 剪切帶厚度因子ξd/2對參數r0的依賴關系Fig. 3 Shear band thickness factor ξd/2 vs. parameter r0

(25)

(26)

(27)

此式為假定的參數β值提供了一個后驗條件．正確的結果應使β的后驗值與通過設定r0得到的設定值相等,據此確定參數從而獲得問題的完全解．

圖4 對應r0設定值參數β的先驗和后驗值Fig. 4 Priori and posteriori values of β vs. parameter r0

圖5 半帶厚內y(ξ)和z(ξ)Fig. 5 Variations of y(ξ)and z(ξ)in the half shear band

圖6 半帶厚內和Fig. 6 Variations of z(ξ) in the half shear band

圖7 剪切帶內正則化應力、偶應力分布Fig. 7 Distributions of the normalized stress and the couple stress components in the shear band

圖8 剪切帶內正則化應變率、微曲率率和剪切速率分布Fig. 8 Distributions of the normalized strain rate, the micro-curvature rate and the shear velocity in the shear band

3.4 剪切帶厚度因子的作用

(28)

4 結 論

剪切帶的厚度具有土體顆粒粒徑尺度相近的量級,帶內變形量的劇烈變化只有采用高階連續介質力學模型才能恰當模擬．采用筆者[11]先期提出的微極亞塑性模型模擬的理想顆粒土,剪切帶的臨界狀態(無剪脹流動狀態)可以通過一組流動法則、強度條件和關于Cosserat轉動角速度的一個非線性常微分方程描述．該方程的解析求解具有難度．筆者在本文中結合剪切帶的基本特征,討論了該微分方程的一些主要特點和關鍵參數的取值范圍,通過應用能量守恒關系補充完備了求解條件,從而能夠利用數值積分的方法求出剪切帶厚度因子及帶內應力、變形率各分量的分布規律等完全解．這一解答揭示了顆粒土失效時凝聚在小尺度范圍內的材料剪切流動機制,而剪切帶厚度因子與細觀特征長度對于本構模型中細觀強度參數的確定具有重要作用．此外,形如式(11)這樣有趣的非線性常微分方程也值得應用數學研究者們的關注．

附錄 方程(12)的推導

(A1)

(A2)

其中

(A3)

上式兩邊乘方后可求得R的下列表達式:

(A4)

代入式(A2)后得到關于θ的常微分方程(11)．