非飽和重塑黃土應(yīng)變控制等應(yīng)力比試驗三維離散元分析

李 濤，蔣明鏡，李立青

(1.中國民航大學 機場學院， 天津 300300； 2.天津大學 建筑工程學院， 天津 300072)

黃土主要分布在干旱和半干旱地區(qū)，黃土工程中遇到的土體多處于非飽和狀態(tài)，非飽和黃土試驗和理論研究一直是巖土工程界的熱門課題之一[1-3]。非飽和黃土的試驗研究主要針對等向壓縮、側(cè)限壓縮、常規(guī)三軸和真三軸應(yīng)力路徑等[4-7]，但實際工程中如黃土路基和填方工程[8]可能出現(xiàn)小主應(yīng)力和大主應(yīng)力共同增加的等應(yīng)力比路徑[9]。在等應(yīng)力比路徑下，隨著應(yīng)力比的變化，土體可能發(fā)生壓縮控制和剪切控制的轉(zhuǎn)換，難以通過常見應(yīng)力路徑試驗反映，因此，有學者通過等應(yīng)力比路徑研究黃土的力學性質(zhì)并應(yīng)用到本構(gòu)模型研究[9-11]。

黃土宏觀力學性質(zhì)與其微觀特性密不可分。為了研究黃土宏觀力學特性的微觀機制，不少學者開展掃描電鏡觀測、X射線計算機斷層掃描[12]等微觀試驗。為了分析微觀試驗難以全面獲取的接觸方向、接觸力等統(tǒng)計信息，Jiang等[13]采用二維離散單元法(DEM)，使用簡單膠結(jié)接觸模型模擬黃土粒間膠結(jié)，研究了等應(yīng)力比路徑下黃土的力學特性。由于應(yīng)力加載進行高應(yīng)力比壓縮試驗難以反映理想塑性和軟化，模擬采用完備的非飽和重塑黃土三維接觸模型[14]，對非飽和重塑黃土離散元試樣開展了不同應(yīng)力比的應(yīng)變加載等吸力等應(yīng)力比試驗，研究了黃土數(shù)值試樣的應(yīng)力、變形及其對應(yīng)的微觀力學性質(zhì)(平均配位情況、接觸法向和接觸力組構(gòu)等)。

1 黃土接觸模型及其微觀參數(shù)

黃土顆粒/團粒形狀不規(guī)則，表面附著細小黏土顆粒，在粒間接觸上表現(xiàn)出抗彎和抗扭能力，研究采用的黃土接觸模型[14]，基于Jiang等[15]的完整接觸模型，同時采用粒間引力考慮顆粒間的范德華力和毛細力作用。接觸模型已應(yīng)用于非飽和黃土三軸試驗[16]和一維壓縮試驗[17]離散元數(shù)值分析。

假設(shè)顆粒間通過圓形截面接觸，接觸半徑為[15]：

Rc=βR

(1)

式中：β為接觸半徑系數(shù)，R=2R1R2/(R1+R2)，R1和R2為兩顆粒半徑?？紤]粒間引力時，顆粒間力和彎矩的計算示意圖如圖1所示。接觸剛度可由顆粒等效模量和顆粒法切向剛度比計算。

圖1 接觸模型力學響應(yīng)示意圖

(2)

(3)

(4)

式中:Fn為法向接觸力;Fa為粒間引力;μ為粒間摩擦系數(shù);ζc為顆粒局部壓碎系數(shù)，可取ζc=2.1。

接觸模型參數(shù)如表1所示。粒間引力是影響非飽和重塑黃土宏觀力學性質(zhì)的主要微觀變量之一，對于非飽和黃土，范德華力與毛細力相比很小，粒間引力主要是由毛細水引起。粒間引力為[14]：

(5)

式中:d50為中值粒徑;s為基質(zhì)吸力;σv為范德華力系數(shù);ξa，ast和bst為擬合參數(shù)，ast=ca1exp(ca2e0),e0為試樣的初始孔隙比。

表1 離散元模擬微觀參數(shù)[14]

2 數(shù)值試樣制備與離散元實施

采用分層欠壓法分五層制備黃土離散元試樣，顆粒數(shù)目為42 180個，粒徑級配曲線如圖2所示，中值粒徑d50為20 μm。為了制備大孔隙比黃土數(shù)值試樣，在制樣過程中施加范德華力，大小按下式計算：

(6)

式中:σv取4 kPa。試樣生成后在12.5 kPa豎向壓力下側(cè)限壓縮穩(wěn)定(預(yù)壓)，預(yù)壓完成后孔隙比為0.934。

圖2 離散元試樣粒徑級配曲線

預(yù)壓完成后，對黃土離散元試樣上下墻體分別按10/s應(yīng)變率施加荷載，對四個側(cè)墻采用伺服控制，使試樣的應(yīng)力比(η=q/p，p為平均應(yīng)力，q為偏應(yīng)力)保持不變，以實現(xiàn)試樣的應(yīng)變控制等應(yīng)力比加載。應(yīng)變控制對比分級荷載施加方式容易實現(xiàn)較高應(yīng)力比下的等應(yīng)力比加載，因此為本文所采用。試驗過程中保持試樣的基質(zhì)吸力不變(100 kPa)。在數(shù)值分析過程中，隨著試樣調(diào)整，接觸可能消失和生成，新生成接觸自動施加接觸模型并賦予模型參數(shù)。

3 離散元試樣的宏觀性質(zhì)

圖3給出了非飽和重塑黃土離散元試樣壓縮曲線，通過與室內(nèi)試驗結(jié)果的對比可知，在半對數(shù)坐標中，非飽和黃土孔隙比先緩慢減小，在屈服應(yīng)力之后快速減小，離散元試樣的壓縮線趨勢與室內(nèi)試驗[5]相同。離散元孔隙比變化數(shù)值上與室內(nèi)試驗差別較大，主要是因為離散元在顆粒形狀、顆粒團簇破碎等方面與實際黃土顆粒不同造成的。關(guān)于應(yīng)力比的影響，離散元結(jié)果表明在一定的應(yīng)力比范圍內(nèi)(ηM的試樣在壓縮流塑階段出現(xiàn)了剪脹現(xiàn)象，且應(yīng)力比越高，試樣剪脹性越強。

圖3 等應(yīng)力比壓縮曲線

圖4給出了非飽和重塑黃土離散元試樣等應(yīng)力比試驗的偏應(yīng)力-偏應(yīng)變關(guān)系曲線，通過與室內(nèi)試驗結(jié)果對比可知，離散元試樣的偏應(yīng)力-偏應(yīng)變關(guān)系與室內(nèi)試驗[10]結(jié)果趨勢相同。在低應(yīng)力比下(ηM)，試樣主要表現(xiàn)出剪切特性，偏應(yīng)力-偏應(yīng)變關(guān)系類似常規(guī)三軸試驗，隨著偏應(yīng)變的增加，試樣發(fā)生剪切破壞，并伴隨剪縮或剪脹。

圖4 應(yīng)力-應(yīng)變曲線

4 離散元試樣的微觀性質(zhì)

為了探求非飽和重塑黃土宏觀力學性質(zhì)的微觀機制，對離散元試樣的配位情況和組構(gòu)情況進行了分析。

平均配位數(shù)(簡稱“配位數(shù)”)和平均力學配位數(shù)(簡稱“力學配位數(shù)”)可以用來分析土體顆粒的平均接觸情況，采用下式計算[18]：

Z=2Nc/Np

(7)

Zm=(2Nc-Np1)/(Np-Np1-Np0)

(8)

式中：Z為配位數(shù);Zm為力學配位數(shù);Nc為接觸數(shù)量;Np為顆粒數(shù)量;Np0和Np1分別是接觸數(shù)量為0和1的顆粒數(shù)量;Np-Np1-Np0為有效顆粒數(shù)量。

圖5給出了非飽和重塑黃土離散元試樣的平均配位情況，試樣的力學配位數(shù)在半對數(shù)坐標系中隨平均應(yīng)力的增加先緩慢增加后快速增加。通過與壓縮曲線的對比可知，在200 kPa平均應(yīng)力以下，隨著應(yīng)力的增加，配位數(shù)、力學配位數(shù)、有效顆粒數(shù)量發(fā)生了較明顯增加，但試樣孔隙比變化較小。即試樣的顆粒重排列和變形僅發(fā)生微調(diào)整，滯后于配位情況的變化。

圖5 配位性質(zhì)

不同應(yīng)力比離散元試樣的力學配位數(shù)隨平均應(yīng)力變化曲線保持很好的歸一性，與壓縮曲線表現(xiàn)相對應(yīng)。但是配位數(shù)曲線有所不同，除等向壓縮試樣(η=0)外，隨應(yīng)力比的增加，配位數(shù)隨平均應(yīng)力增加的速率變大，相同平均應(yīng)力下的有效顆粒數(shù)量更大?？梢姡m然某平均應(yīng)力下不同應(yīng)力比試樣的力學配位數(shù)相同，但配位情況卻有一定不同，高應(yīng)力比試樣中有更多的顆粒參與了傳力(顆粒排列調(diào)整更充分)。高應(yīng)力比(η>M)試樣在進入流塑階段尾段后，有效顆粒數(shù)量接近最大，試樣力學配位數(shù)基本保持不變，與之對應(yīng)的，試樣的孔隙比、平均應(yīng)力和偏應(yīng)力都基本維持不變，試樣達到臨界狀態(tài)。可見，臨界狀態(tài)對應(yīng)的宏微觀變量都基本保持不變，且有效顆粒數(shù)量接近總顆粒數(shù)量(顆粒均參與傳力)，與此對應(yīng)，顆粒排列經(jīng)歷了充分調(diào)整。

組構(gòu)張量可以用來反映矢量集合方向(接觸法向、接觸力方向)的統(tǒng)計規(guī)律，組構(gòu)張量計算如下[19]：

(9)

式中:ni為矢量的方向余弦。偏組構(gòu)Δ可用來表示矢量集合方向的各向異性程度：

(10)

式中:F1、F2和F3是組構(gòu)張量的三個不變量。

圖6給出了非飽和重塑黃土離散元試樣的接觸方向偏組構(gòu)隨平均應(yīng)力的變化，試樣初始接觸方向偏組構(gòu)是原位應(yīng)力狀態(tài)下的偏組構(gòu)(預(yù)壓后)，接觸方向偏組構(gòu)的變化明顯受試驗應(yīng)力比的影響。隨著平均應(yīng)力的增加，側(cè)限壓縮試驗的應(yīng)力比不斷下降(如圖6(b))，但側(cè)限壓縮的接觸方向偏組構(gòu)基本保持不變(高平均應(yīng)力下輕微下降)。在等應(yīng)力比壓縮中，η=0.00和η=0.75試樣的應(yīng)力比始終小于側(cè)限壓縮試驗的應(yīng)力比，接觸方向偏組構(gòu)隨平均應(yīng)力的增加而不斷減小，應(yīng)力比越小，接觸方向偏組構(gòu)隨平均應(yīng)力減小越快； 試樣的應(yīng)力比先是小于而后大于側(cè)限壓縮試驗的應(yīng)力比，接觸方向偏組構(gòu)隨平均應(yīng)力的增加先減小后增加，應(yīng)力比越大，接觸方向偏組構(gòu)隨平均應(yīng)力增加越快?？梢?，試樣接觸方向偏組構(gòu)的變化與其初始狀態(tài)和加載應(yīng)力路徑關(guān)系密切。

圖7給出了η=0.00和η=1.25試樣接觸方向組構(gòu)玫瑰圖(軸向為各個1°徑向范圍內(nèi)接觸的數(shù)量)的變化。隨著加載的進行，η=0.00試樣各個方向接觸數(shù)量均增加，但水平傾向接觸數(shù)量增加的更快，試樣接觸方向各向異性程度減小。η=1.25試樣總接觸數(shù)量有明顯增加，但水平傾向接觸數(shù)量基本不變，試樣接觸方向各向異性程度增加。

圖8給出了非飽和重塑黃土離散元試樣接觸力偏組構(gòu)隨平均應(yīng)力的變化。施加毛細力后，試樣初始接觸力偏組構(gòu)(0.239)遠大于初始接觸方向偏組構(gòu)(0.074)，說明接觸力的初始各向異性程度大于接觸方向。η≤1.25試樣的總接觸力、強接觸力(大于平均接觸力)和弱接觸力(小于平均接觸力)偏組構(gòu)均表現(xiàn)為隨平均應(yīng)力的增加先下降再增加；η≥1.50試樣接觸力偏組構(gòu)(總、強、弱)隨平均應(yīng)力增加較快，最終達到較大的各向異性程度；側(cè)限壓縮試驗接觸力偏組構(gòu)(總、強、弱)經(jīng)歷了“緩慢下降-快速下降-緩慢變化”的過程。在整個等應(yīng)力比加載試驗(不含等向壓縮試樣)過程中，接觸力偏組構(gòu)始終大于接觸方向偏組構(gòu)；強接觸力偏組構(gòu)大于弱接觸力偏組構(gòu)，即強接觸力各向異性程度很高，對承擔豎向荷載貢獻大。

圖6 接觸方向組構(gòu)特性

圖7 接觸方向組構(gòu)玫瑰圖

圖9給出了接觸力組構(gòu)玫瑰圖，在加載開始時接觸力具有明顯的各向異性；隨著加載的進行，η=0.00試樣水平傾向接觸力數(shù)量明顯增加，豎直傾向接觸力數(shù)量略微減小，各向異性程度明顯降低；η=1.25試樣水平傾向接觸力數(shù)量基本不變，豎直傾向接觸力數(shù)量有一定的增加，各向異性程度增加。

圖8 接觸力組構(gòu)特性

圖9 接觸力組構(gòu)玫瑰圖

離散元試樣反映了實際黃土的范德華力、毛細力、顆粒間摩擦和轉(zhuǎn)動抗力等特征，但對顆粒級配、顆粒形狀進行了簡化，忽略了顆粒破碎。因此，離散元模擬僅能定性分析非飽和重塑黃土配位數(shù)和組構(gòu)的變化規(guī)律，配位數(shù)和偏組構(gòu)數(shù)值會與實際黃土有一定的不同，模擬的定性結(jié)果與微觀觀測結(jié)果[20]是協(xié)調(diào)一致的。

5 結(jié) 論

(1) 離散元模擬結(jié)果能定性再現(xiàn)室內(nèi)等應(yīng)力比試驗的壓縮曲線和剪切特性。在低應(yīng)力比下(ηM)試樣主要表現(xiàn)出剪切特性，可能表現(xiàn)出剪脹性。

(2) 不同應(yīng)力比試樣的力學配位數(shù)在半對數(shù)坐標系中隨平均應(yīng)力的增加先緩慢增加后快速增加，不同應(yīng)力比曲線具有很好的歸一性；高應(yīng)力比試樣比低應(yīng)力比試樣有更多的顆粒參與了傳力；臨界狀態(tài)對應(yīng)的宏微觀變量都基本保持不變。

(3) 試樣接觸方向偏組構(gòu)的變化與其初始狀態(tài)和加載應(yīng)力比關(guān)系密切；側(cè)限壓縮的接觸方向偏組構(gòu)基本保持不變。在加載過程中(等向壓縮試樣除外)，接觸力偏組構(gòu)始終大于接觸方向偏組構(gòu)；強接觸力偏組構(gòu)大于弱接觸力偏組構(gòu)。