朱媛媛， 胡育佳， 程昌鈞， 鄭曉妹

(1.上海師范大學(xué) 信息與機(jī)電工程學(xué)院，上海 200234； 2. 上海理工大學(xué) 機(jī)械工程學(xué)院，上海 200093； 3. 上海大學(xué) 力學(xué)系，上海 200093)

基于DQM的空間軸對(duì)稱流體飽和多孔熱彈性柱體動(dòng)力學(xué)特性研究

朱媛媛1， 胡育佳2， 程昌鈞3， 鄭曉妹1

(1.上海師范大學(xué) 信息與機(jī)電工程學(xué)院，上海 200234； 2. 上海理工大學(xué) 機(jī)械工程學(xué)院，上海 200093； 3. 上海大學(xué) 力學(xué)系，上海 200093)

在熱局部平衡條件下研究了空間軸對(duì)稱不可壓流體飽和多孔熱彈性柱體在表面溫度載荷作用下的動(dòng)力學(xué)特性?；赿e Boer多孔介質(zhì)混合物理論，給出了問題的數(shù)學(xué)模型；其次綜合采用微分求積方法-二階向后差分法-Newton-Raphson迭代法求解了數(shù)學(xué)模型，得到柱體各離散點(diǎn)處未知物理量的數(shù)值結(jié)果，進(jìn)而可分析柱體的動(dòng)力學(xué)特性。為了驗(yàn)證該方法的正確性，計(jì)算了不可壓流體飽和多孔彈性柱體的動(dòng)力固結(jié)問題，并與de Boer等的解析結(jié)果進(jìn)行了比較，兩者吻合良好，也證明DQM有較小的計(jì)算量和較高的精度。最后分別研究和比較了柱體只受到表面外載荷作用下的動(dòng)力學(xué)特性和受到兩種表面溫度載荷與外載荷聯(lián)合作用下的動(dòng)力學(xué)特性，考察了材料的某些參數(shù)對(duì)柱體動(dòng)力學(xué)特性的影響。

流體飽和多孔熱彈性柱體；多孔介質(zhì)混合物理論(PMT)；空間軸對(duì)稱；微分求積法(DQM)；熱動(dòng)力學(xué)特性

近年來，飽和多孔介質(zhì)熱-流-固耦合系統(tǒng)的動(dòng)力學(xué)行為已得到國內(nèi)外學(xué)者的廣泛重視。由于飽和多孔介質(zhì)中固相骨架和孔隙流體之間具有不同的物理和力學(xué)性能，因此固相和液相之間存在著與孔隙結(jié)構(gòu)密不可分的相互作用；同時(shí)，孔隙流體的運(yùn)動(dòng)不僅有擴(kuò)散效應(yīng)，而且常伴隨熱反應(yīng)的發(fā)生，如：核廢料污染物在地下低滲透巖土中的傳熱-傳質(zhì)問題、石油天然氣開發(fā)過程中的流固耦合問題、地?zé)醿?chǔ)層的合理開發(fā)和利用問題中的流固耦合問題等，這些問題都可看作是飽和多孔介質(zhì)熱-流-固耦合作用的結(jié)果。所以，對(duì)飽和多孔介質(zhì)力學(xué)熱力學(xué)性能與行為的研究是典型的多場(chǎng)耦合問題，不僅在土力學(xué)、水文學(xué)等經(jīng)典應(yīng)用領(lǐng)域發(fā)揮重要作用，也已成為許多新興學(xué)科和應(yīng)用技術(shù)發(fā)展的關(guān)鍵，相關(guān)理論和數(shù)值方法的研究具有重要的理論意義和廣泛的應(yīng)用背景。

目前，研究多孔介質(zhì)宏觀力學(xué)行為的理論主要有：Biot理論、多孔介質(zhì)混合物理論和雜交混合物理論等。Biot[1]建立了飽和多孔介質(zhì)熱彈性和熱動(dòng)力理論，成功地應(yīng)用于諸多工程領(lǐng)域，取得了眾多成果。近年來，基于Biot理論，Cui等[2]引入熱孔隙度狀態(tài)表面的概念，提出了飽和土熱-水力-力學(xué)耦合分析的一個(gè)較完整的理論模型。Wu等[3]提出了非飽和土熱-水力-力學(xué)本構(gòu)模型。劉干斌等[4]利用Fourier變換研究了荷載作用下地基土體的熱-水-力耦合動(dòng)力響應(yīng)問題。白冰[5]分析了各向同性飽和多孔介質(zhì)熱-水-力一維耦合系統(tǒng)在循環(huán)溫度荷載下的響應(yīng)。

但是Biot理論是建立在力學(xué)直觀基礎(chǔ)上的，其中耦合質(zhì)量項(xiàng)的引入與連續(xù)介質(zhì)力學(xué)相違背。因此，需要建立更加精確的宏觀數(shù)學(xué)理論來描述飽和多孔介質(zhì)的力學(xué)行為。近年來由de Boer等[6-7]發(fā)展起來的多孔介質(zhì)混合物理論, 采用公理化方法和體積分?jǐn)?shù)概念研究多孔介質(zhì)的力學(xué)問題, 將若干微觀性質(zhì)直接用來描述宏觀性質(zhì)，避免了Biot理論的非相容性和雜交混和物理論[8-9]中的繁雜公式，同時(shí)較容易反映某些動(dòng)力學(xué)和非線性效應(yīng)，為研究飽和土的力學(xué)行為提出了另一條新途徑?；诙嗫捉橘|(zhì)混合物理論，de Boer等[10]利用Laplace變換給出了一維動(dòng)力學(xué)響應(yīng)的解析解。Heider等[11]研究了在飽和多孔半空間介質(zhì)中波的傳播。Hu等[12]研究了黏彈性流體飽和多孔介質(zhì)的動(dòng)力學(xué)特性。de Boer等[13]給出了多孔介質(zhì)熱動(dòng)力學(xué)行為的本構(gòu)關(guān)系。He等[14]建立了熱局部非平衡條件下不可壓飽和多孔彈性介質(zhì)熱-流-固耦合模型，Yang[15]建立了相應(yīng)的Gurtin型廣義變分原理。Qin等[16]建立了非飽和土的熱-水-力耦合分析模型，可反映熱膨脹、熱滲流、水的相變、氣體溶解及土骨架變形等現(xiàn)象的耦合過程。Ehlers等[17]考慮了多孔介質(zhì)中傳質(zhì)傳熱問題。

本文基于de Boer多孔介質(zhì)混合物理論，在熱局部平衡條件下，研究了空間軸對(duì)稱流體飽和多孔熱彈性柱體在表面溫度載荷作用下的動(dòng)力學(xué)特性。首先給出了問題的數(shù)學(xué)模型，然后分別在空間域內(nèi)采用微分求積法(DQM)，時(shí)間域內(nèi)采用二階向后差分格式離散控制方程、邊界條件和初始條件，最后運(yùn)用Newton-Raphson法進(jìn)行迭代求解，可得到柱體各離散點(diǎn)處未知物理量的數(shù)值結(jié)果。為了驗(yàn)證本文方法的正確性，文中計(jì)算了無熱效應(yīng)時(shí)不可壓流體飽和多孔彈性柱體的動(dòng)力學(xué)固結(jié)問題，并與1993年de Boer等利用Laplace變換得到的解析解進(jìn)行了比較，兩者吻合良好，證明DQM具有計(jì)算量小，精度高，數(shù)值穩(wěn)定等優(yōu)點(diǎn)。最后在熱局部平衡條件下，利用本文方法研究了空間軸對(duì)稱流體飽和多孔熱彈性柱體只受到表面外載荷作用下的動(dòng)力學(xué)特性和受到兩種表面溫度載荷與外載荷聯(lián)合作用下的動(dòng)力學(xué)特性，考察了材料的某些參數(shù)對(duì)柱體動(dòng)力學(xué)特性的影響，得到了一些有益的結(jié)論。

1 問題的數(shù)學(xué)描述

1.1 控制微分方程

考察圖1所示的空間軸對(duì)稱不可壓流體飽和多孔熱彈性柱體。設(shè)柱體固相和液相的溫度變化相同θS=θF=θ，即柱體處于熱局部平衡狀態(tài)，基于de Boer多孔介質(zhì)混合物理論，對(duì)各向同性線性熱彈性多孔固相骨架和理想流相，在兩相不可壓、小變形、小流速的假設(shè)下，任意rz平面上柱體的基本未知量滿足如下控制微分方程。

圖1 rz-平面上柱體的物理模型Fig.1 The physical model of cylinder on the rz-plane

質(zhì)量守恒方程為

(1)

忽略體積力，動(dòng)量和動(dòng)量矩守恒方程為

Svwr+βθ,r=0

(2)

能量守恒方程為

(3)

幾何關(guān)系為

(4)

假設(shè)固相材料是線性彈性熱各向同性的，本構(gòu)關(guān)系為

(5)

(6)

1.2 邊界條件

假定柱體上表面為理想排水，給定溫度變化φ(r,t)并承受時(shí)間相關(guān)垂直載荷q(r,t)作用，底部為剛性不排水且絕熱；側(cè)表面的徑向位移受到約束且絕熱非排水(圖1)。有如下邊界條件。

上表面(0≤r≤r1，z=0)的邊界條件

(7a)

底部(0≤r≤r1，z=h)的邊界條件

ur=0,uz=0,p,z=0,θ,z=0

(7b)

側(cè)表面(r=r1,0≤z≤h)的邊界條件

(7c)

對(duì)稱軸(r=0,0≤z≤h)上的對(duì)稱性條件

(7d)

1.3 初始條件

假設(shè)t≤0時(shí)柱體處于靜止，初始條件為

(8)

這樣，式(1)～(6)，邊界條件(7)和初始條件(8)構(gòu)成了不可壓空間軸對(duì)稱流體飽和多孔熱彈性柱體6個(gè)基本未知量(即固相位移ur,uz，流相相對(duì)速度wr,wz，孔隙壓p和變溫θ)滿足的控制微分方程。

2 求解控制方程的DQM

由初邊值問題式(1)～(8)可見，獲得問題的解析解或半解析解是困難的。這里將采用微分求積方法(DQM)在空間域內(nèi)離散控制方程和邊界條件，再采用二階向后差分格式來處理時(shí)間導(dǎo)數(shù)，最后，在離散化的初始條件(8)下運(yùn)用Newton-Raphson迭代方法求離散化代數(shù)方程組的解，從而可得各離散點(diǎn)未知量的數(shù)值結(jié)果。

2.1 微分求積方法簡(jiǎn)介

20世紀(jì)70年代，Bellman等[18-19]提出了一種新的數(shù)值方法，即微分求積法(DQM)，由于DQM不依賴泛函和變分原理，具有公式簡(jiǎn)單、使用方便、精度高、計(jì)算量少等優(yōu)點(diǎn)，引起人們廣泛關(guān)注，并成功應(yīng)用于許多領(lǐng)域[20-22]。DQM思想是將未知函數(shù)對(duì)某方向自變量的偏導(dǎo)數(shù)近似表達(dá)為其各離散點(diǎn)(節(jié)點(diǎn))處相應(yīng)函數(shù)值的加權(quán)和，其中，權(quán)系數(shù)與具體問題無關(guān)，而只與解區(qū)域中所選擇的試函數(shù)和離散點(diǎn)有關(guān)。因此，利用DQM，任何微分方程都能轉(zhuǎn)化為相應(yīng)的代數(shù)方程。本文亦將用DQM來對(duì)基本場(chǎng)方程進(jìn)行空間離散，關(guān)于DQM的基本原理和公式，可以參考文獻(xiàn)[18-19]。

考慮在區(qū)域Ω={(r,z)|0≤r≤r1,0≤z≤h}內(nèi)的函數(shù)Ψ(r,z)，采用Chebyshev-Lobatto多項(xiàng)式零點(diǎn)作為節(jié)點(diǎn)的坐標(biāo)，分別沿r,z方向布置Nr×Nz個(gè)節(jié)點(diǎn)(圖2)，根據(jù)DQM，函數(shù)Ψ(r,z)在節(jié)點(diǎn)r=rξ，z=zη處對(duì)自變量r的n階偏導(dǎo)數(shù)可近似表示為

(9)

圖2 rz-平面上的布點(diǎn)Fig.2 Nodes collocated on the rz-plane

2.2 空間域內(nèi)控制方程的DQ離散化

利用公式(9)在任意rz-平面內(nèi)離散控制微分方程(1)～(3)，可得rz-平面上空間域內(nèi)DQ離散化方程為

(10a)

(10b)

(10c)

(10d)

(10e)

(10f)

其中，ξ=2,3…,Nr-1,η=2,3…,Nz-1。

由幾關(guān)何系(4)，應(yīng)變分量的DQ離散化形式為

(11)

其中，ξ=1,2,3…,Nr,η=1,2,3…,Nz。

本構(gòu)關(guān)系(5)的離散化形式為

(12)

式中，ξ=1,2,3…,Nr,η=1,2,3…,Nz。

同時(shí)，固相有效應(yīng)力分量與總應(yīng)力分量關(guān)系(6)的離散化形式為

(13)

式中，ξ=1,2,3…,Nr,η=1,2,3…,Nz。

式(10)～(13)構(gòu)成基本場(chǎng)方程(1)～(6)空間域內(nèi)DQ離散化形式，它是關(guān)于時(shí)間t的代數(shù)-微分方程組。

2.3 對(duì)稱軸上奇異性條件的處理

當(dāng)r→0

(14)

當(dāng)r→0

(15)

式中，η=1,2，…，Nz。

2.4 邊界條件的DQ離散化

柱體上表面(ξ=1,2,…,Nr,η=1)邊界條件的DQ離散化形式為

(p)ξ1=0, (θ)ξ1=φξ1

(16a)

式中，qξ1,φξ1是ξ=1,2,…,Nr,η=1處外載荷q(r,t)和溫度載荷φ(r,t)的值。

底部(ξ=1,2,…,Nr,η=Nz)邊界條件的DQ離散化形式為

(ur)ξNz=0,(uz)ξNz=0,

(16b)

側(cè)表面(ξ=Nr,η=1,2,…,Nz)邊界條件的DQ離散化形式為

(16c)

對(duì)稱軸(ξ=1,η=1,2,…,Nz)上對(duì)稱性條件的DQ離散化形式為

(16d)

2.5 DQ離散化方程的簡(jiǎn)化

為了提高求解速度，可由公式(10b)～(10c)解出未知量(wr)ξη和(wz)ξη代入 (10a),(10d)～(10g)中，離散化方程(10)轉(zhuǎn)化為在離散點(diǎn)(ξ,η)處只有4個(gè)未知量 (ur)ξη,(uz)ξη,(p)ξη，(θ)ξη的方程組。

2.6 時(shí)間導(dǎo)數(shù)的二階向后差分格式

離散化方程(10)～(16)是一組關(guān)于時(shí)間t的代數(shù)-微分方程，將采用二階向后差分格式來離散函數(shù)φ(t)關(guān)于時(shí)間的導(dǎo)數(shù)

(17)

式中：Δt為時(shí)間步長(zhǎng)；φ(tn,x)為在t=tn時(shí)刻的函數(shù)值。

利用式(17), 對(duì)DQ離散化方程組式(10)～式(16)進(jìn)行時(shí)間離散，并在離散化的初始條件式(8)下，利用Newton-Raphson迭代方法求解離散化代數(shù)方程組，從而可得各離散點(diǎn)各未知量的數(shù)值結(jié)果。

3 流體飽和多孔熱彈性柱的動(dòng)力學(xué)特性

3.1 數(shù)值結(jié)果的驗(yàn)證

de Boer等利用Laplace變換得到了不可壓一維無限長(zhǎng)圓柱體在表面荷載作用下動(dòng)力固結(jié)問題的解析解，本文運(yùn)用DQM求解該問題，并與de Boer等的解析解進(jìn)行比較，以說明本文方法的正確性，結(jié)果的可靠性。考慮高h(yuǎn)=10 m，半徑r1=1 m的不可壓流體飽和多孔彈性體，其上邊界為理想排水，并承受時(shí)間相關(guān)垂直載荷q(t)的作用；底部邊界為剛性不排水；側(cè)表面為不排水邊界，流相和固相的材料參數(shù)見表1。

考慮兩種垂直載荷q(t)：① 階梯載荷q(t)=q0h(t)；② 周期載荷q(t)=q0[1-cos(ωt)]，其中，q0=3 kN/m2，ω=75 s-1，h(t)是Heaviside函數(shù)。

表1 流體飽和土的材料參數(shù)Tab.1 Physical parameters of fluid-saturated soil

圖3和圖4給出了在對(duì)稱軸上不同深度處的沉降uz，其中，點(diǎn)為解析解，實(shí)線和虛線為利用DQM得到的數(shù)值解，分別對(duì)應(yīng)Δt=0.01 s和Δt=0.02 s，計(jì)算中取Nr=Nz=7。可以看到，用本文方法來計(jì)算流體飽和多孔彈性介質(zhì)的動(dòng)力學(xué)特性是有效的，也證明本方法具有很高的精度和收斂性。

3.2 空間軸對(duì)稱流體飽和多孔熱彈性柱的動(dòng)力學(xué)特性

考察高h(yuǎn)=10 m，半徑r1=1 m的體內(nèi)無熱源不可壓空間軸對(duì)稱流體飽和多孔熱彈性柱體，邊界條件如圖1所示，材料參數(shù)見表2。計(jì)算中取Nr=Nz=7，Δt=1 s。

圖3 階梯載荷作用下沉降的時(shí)-程曲線Fig.3 Time-history curves of the settlement under the step loading

圖4 周期載荷作用下沉降的時(shí)-程曲線Fig.4 Time-history curves of the settlement under the cycle loading表2 熱-流-固耦合系統(tǒng)材料參數(shù)Tab.2 Physical parameters of thermo-fluid-solid system

參數(shù)數(shù)值參數(shù)數(shù)值參數(shù)數(shù)值μS0.231×106kN/m2γFR10000N/m3β13.35Pa/℃λS0.346×106kN/m2κF0.0173m/sΘ010℃αS8.9×10-5CcS937J/kg℃nSOS0.6ρSR2610kg/m3cF4186J/kg℃nFOF0.4ρFR1000kg/m3k3.142×106W/m℃

3.2.1 柱體上表面僅受外載荷作用時(shí)的動(dòng)力學(xué)特性

Hu等[12]利用DQM得到了空間軸對(duì)稱流體飽和多孔彈性柱體的動(dòng)力學(xué)特性，其模型中并沒有考慮熱效應(yīng)的影響。假設(shè)在本文給出的流體飽和多孔熱彈性柱體模型中，柱體上表面無溫度變化，即φ(r,t)=0，考慮兩種模型下柱體的動(dòng)力學(xué)特性。

考慮兩種垂直載荷q(t)：① 階梯載荷q(t)=q0h(t)；② 周期載荷q(t)=q0[1-cos(ωt)]，其中，q0=3 kN/m2，ω=0.4 s-1，h(t)是Heaviside函數(shù)。

圖5和圖6為利用DQM求得的柱體對(duì)稱軸不同深度處的沉降uz，流速wz，孔隙壓p和應(yīng)力σzz，圖中實(shí)線為利用本文模型(考慮熱效應(yīng))得到的結(jié)果，圓點(diǎn)為利用文獻(xiàn)[12]中模型(不考慮熱效應(yīng))得到的結(jié)果。由于現(xiàn)在考慮的問題是柱體上表面無溫度改變，可近似認(rèn)為柱體溫度變化θ=0，即常溫狀態(tài)，則利用兩種模型得到的解應(yīng)趨于一致。計(jì)算中也驗(yàn)證了這一點(diǎn)。進(jìn)一步說明本文理論和方法的正確性，結(jié)果的可靠性。此時(shí)柱體動(dòng)力學(xué)特性結(jié)論可參考文獻(xiàn)[12]。

(a) 沉降

(b) 相對(duì)流速

(c) 孔隙壓

(d) 應(yīng)力圖5 階梯載荷作用下柱體的時(shí)-程曲線Fig.5 Time-history curves of the cylinder under the step loading

3.2.2 柱體上表面在溫度載荷和外載荷聯(lián)合作用下的動(dòng)力學(xué)特性

給定柱體上表面外載荷q(r,t)=3 kN/m2。同時(shí)分別考察兩種表面溫度φ(t)對(duì)柱體動(dòng)力學(xué)特性的影響：① 階梯溫度載荷φ(t)=φ0h(t)；② 周期溫度載荷φ(t)=φ0[1-cos(ωt)]，其中，φ0=10 ℃，ω=0.4 s-1，h(t)是Heaviside函數(shù)。

本文考察和分析了達(dá)西滲透系數(shù)κF，熱交換系數(shù)β，熱傳導(dǎo)系數(shù)k等對(duì)柱體動(dòng)力學(xué)特征的影響，為了節(jié)省篇幅，僅列出了達(dá)西滲透系數(shù)κF影響的計(jì)算結(jié)果，并進(jìn)行了討論，對(duì)于β，k等參數(shù)的影響只列出了結(jié)論。

(1) 在階梯溫度載荷和外載荷情況下的結(jié)果

圖7給出了在階梯溫度載荷和外載荷作用下，不同的達(dá)西滲透系數(shù)κF對(duì)對(duì)稱軸不同深度處的沉降uz，相對(duì)流速wz，孔隙壓p及溫度變化θ的影響。

從圖7(a)中可以看到，在不同的κF下，柱體沉降隨著時(shí)間的增大趨于相同穩(wěn)定值。當(dāng)κF較大時(shí)，固結(jié)作用過程快于熱傳導(dǎo)過程，初始階段固相的熱體積膨脹效應(yīng)未得到明顯體現(xiàn)，柱體沉降大。當(dāng)κF較小時(shí)，熱傳導(dǎo)過程快于機(jī)械載荷下的固結(jié)作用，初始階段首先表現(xiàn)為固相的熱體積膨脹效應(yīng)，而后固結(jié)作用才逐漸發(fā)揮作用。

圖7(b)指出，柱體流速隨時(shí)間增大逐漸趨于零。在初始階段，柱體上表面附近的流速大于內(nèi)部的流速。當(dāng)κF增大時(shí)，初始階段流速的峰值增大。

圖7(c)指出，柱體的孔隙壓由初始值逐漸消散至零，且上表面附近孔隙壓的消散速度快于內(nèi)部的消散速度；當(dāng)κF增大時(shí)，孔隙壓消散速度增大。

從圖7(d)中可以看到，柱體溫度隨時(shí)間的增加逐漸上升并由上表面向縱深處傳導(dǎo)和擴(kuò)散，最后達(dá)到等溫狀態(tài)，但κF對(duì)溫度影響較小。

(a) 沉降

(b) 相對(duì)流速

(c) 孔隙壓

(d) 應(yīng)力圖6 周期載荷作用下柱體的時(shí)-程曲線Fig.6 Time-history curves of the cylinder under the cycle loading

(a) 沉降

(b) 相對(duì)流速

(c) 孔隙壓

(d) 溫度圖7 階梯溫度載荷作用下柱體的時(shí)-程曲線Fig.7 Time-history curves of the cylinder under the step temperature loading

同時(shí)，本文考慮了熱交換系數(shù)β，熱傳導(dǎo)系數(shù)k等對(duì)柱體的影響。結(jié)果表明，熱交換系數(shù)β越小，流固兩相之間相互作用力越小，初始階段固相的熱體積膨脹效應(yīng)被抑制，柱體沉降較大。而當(dāng)熱傳導(dǎo)系數(shù)k較大時(shí)，溫度場(chǎng)達(dá)到穩(wěn)態(tài)所需時(shí)間短，初始階段熱傳導(dǎo)過程快于固結(jié)作用過程，柱體沉降較小。

(2) 在循環(huán)溫度載荷和外載荷情況下的結(jié)果

圖8中給出了在循環(huán)溫度載荷和外載荷作用下，達(dá)西滲透系數(shù)κF對(duì)對(duì)稱軸不同深度處的沉降uz，相對(duì)流速wz，孔隙壓p及溫度變化θ的影響。計(jì)算表明，在循環(huán)溫度載荷下，柱體的熱動(dòng)力學(xué)整體性質(zhì)和在階梯溫度載荷下相似，除了時(shí)-程曲線呈周期性變化外。達(dá)西滲透系數(shù)κF較小時(shí)，熱傳導(dǎo)過程快于機(jī)械載荷下的固結(jié)作用過程，初始階段表現(xiàn)為柱體固相的熱體積膨脹效應(yīng)，沉降的振幅小。

同時(shí)，熱交換系數(shù)β越小，流固兩相之間相互作用力越小，初始階段固相的熱體積膨脹效應(yīng)被抑制，柱體沉降大。而當(dāng)熱傳導(dǎo)系數(shù)k較大時(shí)，溫度場(chǎng)達(dá)到穩(wěn)態(tài)所需時(shí)間短，初始階段熱傳導(dǎo)過程快于固結(jié)作用過程，柱體沉降小。

在所考慮的幾種情況下，柱體溫度隨時(shí)間的增大周期性變化，并由上表面向縱深處傳導(dǎo)和擴(kuò)散，各物理量沿縱深方向存在相位差，但周期相同。

(a) 沉降

(b) 相對(duì)流速

(c) 孔隙壓

(d) 溫度圖8 循環(huán)溫度載荷作用下柱體的時(shí)-程曲線Fig.8 Time-history curves of the cylinder under the cycle temperature loading

4 結(jié) 論

本文基于de Boer多孔介質(zhì)混合物理論，綜合采用DQM-二階向后差分法-Newton-Raphson迭代法，研究了空間軸對(duì)稱流體飽和多孔熱彈性柱體在表面溫度載荷作用下的動(dòng)力學(xué)特性。

為了驗(yàn)證本文方法的正確性，計(jì)算了無熱效應(yīng)時(shí)不可壓流體飽和多孔彈性柱體的動(dòng)力固結(jié)問題，并與現(xiàn)有解析結(jié)果進(jìn)行了比較，二者吻合良好，證明DQM具有精度高，計(jì)算量小，數(shù)值穩(wěn)定等優(yōu)點(diǎn)。

本文分別研究了空間軸對(duì)稱流體飽和多孔熱彈性柱體在熱局部平衡條件下，只受到表面外載荷作用下的動(dòng)力學(xué)特性和受到表面溫度載荷與外載荷聯(lián)合作用下的動(dòng)力學(xué)特性，并比較了在柱體上表面給定兩種不同溫度載荷(即階梯溫度載荷和周期溫度載荷)的結(jié)果，考察了材料的某些參數(shù)對(duì)柱體動(dòng)力學(xué)特性的影響。結(jié)果表明在兩種溫度載荷的作用下，除了在周期溫度載荷條件下柱體物理量的時(shí)-程曲線呈周期性變化外，柱體的熱動(dòng)力學(xué)整體特性是相似的：隨著時(shí)間的增加，柱體沉降趨于穩(wěn)定、流速逐漸趨于零、孔隙壓由初始值逐漸消散至零、溫度逐漸上升并由上表面向縱深處傳導(dǎo)和擴(kuò)散；同時(shí)，當(dāng)達(dá)西滲透系數(shù)κF較大或熱交換系數(shù)β較小時(shí)，初始階段柱體固結(jié)作用過程快于熱傳導(dǎo)過程，柱體沉降大；當(dāng)熱傳導(dǎo)系數(shù)k較大時(shí)，溫度場(chǎng)達(dá)到穩(wěn)態(tài)所需時(shí)間短，初始階段熱傳導(dǎo)過程快于固結(jié)作用過程，柱體沉降較小。另外，在周期溫度載荷作用下，柱體溫度隨時(shí)間的增大周期性變化，并由上表面向縱深處傳導(dǎo)和擴(kuò)散，各物理量沿縱深方向存在相位差，但周期相同。

[1] BIOT M A. Theory of elasticity and consolidation for a porous anisotropic solid[J]. Journal of Applied Physics, 1955, 26(2): 182-185.

[2] CUI Y J, SULTAN N, DELAGE P. A thermomechanical model for saturated clays[J]. Canadian Geotechnical Journal, 2000, 37(3): 607-620.

[3] WU W H,LI X K,CHARLIER R, et al. A thermo-hydro-mechanical constitutive model and its numerical modeling for unsaturated soils[J]. Computers & Geotechnics, 2004, 31(2): 155-167.

[4] 劉干斌, 姚海林, 楊洋, 等. 考慮熱-水-力耦合效應(yīng)多孔彈性地基的動(dòng)力響應(yīng)[J]. 巖土力學(xué), 2007, 28(9): 1784-1795.

LIU Ganbin, YAO Hailin, YANG yang, et al. Coupling thermo-hydro-mechanical dynamic response of a porous elastic medium[J]. Rock and Soil Mechanics, 2007, 28(9): 1784-1795.

[5] 白冰. 循環(huán)溫度荷載作用下飽和多孔介質(zhì)熱-水-力耦合響應(yīng)[J]. 工程力學(xué), 2007, 24(5): 87-92.

BAI Bing. Thermo-hydro-mechanical responses of saturated porous media under cycle thermal loading[J]. Engineering Mechanics, 2007, 24(5): 87-92.

[6] DE BOER R. Theoretical poroelasticity-a new approach[J]. Chaos Solitons & Fractals, 2005, 25(4): 861-878.

[7] BOWEN R M. Compressible porous media models by use of the theory of mixtures[J]. International Journal of Engineering Science, 1982, 20(6): 693-735.

[8] SCHREFLER B A. Mechanics and thermodynamics of saturated /unsaturated porous materials and quantitative solution[J]. Applied Mechanics Reviews, 2002, 55(4): 351-388.

[9] PLACIDI L, DELL’LSOLA F, IANIRO N, et al. Variational formulation of pre-stressed solid-fluid mixture theory, with an application to wave phenomena[J]. European Journal of Mechanics A/Solids, 2007, 27(4): 582-606.

[10] DE BOER R, EHLERS W, LIU Z. One-dimensional transient wave propagation in fluid-saturated incompressible porous media[J]. Archive of Applied Mechanics, 1993, 63(1): 59-72.

[11] HEIDER Y, MARKERT B, EHLERS W. Dynamic wave propagation in infinite saturated porous media half spaces[J]. Computational Mechanics, 2012, 49(3): 319-336.

[12] HU Y J, ZHU Y Y, CHENG C J. DQM for dynamic response of fluid-saturated visco-elastic porous media[J]. International Journal of Solids & Structures, 2009, 46(7): 1667-1675.

[13] DE BOER R, KOWALSKI S J. Thermodynamics of fluid-saturated porous media with a phase change[J]. Acta Mechanica, 1995, 109(1): 167-189.

[14] HE L W, JIN Z H. A local thermal non-equilibrium poroelastic theory for fluid saturated[J]. Journal of Thermal Stresses, 2009, 33(8): 799-813.

[15] YANG X. Gurtin-type variational principles for dynamics of a non-local thermal equilibrium saturated porous medium[J]. Acta Mechanica Solida Sinica, 2005, 18(1): 37-45.

[16] QIN B, CHEN Z H, FANG Z D, et al. Analysis of coupled thermo-hydro-mechanical behavior of unsaturated soils based on theory of mixtures I[J]. Applied Mathematics and Mechanics (English Edition), 2010, 31(12): 1561-1576.

[17] EHLERS W, HBERLE K. Interfacial mass transfer during gas-liquid phase change in deformable porous media with heat transfer[J]. Transport in Porous Media, 2016, 114(2): 1-32.

[18] BELLMAN R E, CASTI J. Differential quadrature and long term integration[J]. Journal of Mathematical Analysis & Applications, 1970, 34(2): 235-238.

[19] BELLMAM R E, KASHEF B G, CASTI J. Differential quadrature: a technique for the rapid solution of nonlinear partial differential equations[J]. Journal of Computational Physics, 1972, 10(1): 40-52.

[20] BERT C W, MALIK M. Differential quadrature method in computational mechanics: a review[J]. Applied Mechanics Reviews, 1996, 49(1): 1-28.

[21] CHEN C N. The two-dimensional frames model of the dfferential quadrature element method[J]. Computers & Structures, 1997, 62(3):555-571.

[22] ZHU Y Y, LI Y, CHENG C J. Analysis of nonlinear characteristics for thermoelastic half-plane with voids[J]. Journal of Thermal Stresses, 2014, 37(7): 794-816.

Dynamiccharacteristicsforaspatial-axisymmetricfluid-saturatedporousthermo-elasticcylinderbasedonDQM

ZHU Yuanyuan1, HU Yujia2, CHENG Changjun3, ZHENG Xiaomei1

(1. College of Information, Mechanical and Electrical Engineering, Shanghai Normal University, Shanghai 200234, China；2. College of Mechanical Engineering, University of Shanghai for Science and Technology, Shanghai 200093, China;3. Department of Mechanics, Shanghai University, Shanghai 200072, China)

Dynamic characteristics of an incompressible spatial-axisymmetric fluid-saturated porous thermo-elastic cylinder subjected to a surface temperature loading were studied in case of local thermal equilibrium. Firstly, the mathematical model of the problem was established based on de Boer porous media theory. Then, the differential quadrature method, the second-order backward difference scheme and Newton-Raphson iterative method were synthetically used to solve the mathematical model and obtain the numerical results of the unknown quantities at every discretized point, and the dynamic characteristics of the cylinder were further studied. In order to verify the validity of the proposed method, the dynamic consolidation problem of an incompressible fluid-saturated porous elastic cylinder was computed with this method. The obtained numerical results agreed well with the analytical ones published by de Boer et.al, it was shown that the proposed method has two advantages: smaller amount of computation and higher accuracy. Finally, the dynamic characteristics of a spatial-axisymmetric fluid-saturated porous thermo-elastic cylinder subjected to mechanical load or both mechanical load and temperature load were studied and compared, the effects of material’s some parameters on the dynamic characteristics of the cylinder were investigated.

fluid-saturated porous thermo-elastic cylinder； porous media theory(PMT)； spatial-axisymmetric problem； differential quadrature method (DQM)； thermo-dynamic

上海市自然科學(xué)基金(15ZR1431600)

2017-01-13 修改稿收到日期：2017-04-06

朱媛媛 女，博士，副教授，1971年3月生

胡育佳 男，博士，副教授，1979年3月生

TU311

A

10.13465/j.cnki.jvs.2017.23.013