衰變熱源作用下飽和多孔介質(zhì)熱固結問題的擴展精細積分法1)

王路君艾智勇

?(同濟大學地下建筑與工程系,巖土及地下工程教育部重點實驗室，上海200092)

?(浙江大學巖土工程研究所,軟弱土與環(huán)境土工教育部重點實驗室，杭州310058)

衰變熱源作用下飽和多孔介質(zhì)熱固結問題的擴展精細積分法1)

王路君?,?,2)艾智勇?,3)

?(同濟大學地下建筑與工程系,巖土及地下工程教育部重點實驗室，上海200092)

?(浙江大學巖土工程研究所,軟弱土與環(huán)境土工教育部重點實驗室，杭州310058)

熱源作用下飽和多孔介質(zhì)熱固結效應是土木及能源工程領域的一個重要課題.由于問題的復雜性，已有的研究大多將介質(zhì)假定為均勻各向同性，且將熱源假定為恒定強度.實際工程中，天然飽和多孔介質(zhì)常表現(xiàn)出明顯的分層特性，熱源強度也存在衰變性，為此本工作采用擴展精細積分法對衰變熱源作用下層狀飽和多孔介質(zhì)的熱固結問題進行研究.借助于積分變換，將飽和多孔介質(zhì)熱固結問題的偏微分方程轉化為變換域內(nèi)的常微分方程；然后對飽和多孔介質(zhì)微層元進行合并消元，并結合邊界條件，推導出衰變熱源作用下層狀飽和多孔介質(zhì)熱固結問題在積分變換域內(nèi)的擴展精細積分解；對所得解答進行相應的數(shù)值積分逆變換，可獲得所求溫度、超靜孔壓及豎向位移在物理域內(nèi)的解答.基于上述求解過程，編制相應的計算程序進行數(shù)值計算，通過與已有文獻對比，驗證本文擴展精細積分法在求解層狀飽和多孔介質(zhì)熱固結問題中的適應性和正確性；最后通過幾組算例，分析熱源衰變周期、熱源埋深及介質(zhì)的成層性對熱固結效應的影響.結果表明：熱源衰變周期對溫度和超靜孔壓的峰值、以及達到峰值的時間均有明顯影響，衰變周期越長，二者峰值均越大，且達到峰值所需時間越長；熱源埋深對超靜孔壓及豎向位移變化影響顯著，深埋熱源作用時熱源兩側豎向位移呈對稱分布，而淺埋熱源兩側則無此現(xiàn)象；飽和多孔介質(zhì)的分層特性對熱固結效應影響明顯.

熱源，飽和多孔介質(zhì)，熱固結，擴展精細積分法

引言

熱源作用下飽和多孔介質(zhì)熱固結問題是土木、能源、環(huán)境等領域的熱點課題之一.飽和多孔介質(zhì)中，熱源作用會引起孔隙水壓力的產(chǎn)生及消散，而孔隙水壓力的變化會對固體骨架變形及強度產(chǎn)生重要影響.因此，熱固結問題是滲流場、溫度場和應力應變場相互耦合作用的結果，該課題在放射性核廢料地下處置、地熱資源開發(fā)、能量樁等工程中具有重要的研究和應用價值[1-5].在Biot固結理論[6]的基礎上，Biot[7]建立了飽和多孔介質(zhì)熱–水–力耦合問題的基本控制方程，并對方程中的耦合系數(shù)進行了明確的定義和完整的物理解釋，較早開展了飽和多孔介質(zhì)熱固結響應的研究.由于Biot[7]所建立的熱–水–力耦合理論是位移、孔壓及溫度變化耦合的復雜偏微分方程，故在求解上遇到了一定的困難.因此，國內(nèi)外諸多學者對該理論進行簡化和修正，并在求解方法上進行了諸多探索.Booker和Savvidou[8]，Savvidou和Booker[9]分析了飽和巖土介質(zhì)熱固結問題的作用機理，并給出了球形和點熱源作用下熱固結問題的解析解.Mctigue[10]對熱源作用于半無限多孔飽和介質(zhì)表面的熱–水–力耦合問題進行了研究，并得到了其解析解.Bai和Abousleiman[11]基于熱–水–力耦合理論，探討了完全耦合、部分耦合及完全非耦合三種理論的適用條件，并給出了三種情況下的解析解答.白冰等[4,12]基于飽和多孔介質(zhì)熱–水–力耦合控制方程，借助于積分變換對熱彈性固結模型進行求解，并給出溫度、超靜孔壓和位移演化過程的解析式.鄭榮躍等[13]研究了半無限地基在內(nèi)置點熱力源作用下的響應問題，并給出熱力源作用下應力、位移、孔隙水壓力的解.吳瑞潛等[14]建立了變載荷作用下飽和土體一維熱固結問題的解析解.Lu和Lin[15]對衰變熱源作用下多孔彈性半空間的熱固結問題進行了研究，并給出了解析解.Selvadurai和Suvorov[16]對固體骨架為Hooke彈性或彈塑性體的飽和介質(zhì)進行研究，并分析了不同條件下的熱–水–力耦合響應.

天然介質(zhì)經(jīng)過長期的沉積過程，往往表現(xiàn)出明顯的分層特性.為更真實地描述熱固結過程，部分學者開展了層狀介質(zhì)熱固結問題的研究.Giraud等[17]推導了兩層介質(zhì)一維熱固結問題的半解析解答.白冰[18]采用解析方法，對變溫載荷作用時的雙層半無限飽和介質(zhì)的熱固結問題進行了研究.Ai和Wang[5]采用解析層元法求解了層狀飽和介質(zhì)軸對稱熱固結問題.目前針對層狀彈性體系，常用的求解方法有傳遞矩陣法[19-21]、有限層法[22-23]、剛度矩陣法或解析層元法[5,24-26]等.對于層狀彈性體及飽和地基等本構方程較為簡單的情況，采用上述方法往往可推導出其顯式解析解.但對于較復雜的求解模型，如橫觀各向同性體及熱–水–力耦合問題，其控制方程較復雜，極難采用上述方法推導出其顯式解析解.另一方面，目前主流的數(shù)值計算方法，如有限元法等，在求解熱固結問題時需耗費大量時間和內(nèi)存，且需較強的經(jīng)驗性.鑒于此，本文擬基于具有穩(wěn)定性好、計算效率與精度高特點的精細積分法[27-29]，并結合積分變換，建立衰變熱源作用下層狀飽和多孔介質(zhì)熱固結問題的擴展精細積分解.與文獻[28]所給出的精細積分法相比，本文方法不僅可以求解表面受荷問題，而且還可以求解層狀介質(zhì)內(nèi)部任意位置的受荷問題；另外，本文方法的整個求解過程均在積分變換域內(nèi)進行，是精細積分法[27-28[30-31].由于本文方法是精細積分法[28]應用的一個擴展，為此本文將其稱之為“擴展精細積分法”.首先推導出熱固結問題在Laplace-Hankel變換域內(nèi)的常微分矩陣方程；然后對介質(zhì)微層元進行合并消元，得到該問題在變換域內(nèi)的擴展精細積分解；最后對該解進行相應的積分逆變換，可得其在物理域內(nèi)的解答.

1 熱固結問題的常微分矩陣方程

軸對稱條件下，不考慮體力時，飽和多孔介質(zhì)熱固結問題的平衡方程為

式中，σr,σ?,σz分別為r,?,z方向的正應力，σrz為r-z面上的剪應力.

結合廣義熱彈性Hooke定律和有效應力原理，可得到用位移和溫度增量表示的本構關系

熱量在介質(zhì)內(nèi)的傳導符合Fourier導熱定律，則0至t時間內(nèi)，z方向的熱流量Qθ可定義為

式中，K為熱傳導系數(shù).

根據(jù)Fourier定律與能量守恒方程，熱傳導方程可表示為[9]

式中，κ=K/m，m=ρC，ρ和C分別表示密度和比熱容，

據(jù)Darcy定律，初始時刻至t時刻豎向流量為

式中，c為滲透系數(shù).

結合Darcy定律和滲流連續(xù)條件，滲流連續(xù)方程可表示為

式中，αu=3αs(n-1)-3nαf，αs和αf分別為固體顆粒和孔隙水的線膨脹系數(shù)，n表示孔隙率.

由于直接求解上述偏微分方程較困難，本文借助于積分變換方法，將偏微分方程轉化為易于求解的常微分方程.函數(shù)f(r,z,t)關于時間t的Laplace變換及其逆變換為[32]

式中，s表示關于t的Laplace變換參數(shù)，

式中，ξ是關于r的Hankel變換參數(shù)，Jm(ξr)為m階Bessel函數(shù).

對式 (1)～式 (6)進行關于時間t和坐標r的Laplace及Hankel變換，所得常微分方程表示為矩陣形式

2 層狀體系的精細積分法

對于兩點邊值問題的精細積分法，文獻[28]給出了嚴格的理論推導.對于z向任意微層元[za,zb]，其上下界面狀態(tài)向量如圖1.

圖1 微層元示意圖Fig.1 Diagram of an arbitrary layer element

對于式(9)所描述的兩端邊值問題，其兩端狀態(tài)向量存在如下關系[28]

式中，F(xiàn),G,Q,E為待求關系矩陣，它們建立了微層元兩端位移和應力向量間的關系，均為關于za，zb的函數(shù).若層元厚度?=?z=zb-za非常小，則關系矩陣F,G,Q,E可進行Taylor級數(shù)展開，此時其表達式中僅包含層元厚度?及矩陣Φi(i=1,2,3,4)，而矩陣Φi只與介質(zhì)的材料參數(shù)有關，即層元厚度?和介質(zhì)材料參數(shù)一經(jīng)確定，則關系矩陣F,G,Q,E可推導求出，具體可參考文獻[34].

式(10)建立了任意微層元上下界面狀態(tài)變量之間的關系.對于任意兩相鄰微層元1和2，如圖2所示.

圖2 相鄰微層元的合并Fig.2 Combination of adjacent layer elements

根據(jù)交界面處變量的協(xié)調(diào)條件，可以推導出za和zc界面上狀態(tài)變量間的關系式[34]

式中，F(xiàn)3,G3,Q3,E3為合并后形成的微層元3的關系矩陣，其具體表達式為

由層元合并消元次序無關定理[28]知，式 (12)可看作微層元消元合并的一個遞歸表達式.若將微層元3與其相鄰的微層元4繼續(xù)合并形成微層元5，并將微層元序號3,4,5依次用1,2,3替換，則微層元5上下界面處狀態(tài)向量間關系仍可用式(11)和式(12)表述.可見，對于多層介質(zhì)中的任一層元，將其劃分為多個微層元，再進行微層元間的合并操作，最終可建立該層上下表面狀態(tài)變量間的關系.

由以上推導可知，只要?足夠小，則數(shù)值計算結果的誤差是微小的.但當?過度小時，計算時可能會因計算機存儲精度問題而導致有效位數(shù)丟失.為避免該類問題，遞歸式(12)需采用如下表達式

式中，J=(I+G1Q2)-1，K=(I+Q2G1)-1.式(12)與式(13)實際上是相同的，只是將原來的F,E分別表示為I與F#,E#之和的形式，這是由于?很小時，F(xiàn)#,E#與I相比為極小的矩陣.若直接用F和E進行計算，F(xiàn)#和E#可能會因計算機存儲精度問題被消去而得不到精確結果.故遞歸操作時僅涉及F#和E#運算，以避免發(fā)生精度丟失的問題.

3 內(nèi)部受荷問題的擴展精細積分法

如圖3所示，對于內(nèi)部作用有熱源的多層介質(zhì)體系，根據(jù)熱源作用深度HF和計算點深度HC，可將多層體系分成三個部分，為便于表述，將其稱之為“層塊”.三個層塊分別為：層塊1[a,b]，介質(zhì)體系表面與熱源作用(或計算點)深度間的層塊；層塊2[b,c]，熱源作用深度與計算點間的層塊；層塊3[c,d]，計算點(或熱源作用)深度與介質(zhì)體系底面間的層塊.各層塊可能是單個自然層，也可以是某個自然層的一部分，或是由多個自然層組合而成.如圖中虛線所示，各層塊內(nèi)部微層元采用式(12)或式(13)進行消元凝聚，得到各自的關系矩陣，分別記為根據(jù)HF和HC的關系，可分不同情況進行討論.

圖3 層狀體系的三個層塊Fig.3 Three blocks of a layered system

當HF＞HC時，各層塊的關系矩陣方程如下：

層塊1

界面c處可能作用有廣義載荷P?；驈V義位錯PΛ，此時有如下關系

結合式(15)～式(17)，消去c界面的狀態(tài)向量可得

由式(14)和式(18)，可得計算點界面b處廣義應力和位移的表達式

當HF＜HC時，采用類似的推導方法，可建立計算點界面c處應力和位移的表達式

上述推導建立了廣義載荷作用于層狀介質(zhì)體系內(nèi)部時的一般性解答，對于HF=HC,HC=0或HF=0等工況，屬于以上解答的特例.計算時只需根據(jù)相應工況，將三層塊中不存在層塊的關系矩陣分別用I,0,0,I代替，而后采用上述方法，可實現(xiàn)層狀體系特殊受荷時的解答.

當層狀體系內(nèi)部分別作用有恒定總強度為Q0的點熱源和半徑為的圓形熱源時，有

若內(nèi)部作用的熱源為衰變熱源時，假定其初始強度為Q0，其強度隨時間衰減Q=Q0e-γt，其中γ=ln2/t′，t′為熱源的半衰期，此時有

結合式(21)～式(24)，可得到所求問題在積分變換域內(nèi)的解答.

4 數(shù)值計算與分析

第3節(jié)得到了問題在變換域內(nèi)的解答，對其進行相應的積分逆變換，可求得其物理域內(nèi)的最終解.本文Laplace數(shù)值逆變換采用FT[35]法

由式(25)可知，g(t)的求解表達式中僅含有一個自由的參數(shù)M，即累加求和項數(shù).為控制g(t)求解計算中的舍入誤差，參考文獻[36]的建議，對精度要求需定義M，即所需計算精度的有效位數(shù).

因此，F(xiàn)T法的求解可概述為：確定積分變換式，并指定計算參數(shù)t和M的值.首先設定與計算精度有關的參數(shù)M(參考文獻[30,36]的建議，本文M值取為10)，隨后根據(jù)式(25)和式(26)計算g(t,M)的值.參數(shù)M對計算精度及效率的影響分析可參考文獻[30,36].

本文Hankel數(shù)值逆變換借鑒文獻[37]的方法實現(xiàn).考慮到Bessel函數(shù)是震蕩衰減函數(shù)，Gauss積分點或正或負，為保證數(shù)值計算的穩(wěn)定性，計算時采用零點分段的方法.首先，選取Bessel函數(shù)兩相鄰零點為一積分區(qū)段，將原Hankel逆變換的半無限積分區(qū)間劃分為多個積分區(qū)段；然后，每個區(qū)段采用Gauss-Legendre法進行數(shù)值積分計算；最后，將每個區(qū)段的積分結果進行疊加可得到最終的積分值.此法在理論上可表達為：

將N個子積分區(qū)段的結果進行疊加，可得整個區(qū)間的積分值

4.1 驗證

為驗證本文方法及程序的正確性，本文對點熱源作用下飽和多孔介質(zhì)熱固結問題進行計算，并將結果與文獻[9]進行對比，如圖4所示.點熱源強度為Q0且保持不變，熱源埋置于介質(zhì)的深部.無量綱因子定義見圖4，由圖可知，超靜孔壓隨著時間的發(fā)展先逐漸增大，達到一峰值后而逐漸消散.另外，對比圖中結果可知，本文結果與文獻[9]的成果在各時刻均顯示出較高的吻合度，從而證明本文方法在求解該類問題中的適用性.

圖4 點熱源作用時超靜孔壓隨時間變化Fig.4 Variation of excess pore pressure with time due to a point heat source

4.2 熱源衰變周期的影響分析

圖5(a)和圖5(b)展示了單層介質(zhì)熱固結問題中熱源衰變周期對溫度和超靜孔壓的影響.計算點位于z軸，與點熱源間距為b.泊松比μ=0.3保持不變，參數(shù)設置如圖5(a)所示，無量綱因子與4.1節(jié)設置相同.由圖5(a)可知，衰變周期對溫度峰值及達到峰值所需時間影響明顯；衰變周期越大，溫度峰值越大，且達到峰值所需時間越長；當熱源為無衰變熱源時(t0=∞)，溫度隨時間的發(fā)展而增大，最終趨于穩(wěn)定.由圖5(b)可知，衰變周期越長，超靜孔壓峰值越大且達到峰值所需時間越長；衰變周期較短時會有負超靜孔壓出現(xiàn)，這與熱源的迅速冷卻有關.

圖5 衰變周期對熱固結效應的影響Fig.5 E ff ect of the heat source’s half-life on thermal consolidation

4.3 熱源埋深的影響分析

本節(jié)通過兩組算例研究熱源埋深對多孔飽和介質(zhì)熱固結響應的影響.圖6(a)和圖6(b)分別展示了點熱源埋深h對超靜孔壓及豎向位移的影響，假定介質(zhì)厚度為H，且H=10h0，h0為常數(shù)，泊松比μ=0.25，無量綱因子設置仍與4.2節(jié)相同.熱源埋深考慮三種情況，分別為h=h0，h=2h0和h=5h0.由圖6(a)可知，隨著計算點與熱源點間距離的增大，超靜孔壓峰值逐漸減小，且達到峰值時間逐漸延長.從圖6(b)可以發(fā)現(xiàn)，熱源埋置較深時，表面豎向位移峰值較大，但熱源埋置深度對豎向位移的穩(wěn)定值影響不大.

圖6 熱源埋深對超靜孔壓及豎向位移的影響Fig.6 E ff ect of the heat source’s buried depth on the variation of excess pore pressure and vertical displacement

接下來分析深埋和淺埋熱源作用下豎向位移隨時間的變化規(guī)律，具體如圖7.熱源為圓形熱源，其半徑和總強度分別為r0和Q0，深埋和淺埋深度分別為1 000r0和10r0，無量綱因子見圖7(a).由圖知，深埋熱源上下兩側等距離處的位移變化曲線呈對稱分布，并分別向兩側發(fā)生膨脹變形.當τ＜1時，淺埋熱源兩側位移變化基本呈對稱分布，這是由于在較短時間內(nèi)，熱量僅傳導至周圍較小的區(qū)域；隨著時間的發(fā)展，熱量繼續(xù)向更深處傳遞，由于為淺埋熱源，熱量傳遞至介質(zhì)表面而使其上部介質(zhì)膨脹量不再產(chǎn)生變化，熱源下部由于為半空間而導致豎向位移繼續(xù)變化，下側介質(zhì)不斷膨脹使熱源及計算點位移逐漸向上發(fā)展，從而引起圖7(b)中現(xiàn)象的產(chǎn)生.

圖7 不同埋深條件下豎向位移隨時間的變化曲線Fig.7 Evolutions of vertical displacement against time with di ff erent conditions of buried depths

4.4 介質(zhì)成層性的影響分析

實際工程中的天然介質(zhì)往往呈層狀分布，本節(jié)以四層飽和多孔半空間為例分析層狀特性對熱固結響應的影響.如圖8(a)所示，1～3層為有限深度層，第4層為半空間，選5種不同工況進行分析，如表1所示.

圖8 不同時刻多層介質(zhì)超靜孔壓沿深度分布Fig.8 Distribution of excess pore pressure along depth of multilayered medium at di ff erent time

圖8 不同時刻多層介質(zhì)超靜孔壓沿深度分布(續(xù))Fig.8 Distribution of excess pore pressure along depth of multilayered medium at di ff erent time(continued)

表1 多層飽和半空間計算參數(shù)Table 1 Parameters of the multilayered half-space

其他參數(shù)有如下關系：μi=0.25(i=1,2,3,4)，?h1:?h2:?h3:h=5:3:2:10，G1:G2:G3:G4=8:5: 2:1，h為點熱源的埋置深度，其強度為Q0，其他無量綱參數(shù)設置見圖8(d).

圖8展示了不同時刻z軸上不同深度計算點的超靜孔壓分布.由圖可知，在較早時刻 (τ=0.02,0.05)，超靜孔壓主要產(chǎn)生于熱源周圍區(qū)域，且距離熱源越近，其值越大；由于在不同介質(zhì)層分界面處，介質(zhì)參數(shù)出現(xiàn)突變，導致分布曲線在分界面處出現(xiàn)明顯的折點；隨著時間的延長，超靜孔壓的區(qū)域范圍逐漸擴大，且分界面處的曲線逐漸趨于平滑.

5 結論

本文在Laplace-Hankel變換域內(nèi)推導出飽和多孔介質(zhì)熱固結問題的常微分控制方程，并借助于層狀體系任意深度作用載荷時的精細積分法和數(shù)值逆變換技術，得到了衰變熱源作用下層狀飽和多孔介質(zhì)熱固結問題的解答.通過與已有文獻對比，驗證了本文方法在求解熱固結問題時的正確性.最后，通過算例分析了熱源衰變周期、熱源埋深及介質(zhì)成層性對熱固結效應的影響，結果表明：

(1)熱源衰變周期對溫度和超靜孔壓的峰值、達到峰值的時間均有明顯影響；衰變周期越大，二者峰值均越大，且達到峰值所需時間越長.

(2)熱源埋深對超靜孔壓及豎向位移變化影響顯著；深埋熱源作用時熱源兩側豎向位移呈對稱分布，淺埋熱源兩側位移曲線則無此現(xiàn)象.

(3)本文方法能有效求解多層介質(zhì)的熱固結問題，介質(zhì)的分層特性對熱固結效應影響明顯.

1 Delage P,Sultan N,Cui YJ.On the thermal consolidation of Boom clay.Canadian Geotechnical Journal,2000,37(2):343-354

2 王鐵行,李寧,謝定義.土體水熱力耦合問題研究意義、現(xiàn)狀及建議.巖土力學,2005,26(3):488-493(Wang Tiehang,Li Ning, Xie Dingyi.Necessity and means in research on soil coupled heatmoisture-stress issues.Rock and Soil Mechanics,2005,26(3):488-493(in Chinese))

3 蔣中明,Dashnor H.核廢料貯存庫圍巖體熱響應耦合場研究.巖土工程學報,2006,28(8):953-956(Jiang Zhongming,Dashnor H. Studies on coupled fielof thermal response in rock mass of nuclear waste repository.Chinese Journal of Geotechnical Engineering,2006,28(8):953-956(in Chinese))

4 Bai B,Guo LJ,Han S.Pore pressure and consolidation of saturated silty clay induced by progressively heating/cooling.Mechanics of Materials,2014,75:84-94

5 Ai ZY,Wang LJ.Axisymmetric thermal consolidation of multilayered porous thermoelastic media due to a heat source.International Journal for Numerical and Analytical Methods in Geomechanics, 2015,39(17):1912-1931

6 Biot MA.General theory of three-dimensional consolidation.Journal of Applied Physics,1941,12(2):155-164

7 Biot MA.Thermoelasticity and irreversible thermodynamics.Journal of Applied Physics,1956,27(3):240-253

8 Booker JR,Savvidou C.Consolidation around a spherical heat source.International Journal of Solids and Structures,1984,20(11-12):1079-1090

9 Savvidou C,Booker JR.Consolidation around a heat source buried deep in a porous thermoelastic medium with anisotropic fl w properties.International Journal for Numerical and Analytical Methods in Geomechanics,1989,13(1):75-90

10 Mctigue DF.Thermoelastic response of fluid-saturate porous rock.JournalofGeophysicalResearchAtmospheres,1986,91(B9):9533-9542

11 Bai M,Abousleiman Y.Thermoporoelastic coupling with application to consolidation.International Journal for Numerical and Analytical Methods in Geomechanics,1997,21(2):121-132

12 白冰.巖土介質(zhì)非穩(wěn)態(tài)熱固結耦合問題的熱源函數(shù)法.力學學報, 2004,36(4):427-434(Bai Bing.Heat source function method for coupling analyses of thermal consolidation in saturated soil.Acta Mechanica Sinica,2004,36(4):427-434(in Chinese))

13 鄭榮躍,劉干斌,梧松.半空間飽和土內(nèi)置點載荷作用下的熱彈性波動.力學學報,2008,40(3):413-420(Zhen Rongyue,Liu Ganbin, Wu Song.Coupling thermo-hydro-mechanical dynamic response of saturated soil subjected to internal excitation.Chinese Journal of Theoretical and Applied Mechanics,2008,40(3):413-420(in Chinese))

14 吳瑞潛,謝康和,程永鋒.變載荷下飽和土一維熱固結解析理論.浙江大學學報(工學版),2009,43(8):1532-1537(Wu Ruiqian,Xie Kanghe,Cheng Yongfeng.Analytical theory for one-dimensional thermal consolidation of saturated soil under time-dependent loading.Journal of Zhejiang University(Engineering Science),2009, 43(8):1532-1537(in Chinese))

15 Lu JCC,Lin F.Thermal consolidation of a poroelastic full space subjected to a decaying point heat source//Proceedings of the 2nd International ISCM Symposium and the 12nd International EPMESC Conference,2010:407-412

16 Selvadurai APS,Suvorov AP.Thermo-poromechanics of a fluid fille cavity in a fluid-saturate geomaterial//Proceedings of the royal society a mathematical physical and engineering sciences, 2014,470(2163):20130634

17 Giraud A,Homand F,Rousset G.Thermoelastic and thermoplastic response of a double-layer porous space containing a decaying heat source.International Journal for Numerical and Analytical Methods in Geomechanics,1998,22(2):133-149

18 白冰.變溫度載荷作用下半無限成層飽和介質(zhì)的熱固結分析.應用數(shù)學和力學,2006,27(11):1341-1348(Bai Bing.Thermal consolidation of layered porous half-space to variable thermal loading.Applied Mathematics and Mechanics,2006,27(11):1341-1348(in Chinese))

19 Yue ZQ,Yin JH.Backward transfer-matrix method for elastic analysis of layered solids with imperfect bonding.Journal of Elasticity, 1998,50(2):109-128

20 艾智勇,吳超.滲透各向異性可壓縮地基固結的平面應變分析.力學學報,2009,41(5):801-807(Ai Zhiyong,Wu Chao.Analysis on plane strain consolidation of a multi-layered soil with anisotropic permeability and compressbility constituents.Chinese Journal of Theoretical and Applied Mechanics,2009,41(5):801-807(in Chinese))

21 趙宇昕,陳少林.關于傳遞矩陣法分析飽和成層介質(zhì)響應問題的討論.力學學報,2016,48(5):1145-1158(Zhao Yunxin,Chen Shaolin.Discussion on the matrix propagator method to analyze the response of saturated layered media.Chinese Journal of Theoretical and Applied Mechanics,2016,48(5):1145-1158(in Chinese))

22 Booker JR,Small JC.Finite layer analysis of consolidation.I.International Journal for Numerical and Analytical Methods in Geomechanics,1982,6(2):151-171

23 宰金珉,梅國雄.有限層法求解三維比奧固結問題.巖土工程學報,2002,24(1):31-33(Zai Jinmin,Mei Guoxiong.Finite layer analysis of three dimensional Biot consolidation.Chinese Journal of Geotechnical Engineering,2002,24(1):31-33(in Chinese))

24 鐘陽,耿立濤.多層彈性平面問題解的精確剛度矩陣法.巖土力學,2008,29(10):2829-2832(Zhong Yang,Geng Litao.Explicit solution of multiplayer elastic plane by exact sti ff ness matrix method.Rock and Soil Mechanics,2008,29(10):2829-2832(in Chinese))

25 艾智勇,曹國軍,成怡沖.平面應變Biot固結的解析層元.力學學報,2012,44(2):401-407(Ai Zhiyong,Cao Guojun,Cheng Yichong.Analytical layer-element of plane strain Biot’s consolidation.ChineseJournalofTheoreticalandAppliedMechanics,2012,44(2): 401-407(in Chinese))

26 艾智勇,王路君,曾凱.穩(wěn)定溫度場下層狀路面體系的解析層元解.同濟大學學報 (自然科學版),2014,42(11):1665-1669(Ai Zhiyong,Wang Lujun,Zeng Kai.Analytical layer-element solution for layered pavement in stable temperature fieldJournal of Tongji University(Natural Science),2014,42(11):1665-1669(in Chinese))

27 鐘萬勰.結構動力方程的精細時程積分法.大連理工大學學報, 1994,34(2):131-136(Zhong,Wanxie.On precise time-integration method for structural dynamics.Journal of Dalian University of Technology,1994,34(2):131-136(in Chinese))

28 鐘萬勰.彈性力學求解新體系.大連：大連理工大學出版社，1995 (Zhong Wanxie.A New Systematic Methodology for Theory ofElasticity.Dalian:Dalian University of Technology Press,1995(in Chinese))

29 韓澤軍,林皋,李建波.二維層狀地基格林函數(shù)的求解.土木工程學報,2015,48(10):99-107(Han Zejun,Lin Gao,Li Jianbo.The solution of Green’s functions for two-dimensional layered ground.China Civil Engineering Journal,2015,48(10):99-107(in Chinese))

30 Ai ZY,Cheng YC.Extended precise integration method for consolidation of transversely isotropic poroelastic layered media.Computers&Mathematics with Applications,2014,68(12):1806-1818

31 Wang LJ,Ai ZY.Plane strain and three-dimensional analyses for thermo-mechanical behavior of multilayered transversely isotropic materials.International Journal of Mechanical Sciences,2015,103: 199-211

32 Talbot A.The accurate numerical inversion of Laplace transforms.Journal of Institute of Mathematics and Its Application,1979,23(1): 97-120

33 Sneddon IN.The Use of Integral Transform.New York:McGraw-Hill,1972

34 Zhong WX,Lin JH,Gao Q.The precise computation for wave propagation in stratifie materials.International Journal for Numerical Methods in Engineering,2004,60(1):11-25

35 Bailey DH,Swarztrauber PN.A fast method for the numerical evaluation of continuous Fourier and Laplace transforms.SIAM Journal on Scientifi Computing,1994,15(5):1105-1110

36 Abate J,Valko PP.Multi-precision Laplace transform inversion.International Journal for Numerical Methods in Engineering,2004, 60(5):979-993

37 Ai ZY,Yue ZQ,Tham LG,et al.Extended Sneddon and Muki solutions for multilayered elastic materials.International Journal of Engineering Science,2002,40(13):1453-1483

EPIM FOR THERMAL CONSOLIDATION PROBLEMS OF SATURATED POROUS MEDIA SUBJECTED TO A DECAYING HEAT SOURCE1)

Wang Lujun?,?,2)Ai Zhiyong?,3)

?(Department of Geotechnical Engineering,Key Laboratory of Geotechnical and Underground Engineering of Ministry of Education,Tongji University,Shanghai200092,China)

?(Institute of Geotechnical Engineering,Key Laboratory of Soft Soils and Geoenvironmental Engineering of Ministry of Education,Zhejiang University,Hangzhou310058,China)

The thermal consolidation of saturated porous media subjected to a heat source is an important subject in civil engineering and energy engineering.For the complexity of the problem,the porous media are usually treated as homogeneous isotropic media and the heat source is assumed to be a heat source with constant strength in the existingstudies.In engineering practice,natural saturated porous media usually show obvious layered characteristics and the heat source is decaying with time.In this case,the extended precise integration method(EPIM)is presented in this study to investigate the thermal consolidation problems of layered saturated porous media subjected to a decaying heat source. The partial di ff erential equations are reduced to ordinary ones by means of the integral transform techniques.Combining the adjacent layer elements and considering the boundary conditions,the EPIM solutions in the transformed domain of the problems are deduced.With the aid of corresponding numerical integral inversion,the temperatures,excess pore pressures and vertical displacements in the physical domain are obtained.A numerical example with the corresponding calculation program is performed to compare with the existing results,which confir the applicability and validity of the presented method in dealing with the thermal consolidation problems of layered saturated porous media.Finally, numerical examples are carried out to analyse the influenc of the heat source’s half-life and buried depth,as well as the stratificatio of medium on the thermal consolidation behaviour.Numerical results show that:the decay period of heat sources has significan influenc on the peak values and peak time of temperature and excess pore pressure,the longer the decay period,the greater the peak values and the longer the peak time of temperature and excess pore pressure; burial depths have obvious influenc on the variations of excess pore pressure and vertical displacement,the evolutions of vertical displacements against time on both side of the deeply buried heat source are symmetrical,while there is no such phenomenon for the shallow heat source;stratificatio characteristics of the saturated porous media shows prominent e ff ects on the thermal consolidation.

heat source,saturated porous media,thermal consolidation,extended precise integration method

O302

A

10.6052/0459-1879-16-272

2016–09–28收稿，2017–01–04錄用,2017–01–09網(wǎng)絡版發(fā)表.

1)國家自然科學基金資助項目(50578121,41672275).

2)王路君，博士，主要研究方向：巖土工程及地下工程.E-mail：wanglujun007@163.com

3)艾智勇，教授，主要研究方向：巖土工程及地下工程.E-mail：zhiyongai@#edu.cn

王路君，艾智勇.衰變熱源作用下飽和多孔介質(zhì)熱固結問題的擴展精細積分法.力學學報,2017,49(2):324-334

Wang Lujun,Ai Zhiyong.EPIM for thermal consolidation problems of saturated porous media subjected to a decaying heat source.Chinese Journal of Theoretical and Applied Mechanics,2017,49(2):324-334