基于碰撞模型的斜坡滚石颗粒速度预测1)

魏新容 段绍臻,2) 孙金龙 王文达

∗(兰州理工大学土木工程学院，兰州 730050)

† (兰州理工大学甘肃省土木工程防灾减灾重点实验室，兰州 730050)

引言

随着西部大发展的稳步推进，西部山区的工程结构密度与规模逐渐增大.鉴于山区地质条件较差，滑坡落石灾害日益突出[1-2]，滚石灾害精确预测显得尤为重要.

滑坡滚石的运动形式主要分为4 种[3]：(1) 滚石的滑动；(2) 滚石的自由飞落；(3) 滚石的碰撞；(4)滚石的滚动.其中滚石碰撞是滑坡落石撞击过程中最复杂、最不确定的过程，也是最基础的工作和内容[4-5]，滚石碰撞弹跳的最大高度、最远距离、碰后速度及冲量大小等是预估滚石灾害范围、被动防护装置选址和优化设计的重要参数[6].

由于滚石碰撞后对周围物体会造成巨大的冲击从而导致事故的发生，目前对滚石--壁面冲击力、冲击/碰撞后滚石速度以及运动距离受到较多学者的关注[4-5,7-10].作为滚石碰撞过程中最直接的研究方法——试验研究来说，研究者重点分析滚石碰撞过程中影响因素对滚石--壁面碰撞后速度[11-12]或对碰后与速度相关物理量的影响[13-17].然而，由于滚石尺寸较大，开展全模型试验不仅成本高，滚石碰撞过程中不可控因素较多(如：坡角[18]、入射速度、入射角度、被撞颗粒位置[19]以及运动速度等随机因素)等，造成试验过程中捕捉有效信息困难.因此研究者提出利用理论模型来弥补这一不足.

滚石灾害多发于斜坡，坡体既包括整体稳定的斜坡，也包括正在运动中的坡体和崩塌体等.但是已有的碰撞理论模型中大多将滚石碰撞模型简化为质点的碰撞[20-22]或圆盘/球与静止平板的碰撞[7,23-27]；而滚石与斜板碰撞理论模型比较少见[3].同时Pei 等[28]在研究中发现坡体表层松散滑动易诱发小型基地岩块失稳，更容易发生滚石灾害.张正雄等[29]也发现考虑滚石碰撞到运动滚石时可大大提高滚石落下的个数，且相较于单块岩体沿静止坡面滚下具有更大的冲击力.也就是说运动中的坡体和崩塌体更易诱发和扩大滚石灾害，其坡体和崩塌体的运动情况对滚石碰撞过程的影响不容忽视.但已有研究中关于滚石冲击到运动中坡体的碰撞理论模型鲜见报道.因此导致了滚石碰撞理论模型预测结果与实际结果出入较大，滚石碰撞模型精准预测仍面临着巨大挑战.

本文针对崩塌的滚石冲击到运动的坡体上这一自然现象，基于颗粒接触理论[30-33]，建立滚石与斜坡碰撞过程的理论模型，得到碰撞后滚石颗粒速度的解析解，重点分析了斜坡上被碰颗粒速度和坡角对碰撞后滚石颗粒运动特性(颗粒反弹速度和角度以及反弹旋转角速度)的概率分布的影响，给出碰后滚石颗粒反弹速度预测.

1 滚石颗粒随机撞击模型的建立

1.1 理论模型的选取

自然界中，崩塌的滚石颗粒带有一定的速度冲击到复杂地形的斜坡或者斜坡上的其他运动滚石颗粒上，入射速度和角度、碰撞位置等都是随机的，并且斜坡滚石颗粒碰撞过程时间短、受力复杂，滚石颗粒形状各异，使得直接建模困难较大.为描述碰撞的复杂性以及普遍性，将滚石与斜坡碰撞进行如下处理，建立等效模型.虽然滚石运动为三维运动，但滚石运动的主方向仍可近似认为是平行于入射速度平面内的二维运动[32].为简化模型，采用已有文献中常用的处理方法，将非规则形态的滚石颗粒碰撞简化为二维均质圆盘之间的碰撞，且碰撞为非弹性的[32]，碰撞过程如图1 所示，即：等效半径RA的崩塌滚石颗粒A 以某一速度直接与斜坡(RB→∞)或斜坡上某一运动的等效半径为RB的滚石颗粒发生碰撞，碰撞的位置为J1，其中崩塌滚石A 的入射速度大小记为VA，入射角度为θ(将与水平方向x轴所成的角度记为θ)；为了反映运动坡体对滚石碰撞过程的作用，可由坡体上被碰颗粒B 的运动情况来承载，将坡体上被碰颗粒B 在碰撞前的速度作为运动坡体的速度记为VB，坡体上被碰滚石颗粒B 与斜坡的接触位置为J2，且VB与水平方向的夹角记为启动角α；将颗粒A，B 质心连线与水平方向的夹角记为撞击角β，由于斜坡地形复杂，撞击角β 随着地形的变化而变化；斜坡坡角设为γ.

在图1 所示的滚石颗粒撞击模型中，斜坡角γ、滚石A 的入射速度VA、入射角θ 和斜坡上滚石B 的运动速度以及启动角α 均受滚石尺寸以及斜坡等因素决定且是随机变化的，滚石颗粒A 和B 的碰撞点J1，颗粒B 与斜坡的接触点J2也是随机变化的，即碰撞角度β 是随机变化的.

图1 斜坡滚石颗粒碰撞模型Fig.1 The model of rock impacting on the slope

1.2 基本方程的建立

由上述滚石颗粒撞击模型[32-33]以及模型假设，并且考虑滚石颗粒A 和B 撞击的两个阶段(压缩阶段、恢复阶段).

(1)对于滚石A 和B 的碰撞压缩阶段，即：从碰撞开始到碰撞压缩量达到最大时，所遵循的力学方程如下：

颗粒A 切向方向的碰撞冲量表示为

颗粒A 法向方向的碰撞冲量表示为

颗粒A 依据冲量矩定理得到

颗粒B 的x方向的冲量为

颗粒B 的y方向的冲量为

颗粒B 依据冲量矩定理得

(2) 颗粒A 和B 碰撞完成后整个碰撞过程所遵循的力学方程如下：

颗粒A 的碰撞后切向冲量表示为

颗粒A 的碰撞后的法向冲量表示为

颗粒A 依据冲量矩定理得

颗粒B 的x方向的冲量为

颗粒B 的y方向的冲量为

颗粒B 依据冲量矩定理得

1.3 连接条件

上述模型的基本方程共12 个，但是有20 个未知量，显然方程组不封闭，因此需要补充连接条件.

根据两颗粒碰撞时的接触点处的连接条件得到方程，对于两颗粒的接触点J1得：x方向的平衡方程

y方向的平衡方程

对于B 颗粒与地面接触点J2得：

x方向的平衡方程

y方向的平衡方程

各个冲量之间的关系为

由于碰撞为不完全非弹性碰撞，所以本文考虑的法向、切向碰撞恢复系数，即k1,k2采用的是Xie 等[34]恢复系数公式，在此不再赘述.

2 结果分析

将方程(1)∼(20)进行联立求解，得到含有坡角γ，入射速度VA，入射角度θ，碰撞角度β，颗粒B 的速度VB以及启动角度α 等参数的反弹颗粒速度(UAx，UAy，ωA)的解析解.

在斜坡滚石运动过程中，这些变量的取值是有物理意义以及一定范围的，对于岩质落石来说，密度约为2500 kg/m3；根据文献[35]的试验选取滚石颗粒入射角度θ 范围(20◦，90◦)，碰撞角度β 范围与θ 有关[36]且β 取值范围为(θ，θ+90◦).

我们利用等间距直方图理论[37]计算滚石颗粒反弹速度以及角度的频率.假设样本空间包括所有的反弹速度和角度，将其作为随机变量X，其中随机变量个数为n.假设a=min{X}，b=max{X}，将区间[a−ε,b+ε] 等分为m个小区间，其中ε 是小的正数，统计出落在小区间i中(i=1,2,···,m) 样本的频数fi，得到频率P(xi)=fi/n.

计算过程主要分为如下几个部分：(1)划分入射角度θ、碰撞角度β 所在的区间，且假设上述随机变量在区间内为等概率分布.(2)由于碰撞角度β 范围与入射角度θ 有关，一个入射角θ 对应一个β 的区间，选取β 每个子区间中的均值代入反弹颗粒解析解中.(3) 重复过程(2)得到不同入射角θ 下每个子区间下的值作为样本.(4)确定固定坡角下所有的样本空间并输出频率.

2.1 滚石--壁面碰撞模型验证与分析

将本模型中解析解与章广成等[35]滚石坡面的试验结果进行对分析.理论模型得到的反弹颗粒速度与其试验值对比如图2 所示，其中验证模型中所采用的滚石试样的物理力学性质如表1 中所示.图2 中离散的方块代表反弹颗粒速度的均值(实心代表试验值，空心代表理论值)，I 代表颗粒反弹速度的方差值.从上图中可以看出该模型的理论值均值以及方差可以涵盖了试验值结果，且吻合较好.

图2 理论结果与试验结果对比Fig.2 Comparison of theoretical results with experimental results

表1 模型计算数据与试验数据对比Table 1 Comparison of test data and model calculation data

从图2 可以看出，由于滚石碰撞的位置不同，即碰撞角不同，在相同的碰撞条件下，其反弹速度不是一个固定值，而是在一定范围内变化，再次说明了滚石碰撞位置影响颗粒反弹速度，颗粒反弹速度随撞击位置而变化，应有一个分布范围.因此，上述对比分析说明了该理论模型具有适用性，同时表明影响滚石--壁面碰撞过程的因素都是随机的，且对滚石运动产生明显影响.

2.2 坡体上被碰颗粒B 的速度VB 对反弹滚石运动特性概率分布的影响

斜坡上被碰滚石颗粒B 的速度VB反映了坡体运动的情况，已有文献[28-29]指出坡体的运动对滚石碰撞过程的影响不容忽视.因此，下面分析VB对反弹滚石颗粒运动特性概率分布的影响.

在此，我们选取了碰撞恢复系数k1=0.63，k2=0.21[35]；据文献[38]选取滚石颗粒的入射速度为(5 m/s，10 m/s，15 m/s，20 m/s)，同样假设滚石入射速度在该范围内为等概率分布，并重复上述频率的计算过程，得到坡角为γ=55◦时滚石颗粒撞击到静止坡面和运动坡面后，反弹速度和角度以及反弹角速度的频率(如图3 所示).

由图3 可知，反弹滚石颗粒的速度和角度以及反弹旋转角速度均服从高斯分布(y=y0+(A/e[−2(x−xc)/w]2)，坡体上被碰颗粒B 的速度对反弹颗粒线速度和角度以及反弹角速度概率分布形式定性上影响不大，但对其概率分布定量上有显著影响，例如随着坡面运动速度逐渐增大，反弹速度和反弹旋转角速度的高斯分布的中心值xc呈线性增大如图4(a)所示，其中反弹角速度的中心值影响显著；如图4(b)所示，反弹角度的高斯分布的中心值xc呈线性减小.即：当崩塌滚石碰撞到斜坡上其它运动滚石上时(VB0)，反弹速度和旋转角速度均值均要高于碰撞到斜坡上其它静止滚石上(VB=0)时的情况，坡体上被碰滚石颗粒的运动更利于滚石颗粒的反弹，且反弹后速度和旋转角速度增大，这一结论与文献[29]中的结论一致.也就是说在滚石灾害预测中，进入下一碰撞周期时的能量增大，灾害的破坏能力和范围也更大.

图3 VB 不同时，碰后反弹颗粒速度以及角度概率分布(VB=0，1，2，3，4，5 m/s，RA=RB=0.2 m)Fig.3 Rebound particle velocity and angular probability distribution after collision when the particle B’velocity is different(VB=0，1，2，3，4，5 m/s，RA=RB=0.2 m)

图4 反弹颗粒速度和旋转角速度的概率分布的中心值xc 随坡体上颗粒B 运动情况的趋势图(VB=0，1，2，3，4，5 m/s，RA=RB=0.2 m)Fig.4 The variation of the center value xc of the probability distribution of the rebound velocity and angular velocity with the particles B velocity(VB=0，1，2，3，4，5 m/s，RA=RB=0.2 m)

2.3 坡角对反弹滚石运动特性概率分布的影响

何思明指出在大于45◦的坡角下，更容易发生滚石灾害，且调查发现，最容易发生滚石灾害的坡度为55◦∼70◦[39].因此，坡角对滚石灾害的影响不容忽视.接下来进一步利用等间距直方图理论，重复上述频率计算过程得到了滚石与静止床面碰撞时，不同的坡角下反弹速度、角度以及反弹角速度的频率(如图5 所示).

图5 静止坡面，不同坡角下碰后反弹颗粒速度以及旋转角速度的概率分布(γ=10◦，20◦，30◦，45◦，55◦，65◦，RA=RB=0.2 m)Fig.5 Probability distribution of rebound particle velocity and rotational angular velocity on a stationary slope with different slope angles(γ=10◦，20◦，30◦，45◦，55◦，65◦，RA=RB=0.2 m)

图5 静止坡面，不同坡角下碰后反弹颗粒速度以及旋转角速度的概率分布(续)(γ=10◦，20◦，30◦，45◦，55◦，65◦，RA=RB=0.2 m)Fig.5 Probability distribution of rebound particle velocity and rotational angular velocity on a stationary slope with different slope angles(γ=10◦，20◦，30◦，45◦，55◦，65◦，RA=RB=0.2 m)(continued)

由图5 可知，反弹滚石颗粒速度和角度以及反弹旋转角速度均服从高斯分布形式(y=y0+(A/())e[−2(x−xc)/w]2)，从图5(b)中，可以看出坡角为45◦时，反弹角的概率密度曲线虽然仍服从高斯分布，但是与其他角度的曲线图趋势有明显不同，由图5(b)可知，γ=45◦的反弹角度概率密度相比γ=10◦，γ=20◦和γ=30◦的反弹角度概率密度明显增大；而当坡角大于45◦时，反弹角度概率密度明显增高，且概率密度峰值发生左移.即针对滚石反弹角度来说，坡角45◦是一个临界坡度值，这一结论与文献[39]定性上是一致的.

坡角对反弹线速度和角度以及反弹旋转角速度概率分布形式定性上影响不大，但对其概率分布定量上有明显影响，由图6 所示，随着坡角的逐渐增大，反弹速度和角度的高斯分布的中心值xc逐渐减小，反弹旋转角速度的高斯分布的中心值xc先增大后逐渐减小.

3 结论与讨论

3.1 结论

基于上述理论分析，得出以下的初步结论：

(1)本文建立了普适的斜坡滚石碰撞模型，该碰撞模型可描述滚石与斜坡、滚石与滚石、滚石与平板，滚石与运动床面的碰撞过程.

(2) 得到了滚石碰撞后反弹滚石颗粒速度解析解，这一解析解包含了影响滚石碰撞过程的6 个随机因素，即：坡角、坡体上被碰颗粒B 的速度和角度、入射速度和角度、撞击角度等随机因素.

图6 反弹颗粒速度和旋转角速度的概率分布的中心值xc 随坡角变化趋势图(γ=10◦，20◦，30◦，45◦，55◦，65◦，RA=RB=0.2 m)Fig.6 The variation of the center value xc of the probability distribution of the rebound velocity and angular velocity with slope angle(γ=10◦，20◦，30◦，45◦，55◦，65◦，RA=RB=0.2 m)

(3)该理论模型结果显示，反弹滚石速度和角度以及反弹旋转角速度的概率分布均服从高斯分布形式，坡体上被碰滚石颗粒B 的速度和坡角对反弹滚石速度和角度以及反弹旋转角速度的概率分布定性上没有影响，但对高斯分布的中心值参数有显著影响.

3.2 讨论

该理论碰撞模型反映的是滚石--斜坡单次碰撞过程，对滚石运动演化全过程的预测存在一定的局限性，不能直观表述滚石多次、连续的碰撞过程.但是本模型中得到的碰后速度与滚石进入下次碰撞的入射速度和能量有密切关系，且碰后速度的变化是计算滚石冲击力、冲击能量以及运动距离的关键输入参数，可为滚石灾害范围预测以及被动防护装置结构设计提供一定的理论依据.

另外，实际落石运动演化过程是不规则形状滚石三维运动过程，简化后的二维圆盘碰撞理论模型在主要运动方向上可以体现碰撞过程的运动特性，反映碰后速度变化，但是仍然会忽视了垂直于圆盘方向上的速度变化和能量耗散，并且简化后的圆盘忽略了不规则滚石形心和质心不重合带来的接触点和力的作用点变化的影响，使得预测结果与实际落石速度以及运动轨迹存在一定的误差.因此，在进一步的研究中，试图将二维圆盘碰撞拓展为三维球体，并将滚石的不规则形状考虑进去，使得滚石-斜坡碰撞模型更接近实际情况，进一步提高滚石-斜坡碰撞理论模型的预测精度.

附录

RA/m:崩塌颗粒半径

RB/m:坡体上被碰颗粒半径

J1:颗粒A 和颗粒B 碰撞时的接触位置

J2:坡体上被碰滚石颗粒B 与斜坡的接触位置

VA/(m·s−1):崩塌颗粒A 的入射速度

VB/(m·s−1):颗粒B 在碰撞前的速度

α/(◦):VB与水平方向的夹角

β/(◦):颗粒A，B 质心连线与水平方向的夹角

γ/(◦):斜坡坡角

θ/(◦):崩塌颗粒A 入射时与水平方向x轴所成的角度

mA/(kg·m−3):崩塌颗粒A 的质量

mB/(kg·m−3):崩塌颗粒B 的质量

/(kg·m·s−1):颗粒A 与B 之间的切向碰撞冲量

/(kg·m·s−1):颗粒A 与B 之间的法向碰撞冲量

UAx/(m·s−1):颗粒A 在碰撞结束时质心速度的水平分量

UAy/(m·s−1):颗粒A 在碰撞结束时质心速度的垂向分量

ωA/(rad·s−1):颗粒A 在碰撞结束时角速度

UBx/(m·s−1):颗粒B 在碰撞结束时质心速度的水平分量

UBy/(m·s−1):颗粒B 在碰撞结束时质心速度的垂向分量

ωB/(rad·s−1):颗粒B 在碰撞结束时角速度

Sτ/(kg·m·s−1):颗粒A,B 接触点处的切向碰撞冲量

Sn/(kg·m·s−1):颗粒A,B 接触点处的法向碰撞冲量

/(kg·m·s−1):碰撞过程中压缩变形达到最大时，斜坡作用在颗粒B 接触点处的切向冲量

/(kg·m·s−1):碰撞过程中压缩变形达到最大时，斜坡作用在颗粒B 接触点处的法向冲量

S1/(kg·m·s−1):碰撞结束时B 受到斜坡的切向冲量

S2/(kg·m·s−1):碰撞结束时B 受到斜坡的法向冲量

k1:法向碰撞恢复系数

k2:切向碰撞恢复系数