泥石流是我国最主要的地质灾害之一，修建在泥石流高发区的铁路和公路桥梁，在建设运营过程中都面临艰巨的泥石流防治任务^[1-2].泥石流对桥梁的危害方式主要包括冲击、冲刷和淤埋，尤以大块石对桥梁的冲击力最为显著^[3-4].实践表明，许多桥梁被泥石流毁坏多与泥石流大块石对桥墩的冲击破坏有关^[5].因此，有关泥石流大块石的冲击力计算成为泥石流灾害防治研究的热点.目前，国内外计算泥石流大块石冲击力的主要方法较多.如章书成^[6]以结构力学为基础，将柱型构筑物简化为悬臂梁、简支梁模型来计算泥石流大块石作用于梳子坝支墩、桥梁等柱型结构上的冲击力.工程实践表明，此法计算结果比估计的承载力大^[7].Thornton^[8]视泥石流块石与构筑物分别为刚性体和塑性体来计算泥石流块石对构筑物的冲击力，其适用于强度较低的构筑物被冲击的情况^[3].Yamaguchi^[9]、Mizuyama^[10]、Braccesi等^[11]和黄宏斌等^[12]以Hertz弹性碰撞理论为基础，采用修正系数的方法来综合考虑系统材料、冲击过程中摩擦、断裂等因素对冲击力的综合影响.这种基于Hertz理论进行修正的方法是目前国内外计算泥石流大块石冲击力的主流方法之一.但其最大的不足是修正系数力学意义模糊，一般通过工程经验或实验确定，具有主观性强、成本高昂或耗时周期长等缺点.

本文针对泥石流中大块石对桥墩的冲击这一具体形式进行研究.以Hertz碰撞理论为基础，根据冲击系统中对象的材料特性、冲击物体的相对大小及结构变位对冲击力进行修正，其考虑因素全面，理论依据充分，参数明确，使用方便，可为泥石流多发区的桥梁设计及防护提供可信的计算方法.

1 Hertz碰撞理论

Hertz在假设两球接触面半径为a的圆的前提下，给出了两个球体在压力P作用下的完备解，并给出了碰撞时的冲击力公式(1)、(2) (图 1)^[13].

接触压力P(r)分布公式为

$P\left( r \right) = \frac{{2E{R^{\frac{1}{2}}}{\delta ^{\frac{3}{2}}}}}{{{\rm{ \mathsf{ π} }}{a^2}}}{\left[ {1 - {{\left( {\frac{r}{a}} \right)}^2}} \right]^{\frac{1}{2}}}.$

(1)

式中：$E = \frac{{{E_1}{E_2}}}{{\left( {1 - {\mu ^2}_1} \right){E_2} + \left( {1 - {\mu ^2}_2} \right){E_1}}}$为等效模量，$R = \frac{{{R_1}{R_2}}}{{{R_1} + {R_2}}}$为等效半径，a为接触半径，δ=δ₁+δ₂为总压缩量.其中μ₁、μ₂分别为两小球的泊松比，E₁、E₂分别为两小球的弹性模量，R₁、R₂分别为两个球的半径，δ₁、δ₂分别为两个球的压缩量.

两球碰撞时的冲击力

$F = \frac{{4\sqrt R }}{{3{\rm{ \mathsf{ π} }}({k_1} + {k_2})}}{\delta ^{\frac{3}{2}}}.$

(2)

式中${k_1} = \frac{{1 - \mu _1^2}}{{{\rm{ \mathsf{ π} }}{E_1}}}$，${k_2} = \frac{{1 - \mu _2^2}}{{{\rm{ \mathsf{ π} }}{E_2}}}$分别为两球的刚度.

根据牛顿第二定律，则有

$\frac{{{{\rm{d}}^2}\delta }}{{{\rm{d}}{t^2}}} = \frac{{{{\rm{d}}^2}{\delta _1}}}{{{\rm{d}}{t^2}}} + \frac{{{{\rm{d}}^2}{\delta _2}}}{{{\rm{d}}{t^2}}} = (\frac{1}{{{m_1}}} + \frac{1}{{{m_2}}})\left( { - F} \right),$

(3)

由式(2)、(3)得

$m\frac{{{{\rm{d}}^2}\delta }}{{{\rm{d}}{t^2}}} + \frac{{4\sqrt R }}{{3{\rm{ \mathsf{ π} }}({k_1} + {k_2})}}{\delta ^{\frac{3}{2}}} = 0.$

(4)

对式(4)进行积分，则弹性碰撞力

${F_{\rm{e}}} = {\left( {\frac{{5m{V^2}}}{4}} \right)^{\frac{3}{5}}}{\left( {\frac{{4\sqrt R }}{{3{\rm{ \mathsf{ π} }}({k_1} + {k_2})}}} \right)^{\frac{2}{5}}}.$

(5)

式中$m = \frac{{{m_1}{m_2}}}{{{m_1} + {m_2}}}$为等效质量，其中m₁、m₂分别为两小球的质量.

2 模型的建立

桥墩截面形状一般为圆形、方形及其复合形状，墩底置于足够深度的位置，目前大多将其简化为悬臂梁结构和简支梁结构.实际上，桥墩结构不是纯粹的悬臂结构，更不是简支结构.工程实践中应合理地考虑上部桥梁对墩顶的约束作用.因而桥墩结构应为墩底为固定支座，墩顶为定向支座的一种复合约束结构.在块石冲击力作用下，由于桥梁与墩顶之间的摩擦力对墩顶的约束从而保证墩顶不发生位移时，按照桥墩两端为固定支座考虑；当冲击力足够大时，桥墩可往一个方向产生一定的位移，此时桥墩结构为墩底为固定支座、墩顶为定向支座外加恒定摩擦反力G·f作用的复合结构，其中G为桥梁上部结构作用在墩顶的重力，f为支座的动摩擦系数.

将大块石简化为均匀球体，以速度V撞击桥墩，撞击位置位于l_i处；桥墩截面为圆形，总长l.可以建立如图 2所示的泥石流大块石冲击计算模型.

3 冲击力公式的修正

目前基于Hertz碰撞理论的冲击力计算公式中，推导出的绝大多数冲击力公式是基于球体的应力应变关系为弹性碰撞这一假设条件得来，而实践中大多数工程材料为应力应变关系较复杂的弹塑性体，直接使用Hertz弹性碰撞的冲击力公式所计算出的结果必然偏大；大多数公式只考虑冲击物体的大小对冲击力的影响，事实上，被冲击物的大小也是影响冲击力大小的关键所在；物体碰撞过程中，被撞物体均有不同程度的变形，从能量守恒的观点来看，物体的被撞击位置不同，其变形程度也不同，从而相应的冲击力也一定不同.本文基于以上3个方面，量化冲击系统中对象的材料性质、相对尺寸大小及构筑物的结构变位对冲击力的影响，最终导出泥石流中大块石对桥梁冲击力大小的计算公式.

当被冲击物的半径R₂和质量m₂均趋于无穷大，即R₂=∞，m₂=∞时，由式(1)和(5)中等效质量和等效半径表达式可得R=R₁，m=m₁，由此导出冲击力计算公式为

${F_{\rm{e}}} = \left[ {\frac{4}{{3{\rm{ \mathsf{ π} }}{{({k_1} + {k_2})}^{0.4}}}}} \right]{\left( {\frac{{5\rho {{\rm{ \mathsf{ π} }}^2}}}{4}} \right)^{0.6}}{V^{1.2}} \times {R^2},$

(6)

式中ρ为冲击物体的密度.

式(6)仅适用于碰撞系统材料性质为弹性、被撞击物体尺寸为无穷大及忽略被撞击物结构变位的情况，因此，在工程实践中，直接应用式(6)计算出的冲击力势必会与实际值有较大差别.为此，引入相应的修正系数来考虑碰撞过程中系统对象的材料性质、相对尺寸大小和结构变位对冲击力的影响，则实际冲击力计算公式如下：

$F = {k_{\rm{e}}}\left[ {\frac{4}{{3{\rm{ \mathsf{ π} }}{{({k_1} + {k_2})}^{0.4}}}}} \right]{\left( {\frac{{5\rho {{\rm{ \mathsf{ π} }}^2}}}{4}} \right)^{0.6}}{V^{0.8}} \times {R^2}.$

(7)

式中k_c=k_c1k_c2k_c3为综合修正系数，其中k_c1，k_c2，k_c3分别为冲击系统的材料修正系数、球径比(R₂/R₁)修正系数和结构变位修正系数.

可以看出，式(6)与式(7)中冲击物速度的幂次方不同，这主要是本文将名义材料修正系数中的速度分项合并到冲击力计算公式，方便单独量化材料性质对冲击力的影响所致，具体见修正系数k_c1，k_c2，k_c3的推导过程.

3.1 材料修正系数k_c1

国内外学者发现，基于Hertz接触力学所提出的冲击力计算公式中，其计算结果往往比实际观测结果大1~2个数量级^[14]，这是由于实际工程材料与Hertz接触力学假定的弹性材料在力学效应上存在显著区别.何思明等以Hertz碰撞理论、Thornton接触力学模型和牛顿第二定律，推导出被撞结构为无穷大时的最大压缩量和冲击力计算公式^[15]为

${\delta _{{\rm{max}}}} = {\left[ {\frac{{15m{{(V\cos \theta )}^2}}}{{16{R^{\frac{1}{2}}}E}}} \right]^{\frac{2}{5}}}.$

(8)

${F_{{\rm{ep}}}} = \frac{4}{3}E{R^{\frac{1}{2}}}\delta _{\max }^{\frac{3}{2}} - 2{\rm{ \mathsf{ π} }}\int_0^{{a_{\rm{p}}}} {\left[ {P\left( r \right) - {P_y}} \right]} rdr.$

(9)

式中：θ为冲击角度，a_p为构筑物塑性区半径，P_y为接触屈服压应力.

可以看出，基于Hertz接触力学推导的冲击力计算公式(6)和(9)仅在碰撞系统的材料上有区别.因此，可根据式(6)、(9)推导出材料修正系数k_c1.

根据Hertz接触力学理论和Thornton接触力学模型，存在如式(10)、(11)^[16]的关系：

$R\delta = a_{\rm{p}}^2 + a_{\rm{y}}^2 = {a^2},$

(10)

${\rm{ \mathsf{ π} }}R{P_{\rm{y}}} = 2E{a_{\rm{y}}}.$

(11)

式中a_y为初始屈服对应的接触面半径.

假设冲击过程中被冲击构筑物材料满足Mises屈服准则^[17]，则初始屈服应力服从式(12)的关系：

${P_{\rm{y}}} = \left( {1.234 + 1.256\mu } \right)Y.$

(12)

式中Y为接触材料的屈服强度.

综合式(1)、(8)~(12)，可推导出考虑冲击作用为正碰，被冲击结构材料为弹塑性材料时的冲击力公式：

${F_{{\rm{ep}}}} = \frac{4}{3}E{R^{\frac{1}{2}}}{\left( {\frac{{15m{V^2}}}{{16{R^{\frac{1}{2}}}E}}} \right)^{\frac{3}{5}}}{\left( {\frac{{{a_{\rm{y}}}}}{a}} \right)^3} + {\rm{ \mathsf{ π} }}{P_{\rm{y}}}a_{\rm{p}}^2.$

(13)

文献[15]的研究结果表明，块石与混凝土组成的冲击系统的临塑冲击速度很小，实践中泥石流速度一般都会超过临界速度，其势必在接触面上形成大范围的塑性变形区，此时材料初始屈服半径较塑性区半径至少小一个数量级，因此，名义材料修正系数为

$k_{{\rm{c1}}}^* = \frac{{{F_{{\rm{ep}}}}}}{{{F_{\rm{e}}}}} = \frac{{3{\rm{ \mathsf{ π} }}\left( {1.234 + 1.256\mu } \right)Y}}{{{{(1{\rm{ }}280{\rm{ \mathsf{ π} }}{E^4})}^{0.2}}}}{V^{ - 0.4}}.$

(14)

实际材料修正系数只考虑材料对冲击力的影响，因此，将名义材料修正系数中的速度分项合并到式(6)中的速度分项，得出实际冲击力计算式(7)中最终的材料修正系数为

${k_{{\rm{c}}1}} = \frac{{3{\rm{ \mathsf{ π} }}\left( {1.234 + 1.256\mu } \right)Y}}{{{{(1{\rm{ }}280{\rm{ \mathsf{ π} }}{E^4})}^{0.2}}}}.$

(15)

3.2 球径比修正系数k_c2

由上述可知，冲击系统最大压缩量、塑性区半径与等效半径R分别呈反相关、正相关关系，而等效半径R随被冲击物半径的增长而增加，进而表明了冲击力大小与被冲击物尺寸大小有关.

为量化被冲击物尺寸对冲击力的影响，定义球径比为n₁=R₂/R₁，其中R₁、R₂是冲击物和被冲击物的代表半径.以被冲击物代表半径为无穷大作基准，则根据式(5)，可求出不同球径比下冲击力与被冲击物为半无限空间体时冲击力的比值F_en/F_e∞，其中F_en和F_e∞分别为不同球径比下的冲击力和被撞结构代表半径为无穷大的冲击力.则球径比修正系数k_c2最终表达式如下：

${k_{{\rm{c}}2}} = \frac{{{F_{{\rm{en}}}}}}{{{F_{{\rm{e}}\infty }}}} = \frac{{n_1^2}}{{{{(1 + n_1^3)}^{0.6}}{{(1 + {n_1})}^{0.2}}}}.$

(16)

3.3 结构变位修正系数k_c3

混凝土桥墩在不同大块石冲击力作用下，墩顶将会出现不同的变形响应.当冲击力在墩顶产生的附加力小于桥梁上部结构与墩顶的最大摩擦力时，墩顶将保持不动，此时桥梁结构视为两端固定的力学模型；当冲击力足够大时，墩顶产生一定量沿冲击方向的位移时，桥梁结构视为墩底为固定支座、墩顶为定向支座附加恒定动摩擦力的力学模型.

根据机械能守恒定律，冲击过程中将动能转化成桥墩的应变能，则由结构力学可推导出桥墩应变能公式：

墩顶无位移时

${V_{\varepsilon 1}} = \frac{{{F_i}^2l_i^3{{(l - {l_i})}^3}}}{{6EI{l^3}}}.$

(17)

墩顶有位移时，令Gf=n₂F，考虑到上部桥梁的自重、支座设计摩擦系数及大块石冲击荷载，一般情况下n₂≤0.2^[18]，则有

${V_{\varepsilon 2}} = \frac{{{F_i}^2\left( {4{n^2}_2{l^4} - 15{n_2}{l^2}_i{l^2} + 4{l^3}_il(1 + {n_2}) + 3{l^4}_i} \right)}}{{24EIl}}.$

(18)

冲击过程中，冲击在离固定端越近的位置其相应的冲击力越大.泥石流冲击桥墩时，其离墩底最近的冲击位置为R处，以此位置冲击力为基准，令冲击力结构变位修正系数k_c3=F_i/F_R，则有：

墩顶无位移时

${k_{{\rm{c}}3}} = {\left( {\frac{{{R^3}{{\left( {l - R} \right)}^3}}}{{l_i^3{{(l - {l_i})}^3}}}} \right)^{\frac{1}{2}}}.$

(19)

墩顶有位移时

${k_{{\rm{c}}3}} = {\left( {\frac{{4n_2^2{l^4} - 15{n_2}{R^2}{l^2} + 4{R^3}l(1 + {n_2}) + 3{R^4}}}{{4n_2^2{l^4} - 15{n_2}l_i^2{l^2} + 4l_i^3l(1 + {n_2}) + 3l_i^4}}} \right)^{\frac{1}{2}}}.$

(20)

冲击修正系数k_c3表达式的选择，需先经试算后再确定.首先假设在墩顶无位移的情况下求得相应的冲击力，此时桥墩两端为固定支座，可求出桥梁上部固定端的剪力.其次将桥梁上部固定端的剪力和摩擦反力比较，当上部固定端的剪力小于摩擦反力，即F(l-l_i)²(l+2l_i)/l³ < Gf，则假设成立，此时冲击修正系数k_c3表达式为式(19)；反之，假设不成立，此时墩顶会发生位移，冲击修正系数k_c3表达式为(20).

4 算例与分析

为验证本文方法的合理性，以汶川县登基沟处桥梁受泥石流大块石冲击力计算予以说明.汶川县登基沟泥石流沟所携带大块石为花岗岩，沟口桥墩采用C35钢筋混凝土浇筑而成.根据登基沟泥石流观测数据，取泥石流大块石速度为4.6 m/s研究泥石流大块石对桥墩的冲击力，其余影响块石与桥墩冲击力的相关计算参数见表 1、2.

代入计算参数，可以求得k_c1=0.262，k_c2=0.633，k_c3=0.401，k_c=0.067，最大冲击力F=2 280.230 kN.如直接采用Hertz弹性理论进行计算，其理论冲击力F′=34 033.284 kN，远大于目前国内外通过现场观察所得的经验估计值^[11-12].本例计算结果表明，在大块石冲击墩柱一类的冲击力计算中，系统材料和刚度对冲击力计算结果影响显著，综合考虑冲击系统特性计算出的冲击力仅为Hertz弹性理论计算结果的6.7%，所得计算冲击力更贴近实际情况，表 3为典型泥石流块石冲击灾害过程中冲击力预估值与Hertz理论计算值比较.

此外，以大块石速度为变量绘制了登基沟泥石流大块石冲击力随大块石冲击速度关系曲线，由图 3可以看出, 基于弹性材料的Hertz理论计算出的冲击力随冲击速度的增加而急剧增大，而本文的冲击力计算随速度的增加而增长较为平缓，这主要是由于计算采用的冲击系统材料为典型的弹塑性材料；与Mizuyama^[10]和Yamaguchi^[9]公式计算结果相比，本文冲击力计算结果偏小，主要是由于本文提出的计算公式综合考虑了冲击系统的材料及刚度特性.

由式(15)可知，冲击系统材料确定，冲击力计算中的材料修正系数随之确定.由于泥石流中常见块石弹性模量为10.5~135.5 GPa，泊松比为0.10~0.35；工程中常用混凝土弹性模量为22.0~32.5 GPa，泊松比为0.18~0.20.根据式(1)中等效模量表达式，可得等效模量为8.3~30.8 GPa.结合式(15)进行计算，绘制材料修正系数与等效模量的关系曲线.由图 4可以看出，材料修正系数随等效模量的增加而增加，其值在0.20~0.45.表明随系统材料等效模量的增加，系统弹性冲击特性愈加明显.

由式(16)及图 5可知，球径比修正系数随球径比的增大而增加.计算表明n₁=2时，k_c2=0.859；n₁=5时，k_c2=0.960；n₁=10时，k_c2=0.981.可以看出，当0≤n₁ < 2时，球径比修正系数增长显著；当2≤n₁ < 5时，球径比修正系数增长相对缓慢；当5≤n₁ < 10时，球径比修正系数增长变化不大；当n₁≥10时，球径比修正系数基本稳定.表明当被冲击结构代表尺寸较小时，应重视构筑物尺寸对冲击力的影响.当其代表尺寸增加到一定程度时，可不再考虑构筑物尺寸大小对冲击力的折减效应.

由式(19)、(20)及图 6可知，桥梁墩顶不同的位移方式对冲击力的影响显著不同.对于墩顶不发生位移的情况下，从桥墩底部至桥墩中部其修正系数逐渐减小，衰减速率逐步降低，在桥墩中部位置时修正系数达最小值；对于墩顶沿撞击方向发生了一定位移的情况下，结构变位修正系数不仅与冲击位置有关，也与桥梁上部结构作用于墩顶的摩擦力有关，上部结构提供的摩擦力越大，相应的结构变位修正系数越大，增长速率逐步降低.在上部结构一定的情况下，结构变位修正系数自墩底到墩顶的方向先增大后减小.

5 结论

1)冲击系统的材料性质、相对尺寸大小及冲击位置对冲击力影响显著，工程设计中应加以考虑.

2)综合考虑冲击系统的材料性质、相对尺寸大小及构筑物结构变位情况所计算出的最大冲击力远小于Hertz弹性碰撞理论所得的冲击力，更贴近实际情况.

3)本文提出的修正冲击力计算公式可弥补工程经验和模型试验方法确定修正参数主观性强、离散性大的不足.其严格按照力学原理推导，理论依据充分，参数明确，修正计算采用显示解，使用方便，计算结果更为贴近工程预估值，可为泥石流区地质灾害设计提供有效参考.