0 引言

Terzaghi^[1]一维固结理论为现代土力学的发展奠定了基础，在工程中应用也最为广泛。但随着工程地质条件越来越复杂，对一些工程性质较差的特殊土层，Terzaghi一维固结理论计算与实测值间存在较大偏差，造成偏差的主要原因有未考虑大变形和次固结问题。在Terzaghi一维固结理论的基本假设中，渗透系数k和压缩系数a假设为常数值，而在工程实际中，渗透系数k、压缩系数a和孔隙比e是随有效应力增大而变化的，且变化幅值较大^[2-6]；特别是高压缩性软土除具有明显的大变形特性外^[7-8]，还具有较强的流变特性，即在主固结完成后，次固结沉降显著^[9-10]。故Davis等^[11-12]采用e-lg p曲线改进了压缩变形的线性关系，获得一维非线性固结方程的解析解。Gibson等^[13]在固结理论中引入了大应变固结理论，推导了大变形固结基本控制方程；Poskitt^[14]、谢康和等^[15]均考虑了k，a的变化，得到相近的一维大变形固结方程解析解。李传勋等^[16-17]推导出基于非达西渗流软土一维非线性固结控制方程，进而分析了非达西渗流下的固结速率。总之在一维非线性固结理论的研究中，主要难点是提出的非线性固结微分方程只有在特殊情况下才可求得解析解，面对这一问题，引入了半解析法和有限差分法，李冰河等^[18]采用比较公认的e-lg σ和e-lg k_v关系式，用有限差分法得到了非线性固结理论的半解析解。Mesri^[19]，Chen^[20]，Abbasi等^[21]运用这两种方法求解一维非线性固结方程取得了很多成果。近年来，有关学者针对土的竖向固结系数也进行了修正研究，高俊等^[22]提出初始状态下的渗透系数预测公式，对固结系数表c_v计算式进行修正，从而提高Terzaghi一维固结理论的精确度。

半解析法是把解函数一部分用精确解求解，一部分数值逼近，所以能得到一个半解析解，半解析解精度低于解析解，但是高于数值解，也就说它是一个精度更高的数值解，而有限差分法属于数值分析法，解为数值解。对于多元且复杂的偏微分方程，数值解不仅计算精度低误差大，而且适用范围小，故研究非线性固结微分方程的解析解意义重大。而针对于土固结系数c_v的修正方面，没有将孔隙比e和渗透系数k随时间变化的规律作为非线性固结理论的先天条件，仅修正c_v这个参数，与实际还是存在很大误差。

现以泥炭质土为研究对象，选取滇池附近2个场地的泥炭质土进行固结-渗透联合试验，对其e-σ′和k-σ′关系进行研究，提出适合于泥炭质土的非线性固结理论模型，并获得一定条件下的解析解，从而对Terzaghi一维固结理论进行进一步的修正，最后结合试验分析探讨修正固结方程的合理性。

1 土样基本性质及分类

试验土样取自昆明滇池周围，属第四系全新统冲湖积层，颜色灰黑色，高压缩性，孔隙比e, 含水率和有机质含量高。按照ASTM^[23]等对泥炭质土分类方法，滇池湖相沉积泥炭质土属高分解度土。利用所选场地正在进行基坑开挖的便利，采用人工取原状块样到实验室，取样时尽量保证每组土样性质差异性最低。3^#1组取自昆明市广福路和滇池路交叉口附近的锦慧金融中心，距离滇池水域1.5 km左右，3^#2，3^#3，3^#4组取自昆明市官渡区渔村社区，距离滇池水域1 km左右，试验土样的物理力学指标如表 1所示。

2 试验方案

研究高应力范围内孔隙比和有效应力之间的相关关系时，因改进STQY1-3型渗压仪施加的固压范围有限，故采用常规固结试验仪，按0→12.5 →25→50→100→200→300→400→600→800→1 200→1 600 kPa的顺序加载，每级荷载下持续24 h。试验步骤按试验标准^[24]进行。

固结-渗透联合试验采用集压缩渗透于一体的改进STQY1-3型渗压仪，能施加的最大固结应力1 000 kPa，最大渗透压力为250 kPa，出力灵敏度为1 kPa，试验采用专用环刀。土样面积30 cm²，高2 cm。试验前土样进行抽气饱和，饱和完成后测饱和度是否满足要求，平行试验3组。固结压力加载级数为25，50，100，200，300，400，500，600 kPa，每级固压下试样变形稳定标准为：变形速率小于0.01 mm/d，双向排水。打开渗水阀门，通过空气压缩机施加远小于固结压力的渗透压力，进行该级压力下的变水头固结渗透试验，渗透系数k按式(1)计算:

(1)

式中，a为变水头管的内横截面积；A为试样面积；t₁，t₂为前后两次的测试时间；h₁，h₂分别为t₁，t₂对应的初始及终止水头和试样高度。

渗透试验具体步骤如下：

(1) 采用专用环刀取样后测定试样密度，抽气饱和，饱和完成后测试饱和度是否满足要求。

(2) 清洗缸体，对止水圈、螺母等密封设备进行密封性检查。

(3) 将饱和后的试样装进缸体，密封试样，安装完成后，调节百分表，使量程为0~10 mm的指针在0~9之间，精度0.01 mm的指针调节到0。

(4) 测定原状土渗透系数k₀后，施加相应的固结应力。固结标准为：每级加压时长必须超过24 h，且变形速率不大于0.01 mm/h。固结稳定后，向缸体内充满水，根据固压和不同土质施加相应的渗透压差，打开A，B管水头控制阀，为保证渗流稳定，等5 min开始记录玻璃管液面下降高度与时间的关系，并整理数据。

(5) 卸压、拆除仪器，计算k及数据分析。

3 分析模型

泥炭质土随着施加固压的不同会呈现不同孔隙比e，e随施加固压的增大呈减小趋势，e和渗透系数k间存在对应关系，国内外学者在建立k与e之间的非线性渗透模型上取得了大量研究成果^[26-29]。同理在k和有效应力σ′之间也存在一个确定的对应关系，本研究以试验数据为依托对泥炭质土的σ′和k之间对应关系进行研究，并建立σ′和k之间的模型。

3.1 孔隙比变化分析

在Terzaghi一维固结理论中，固结推导过程中将e和σ′间的数值关系视为线性，即de=-adσ′，a保持不变。根据泥炭质土在高应力条件下的固结试验结果，分别作出e-σ′及e-lg σ′曲线，如图 1所示。

从图 1中可看出，泥炭质土e-σ′曲线的线性拟合度较低，采用e-σ′线性关系进行计算会产生很大误差。无论是重塑土还是原状土压缩系数a均随σ′的增大而减小，即Terzaghi固结理论存在假设的a值前期偏小后期偏大的情况。对比分析e-σ′和e-lg σ′曲线发现，用线性拟合e-lg σ′曲线的相关系数R²比e-σ′曲线的高，R²都大于0.8，重塑土的线性相关度比原状土的高。

低应力(25~600 kPa)下泥炭质土的e-lg σ′曲线如图 2所示。

从图 2可看出，其e随σ′变化规律与在高应力条件下相同。高应力和低应力条件下泥炭质土σ′和e间关系式可表示为：

(2)

式中，e₀为初始孔隙；C_c为土的压缩指数。

3.2 渗透系数变化分析

泥炭质土渗透系数k和有效应力σ′关系曲线如图 3所示。从图 3中可看出，泥炭质土的k随σ′增大而逐渐减小，泥炭质土的初始渗透系数受初始孔隙比e₀的影响，即e₀越大渗透系数初始值也越大。k的变化幅度受压缩系数的影响，3^#1组土样的压缩系数最大，3^#2组土样的压缩系数最小，从图中可以明显看出，在压缩过程中3^#1组土样渗透系数的变化幅度比3^#2组土样大得多。泥炭质土lg k-lg σ′曲线上切线斜率定义为渗压指数C_k，C_k随有效应力的增大逐渐增大。且C_k和土压缩系数相关，泥炭质土压缩系数越大，在固压作用下越容易被压缩，渗透系数的变化幅度也越大，从而导致高压缩性土的C_k值比低压缩性的大。同样可从图 3中看出，4组土样的lg k-lg σ′曲线随有效压应力增大呈集中趋势，即随着有效压应力的增大，各组泥炭质土的lg k值越来越接近。

对比k-σ′和lg k-lg σ′曲线可明显看出，lg k-lg σ′曲线线性相关度更高，对泥炭质土lg k-lg σ′间的关系进行线性拟合，各组土样线性R²均大于0.85，说明可用线性关系反应k和σ′间关系。从数学的角度可知对于函数y=log_ax，自变量x不可取0，结合实际试验情况，分析lg k-lg σ′之间关系，可将p=0时的情况用p=1 kPa近似代替，故lg k-lg σ′间关系式为：

(3)

式中，k₀为初始渗透系数; C_k为lg k-lg σ′曲线图拟合的斜率，本研究定义为土的压缩渗透指数。

有关学者对黏性土渗透特性的研究是基于k-e间规律，建立了适用于不同地区不同性质软黏土k-e间的模型^[25]。本研究则基于k-σ′间规律建立适用于滇池泥炭质土k-σ′间模型。根据图 3可得到泥炭质土C_k的取值范围为1.0~2.0。

3.3 非线性固结过程分析

图 4为均质土层某时刻饱和土层中孔隙水压力(或有效应力)的分布图，从图可看出，当其顶部作用有连续均布荷载p₀，在Δt时间段内单元体渗出水的体积量等于孔隙体积的减少量即：

(4)

对式(4)化简得：

(5)

根据式e=e₀-C_clg σ′，可得

(6)

根据渗透系数和有效固结应力的关系式lg k=k₀-C_klg σ′可得：

(7)

以式 ，式(7)和式(6)代入式(5)得：

(8)

化简得：

(9)

令 得

(10)

总应力p等于有效应力σ′和超静孔隙水压力u的和，故式(10)可改写成：

(11)

式(11)即为修正饱和土一维非线性固结微分方程，其中c′ _v，p，C_k均为常数。

4 一维非线性固结微分方程求解 4.1 单面排水情况 4.1.1 边界条件

假设土层为单面排水，起始超静孔隙水压力为线性分布，如图 5所示。

设土层排水面的起始总应力为p₁，不透水面的起始总应力为p₂，两者的比值为：

(12)

深度z处的初始超孔隙水压力和总应力相等即：

(13)

当0≤t≤∞，z=H时，=0，对式u+σ=p_z对z求偏导为：，此时=0，故 。

使用式(10)为基础，初始及边界条件为：

边界条件1：< t < ∞，z=0时，σ′=p₁。

0 < t < ∞，z=H时，,

t=∞，0≤z≤H时，。

初始条件1：t=0, 0≤z≤H时，σ′=0。

如若使用式(11)为基础，初始及边界条件为：

边界条件2：0 < t < ∞，z=0时，u=0。

0 < t < ∞，z=H时，=0，

t=∞，0≤z≤H时，u=0。

初始条件2：t=0，0≤z≤H时，u= 。

4.1.2 方程求解

根据试验结果可知滇池泥炭质土的C_k取值范围为1.0~2.0。

(1) C_k取值为1时：

将C_k=1代入式(11)修正饱和土一维非线性固结微分方程，可化简为：

(14)

式(14)与经典饱和土的一维固结微分方程结构相同，可采用相同的方法求解，在该边界条件下的求解结果为：

(15)

式中，m为正奇数，m=1, 3, 5, 7, …。

参考经典饱和土的一维固结方程，令，可把方程改写为：

(16)

在实际中常取第一项的值，即取m=1得：

(17)

(2) 当C_k取2时：

将C_k=2代入式(11)修正饱和土一维非线性固结微分方程，可化简为：

(18)

引入激活函数:

(19)

根据激活函数tanh x的性质，

(20)

(21)

可得到式(18)的一个特解为：

(22)

4.2 双面排水情况

土层为双面排水时，令α=p₁/p₂，土层厚度为2H，初始和边界条件：

边界条件3：0 < t < ∞，z=0时，u=0;0 < t < ∞，z=H时，u=0;t=∞，0≤z≤H时，u=0。

初始条件3：t=0，0≤z≤H时，u= 。

C_k取值为1时：将C_k=1代入式(11)修正饱和土一维非线性固结微分方程，可求得特解：

(23)

同理在实例中常取第1项值，即m=1可得:

(24)

式中，m为正奇数(1, 3, 5, 7, …)；H为压缩土层最远的排水距离，当土层为单项排水时，H取土层的厚度。

竖向固结时间因数T′ _v按式(25)计算：

(25)

5 试验验证和算例分析

根据经典的太沙基固结理论可知，单面排水情况下任意时刻土层固结度可表示为：

(26)

本研究对土层的平均固结度采用同样的概念，即某一时刻的有效应力图面积和最终有效应力图面积之比值，称为竖向平均固结度 ，用公式表示为：

(27)

将式(13)、式(17)代入式(27)可得到相应的任一时刻的土层平均固结度：

(28)

任一时刻土层的沉降量同Terzaghi固结理论为：

(29)

土层为双面排水时，将式(13), (24)代入式(28)得任意时刻土层平均固结度 ：

(30)

5.1 试验验证

为验证本研究推导公式的准确性，在渔村社区再次取土样做固结试验以验证本研究公式的可靠性(所取土样C_k=1)，取样深度―7 m，所取土样参数如表 2所示。试样所用环刀尺寸同上，试样做抽气饱和处理，饱和后快速装入缸体，分别施加50，200，500 kPa荷载，并记录试验数据。

图 6为50，200，500 kPa荷载下Terzaghi一维固结理论，本研究非线性固结理论预测的和试验得出的S_t-t关系曲线。

试样固结初期，初始沉降S_t=0=0，随着时间发展，沉降S_t逐渐增大，并最终趋于某个值，Terzaghi固结理论和本研究非线性固结理论都随着固结过程的发展最终稳定值S_∞=aσH/ (1+e₀)。Terzaghi固结理论与实际误差很大，2 cm高试样用Terzaghi固结理论预测发现固结完成所用的时间很短，而本研究非线性固结理论更接近实际。特别是当施加的固结荷载越大，预测的S_t-t关系曲线越接近实际。但泥炭质土结构复杂，e-lg σ′和lg k-lg σ′两个模型和实际存在一定误差，且非线性固结理论推导过程中未考虑到土瞬时沉降等影响因素，故理论值与试验值还存在一定的误差。

5.2 算例分析

假设存在5 m厚的饱和泥炭质土层，土的物理力学指标同3^#3组，余下物理指标如表 3所示。采用堆载压实固结法，上覆均布荷载，在单项排水情况下进行一维非线性应变固结，如图 5所示。

图 7为本研究公式预测泥炭质土-T′ _v关系曲线。

从图 7可看出，预测的初始值为16.8%，存在瞬时固结。随着T′ _v的增大，逐渐增大，并最终趋于稳定=1。在C_k=1时，本研究非线性固结推导式预测的泥炭质土一维固结 -T′ _v关系曲线与Terzaghi固结理论一致。经典Terzaghi固结理论是本研究非线性固结式在C_k=1时的一种特殊情况，也从另一个角度说明本研究固结式的实用性和概括性。

图 8为泥炭质土地基平均固结度与时间关系曲线。

由图 8可见，随着固结时间的推移，固结度逐渐增大，并最终无限趋近于1，但本研究非线性固结理论预测稳定所需时间远远大于运用经典Terzaghi固结理论预测固结稳定时间，其本质原因是经典Terzaghi一维固结方程忽略了泥炭质土施加竖向应力后压缩系数a和渗透系数k随着σ′的增大逐渐减小的影响，同时泥炭质土e-lg σ′模型中的C_c值要远大于其压缩系数a，从而导致本研究非线性固结理论中的c′ _v要远小于Terzaghi固结理论中的c_v，故本研究一维非线性固结理论预测固结稳定所需时间要远大于经典Terzaghi固结理论预测固结稳定时间。

综上所述，算例分析中得到本研究一维非线性固结理论和Terzaghi一维固结理论的-T′ _v关系曲线一致，充分说明本研究提出的一维非线性固结理论的合理性、实用性和概括性。在-t曲线上存在很大差别，随固结时间推移，固结度均逐渐增大，本研究一维非线性固结理论预测固结稳定所需的时间要远大于Terzaghi饱和一维固结理论预测值，更符合工程实际，具有良好的工程实际意义及应用前景。

6 结论

(1) 以滇池泥炭质土为研究对象，进行常规固结试验和固结-渗透联合试验，得出适用于泥炭质土e-σ′，k-σ′之间关系模型，其中利用k-σ′模型取代k-e模型使得非线性固结微分方程得到了极大的简化。

(2) 考虑文中泥炭质土e-σ′，k-σ′间模型，从而推导出了适用泥炭质土的一维非线性固结微分方程，为解决含泥炭质土的工程问题提供理论支撑。

(3) 通过数学计算，解出了在C_k=1条件下非线性固结微分方程的解析解，并且对C_k=2条件下非线性固结微分方程进行分析，得出了一个特殊解。

(4) 在C_k=1的情况下，得到的饱和一维非线性固结微分方程与Terzaghi一维固结方程结构相同，仅竖向固结系数c_v存在差异；在C_k=2情况下，通过引入激活函数得到非线性固结微分方程的一个特解。

(5) 对C_k=1类型泥炭质土进行固结试验研究发现，在固结过程中Terzaghi一维固结理论预测沉降速率比实测快，到达沉降稳定状态所需时间短；推导的一维非线性固结理论预测的各时间段沉降值和实测值接近，预测到达沉降稳定所需时间的误差小，且上覆荷载越大S_t-t关系曲线与实测值越接近，推导的一维非线性固结理论预测结果越精确。

(6) 通过算例分析得出，使用本研究推导的一维非线性固结理论预测的-t关系曲线与Terzaghi一维固结理论存在很大差异，使用推导的一维非线性固结理论预测固结稳定所需的时间远远要大于Terzaghi一维固结理论，而两者-T′ _v关系曲线为一致，充分说明了本研究提出的一维非线性固结理论的合理性、实用性和概括性。