在岩土工程建设中必不可少会产生大量建筑废弃泥浆，为响应国家环境污染治理要求，经济环保的废弃泥浆处理技术成为社会发展的迫切需要。目前，中国废弃泥浆脱水技术主要有泥浆固化处理和泥浆固液分离^[1]。其中，泥浆固液分离相较于泥浆固化处理不受场地和天气等的影响，因而更具有优势。泥浆固液分离技术分为生物方式、化学方式和物理方式脱水。其中，生物脱水方式采用微生物发酵法^[2]，但仅对有机物含量较高的泥浆有明显的脱水效果，而对于建筑废弃泥浆的脱水效果并不明显。化学脱水方式采用分离设备和化学药剂，工程实践^[3-5]表明，化学方式脱水确实有一定的效果，但该方法操作流程复杂，仪器设备价格昂贵，此外由于添加了化学药剂，导致泥浆清液被分离出来之后还需进一步处理和检测达标后才能排放。物理脱水方式包括滤袋脱水、压滤脱水、真空负压脱水、离心脱水等方法。其中，文献[3]和文献[6]的试验表明，离心脱水方式对泥浆处理有较好的脱水效果，且具有经济、快速等优势。因此，离心脱水技术在处理主要由黏土组成的建筑废弃泥浆过程中更具有优势和前景，先后有不少学者在从事其实践研究。文献[7]采用卧螺离心机设计了污泥脱水工艺流程，并采用正交试验方法得到了该工艺流程的最佳处理参数，从而有效实现了污泥减量化处理。文献[8]提出了一个可以预测卧螺离心机动态行为的物理模型，该模型能够考虑沉降过程、滤饼固结和泥沙输送，进而分析了转速、差速等相关参数对离心脱水过程的影响。文献[9]通过调试离心设备参数对底泥进行脱水处理，从而实现了底泥的安全高效、无二次污染处理。这些文献表明，离心脱水技术在工程实践中应用已较为广泛，但多针对离心设备的相关参数进行调试，从而在特定工况下实现对泥浆的有效脱水，而对于泥浆在离心机内部脱水过程的内在机理研究相对薄弱。

泥浆离心脱水本质上是将大孔隙比、高含水率泥浆中的水从泥浆土骨架中排出的过程，因此，泥浆离心脱水也可视为软土大变形固结过程。文献[11]在其分析软土大变形固结的CS2模型^[10]基础上，针对离心加载情况进行了试验分析，并考虑离心加速度随着试样深度的线性变化。实际上，固结理论分析多假定孔隙水的渗流符合达西定律。但是，越来越多的渗透实验研究^[12-15]表明，黏性土中的渗流偏离达西定律，并分别提出了多种不同形式的非达西渗流拟合式，其中最为常用的当属Hansbo^[15]渗流公式。已有文献[16-17]把该渗流模型引入了固结分析中，认为该渗流延缓了饱和黏性土的渗透固结。因此，本文拟引入Hansbo渗流方程描述离心脱水过程中水的排出行为，建立考虑泥浆非达西渗流特性的离心脱水(centrifuge dewatering of non-Darcian flow，CDN)模型, 该分析模型具有模块化的优点，可同时考虑自重、加卸载等荷载形式。然后，通过设计室内离心模型试验对其进行验证。最后分析了相关参数对泥浆离心脱水过程的影响。

1 泥浆离心脱水的数值分析模型 1.1 CDN模型

在离心状态下，如图 1(a)所示，以远离轴心一端为底部，泥浆试样在底部以角速度ω开始绕转轴旋转，所产生的离心加速度为ω²r，离心加速度沿径向指向转轴方向，所产生的离心加速度不断减小，直至转轴处其大小为0。而在只有重力的情况下，试样内部产生的重力加速度是恒定的, 如图 1(b)所示。

泥浆试样的离心脱水过程如图 2所示，泥浆试样以角速度ω绕转轴旋转，在离心作用下出现不断位移的土水界面，朝向转轴方向的为泥浆上清液，背向转轴方向的为脱水后的泥浆试样。

对于初始时刻，如图 3(a)所示，假定泥浆初始孔隙比e₀处处相同，初始厚度为H₀，静水面高度为H_w且在离心脱水过程中保持不变。以其底面为基准面，沿径向朝向轴心方向为正方向，建立z坐标系。泥浆体自下而上被均匀划分为n个单元，且每个单元都具有相同的横截面面积和初始厚度L₀，每个单元的位置以单元节点至基准面的距离表示，且初始位置为z_0,j。经过时间t后，见图 3(b)，泥浆体厚度变为H^t，泥浆体单元厚度和位置分别变为L_j^t和z_j^t。

1.2 本构关系

泥浆体的脱水过程类似于正常固结土的固结过程，这里引入经典的非线性经验公式(1)来描述其压缩过程。

$ e-e_{\text {ref }}=C_{\mathrm{c}}\left(\lg \sigma_{\text {ref }}^{\prime}-\lg \sigma^{\prime}\right) $

(1)

式中：C_c为泥浆的压缩系数，e、e_ref分别为对应于有效应力σ′和某参考应力σ′_ref的孔隙比。

考虑黏性土中的渗流偏离达西定律，因此，这里引入Hansbo渗流方程^[15]描述泥浆中水的排出。为便于分析，类似文献[17]，将该方程写为

$ v=-k i $

(2)

$ k=\left\{\begin{array}{l} \frac{K|i|^{m-1}}{m i_1^{m-1}}, |i| \leqslant i_1 \\ \frac{K\left[1-i_1(m-1)\right]}{m|i|}, |i|>i_1 \end{array}\right. $

(3)

式中：v为渗流速度，i为水力梯度，i₁和m为非达西渗流参数，k为等效渗透系数，K为渗透系数，随孔隙比e的变化而变化，采用文献[18]中提出的关系式(4)表达。

$ \lg K=\frac{e-e_{\mathrm{ref}}}{C_{\mathrm{K}}}+\lg K_{\mathrm{ref}} $

(4)

式中: C_k为渗透指数，K_ref为对应e_ref的渗透系数。

1.3 单元应力及孔隙压力

CDN模型中，单元j节点处的应力主要由在离心作用下上部单元产生的离心力所造成的，采用文献[11]表征离心机特性参数的离心因子N来描述离心作用的大小。在泥浆脱水的实际工程中，其所受重力可以忽略，底部离心因子N₀可表示为

$ N_0=\frac{\omega^2 r}{g} $

(5)

离心因子大小在z方向上线性变化，则t时刻单元j节点处产生的离心因子N可表示为

$ N_j^t=N_0-\alpha z_j^t $

(6)

式中$ \alpha=\frac{\omega^2}{g}$。

时刻t单元j的饱和重度γ_j^t可表示为

$ \gamma_j^t=N_j^t \gamma_{\mathrm{w}} \frac{G_{\mathrm{s}}+e_j^t}{1+e_j^t} $

(7)

式中：γ_w为水的重度，G_s为泥浆中的土粒比重，e_j^t为在时刻t单元j的孔隙比。

其中泥浆顶面处孔隙静水压为

$ u_1^t=\int_{H^t}^{H_{\mathrm{w}}} \gamma_{\mathrm{w}} N \mathrm{~d} z=\gamma_{\mathrm{w}}\left[N_0\left(H_{\mathrm{w}}-H^t\right)-\frac{\alpha}{2}\left(H_{\mathrm{w}}^2-H^{t2}\right)\right] $

(8)

综合式(6)、(7)、(8)可计算时刻t单元j中点处的应力, 即

$ \begin{aligned} & \sigma_j^t= u_1^t+\gamma_{\mathrm{w}}\left(N_j^t-\alpha \frac{L_j^t}{4}\right) \frac{L_j^t}{2} \frac{G_{\mathrm{s}}+e_j^t}{1+e_j^t}+ \\ & \gamma_{\mathrm{w}} \sum\limits_{i=j+1}^{R_j} L_i^t N_j^t \frac{G_{\mathrm{s}}+e_j^t}{1+e_j^t} \end{aligned} $

(9)

式中$ j=1, 2, \cdots, R_j, L_j^t$为时刻t泥浆体单元j的厚度。

在时刻t，由节点j的孔隙比e_j^t根据式(1)计算相应的有效应力σ′_j^t，相应的孔隙水压力可表示为

$ u_j^t=\sigma_j^t-\sigma_j^{\prime t} $

(10)

1.4 渗流速度及土水界面位移

在离心力场中，由于z方向上离心加速度竖直向上不断减小，因此相应水力梯度有所不同，其中任意边界下水的总势能为

$ \varphi_{\mathrm{b}}=g N_0 H_{\mathrm{w}}-\frac{g \alpha}{2} H_{\mathrm{w}}^2 $

(11)

在节点j处时刻t的水的势能表达式为

$ \varphi_j^t=g \int_0^{z_j^t} N \mathrm{~d} z+\frac{u_j^t}{\rho_{\mathrm{w}}}=g\left(N_0 z_j^t-\alpha \frac{z_j^{t 2}}{2}\right)+\frac{u_j^t}{\rho_{\mathrm{w}}} $

(12)

这样，节点j处时刻t的等效水力梯度i_j^t为

$ i_j^t=\frac{\varphi_{j+1}^t-\varphi_j^t}{g\left(z_{j+1}^t-z_j^t\right)} $

(13)

相邻单元的渗流示于图 4，规定沿半径指向转轴的渗流方向为正，以单元j和j-1为例，相邻单元间的相对渗流速度根据相应的等效水力梯度i_j^t计算。

$ v_{\mathrm{rf}, j}^t=-k_{s, j}^t i_j^t $

(14)

式中k^t_{s, j}为时刻t单元j中点处的等效渗透系数，当半径R比土体厚度大得多时，k^t_{s, j}可表示为

$ k_{s, j}^t=\frac{k_{j+1}^t k_j^t\left(L_{j+1}^t+L_j^t\right)}{L_{j+1}^t k_j^t+L_j^t k_{j+1}^t} $

(15)

其中k_j^t为时刻t单元j中点处的渗透系数，按式(3)求得。

在泥浆脱水过程中，假定泥浆体靠近转轴端的边界为完全透水边界，而远离转轴端的边界为完全不透水边界，根据式(11)和式(13)可得这两个边界处的等效水力梯度分别为

$ \begin{gathered} i_1^t=\frac{g N_0 H_{\mathrm{w}}-\frac{g \alpha}{2} H_{\mathrm{w}}^2-\varphi_1^t}{g\left(H^t-z_1^t\right)} \end{gathered} $

(16)

$ i_{\mathrm{R}}^t=0 $

(17)

这样，两个边界处的渗流速度可分别表示为

$ v_{\mathrm{r}, 1}^t=-k_1^t \frac{g N_0 H_{\mathrm{w}}-\frac{g \alpha}{2} H_{\mathrm{w}}^2-\varphi_1^t}{g\left(H^t-z_1^t\right)} $

(18)

$ v_{\mathrm{rf}, \mathrm{R}}^t=0 $

(19)

当时刻t相应位置处的相对速度求出后，则相继可获得t+Δt时刻不同单元厚度L_j^t+Δt和孔隙比e_j^t+Δt，进而获得泥浆体厚度H^t+Δt、土水界面位移S^t+Δt和平均固结度U_t。相应的计算式分别为

$ L_j^{t+\Delta t}=L_j^t-\left(v_{\mathrm{r} f, j}^t-v_{\mathrm{r} f, j-1}^t\right) \Delta t $

(20)

$ e_j^{t+\Delta t}=\frac{L_j^{t+\Delta t}\left(1+e_{0, j}\right)}{L_0-1} $

(21)

$ H^{t+\Delta t}=\sum\limits_{j=1}^R L_j^{t+\Delta t} $

(22)

$ S^{t+\Delta t}=H_0-H^{t+\Delta t} $

(23)

$ U_t=\frac{S^{t+\Delta t}}{S_{\mathrm{f}}} $

(24)

式中：e_{0, j}为泥浆层中j单元的初始孔隙比，S_f为土水界面最终位移量。

1.5 时间增量

前述迭代计算过程类似于显式有限差分法，是条件收敛的。为保证迭代计算的收敛以及结果的精度，时间增量Δt采用文献[11]提出的方法确定。具体计算公式为

$ \Delta t=\min \left[\frac{0.4 \gamma_{\mathrm{w}} \alpha_{v, j}^t\left(L_j^t\right)^2}{k_j^t\left(1+e_j^t\right)}, \left|\frac{0.01 L_0\left(e_{0, j}-e_{\mathrm{f}, j}\right)}{\left(1+e_{0, j}\right)\left(v_{\mathrm{rf}, j}^t-v_{\mathrm{rf}, j+1}^t\right)}\right|\right] $

(25)

式中α^t_{v, j}为单元j在t时刻的压缩系数，可根据压缩指数C_c与有效应力的关系获得，e_{f, j}为单元j的最终孔隙比，可由最终应力按式(1)计算而得。

2 模型试验 2.1 离心作用下泥浆脱水过程对比

文献[11]曾对文献[20]中给出的一个新加坡软土离心试验的案例进行数值模拟。相关参数为：转速ω=23.87 rad/s，土体初始厚度H₀=260 mm，水位H_w=260 mm，初始孔隙比e₀=3.26，初始有效应力σ′=0.189 kPa，初始渗透系数k=7.54×10^-9 cm/s，土粒比重G_s=2.65，且孔隙比与有效应力和渗透系数的关系分别为

$ e=2.76-0.69 \lg \sigma^{\prime} $

(26)

$ e=12.52+1.14 \lg k $

(27)

这里模拟时分为两种情况：1)将Hansbo渗流方程中的参数取值为m=1，i₁=0，退化到Darcy渗流情况；2)假定新加坡软土有非达西渗流特征，将Hansbo渗流方程中的参数取值为m=2，i₁=10。模拟结果与文献[11]的结果均示于图 5。可以看出，Darcy渗流情况和文献[11]中的数值模拟是一致的，说明本文的计算程序是正确的。但图 5表明按Hansbo渗流分析可以更好地模拟新加坡软土的沉降情况。

2.2 与室内试验结果的对比

本文试验泥浆取自河南周口某工地。经室内试验测得其主要物理参数如下：土粒比重G_s=2.71，液限w_L=39.4%，塑限w_P=22.3%，塑性指数I_P=17.1。通过压缩试验获得其压缩指数C_c=0.39。在研究其渗透性时，为便于试样成形和考虑渗透试验的要求，利用静压法制作了孔隙比为1.05~1.52的一组重塑试样进行渗透试验，测定该组试样的渗透系数为1.93×10^-6~1.19×10^-5 cm/s，渗透指数C_k=0.62，m=1.30~1.86，i₀=20~30，说明试样具有显著的非达西渗流特征。

2.2.1 试验设备

为进一步进行室内试验的模拟，自行设计了一种离心脱水试验装置，如图 6所示。主要性能指标为：离心机转盘最外侧边缘处半径为180 mm，有效离心半径为170 mm；转盘内可放置2~4个泥浆试筒，泥浆试筒高60 mm、直径为40 mm；电机调速器设置有5个档位，本次试验选用1、3、5三个档位，对应的电机转速分别为275、570、825 r/min。

本试验装置通过电机提供设备动力。试验时可将电机调速器调至所需档位，由电机带动转盘旋转，从而对转盘内试筒中的泥浆产生离心力场，试样筒内开始出现土水分界面。随着脱水时间的延长，土水分界面向外位移，可通过设置在转盘内的无线摄像机随时进行监测。

2.2.2 试验过程

将泥浆缓慢倒入高为60 mm的圆柱形试筒(试样有效高度为57 mm)中，倾倒过程中尽量避免气泡产生。随后用刮土刀刮除试筒上部多余泥浆，然后用保鲜膜封住顶部并用胶带缠绕封口。制备完毕的泥浆试样如图 7所示。

本次试验按孔隙比为2.71制备泥浆。考虑到时间问题，试验时的泥浆孔隙比与之前略有差异，见表 1。泥浆试样分为4组，分别命名为A-1、A-2、A-3、A-4，其中A-1为自重对照组，其他3组为离心试验组，且每组分别制作两个泥浆试样，依次放入到转盘中的两个试样槽中。在进行室内离心试验时，将自重试验组试样A-1竖直静置于水平桌面上，采用无线摄像机记录其土水界面位移过程，泥浆离心脱水的试验方案示于表 1。

2.2.3 室内试验结果及模型验证

根据摄像机所记录的画面，对不同时刻土水界面位移量进行读数，并绘制泥浆试样土水界面位移量与时间关系曲线示于图 8，可以看出，土水界面最终位移量随着电机转速的增大而增大，且随时间的延长而逐渐趋于稳定。

用上述CDN模型对模型试验过程分别进行模拟。参数取值如下: 比重G_s=2.71，压缩指数C_c=0.39，渗透指数C_k=0.62，根据C_k按式(4)计算初始渗透系数K₀=1.69×10^-7 cm/s，Hansbo渗流参数由试算取值为m=1.26，i₁=3.58，其他参数见表 2。

图 8也给出了CDN模型的数值解。可以看出，总体上CDN模型可以较好地拟合泥浆离心试验的结果。

3 参数分析

参照文献[21-23]和室内试验，本节参数取值如下：转鼓半径r=300 mm，电机转速n=1 000 r/min，土粒比重G_s=2.7，压缩指数C_c=0.4，渗透指数C_k=0.6，参考应力σ′_ref=28.0 kPa，参考孔隙比e₀=1.26，初始孔隙比e₀=2.7，初始有效应力σ₀=7.03× 10^-3 kPa，初始含水率w₀=100%，初始渗透系数K₀=1.0×10^-7 m/s，非达西渗流参数m=1.5，i₁=5.0，单元数n=100，泥浆初始高度h₀=h_w=150 mm，下面单因素分析中，其他参数不变。

3.1 电机转速的影响

为考察电机转速对脱水过程的影响，这里取电机转速分别为500 r/min、750 r/min、1 000 r/min、1 250 r/min、1 500 r/min，图 9给出了不同电机转速下的土水界面位移与时间关系曲线，可以看出，随着电机转速的增大，泥浆最终土水界面位移量不断增加，同时，泥浆土水界面位移量线变得接近，说明随着电机转速的增加对土水界面位移量的影响逐渐减小，最终土水界面位移量增加的趋势不断减缓，此外，电机转速越大，泥浆开始脱水的阶段就越早，且到达最终土水界面位移量的所需时间也越短。

图 10给出了不同电机转速下的土水界面移动速度与时间的关系曲线。图 10表明，土水界面移动速度及土水界面移动速度峰值与电机转速呈正相关关系，电机转速越大，土水界面移动速度达到峰值的时间越早，且土水界面移动速度在达到峰值之后，其速度快速下降，不同转速下的泥浆土水界面移动速度曲线很快趋于一致。

图 11给出了不同电机转速下的平均固结度与时间的关系曲线。不难看出，电机转速越快，泥浆脱水固结进程开始就越早，且越早完成固结，因此提高电机转速，将会提高离心脱水固结的效率。

3.2 泥浆初始厚度的影响

泥浆初始厚度，即相对于离心机内部液环层的高度，其高度由卧螺离心机的转鼓直径和溢流板直径共同决定。为考察泥浆初始厚度对脱水过程的影响，这里取泥浆初始厚度分别为100 mm、150 mm、200 mm、250 mm。图 12给出了不同初始厚度的土水界面位移与时间关系曲线，可以看出，随着泥浆初始厚度的增加，最终土水界面位移量不断增加，这与文献[11]中所得出的结论一致。同时，初始厚度越大，达到最终土水界面位移量所需的时间越长。同时发现界面位移与时间关系曲线中出现初期与后期不同的性状，其主要原因在于，离心脱水试验中，泥浆初始厚度较大的泥浆试样，其较多部分更接近转轴中心，那么这部分泥浆所受到的离心力较小，因此在初始时刻泥浆厚度较大的试样脱水速度较慢。随着时间的延长，土水界面将远离转轴，即更多的土颗粒将处于更大的半径处，这样其受到的离心力比初始厚度较小的试样更大，脱水速度自然加快。

图 13给出了不同初始厚度的土水界面移动速度与时间的关系曲线，可以看出，泥浆初始厚度越大，达到土水界面移动速度峰值越大，且达到速度峰值的时间随泥浆初始厚度的增加而延长。但泥浆初始厚度越大，土水界面位移初期速度就越小，且位移速度在达到峰值之后，其速度下降就越快。

3.3 泥浆压缩指数C_c的影响

为考察泥浆压缩指数C_c对脱水过程的影响，这里取泥浆压缩指数C_c分别为0.4、0.6、0.8、1.0。图 14给出了不同压缩指数C_c的泥浆土水界面位移与时间的曲线，可以看出，泥浆最终土水界面位移量随压缩指数C_c的增大而增大，同时，在泥浆脱水的初期，不同压缩指数的泥浆其土水界面位移曲线近乎一致，后期达到最终土水界面位移的时间随压缩指数C_c的增大而延长。

图 15给出了不同压缩指数C_c的泥浆土水界面移动速度与时间的关系曲线，可以发现，在土水界面位移的初期，土水界面移动速度差别非常大，随后速度曲线近乎一致，达到峰值及其对应的时间无差异，可见C_c对于泥浆的脱水速度影响是微弱的。

3.4 泥浆非达西渗流参数的影响

渗流过程采用Hansbo渗流方程，来探讨非达西渗流过程对于泥浆离心脱水过程的影响。其中，对于非达西渗流参数m和i₁的分析，当分析参数m对脱水过程的影响时，对i₁取5.0，参数m取1.2、1.5、1.8；当分析参数i₁对脱水过程的影响时，对m取1.5，i₁取2.0、5.0和8.0。

图 16和图 17给出了在不同电机转速下泥浆与非达西渗流参数的关系曲线。可以看出，非达西渗流参数m和i₁会减缓土水界面位移的过程，且随着m和i₁的增大，减缓的趋势会增加。同时，随着电机转速的增加，非达西渗流参数m和i₁对于脱水过程的影响逐渐减弱，即只有在低转速下，非达西渗流参数m和i₁才会对泥浆的脱水过程产生明显影响，在高转速下，几乎可以忽略非达西渗流参数m和i₁的影响。因此，在对黏土泥浆进行离心脱水过程中选取较高的电机转速，可以削弱黏土泥浆的非达西渗流的影响，从而获得更好更快的脱水效果。

4 结论

基于分段线性化方法，建立了在离心状态下能够考虑泥浆非达西渗流特性的CDN模型，通过与已有研究成果及室内试验结果对比，验证了模型的适用性及算法的准确性，最后分析了相关参数对泥浆脱水过程的影响，得出的主要结论如下：

1) 随着电机转速、泥浆初始厚度和压缩指数C_c的增大，泥浆最终土水界面位移量会增大，而非达西渗流参数对其无影响。

2) 提高电机转速，初期土水界面移动速度及土水界面移动速度峰值会随之增大，达到土水界面移动速度峰值所需时间缩短，而且会缩短离心脱水固结所用的时间；增加泥浆初始厚度，初期土水界面移动速度降低，土水界面移动速度峰值增大，达到土水界面移动速度峰值所需时间延长。

3) 在相同离心状态下，初期土水界面移动速度随压缩指数C_c的增大而减小，但C_c的取值对土水界面移动速度峰值与达到土水界面移动速度峰值所需时间无影响。土水界面移动速度随非达西渗流参数m和i₁的增大而减小，但提高电机转速，可以削弱非达西渗流参数对土水界面移动速度的影响。