半无限空间饱和土一维热固结精确积分解

钮家军; 王秀凯; 凌道盛; 巩师林; 石吉森

doi:10.11779/CJGE202108007

半无限空间饱和土一维热固结精确积分解 English Version

钮家军^{1, 2,},
王秀凯^{1, 2},
凌道盛^{1, 2, ,},
巩师林^{1, 2},
石吉森^{1, 2}

1.
浙江大学岩土工程研究所,浙江　杭州　310058
2.
浙江大学软弱土与环境土工教育部重点实验室,浙江　杭州　310058

基金项目:

国家重点研发计划项目 2016YFC0800200

国家自然科学基金项目 51578502

详细信息

作者简介:
钮家军(1993— ),男,博士研究生,主要从事多孔介质热–水–力耦合现象的研究。E-mail：niu0918@zju.edu.cn

通讯作者:
凌道盛, E-mail：dsling@zju.edu.cn

中图分类号: TU431
计量
- 文章访问数: 241
- HTML全文浏览量: 53
- PDF下载量: 129
出版历程
- 收稿日期: 2020-09-08
- 网络出版日期: 2022-12-02
- 刊出日期: 2021-07-31

Exact integral solutions for one-dimensional thermal consolidation of semi-infinite saturated soils

NIU Jia-jun^{1, 2,},
WANG Xiu-kai^{1, 2},
LING Dao-sheng^{1, 2, ,},
GONG Shi-lin^{1, 2},
SHI Ji-sen^{1, 2}

1.
Institute of Geotechnical Engineering, Zhejiang University, Hangzhou 310058, China
2.
MOE Key Laboratory of Soft Soils and Geoenvironmental Engineering, Zhejiang University, Hangzhou 310058, China

摘要

摘要: 基于Coussy提出的饱和多孔介质热–水–力耦合理论,给出考虑热渗效应和等温热流效应的饱和土热固结耦合控制方程,然后通过正弦和余弦变换方法,结合Hamilton-Cayley定理,给出了两类任意非齐次边界条件下,半无限空间饱和土一维热固结精确积分解的求解方法,避免了Laplace变换处理复杂边界问题时的数值反演。最后,通过算例验证了本文解析解的正确性,并分析了半无限空间表面透水或不透水时,饱和土体在季节性温度荷载作用下的温度、超静孔压和位移响应特性。结果表明：温度、超静孔压和位移受季节性温度荷载影响而呈现周期性波动,表面不透水情形的土体超静孔压和位移响应明显高于相应的表面透水情形的土体响应。
- 半无限饱和土 /
- 一维热固结 /
- 正弦和余弦变换 /
- 精确积分解 /
- 季节性温度荷载
Abstract: Based on the thermo-hydro-mechanical coupled theory for saturated porous media proposed by Coussy, two coupled governing equations for thermal consolidation of saturated soils are presented, taking thermo-osmosis and thermo- filtration into account. The sine and cosine transforms as well as the Hamilton-Cayley law are used to obtain exact integral solutions for one-dimensional thermal consolidation of semi-infinite saturated soils with two types of arbitrary nonhomogeneous boundary conditions. It helps avoid the errors of numerical inverse Laplace transform when dealing with complex boundary conditions. Finally, case studies are provided to validate the proposed analytical solutions and to investigate the thermal responses of saturated soils with seasonal thermal loads on the permeable or impermeable surface. The results show that temperature change, excess pore water pressure and displacement all fluctuate periodically due to seasonal thermal loads. The responses of the excess pore water pressure and displacement of the saturated soils with an impermeable surface are much higher than those with a permeable surface.
- semi-infinite saturated soil /
- one-dimensional thermal consolidation /
- sine and cosine transform /
- exact integral solutions /
- seasonal thermal load

HTML全文

0. 引言

饱和土土颗粒和孔隙流体常具有热胀冷缩的特性,由于孔隙流体的热膨胀系数通常大于土颗粒的热膨胀系数,温度升高会导致土体内产生超静孔压,超静孔压的消散与土体的热胀冷缩都将导致土体的变形,这一现象在岩土工程中常常被忽略。但随着中国不断推进地热资源开发和地下能源存储等工程的建设,温度场影响下的岩土体热–水–力耦合响应分析已成为研究的热点。

许多学者针对温度变化引起的饱和弹性土体热固结响应展开了解析研究。Booker等^[1-2]最早系统进行了相关工作,如通过Laplace变换、Fourier变换先后提出了点热源和球热源周围土体的热固结解析解。Mctigue^[3]将由热传导方程得到的温度场解代入渗流场方程,得到了半无限空间饱和土体表面透水或不透水时,土体内部超静孔压在表面恒定温度荷载作用下的热响应解。Giraud等^[4]建立了内部存在指数衰减热源的半无限空间饱和土体简化计算模型,并通过Laplace变换求解得到了温度、超静孔压和位移随深度的分布。Smith等^[5-6]提出了力场–渗流场–温度场全耦合的饱和土固结理论,并提出了在Laplace变换域上解决热弹性固结问题的边界元法。Bai等^[7]通过有限Fourier变换给出了单侧排水的单层饱和土体热固结解,分析了控制方程非耦合、部分耦合及全耦合情况下温度、超静孔压与位移的变化区别。Zhou等^[8-9]建立了考虑热渗效应与等温热流效应的热固结全耦合控制方程,并用Laplace变换给出了线热源、柱热源与球热源作用下的饱和土体热固结解析解。Blond等^[10]考虑了土颗粒和孔隙流体的可压缩性,给出了半无限空间土体在表面简谐温度作用下超静孔压变化最大值的解,并确定了孔压最大值所在位置。白冰^[11-13]利用Laplace变换和数值逆变换,给出了半无限饱和介质在周期变化温度作用时,土体内部一维热响应的半解析解。Liu等^[14-15]建立了热–水–力耦合的动力形式控制方程,并采用Laplace变换得到了轴对称和球对称情况的土体热动力响应。He等^[16]与Yang等^[17]分别假设土骨架与孔隙流体处于局部非热平衡态,用Laplace变换的方法求解土体内渗流场与温度场分布,并与局部热平衡态的热固结解进行对比。Ai等^[18-20]、Wang等^[21-22]在Laplace变换域引入解析层元法,对层状多孔介质的热–水–力耦合响应进行了大量研究。

从前人的研究了解到,饱和土热固结的解析解可用于揭示饱和弹性土体受温度荷载作用的固结响应规律,也可用于验证相关数值解的正确性,具有十分重要的意义。但学者们广泛应用了Laplace变换对热固结问题进行求解,复杂边界条件下得到的Laplace变换域解无法直接进行逆变换得到其对应的时域表达式,数值逆变换的参数又常难以确定。此外,学者们对半无限空间饱和土体在季节性温度荷载作用下的一维热固结问题研究较少,表面透水或不透水对热固结的影响未进行对比分析。因此,本文借鉴钮家军等^[23]提出的多场耦合方程求解方法,采用正余弦变换与Hamilton–Cayley定理,给出两类非齐次边界条件下,半无限饱和土体一维热固结精确时域显式积分解的求解方法,避免了Laplace数值逆变换可能带来的误差。最后通过算例验证本文解的正确性,并对比研究表面透水和不透水的半无限空间饱和土体在季节性温度荷载作用下的一维热固结响应。

1. 饱和土热–水–力耦合控制方程

Smith等^[5-6]曾提出饱和土热–水–力耦合的平衡方程、质量守恒方程与能量守恒方程；Zhou等^[8]进一步提出考虑了热渗效应、等温热流效应及物性参数变化的饱和多孔介质全耦合理论模型；Coussy^[24]从热力学第一定律和第二定律角度,推导得到了对称形式的饱和多孔介质耦合控制方程；Bai等^[25-26]考虑了渗流过程对热传导的影响,给出了饱和土多场全耦合的控制方程。本文首先给出土体的本构关系及修正的Darcy定律和Fourier定律：

$σ_{i j} = 2 G ε_{i j} + λ ε_{v} δ_{i j} - 3 K_{0} α_{s} T δ_{i j} - α p δ_{i j}$ ,

(1)

$q_{w} = - \frac{k_{w}}{γ_{w}} \nabla p - S_{w} \nabla T$ ,

(2)

$q_{T} = - S_{T} \nabla p - K_{T} \nabla T$ 。

(3)

式中 $ε_{v} = ε_{x} + ε_{y} + ε_{z} ； K_{0} = λ + \frac{2}{3} G ； α = 1 - \frac{K_{0}}{K_{s}} ； S_{T} \approx$ $T_{0} S_{w}$ ； $σ_{i j}$ 为总应力；λ,G为Lame常数； $ε_{v}$ 为体应变；K₀为土骨架体积模量； $α_{s}$ 为土颗粒的线热膨胀系数；p为土体内超静孔压；T为区别于初始温度T₀的土体温度改变量；α为Biot系数； $γ_{w}$ 为孔隙流体的重度；k_w和K_T分别为土体的渗透系数与热传导系数；S_w和S_T分别为土体的热渗效应系数与等温热流效应系数。

根据修正的Darcy定律可知,当式（2）右侧两项的差值不能忽略不计时,即土渗透系数较小或温度梯度较大时,应当考虑热渗效应的影响,这与Carnahan^[27]、Horseman等^[28]的研究结论一致。等温热流效应同理。将考虑热渗效应、等温热流效应的Darcy定律和Fourier定律代入Coussy^[24]提出的饱和多孔介质热–水–力耦合控制方程,可以得到

$G \nabla^{2} u + (λ + G) \nabla ε_{v} - α \nabla p - β \nabla T = 0$ ,

(4)

$\frac{k_{w}}{γ_{_{w}}} \nabla^{2} p + S_{w} \nabla^{2} T = α \frac{\partial ε_{v}}{\partial t} + α_{p} \frac{\partial p}{\partial t} - Y \frac{\partial T}{\partial t}$ ,

(5)

$\frac{S_{T}}{T_{0}} \nabla^{2} p + \frac{K_{T}}{T_{0}} \nabla^{2} T = β \frac{\partial ε_{v}}{\partial t} - Y \frac{\partial p}{\partial t} + \frac{c_{v}}{T_{0}} \frac{\partial T}{\partial t}$ 。

(6)

其中,

$β = 3 α_{s} K_{0}, Y = 3 (α - n) α_{s} + 3 n α_{w}$ ,

$α_{p} = n / K_{w} + (α - n) / K_{s}, c_{v} = (1 - n) ρ_{s} c_{s} + n ρ_{w} c_{w}$ 。

式中 u为土骨架的位移；t为时间；n为土体的初始孔隙率；K_s和K_w分别为土颗粒和孔隙流体的体积模量；ρ_s和ρ_w分别为土颗粒和孔隙流体的密度；α_s和α_w分别为土颗粒和孔隙流体的线热膨胀系数；c_s和c_w分别为土颗粒和孔隙流体的质量比热容。

考虑一维情况,将式（4）左右两侧分别积分,可以得到

$ε_{v} = \nabla \cdot u = \frac{α p + β T + q (t)}{λ + 2 G} = \frac{α p + β T + q (t)}{E_{c}}$ 。

(7)

式中 q(t)为关于时间的函数,由位移边界确定,常与外荷载有关；ε_v为体应变,以拉为正。

将式（7）代入式（5）和（6）,可以得到以超静孔压p与温度改变量T为未知量的饱和土体热–水–力耦合控制方程,其一维形式为

$\begin{array}{l} \frac{k_{w}}{γ_{w}} \frac{\partial^{2} p}{\partial x^{2}} + S_{w} \frac{\partial^{2} T}{\partial x^{2}} \\ = (\frac{α^{2}}{E_{c}} + α_{p}) \frac{\partial p}{\partial t} + (\frac{α β}{E_{c}} - Y) \frac{\partial T}{\partial t} + \frac{α}{E_{c}} \frac{\partial q (t)}{\partial t} \end{array}$ ,

(8)

$\begin{array}{l} \frac{S_{T}}{T_{0}} \frac{\partial^{2} p}{\partial x^{2}} + \frac{K_{T}}{T_{0}} \frac{\partial^{2} T}{\partial x^{2}} \\ = (\frac{α β}{E_{c}} - Y) \frac{\partial p}{\partial t} + (\frac{β^{2}}{E_{c}} + \frac{c_{v}}{T_{0}}) \frac{\partial T}{\partial t} + \frac{β}{E_{c}} \frac{\partial q (t)}{\partial t} \end{array}$ 。

(9)

为便于表示,将式（8）,（9）简写为

$g_{11} \frac{\partial^{2} p}{\partial x^{2}} + g_{12} \frac{\partial^{2} T}{\partial x^{2}} + c_{11} \frac{\partial p}{\partial t} + c_{12} \frac{\partial T}{\partial t} = q_{1} (t)$ ,

(10)

$g_{21} \frac{\partial^{2} p}{\partial x^{2}} + g_{22} \frac{\partial^{2} T}{\partial x^{2}} + c_{21} \frac{\partial p}{\partial t} + c_{22} \frac{\partial T}{\partial t} = q_{2} (t)$ 。

(11)

式中, $g_{11} = \frac{k_{w}}{γ_{w}}, g_{12} = S_{w}, c_{11} = - \frac{α^{2}}{E_{c}} - α_{p}, c_{12} = Y - \frac{α β}{E_{c}}$ , $g_{21} = \frac{S_{T}}{T_{0}}, g_{22} = \frac{K_{T}}{T_{0}}, c_{21} = Y - \frac{α β}{E_{c}}, c_{22} = - \frac{β^{2}}{E_{c}} - \frac{c_{v}}{T_{0}}$ 。

2. 半无限空间饱和土一维热固结

2.1 半无限空间饱和土一维热固结定解问题

为研究季节性温度变化对半无限空间饱和土体固结的影响,假设土体表面作用有大范围的均布温度荷载,饱和土的温度传导、渗流与变形均发生在沿深度方向,简化为一维热固结问题。设土体表面为坐标原点,坐标轴x方向垂直向下,建立如图1所示的半无限空间饱和土一维热固结计算模型。

图 1 半无限空间饱和土一维热固结计算模型

Figure 1. Model for one-dimensional thermal consolidation of semi-infinite saturated soils

下载: 全尺寸图片幻灯片

设土体位移、超静孔压与温度变化量的初始条件与无穷远处的自然边界条件分别为

$u (x, 0) = 0, p (x, 0) = 0, T (x, 0) = 0$ ,

(12)

$u (\infty, t) = 0, p (\infty, t) = 0, T (\infty, t) = 0$ 。

(13)

考虑上表面为自由边界,超静孔压和温度变化量为以下两类边界条件：

$p (0, t) = f_{1} (t), T (0, t) = f_{2} (t) (边界条件 (A))$ ,

(14)

$p_{, x} (0, t) = f_{1} (t), T_{, x} (0, t) = f_{2} (t) (边界条件 (B))$ 。

(15)

式中,f₁(t)和f₂(t)为任意关于时间t的函数。边界条件（A）代表已知半无限空间土体表面超静孔压p和温度改变量T随时间变化的工况；根据修正的Darcy定律和Fourier定律,边界条件（B）代表已知半无限空间土体表面渗流流量q_w和热流流量q_T随时间变化的工况。

式（10）～（15）分别给出了饱和土热固结的控制方程,边界条件和初始条件,构成了半无限空间饱和土一维热固结的定解问题。

2.2 求解过程

采用正弦变换和余弦变换方法对上述定解问题进行求解,正弦变换、余弦变换及其相应逆变换的定义分别为

$\begin{array}{l} s (ω) = \int_{0}^{\infty} f (x) \sin (ω x) d x, \\ f (x) = \frac{2}{π} \int_{0}^{\infty} s (ω) \sin (ω x) d ω 。 \end{array}}$

(16)

$\begin{array}{l} c (ω) = \int_{0}^{\infty} f (x) \cos (ω x) d x, \\ f (x) = \frac{2}{π} \int_{0}^{\infty} c (ω) \cos (ω x) d ω 。 \end{array}}$

(17)

（1）边界条件(A)的求解

对超静孔压p和温度改变量T实行正弦变换,定义为

${\bar{p}}_{s} (ω, t) = \int_{0}^{\infty} p (x, t) \sin (ω x) d x$ ,

(18)

${\bar{T}}_{s} (ω, t) = \int_{0}^{\infty} T (x, t) \sin (ω x) d x$ 。

(19)

根据正弦变换的性质,结合边界条件,可将式（10）,（11）转化为如下形式：

$\begin{matrix} g_{11} [- ω^{2} {\bar{p}}_{s} (ω, t) + ω p (0, t)] + g_{12} [- ω^{2} {\bar{T}}_{s} (ω, t) + \\ ω T (0, t)] + c_{11} \frac{\partial {\bar{p}}_{s} (ω, t)}{\partial t} + c_{12} \frac{\partial {\bar{T}}_{s} (ω, t)}{\partial t} = {\bar{q}}_{s1} \end{matrix}$ ,

(20)

$\begin{array}{l} g_{21} [- ω^{2} {\bar{p}}_{s} (ω, t) + ω p (0, t)] + g_{22} [- ω^{2} {\bar{T}}_{s} (ω, t) + \\ ω T (0, t)] + c_{21} \frac{\partial {\bar{p}}_{s} (ω, t)}{\partial t} + c_{22} \frac{\partial {\bar{T}}_{s} (ω, t)}{\partial t} = {\bar{q}}_{s2} \end{array}$ 。

(21)

式（20）,（21）记为矩阵微分方程形式,即

$C \frac{d X (t)}{d t} = ω^{2} G X (t) + H_{s} (t)$ ,

(22)

式中,

$X (t) = {\begin{matrix} {\bar{p}}_{s} (ω, t) \\ {\bar{T}}_{s} (ω, t) \end{matrix}}, C = [\begin{matrix} c_{11} & c_{12} \\ c_{21} & c_{22} \end{matrix}], G = [\begin{matrix} g_{11} & g_{12} \\ g_{21} & g_{22} \end{matrix}]$

$H_{s} = {\begin{matrix} {\bar{q}}_{s1} - g_{11} ω f_{1} (t) - g_{12} ω f_{2} (t) \\ {\bar{q}}_{s 2} - g_{21} ω f_{1} (t) - g_{22} ω f_{2} (t) \end{matrix}}$ 。

对应的初始条件的正弦变换为

${\bar{p}}_{s} (ω, 0) = 0, {\bar{T}}_{s} (ω, 0) = 0$ 。

(23)

（2）边界条件(B)的求解

与边界条件(A)的求解类似,对超静孔压p和温度改变量T进行余弦变换,记为

${\bar{p}}_{c} (ω, t) = \int_{0}^{\infty} p (x, t) \cos (ω x) d x$ ,

(24)

${\bar{T}}_{c} (ω, t) = \int_{0}^{\infty} T (x, t) \cos (ω x) d x$ 。

(25)

根据余弦变换的性质,结合边界条件,可将式（10）,（11）转化为如下形式：

$\begin{array}{l} g_{11} [- ω^{2} {\bar{p}}_{c} (ω, t) - p_{, x} (0, t)] + g_{12} [- ω^{2} {\bar{T}}_{c} (ω, t) - \\ T_{, x} (0, t)] + c_{11} \frac{\partial {\bar{p}}_{s} (ω, t)}{\partial t} + c_{12} \frac{\partial {\bar{T}}_{s} (ω, t)}{\partial t} = {\bar{q}}_{c1} \end{array}$ ,

(26)

$\begin{array}{l} g_{21} [- ω^{2} {\bar{p}}_{c} (ω, t) - p_{, x} (0, t)] + g_{22} [- ω^{2} {\bar{T}}_{c} (ω, t) - \\ T_{, x} (0, t)] + c_{21} \frac{\partial {\bar{p}}_{s} (ω, t)}{\partial t} + c_{22} \frac{\partial {\bar{T}}_{s} (ω, t)}{\partial t} = {\bar{q}}_{c2} \end{array}$ 。

(27)

将式（26）,（27）记为矩阵微分方程形式,即

$C \frac{d X (t)}{d t} = ω^{2} G X (t) + H_{c} (t)$ ,

(28)

式中,

$X (t) = {\begin{matrix} {\bar{p}}_{c} (ω, t) \\ {\bar{T}}_{c} (ω, t) \end{matrix}}, C = [\begin{matrix} c_{11} & c_{12} \\ c_{21} & c_{22} \end{matrix}], G = [\begin{matrix} g_{11} & g_{12} \\ g_{21} & g_{22} \end{matrix}]$ ,

$H_{c} = {\begin{matrix} {\bar{q}}_{c1} + g_{11} f_{1} (t) + g_{12} f_{2} (t) \\ {\bar{q}}_{c 2} + g_{21} f_{1} (t) + g_{22} f_{2} (t) \end{matrix}}$ 。

对应的初始条件的余弦变换为

${\bar{p}}_{c} (ω, 0) = 0, {\bar{T}}_{c} (ω, 0) = 0$ 。

(29)

（3）常微分方程的解

式（22）,（28）均为关于时间t的常微分方程,其可统一表示为标准形式：

$\frac{d}{d t} X (t) = ω^{2} N X (t) + H (t)$ ,

(30)

式中,N=C^-1G,对边界条件(A),H=C^-1H_s；对边界条件(B),H=C^-1H_c。

方程（30）有以下形式的解：

$X (t) = e^{ω^{2} N t} X (0) + e^{ω^{2} N t} \int_{0}^{t} e^{- ω^{2} N ξ} H (ξ)d ξ$ 。

(31)

根据Hamilton-Cayley定理, $e^{N t}$ 写为如下形式：

$e^{N t} = α_{0} (t) I + α_{1} (t) N$ ,

(32)

其中,

$α_{0} (t) = \frac{r_{1} e^{r_{2} t} - r_{2} e^{r_{1} t}}{r_{1} - r_{2}}, α_{1} (t) = \frac{e^{r_{1} t} - e^{r_{2} t}}{r_{1} - r_{2}}$ ,

(33)

式中,r₁和r₂为矩阵N的特征值。当r₁=r₂时,其表达式应为

$α_{0} (t) = (1 - r_{1} t) e^{r_{1} t}, α_{1} (t) = t e^{r_{1} t}$ 。

(34)

因此,将初始条件（23）,（29）, $e^{N t}$ 与H(t)代入式（31）,逐步计算即可求得超静孔压p和温度改变量T在变换域上的解,再根据相对应的正弦逆变换（16）或余弦逆变换（17）即可得到其对应的时域精确显式积分解。将超静孔压p和温度改变量T的解代入式（7）,即可得到土体体应变 $ε_{v}$ 的精确显式积分解。将体应变沿深度方向积分,即可得到饱和土体各深度处的位移解u。

3. 算例分析

3.1 饱和土体热–水–力耦合物性参数

假定土体的初始温度T₀=300 K,表面外荷载q(t) =0。参考Bai^[11]给出的土体物性参数,本文土体的热–水–力基本参数取值见表1。

表 1 饱和土体的材料参数

Table 1. Material parameters of saturated soils

孔隙率n	泊松比μ	土体弹性模量E/MPa	土颗粒体积模量K_s/GPa	土颗粒密度ρ_s/(kg·m^-3)	土颗粒比热容c_s/(J·kg^-1·K^-1)	土颗粒热膨胀系数α_s/K^-1	孔隙流体体积模量K_w/GPa	孔隙流体密度ρ_w/(kg·m^-3)	孔隙流体比热容c_w/(J·kg^-1·K^-1)	孔隙流体热膨胀系数α_w/K^-1	土体热传导系数K_T/(J·s^-1·m^-1·K^-1)	热渗系数S_w/(m²·s^-1·K^-1)
0.4	0.3	5	20	2600	800	1.5×10^-5	5	1000	4200	2×10^-5	1	5×10^-11

下载: 导出CSV

| 显示表格

因为土体某深度处温度升高越快,渗流越慢,超静孔压响应越大；相反,温度升高越慢,渗流越快,超静孔压响应越小,所以饱和土中的超静孔压响应是渗流和温升快慢的耦合作用结果。Bai^[12]在研究中引入了无量纲固结系数与热扩散系数的比值c/κ来表示两者的相对速度。本文类似引入简化的无量纲相对扩散数 $v_{1} / v_{2}$ ,其中v₁为太沙基一维固结理论定义的固结系数, $v_{1}$ =k_wE_s/ $γ_{w}$ ； $v_{2}$ 为根据热传导方程定义的温度变化系数, $v_{2}$ =K_T/c_v。

3.2 解析解正确性验证

Mctigue^[3]提出直接求解非耦合的热传导方程,得到表面温度突增时,半无限饱和土体内的温度响应,再将温度场的解代入渗流场方程,得到土体内超静孔压的解析解。为验证本文解析解的正确性,将本文解析方法计算得到的土体热固结响应结果与Mctigue^[3]方法的结果进行对比。所以,将本文温度场控制方程简化为热传导方程,土体的初始条件和无穷远处的边界条件不变,表面边界条件设为透水边界,施加突增温度荷载：

$p (0, t) = 0, T (0, t) = 25$ 。

(35)

图2,3给出了分别采用本文的积分解方法和Mctigue^[3]方法得到的不同时刻土体温度改变量T和超静孔压p随深度的分布。由图3可以看出,两种方法计算得到的结果几乎相同,所以本文解析解的正确性可以得到验证。图2,3还反映了突增温度荷载作用下的半无限空间饱和土体的一维热固结响应特性,即温度从表面向深处传播,较深处温度增长比较浅处慢；各深度处温度上升引起的超静孔压峰值几乎相同,并呈波动特征向土体深处传递。

图 2 表面温度突增时不同时刻温度分布图

Figure 2. Distribution of temperature change along depth at different moments subjected to step temperature increase

下载: 全尺寸图片幻灯片

图 3 表面温度突增时不同时刻超静孔压分布图

Figure 3. Distribution of excess pore pressure along depth at different moments subjected to step temperature increase

下载: 全尺寸图片幻灯片

3.3 半无限空间饱和土表面季节性温度荷载

根据半无限空间饱和土体表面是否有覆盖层分为两类边界条件。当土体表面无覆盖层时,土体表面为透水边界,由于土体表面与大气直接接触,受阳光直射,故假定土体表面温度变化与气温相似,为分析季节性周期变化的温度荷载对饱和土体的热固结影响,取土体表面温度周期变化,幅值为25℃,周期为T_p=π×10⁷ s,约为1 a,即边界条件(A)。当土体表面有覆盖层时,如机场跑道、高等级路面等表面有道面板等不透水结构,设定土体表面为不透水边界,又因为覆盖层的存在,外界温度变化无法直接改变土体表面温度,而是通过热量向下传递的形式产生影响,因此根据Fourier定律,假定土体上表面温度荷载以周期变化的热流量表示,即边界条件(B),热流量变化周期与边界条件(A)相同,为T_p=π×10⁷ s,约为1 a。为分析季节性温度荷载作用下,表面渗透条件对饱和土热固结响应的影响,通过计算使热流边界条件(B)情况下得到的土体表面温度变化响应与边界条件(A)表面周期变化温度近似相同：

$p (0, t) = 0, T (0, t) = 25 \sin (2 t / 10^{7}) (边界条件 (A))$ ,

(36)

$p_{x} (0, t) = 0, K_{T} T_{x} (0, t) = - 20 \cos (2 t / 10^{7}) (边界条件 (B)) 。$

(37)

3.4 半无限空间饱和土体一维热固结特性分析

（1）边界条件(A)

图4,5分别给出了半无限饱和土表面透水时,不同深度处温度改变量T和超静孔压p的变化过程。由于表面温度周期性变化,土体不同深度处的温度改变量T和超静孔压p均呈周期性变化,从温度荷载施加后的瞬态响应渐渐过渡到稳定状态,其稳态的变化周期与表面温度周期相同。随着深度增加,温度改变量T变化的幅值减小,且有滞后的趋势。深度0 m处温度变化幅值为25 K,深度5 m处的温度变化不超过±1.98 K。

图 4 边界条件(A)不同深度处温度变化过程

Figure 4. Evolution of temperature change with time at different depths under boundary condition A

下载: 全尺寸图片幻灯片

图 5 边界条件(A)不同深度处超静孔压变化过程

Figure 5. Evolution of excess pore pressure with time at different depths under boundary condition A

下载: 全尺寸图片幻灯片

由于表面为透水边界,温度变化产生的超静孔压可以向两个方向消散,超静孔压p的变化幅值沿深度先增大后减小。深度2 m处,超静孔压稳态变化幅值最大约为6.76 kPa。图6还给出了表面处(x=0)位移随时间的变化,显然位移亦随时间周期性波动,表面最大位移可达4.4 mm。

图 6 边界条件(A)表面位移变化过程

Figure 6. Evolution of surface displacement of soils with time under boundary condition A

下载: 全尺寸图片幻灯片

（2）边界条件(B)

图7,8分别给出了半无限空间饱和土表面为边界条件（B）时,不同深度处温度改变量T和超静孔压p的变化过程。由图7,8可知,土体不同深度处的温度改变量T和超静孔压p受表面周期变化的热流影响而呈现周期变化的特性,稳态变化的周期与表面热荷载周期相同。不同于边界条件（A）的结果,土体温度改变量T和超静孔压p周期变化的幅值都随着深度增大而减小。深度0 m处温度稳态周期变化的幅值最大,约为26.1 K,较边界条件（A）深度0 m处幅值25 K高约4.4%。而深度0 m处超静孔压稳态周期变化的幅值最大,约为16.46 kPa,为边界条件（A）深度2 m处稳态幅值6.76 kPa的2.43倍。原因是深度0 m处温度波动幅值最大,且表面为不透水边界,温度产生的超静孔压只能向无穷远处消散,排水路经长导致了表面的超静孔压相对较大。图9反映了表面处（x=0）位移随时间变化的趋势,表面位移随时间周期性变化,变化最大值约为8.18 mm,比边界条件（A）表面的位移最大值高出85.9%。

图 7 边界条件(B)不同深度处温度变化过程

Figure 7. Evolution of temperature change with time at different depths under boundary condition B

下载: 全尺寸图片幻灯片

图 8 边界条件(B)不同深度处超静孔压变化过程

Figure 8. Evolution of excess pore pressure with time at different depths under boundary condition B

下载: 全尺寸图片幻灯片

图 9 边界条件(B)表面位移变化过程

Figure 9. Evolution of surface displacement of soils with time under boundary condition B

下载: 全尺寸图片幻灯片

综合边界条件（A）和边界条件（B）的结果可知,虽然两种边界条件下的土体温度变化特性类似,但表面不透水边界较透水边界,土体浅层的超静孔压更难以消散,孔压的累积更有可能导致浅层土体发生强度破坏,在机场道路,高等级路面等地基设计与施工过程中应引起重视。

（3）扩散数 $v_{1} / v_{2}$ 的影响

由于在表面不透水情况下,浅层土体更易产生较大的超静孔压和变形,因此本部分保持 $v_{2}$ 不变,采用不同的固结系数 $v_{1}$ ,对边界条件（B）情况下,土体的热固结响应进行进一步研究。图10～12分别给出了不同扩散数 $v_{1} / v_{2}$ 条件下,深度0 m处土体温度改变量T、超静孔压p和位移u的响应。可以发现不同的 $v_{1}$ 对表面温度变化和位移几乎没有影响,但扩散数增大5倍后,超静孔压p的峰值16.46 kPa减小为9.55 kPa,降低了约41.2%；扩散数减小10倍后,超静孔压p的峰值增大至27.48 kPa,提高了约66.9%。不同 $v_{1}$ 对表面位移没有影响的原因在于同样的温度变化产生的总超静孔压是相同的, $v_{1}$ 的不同仅导致了超静孔压沿深度方向更加均匀的分布,如图13。由于假设了土体沿深度方向的压缩模量E_s不变,上表面位移是体应变 $ε_{v}$ 沿深度方向的积分,所以扩散数对土体表面位移几乎无影响。

图 10 边界条件(B)不同扩散数表面温度随时间变化

Figure 10. Evolution of temperature of surface with time with different coefficients under boundary condition B

下载: 全尺寸图片幻灯片

图 11 边界条件(B)不同扩散数表面超静孔压随时间变化

Figure 11. Evolution of excess pore pressure of surface with time with different coefficients under boundary condition B

下载: 全尺寸图片幻灯片

图 12 边界条件(B)不同扩散数表面位移随时间变化

Figure 12. Evolution of displacement of surface with time with different coefficients under boundary condition B

下载: 全尺寸图片幻灯片

图 13 边界条件(B)不同扩散数超静孔压随深度变化

Figure 13. Distribution of excess pore pressure along depth with different coefficients under boundary condition B

下载: 全尺寸图片幻灯片

4. 结论

本文先给出了考虑热渗效应和等温热流效应的饱和土热–水–力耦合固结控制方程,再利用正弦变换和余弦变换,提出了两类任意非齐次边界条件下,半无限空间饱和土体一维热固结响应的精确积分形式解的求解方法,避免了复杂边界Laplace数值逆变换带来的误差。通过算例验证了本文解的正确性,分析了表面透水和不透水的半无限空间饱和土体在季节性周期变化温度荷载作用下的热响应,得到以下4点结论。

（1）土体内部温度从瞬态响应逐渐过渡到稳态周期波动,波动的周期与表面温度荷载相同。温度波动的幅值随深度增加而减小,且深处温度变化滞后于浅处温度变化。

（2）土体内部超静孔压响应呈周期性波动,变化的周期与表面温度荷载相同。表面透水时,孔压波动的幅值随深度增大先增大后减小；表面不透水时,孔压波动的幅值随深度增大而减小。表面不透水情况的最大孔压比表面透水情况高出1倍多。

（3）土体上表面位移随时间周期性变化。表面不透水情况的位移最大值明显高于表面透水情况。

（4）表面不透水且外界温度荷载相同时,土体渗透性越好,深度0 m处的超静孔压波动幅值越低,且超静孔压沿深度方向的分布更趋于均匀化,使得上表面的位移几乎不受影响。

图 1 半无限空间饱和土一维热固结计算模型

Figure 1. Model for one-dimensional thermal consolidation of semi-infinite saturated soils

下载: 全尺寸图片幻灯片

图 2 表面温度突增时不同时刻温度分布图

Figure 2. Distribution of temperature change along depth at different moments subjected to step temperature increase

下载: 全尺寸图片幻灯片

图 3 表面温度突增时不同时刻超静孔压分布图

Figure 3. Distribution of excess pore pressure along depth at different moments subjected to step temperature increase

下载: 全尺寸图片幻灯片

图 4 边界条件(A)不同深度处温度变化过程

Figure 4. Evolution of temperature change with time at different depths under boundary condition A

下载: 全尺寸图片幻灯片

图 5 边界条件(A)不同深度处超静孔压变化过程

Figure 5. Evolution of excess pore pressure with time at different depths under boundary condition A

下载: 全尺寸图片幻灯片

图 6 边界条件(A)表面位移变化过程

Figure 6. Evolution of surface displacement of soils with time under boundary condition A

下载: 全尺寸图片幻灯片

图 7 边界条件(B)不同深度处温度变化过程

Figure 7. Evolution of temperature change with time at different depths under boundary condition B

下载: 全尺寸图片幻灯片

图 8 边界条件(B)不同深度处超静孔压变化过程

Figure 8. Evolution of excess pore pressure with time at different depths under boundary condition B

下载: 全尺寸图片幻灯片

图 9 边界条件(B)表面位移变化过程

Figure 9. Evolution of surface displacement of soils with time under boundary condition B

下载: 全尺寸图片幻灯片

图 10 边界条件(B)不同扩散数表面温度随时间变化

Figure 10. Evolution of temperature of surface with time with different coefficients under boundary condition B

下载: 全尺寸图片幻灯片

图 11 边界条件(B)不同扩散数表面超静孔压随时间变化

Figure 11. Evolution of excess pore pressure of surface with time with different coefficients under boundary condition B

下载: 全尺寸图片幻灯片

图 12 边界条件(B)不同扩散数表面位移随时间变化

Figure 12. Evolution of displacement of surface with time with different coefficients under boundary condition B

下载: 全尺寸图片幻灯片

图 13 边界条件(B)不同扩散数超静孔压随深度变化

Figure 13. Distribution of excess pore pressure along depth with different coefficients under boundary condition B

下载: 全尺寸图片幻灯片

表 1 饱和土体的材料参数

Table 1 Material parameters of saturated soils

孔隙率n	泊松比μ	土体弹性模量E/MPa	土颗粒体积模量K_s/GPa	土颗粒密度ρ_s/(kg·m^-3)	土颗粒比热容c_s/(J·kg^-1·K^-1)	土颗粒热膨胀系数α_s/K^-1	孔隙流体体积模量K_w/GPa	孔隙流体密度ρ_w/(kg·m^-3)	孔隙流体比热容c_w/(J·kg^-1·K^-1)	孔隙流体热膨胀系数α_w/K^-1	土体热传导系数K_T/(J·s^-1·m^-1·K^-1)	热渗系数S_w/(m²·s^-1·K^-1)
0.4	0.3	5	20	2600	800	1.5×10^-5	5	1000	4200	2×10^-5	1	5×10^-11

下载: 导出CSV

参考文献(28)

[1]	BOOKER J R, SAVVIDOU C. Consolidation around a spherical heat source[J]. International Journal of Solids and Structures, 1984, 20(11/12): 1079-1090.
[2]	BOOKER J R, SAVVIDOU C. Consolidation around a point heat source[J]. International Journal for Numerical and Analytical Methods in Geomechanics, 1985, 9(2): 173-184. doi: 10.1002/nag.1610090206
[3]	MCTIGUE D F. Thermoelastic response of fluid-saturated porous rock[J]. Journal of Geophysical Research, 1986, 91(B9): 9533-9542. doi: 10.1029/JB091iB09p09533
[4]	GIRAUD A, ROUSSET G. Thermoelastic and thermoplastic response of a porous space submitted to a decaying heat source[J]. International Journal for Numerical and Analytical Methods in Geomechanics, 1995, 19(7): 475-495. doi: 10.1002/nag.1610190703
[5]	SMITH D W, BOOKER J R. Green’s functions for a fully coupled thermoporoelastic material[J]. International Journal for Numerical and Analytical Methods in Geomechanics, 1993, 17(3): 139-163. doi: 10.1002/nag.1610170302
[6]	SMITH D W, BOOKER J R. Boundary element analysis of linear thermoelastic consolidation[J]. International Journal for Numerical and Analytical Methods in Geomechanics, 1996, 20(7): 457-488. doi: 10.1002/(SICI)1096-9853(199607)20:7<457::AID-NAG805>3.0.CO;2-H
[7]	BAI M, ABOUSLEIMAN Y. Thermoporoelastic coupling with application to consolidation[J]. International Journal for Numerical and Analytical Methods in Geomechanics, 1997, 21(2): 121-132. doi: 10.1002/(SICI)1096-9853(199702)21:2<121::AID-NAG861>3.0.CO;2-W
[8]	ZHOU Y, RAJAPAKSE R K N D, GRAHAM J. A coupled thermoporoelastic model with thermo-osmosis and thermal filtration[J]. International Journal of Solids and Structures, 1998, 35(34/35): 4659-4683.
[9]	ZHOU Y, RAJAPAKSE R K N D, GRAHAM J. Coupled consolidation of a porous medium with a cylindrical or a spherical cavity[J]. International Journal for Numerical and Analytical Methods in Geomechanics, 1998, 22(6): 449-475. doi: 10.1002/(SICI)1096-9853(199806)22:6<449::AID-NAG930>3.0.CO;2-H
[10]	BLOND E, SCHMITT N, HILD F. Response of saturated porous media to cyclic thermal loading[J]. International Journal for Numerical and Analytical Methods in Geomechanics, 2003, 27(11): 883-904. doi: 10.1002/nag.301
[11]	白冰. 循环温度荷载作用下饱和多孔介质热-水-力耦合响应[J]. 工程力学, 2007, 24(5): 87-92. doi: 10.3969/j.issn.1000-4750.2007.05.015 BAI Bing. Thermo-hydro-mechanical responses of saturated porous media under cyclic thermal loading[J]. Engineering Mechanics, 2007, 24(5): 87-92. (in Chinese) doi: 10.3969/j.issn.1000-4750.2007.05.015
[12]	BAI B. Fluctuation responses of saturated porous media subjected to cyclic thermal loading[J]. Computers and Geotechnics, 2006, 33(8): 396-403. doi: 10.1016/j.compgeo.2006.08.005
[13]	BAI B. Thermal consolidation of layered porous half-space to variable thermal loading[J]. Applied Mathematics and Mechanics, 2006, 27(11): 1531-1539. doi: 10.1007/s10483-006-1111-1
[14]	LIU G B, XIE K H, ZHENG R Y. Model of nonlinear coupled thermo-hydro-elastodynamics response for a saturated poroelastic medium[J]. Science in China Series E: Technological Sciences, 2009, 52(8): 2373-2383. doi: 10.1007/s11431-008-0220-8
[15]	LIU G B, XIE K H, ZHENG R Y. Thermo-elastodynamic response of a spherical cavity in saturated poroelastic medium[J]. Applied Mathematical Modelling, 2010, 34(8): 2203-2222. doi: 10.1016/j.apm.2009.10.031
[16]	HE L W, JIN Z H. A local thermal nonequilibrium poroelastic theory for fluid saturated porous media[J]. Journal of Thermal Stresses, 2010, 33(8): 799-813. doi: 10.1080/01495739.2010.482358
[17]	YANG Y, GUERLEBECK K, SCHANZ T. Thermo-osmosis effect in saturated porous medium[J]. Transport in Porous Media, 2014, 104(2): 253-271. doi: 10.1007/s11242-014-0332-5
[18]	AI Z Y, YE Z, ZHAO Z, et al. Time-dependent behavior of axisymmetric thermal consolidation for multilayered transversely isotropic poroelastic material[J]. Applied Mathematical Modelling, 2018, 61: 216-236. doi: 10.1016/j.apm.2018.04.012
[19]	AI Z Y, YE Z, LI Y. Quasti-static response of a multilayered transversely isotropic porothermoelastic material subjected to a cylindrical heat source[J]. Computers and Geotechnics, 2019, 107: 201-213. doi: 10.1016/j.compgeo.2018.11.024
[20]	AI Z Y, WANG L J. Precise solution to 3D coupled thermohydromechanical problems of layered transversely isotropic saturated porous media[J]. International Journal of Geomechanics, 2018, 18(1): 4017121. doi: 10.1061/(ASCE)GM.1943-5622.0001027
[21]	WANG L J, ZHU B, CHEN Y M, CHEN R P, SHI X S. Precise model for predicting excess pore-water pressure of layered soils induced by thermal-mechanical loads[J]. Journal of Engineering Mechanics, 2019, 145(1): 4018114. doi: 10.1061/(ASCE)EM.1943-7889.0001544
[22]	王路君, 艾智勇. 衰变热源作用下饱和多孔介质热固结问题的扩展精细积分法[J]. 力学学报, 2017, 49(2): 324-334. https://www.cnki.com.cn/Article/CJFDTOTAL-LXXB201702009.htm WANG Lu-jun, AI Zhi-yong. EPIM for thermal consolidation problems of saturated porous media subjected to a decaying heat source[J]. Chinese Journal of Theoretical and Applied Mechanics, 2017, 49(2): 324-334. (in Chinese) https://www.cnki.com.cn/Article/CJFDTOTAL-LXXB201702009.htm
[23]	钮家军, 凌道盛, 王秀凯, 等. 饱和单层土体一维热固结精确解[J]. 岩土工程学报, 2019, 41(9): 1715-1723. https://www.cnki.com.cn/Article/CJFDTOTAL-YTGC201909018.htm NIU Jia-jun, LING Dao-sheng, WANG Xiu-kai, et al. Exact solutions for one-dimensional thermal consolidation of single-layer saturated soil[J]. Chinese Journal of Geotechnical Engineering, 2019, 41(9): 1715-1723. (in Chinese) https://www.cnki.com.cn/Article/CJFDTOTAL-YTGC201909018.htm
[24]	COUSSY O. Poromechanics[M]. Chichester: John Wiley & Sons Ltd, 2004.
[25]	BAI B. Thermal response of saturated porous spherical body containing a cavity under several boundary conditions[J]. Journal of Thermal stresses, 2013, 36(11): 1217-1232. doi: 10.1080/01495739.2013.788389
[26]	BAI B, LI T. Irreversible consolidation problem of a saturated porothermoelastic spherical body with a spherical cavity[J]. Applied Mathematical Modelling, 2013, 37(4): 1973-1982. doi: 10.1016/j.apm.2012.05.003
[27]	CARNAHAN C L. Thermodynamic coupling of heat and matter flows in near-field regions of nuclear waste repositories[J]. MRS Proceedings, 1983, 26(1): 1023-1030.
[28]	HORSEMAN S T, MCEWEN T J. Thermal constraints on disposal of heat-emitting waste in argillaceous rocks[J]. Engineering Geology, 1996, 41(1/2/3/4): 5-16.

施引文献(3)

期刊类型引用(2)

1.	闻敏杰，王宗琴，涂园，刘开富. 考虑界面热阻效应的成层饱和土地基热固结分析. 应用基础与工程科学学报. 2025(01): 213-222 . 百度学术
2.	江文豪，李江山，黄啸，程鑫，万勇. 非等温分布条件下考虑半透水边界时饱和黏土的一维固结解析解. 岩土力学. 2022(10): 2744-2756 . 百度学术