一种水库运行期滑坡浸润线和推力优化计算方法

技术领域

本发明涉及一种水库运行期滑坡浸润线和推力优化计算方法，属于地质灾害防治领域。

背景技术

三峡工程于2003年6月正式蓄水发电，2008年以来，库水位每年在145～175m之间周期性波动。库水位周期性涨落是诱发水库库岸滑坡的主要因素之一。目前，库水位波动对滑坡稳定性的影响受到了广泛关注。地下水对滑坡稳定性具有决定性的影响，因此库区涉水滑坡稳定性分析需要首先确定坡体内浸润线的位置。

目前水库滑坡防治工程设计中，通常采用库水位单次下降工况，如在三峡库区，使用175米单次下降至145米水位降落的极端工况，进行涉水滑坡内部的浸润线计算，进而基于剩余推力法计算获取滑坡的剩余推力和稳定性系数。这显然已经不符合近十年来三峡库区实际的运行工况。特别是对于库岸特大水力型滑坡，采用现有方法进行计算，往往使得获取的滑坡剩余推力和稳定性系数无法反映滑坡的真实状态，导致工程的安全储备及造价过高。

基于上述的缺陷，有必要提供一种适用于水库运行期的滑坡浸润线和推力优化计算方法，为库区滑坡防治工程的设计提供技术支撑。

发明内容

本发明的目的在于，提供一种适用于水库运行期的滑坡浸润线和推力优化计算方法。

本发明的技术方案如下：

一种水库运行期的滑坡浸润线和推力优化计算方法，引入了一维非稳定性渗流方程，在单次水位下降时浸润线公式解析解基础上，扩展得到了库水位连续波动时库岸滑坡浸润线的计算公式。将库水位连续波动时库岸浸润线计算公式，引入滑坡稳定性计算的传递系数法中，计算获取滑坡在库水位连续波动时的稳定性和下滑力。

本发明中库水位连续波动时库岸滑坡浸润线的推导过程及计算原理如下：

潜水面的渗流速度由达西定律计算，见公式(1)。

对于一维渗流问题，其中潜水面高度仅随着x改变，与z无关，即H＝H(x，t)。因此，对于单宽流量可由公式(2)计算。

其中ZL是隔水底板的高程，如果将基准面取在隔水底板处，则单宽流量可由公式(3)计算。

单位流入量和流出量之差为该段地下水体积的变化量，有：

△V＝q(x,t)△t-q(x+△x,t)△t (4)

且△V＝[h(x,t+△t)-h(x,t)]△x.1.u

则有：

q(x,t)△t-q(x+△x,t)△t＝[h(x,t+△t)-h(x,t)]△x.1.u

式中u

方程(6)两边都除以△x△t，可得:

令Δt→0，Δx→0，则得：

将裘布衣方程(3)带入公式(8)，则有：

公式(9)即为潜水的一维非稳定性运动微分方程Boussinesq布辛尼斯克方程，整理后写为公式(10)。

这是一个二阶非线性偏微分方程，目前还没有求解析解的方法，通常采用简化方法。简化方法是将括号中的一个含水层厚度h近似地看作常量，用时段始、末潜水流厚度的平均值h

在一维非稳定性渗流方程的基础上，推到获取库水位连续波动作用下岸坡浸润线计算公式。首先做以下几点基本假定：(1)含水层均质、各向同性，侧向无限延伸，具有水平不透水层；(2)初始潜水面水平；(3)潜水流为一维流；(4)库水位以V的速度等速变化；(5)库岸按垂直考虑、库水降幅内的库岸与大地相比小得多，为了简化将其视为垂直库岸。

初始时刻，即t＝0时，各点水位为h(0，0)，设距库岸x处在t时刻的地下水位变幅为：

u(x,t)＝h

断面t＝0时的水位变幅为：

u(x,0)＝h

库水位以v的速度变化时，则在x＝0断面处，有：u(0,t)＝vt；

在x＝∞的断面处，有：u(∞,t)＝0；

由公式(12)可以把上述水位下降的半无限含水层中地下水非稳定渗流归结为下列数学模型：

u(x,0)＝0，0

u(0,t)＝vt，t>0 (15)

u(∞,t)＝0，t>0 (16)

将上述(13)～(16)表述的数学模型利用拉普拉斯积分变换和逆变换，可得到微分方程的解析解，如下：

其中，

令

h

直接使用f(λ)的表达式太复杂，因此前人对其进行了拟合，得到了公式(19)。

则库水位单次等速变化时的浸润线方程如公式(20)所示。

在库水位单次等速变化时的浸润线公式的基础上，可将t时刻库水位的高程写为ht，t-1时刻库水位的高程定义为ht-1，t-1时刻x位置的地下水位定义为h

h

S

上式中，k为渗透系数，单位为m/d；h

由库水位的波动曲线可计算获得不同时步下库水位的波动幅度S

在库水位连续波动的不同时步下的浸润线计算公式的基础上，库水位连续波动的不同时步下滑坡安全系数和剩余下滑力，计算公式如下：

R

T

N

W

上式中，F

根据库水位波动曲线，可由公式(22)确定t时刻的库水位波动值S

本发明的有益效果是：通过解析解的方法计算获取库水位连续波动时的滑坡浸润线曲线，引入库水位连续波动曲线，即可进一步获取库水位连续波动时不同时步下的滑坡剩余推力和稳定性系数，从而使得库区滑坡防治工程设计时，计算获取的滑坡剩余推力和稳定性系数计算更加符合滑坡的真实状态，在保证安全的前提下，能有效节省工程造价，且这种解析解和传递系数法耦合计算方法在工程应用中非常实用。

附图说明

图1为本发明水库运行期滑坡浸润线和推力优化计算流程图。

图2为本发明库岸滑坡浸润线计算模型示意图。

图3为本发明实施例中三峡库区多年库水位波动曲线。

图4为应用实例-塔坪滑坡在三峡库水位连续波动作用下的浸润线位置变化图。

图5为本发明库水位连续波动时滑坡传递系数计算法的计算示意图。

图6为应用实例-塔坪滑坡在三峡库水位连续波动作用下的稳定性变化图。

图7为应用实例-塔坪滑坡使用该方法计算获取的剩余推力和传统方法的对比图。

具体实施方式

本下面将结合附图1-7，对本发明实施例中的技术方案进行清楚、详细的描述，显然，所描述的实施例仅是本发明一部分实施例，而不是全部的实施例，基于本发明中的实施例，本领域普通技术人员在没有做出创造性劳动前提下所获得的所有其他实施例，都属于本发明的保护范围。

实施例：如图1所示，一种水库运行期滑坡浸润线和推力优化计算方法，本发明的实施步骤是：

步骤一、获取库水位波动曲线，如对于三峡库区，可通过多年监测水位线，获取175米蓄水以来，近十年的库水位波动曲线，并进一步根据库水位波动曲线的关键特征点，概化出多年库水位波动概化曲线，如图2所示。

步骤二、通过收集勘察报告，并通过现场调查和试验，可获取滑坡的计算参数。如本实施例中三峡库岸特大水力型滑坡—塔坪滑坡，可通过收集勘察资料，获取滑坡的工程地质参数，如滑坡材料力学参数、结构、地形、水文参数、滑带分布等特征。

步骤三、引入了一维非稳定性渗流方程，推导得到库水位连续波动时库岸滑坡浸润线的计算公式，获取滑坡在库水位连续波动下的浸润线变化曲线。

库水位连续波动时库岸滑坡浸润线的计算公式的推导原理如下：

潜水面的渗流速度由达西定律计算，见公式(1)。

对于一维渗流问题，其中潜水面高度仅随着x改变，与z无关，即H＝H(x，t)。因此，对于单宽流量可由公式(2)计算。

其中ZL是隔水底板的高程，如果将基准面取在隔水底板处，则单宽流量可由公式(3)计算。

单位流入量和流出量之差为该段地下水体积的变化量，有：

△V＝q(x,t)△t-q(x+△x,t)△t (4)

且△V＝[h(x,t+△t)-h(x,t)]△x.1.u

则有：

q(x,t)△t-q(x+△x,t)△t＝[h(x,t+△t)-h(x,t)]△x.1.u

式中u

方程(6)两边都除以△x△t，可得:

令Δt→0，Δx→0，则得：

将裘布衣方程(3)带入公式(8)，则有：

公式(9)即为潜水的一维非稳定性运动微分方程Boussinesq布辛尼斯克方程，整理后写为公式(10)。

这是一个二阶非线性偏微分方程，目前还没有求解析解的方法，通常采用简化方法。简化方法是将括号中的一个含水层厚度h近似地看作常量，用时段始、末潜水流厚度的平均值h

在一维非稳定性渗流方程的基础上，推到获取库水位连续波动作用下岸坡浸润线计算公式。如图3所示，首先做以下几点基本假定：(1)含水层均质、各向同性，侧向无限延伸，具有水平不透水层；(2)初始潜水面水平；(3)潜水流为一维流；(4)库水位以V的速度等速变化；(5)库岸按垂直考虑、库水降幅内的库岸与大地相比小得多，为了简化将其视为垂直库岸。

初始时刻，即t＝0时，各点水位为h(0，0)，设距库岸x处在t时刻的地下水位变幅为：

u(x,t)＝h

断面t＝0时的水位变幅为：

u(x,0)＝h

库水位以v的速度变化时，则在x＝0断面处，有：u(0,t)＝vt；

在x＝∞的断面处，有：u(∞,t)＝0；

由公式(12)可以把上述水位下降的半无限含水层中地下水非稳定渗流归结为下列数学模型：

u(x,0)＝0，0

u(0,t)＝vt，t>0 (15)

u(∞,t)＝0，t>0 (16)

将上述(13)～(16)表述的数学模型利用拉普拉斯积分变换和逆变换，可得到微分方程的解析解，如下：

其中，

令

h

直接使用f(λ)的表达式太复杂，因此前人对其进行了拟合，得到了公式(19)。

则库水位单次等速变化时的浸润线方程如公式(20)所示。

在库水位单次等速变化时的浸润线公式的基础上，可将t时刻库水位的高程写为ht，t-1时刻库水位的高程定义为ht-1，t-1时刻x位置的地下水位定义为h

h

S

上式中，k为渗透系数，单位为m/d；h

由库水位的波动曲线可计算获得不同时步下库水位的波动幅度S

本实施例中，三峡库区特大水力型滑坡—塔坪滑坡在三峡库区多年平均库水位波动下，坡体浸润线曲线计算结果如图4所示。

步骤四、在获取浸润线计算公式的基础上，改进获取库水位连续波动下的传递系数法的计算公式，计算获取不同库水位波动作用下滑坡的稳定性系数和剩余推力，图5为本发明库水位连续波动时滑坡传递系数计算法的计算示意图。

库水位波动作用下滑坡传递系数法的计算公式如下：

R

T

N

W

上式中，F

根据库水位波动曲线，可由公式(22)确定t时刻的库水位波动值S

本实施例中，三峡库区特大水力型滑坡—塔坪滑坡在三峡库区多年平均库水位波动下，使用本发明方法计算获取滑坡的稳定性系数变化曲线如图6所示。

本实施例中，三峡库区特大水力型滑坡—塔坪滑坡在三峡库区多年平均库水位波动下，使用本发明方法计算获取滑坡不同条块剩余推力于传统骤降工况计算获取的剩余推力的对比如图7所示。

塔坪滑坡在多年平均水位运行工况下的最危险安全系数为1.092。计算结果表明，相较现行规范中的175m至145m水位降落工况，最危险安全系数提高了约0.05。

当设计安全系数取值为1.10时，使用本发明方法计算库水位波动工况下的滑坡剩余推力显著小于传统单次下降工况。设桩位置为第8条块时，单次水位下降工况计算获取的桩后设计剩余推力为17081kN/m；本发明方法计算获取的桩后设计剩余推力为15131kN/m。考虑桩前抗力的影响，传统的175m降至145m水位下降工况计算获取的设计剩余推力为5467kN/m；本发明方法计算库水位波动工况下设桩处的剩余推力为3754kN/m。后者相较前者降低了31.3％。

采用直径1m的中口径抗滑桩群进行塔坪滑坡的治理设计。依据现行规范中水位175m降至145m下降工况进行抗滑桩群设计计算时，共需布置1236根，总造价约1.2亿元。而根据本发明提出的方法进行抗滑桩群设计计算时，仅需布置781根，总造价约7900万元，降低约34％。