李保如泥沙起动流速公式修正研究

蔡蓉蓉，张罗号，张红武

（1. 清华大学水沙科学与水利水电工程国家重点实验室，北京 100084；2. 河海大学水利水电学院，江苏南京 210098）

1 现有泥沙起动流速公式的简况

研究泥沙起动问题具有重大的科学意义与实践价值，有助于理解河床冲淤变形规律，对航道运渠等的设计［1］起到决定性作用。研究泥沙起动问题时，往往从泥沙起动流速入手，早在200多年前，法国人杜步华就采用了试验方法观测了泥沙运动与流速的关系［2］，之后国内外科研工作者从各个角度进行了大量相关研究［3］。

沙莫夫于1959年提出的泥沙起动流速公式为［4］：

式中：Vc为泥沙起动流速；γs为泥沙容重；γ 为水流容重；D为泥沙粒径；h为水深。

我国学者曾采用试验的手段研究过泥沙的运动问题［5-6］。例如，李保如等学者早在1953年就建立过泥沙起动流速的如下表达式［7］：

显然，式（2）与式（1）的结构形式是相同的，说明我国学者在这一领域的研究水平已经很高，只是局限于试验资料的粒径范围，没有在当时得出可包含细沙的起动流速公式。

之后李保如［1］进一步研究了粗颗粒及细颗粒泥沙的起动流速问题，引入了当时水流运动研究习用的雷诺数，还考虑了水力学中近壁层流层厚度的概念，进而可以考虑水温及含沙量对起动流速的影响。文献［1］进行受力分析时仅考虑水流推力FD，上举力FL及重力W，建立了泥沙颗粒的受力平衡方程式为：

式中f为“起动系数”，可反映泥沙起动时所受边界的摩擦作用，相当于摩擦系数［1］。

一般天然河道的坡度较小，可近似取sina=0，cosa=1，式（3）可简化为：

文献［1］将所给出的水流推力、上举力及重力的表达式代入上式，即可得出：

式中：C为谢才系数；m在天然情况下可假定为常数，为单颗泥沙承受推力的份数；k为上举力大小与水流推力强度的比值；A1、A2、A3均是常数，取决于泥沙颗粒形状。

将张友龄［8］建立的糙率n与泥沙粒径D的关系式代入曼宁公式，可将谢才系数C表达为如下形式：

式中：K ′为常数；R为水力半径。

李保如指出，摩擦系数f是泥沙粒径D与近壁层流层厚度δ 之比的函数［1］：

由此看出，文献［1］的研究，实际上就相当于将细颗粒的摩擦系数f处理为近壁层流层厚度对细颗粒泥沙“隐蔽”作用的综合影响系数。近壁层流层在现代水力学中又被称为黏性底层，其内黏性切力作用很大，该处流体受到边壁的约束作用，且边壁的绝对粗糙度与该层厚度的比值可反映边壁的摩阻特性［9］。近壁层流层厚度的表达式为［10］：

式中：N为常数，一般取11.6；ν 为水流运动黏滞系数；u*为摩阻流速。

联立式（6）及式（8），可推得D/δ 的表达式为：

联立式（5）及式（6）得：

联立式（7）、式（9）及式（10）得泥沙起动流速的一般形式为：

该式形式与式（1）相同，只不过系数为泥沙起动雷诺数Revd及相对糙率D/R的函数。

李保如利用式（11），结合前人当时的试验资料提出了泥沙起动流速的分区公式：

张红武认为，李保如泥沙起动流速公式适用于黏结力尚不占主导的细颗粒泥沙［11］，对于黏性细颗粒泥沙，使用该式计算的起动流速偏小很多，这是由于当时缺乏典型黏性细颗粒泥沙的试验资料，导致公式在黏性细颗粒泥沙起动流速的走向上表现不佳［12］。

张瑞瑾考虑了细颗粒间的黏结力，认为黏结力与沙粒空隙厚度成反比，与粒径的平方成正比，并与沙粒所受垂直压力有关［3］，借助实测资料率定出4个未知参数，于1961年提出泥沙起动流速统一公式［13］：

该式单位采用kg·m·s。

唐存本［14］认为细颗粒间的黏结力由泥沙颗粒及水膜间的分子力决定，与泥沙干容重有关，进行3个参数的率定后，提出适用于粗颗粒及细颗粒的起动流速公式［3］：

式中：m为指数，一般天然河道可取1/6，平整河床可按m=(1/4.7 )( D/h)0.06计算；γ′为淤积物的干容重；为淤积物的稳定干容重；系数C=0.906×10-4g/cm。该式单位采用g·cm·s。

沙玉清只考虑水流作用力及重力，认为水流作用力系数Cp与泥沙静水沉降时的阻力系数Cd满足指数关系，并使用试验资料于1956年求得泥沙的起动流速公式［2］：

式中ω 为泥沙颗粒沉速，可采用沙玉清公式［15］计算，即滞流区（D=0.0001～0.01 cm）公式：

介流区（D=0.01 ～0.2 cm）公式：

紊流区（D＞0.2 cm）公式：

上述式（19）—（21）单位采用g·cm·s。

沙玉清在1965年对其早期公式（式（18））进行了改进，在“阻力系数”中考虑了散粒体颗粒与黏性沙的区别，引入了分子水膜厚度及孔隙率的影响，提出了能概括粗细颗粒泥沙起动的流速公式［17］：

式中：ε 为孔隙率，稳定值约为0.4；δ 为分子水膜厚度，取0.0001 mm。

窦国仁1958年通过交叉石英丝试验证明附加下压力会对细颗粒泥沙的起动产生影响，考虑了压力水头对黏结力的影响作用，于1960年求得浅水下的泥沙起动流速公式［18］：

式中：Ha是用水柱高度表示的大气压力；δ 是一个水的分子厚度，取3×10-8cm 。

窦国仁1974 年在考虑附加下压力对黏结力影响的基础上又考虑了分子黏结力对泥沙起动的影响，对其早期的公式（式（23））进行修正，提出的泥沙起动流速公式为［19］：

式中：Ks为河床糙度，当D小于等于0.5 mm时，Ks=0.5 mm，当D大于0.5 mm时，Ks=D；δ 为与沙粒缝隙大小相关的特征厚度，δ=0.213×10-4cm，εk=2.56 cm3/s2，系根据交叉石英丝试验成果确定。该式单位采用g·cm·s。

张红武考虑了颗粒间的黏结力及附加下压力，并在黏结力的表达式中引入了水流运动黏滞系数（可反映出近壁层流层厚度，水温及含沙量［20］的影响），在作用于沙粒的有效瞬时流速与断面平均流速的转换中考虑了含沙量及河床摩阻的影响，于2012年提出适用于粗细沙的起动流速统一公式为［20］：

式中：KD为作用于沙粒的有效瞬时流速与断面平均流速的转换系数；δ 参考窦国仁交叉石英丝试验的成果，取为0.213×10-4cm。

1950年代，我国沙玉清及李保如两位学者基于河流动力学原理先后推导出泥沙起动流速的计算公式，其中李保如考虑了水流运动黏滞系数的影响。根据试验观察，同一粒径模型沙的起动流速大小受到水温的影响［21］，说明水流运动黏滞系数是影响起动流速的重要因素，必须将其考虑在内。这一说法可以通过李保如的观点解释，水温升高，近壁层流层厚度减小，对于泥沙颗粒的遮蔽作用减弱，故泥沙的起动流速减小。李保如通过理论与试验研究取得的成果物理意义清晰，代表了当时起动流速研究的较高水平，具有重大意义。

2 李保如泥沙起动流速修正公式

前人推求泥沙起动流速公式时，考虑的作用力及推得的具体表达式各具特色，但都通过实测资料率定出了不止一个系数才得最终的起动流速公式［12］。谢鉴衡［3］指出，“只要适当的选择系数，使这些公式比较符合实测资料还是有可能的”。鉴于李保如公式的物理图形既考虑了泥沙颗粒所处边界的水力学差异，又能体现近壁层流层对细颗粒泥沙的“隐蔽”影响，这在当时的国内外学者中是极为难得的，本文以该公式的物理背景为基础进行修正，期望得到适用于黏性沙及散粒体泥沙的起动流速统一公式。

为此，可将摩擦系数f的函数处理为由近壁层流层厚度对细颗粒泥沙“隐蔽”影响项与可忽略不计该影响的项组成，同时考虑到细颗粒泥沙粒径越细，近壁层流层的“隐蔽”影响越大，参照李保如当时对细颗粒泥沙起动影响的处理思路，进一步将前者表示为“近壁层流层厚度δ 与泥沙粒径D之比”即δ/D 的指数形式，于是将f的函数表示为：

式中k0、m1、C0为待定系数。

直接用曼宁公式表示谢才系数C，引入近壁层流层厚度的表达式（8）推得δ/D 的表达式为：

将上式代入式（26）得：

将式（28）代入式（10），表明在忽略近壁层流层厚度影响时，公式形式与沙莫夫公式相同，故根据能够忽略近壁层流层厚度影响的较粗颗粒泥沙实测起动流速资料，求得能够反映重力作用的系数C0=1.5；利用黏性细颗粒泥沙资料，将近壁层流层厚度对细颗粒泥沙“隐蔽”影响项的系数k0=0.0035，指数m1=2；最终得出如下泥沙起动流速公式为：

该式单位采用kg·m·s。借助试算法，利用式（29），即可求得泥沙的起动流速。

此外，通过上式中的泥沙起动雷诺数Revd，可体现出含沙量对泥沙起动的影响。浑水条件下水流的运动黏滞系数νm［22］为：

式中：Sv为体积含沙量；d50为悬移质泥沙中值粒径，单位为mm。

采用21 组清华大学黄河研究中心水槽试验资料（粒径范围0.001 ～0.193 mm，水深均换算为0.15 m）及文献［3］中窦国仁整理的各家实测资料与从长江实测记录换算而得的资料共计105组泥沙起动资料（粒径范围0.20～78.14 mm，水深均换算为0.15 m），分别求出不同粒径泥沙颗粒近壁层流层影响项及重力作用项所占的比例，绘制成图，如图1所示。可以看出，随着泥沙粒径的逐渐增大，重力作用项所占比例逐渐增加，近壁层流层影响项所占比例逐渐减少。从工程角度讲，泥沙粒径较大时（如大于0.1 mm），可忽略近壁层流层影响项，泥沙粒径较小时（如小于0.01 mm），可忽略重力作用项。

3 李保如泥沙起动流速修正公式的验证

图1 泥沙粒径与近壁层流层影响项及重力作用项的关系

本文除利用上述试验资料外，又利用文献［2］图5中何之泰、俞世煜、姜国干等人整理的共计67组泥沙起动资料（粒径范围0.001～70.63 mm，水深均换算为1 m），取指数为1/6，使用指数流速分布公式将上述水深1 m下的实测起动流速资料转换为水深0.15 m下的起动流速资料，利用式（29）计算得出了水深0.15 m时不同粒径泥沙的计算起动流速（用张友龄公式确定n与D的关系，其系数A按张红武公式［12］计算），并与窦国仁、张瑞瑾等各家著名公式相比较（为便于比较其中张红武公式中谢才系数C取60 m1/2/s），如图2所示。从该图可看出，对于粒径小于0.01 mm的细颗粒泥沙，除唐存本公式计算结果稍大以外，各家公式计算结果相差不大；对于粒径在0.01～0.1 mm之间的泥沙，除沙玉清公式及唐存本公式计算结果略有偏大以外，各家公式计算结果相近；对于粒径在0.1～1 mm之间的泥沙，各家公式计算结果较为接近；对于较粗颗粒泥沙，张红武公式、沙玉清公式，唐存本公式及本文公式与实测点据的吻合程度较高，张瑞瑾公式与窦国仁公式的计算结果偏小。可以看出，与实测资料及各家公式相比，本文公式颇为符合实际，取得了较好的验证效果。

采用式（29）计算得出了粒径为0.013 mm的泥沙在不同水深下（0.15 ～20 m）的计算起动流速，并与各家公式计算结果相比较，如图3所示，本文公式可以较好的反映出泥沙起动流速随水深增加而增大的规律。

图2 各家公式同实测资料的比较结果（水深0.15m）

图3 泥沙起动流速随水深的变化规律（D=0.013mm）

图4 本文公式计算结果与实测资料的比较

为便于比较各家公式的计算结果，使用文献［23］中窦国仁整理的长江资料（粒径范围0.1～0.18 mm）、南科院不同容重模型沙资料（粒径范围0.63～8 mm）、7 组黄河实测资料（粒径范围0.059～0.162 mm）、3 组伏尔加河资料（粒径范围0.383～0.73 mm）、谢鉴衡整理的模型沙资料（粒径范围0.24～0.8 mm），前苏联中亚河资料（粒径范围8.3～13 mm）及3组汾河河津资料（粒径范围0.095～0.169 mm），利用图示法对李保如泥沙起动流速修正公式进行验证，验证结果如图4所示，通过计算求得使用本文公式的计算结果与实测资料的相关系数为0.966，平均相对误差为17.9%。由文献［20］及文献［23］可知各家公式对相同资料的验证结果，平均相对误差分别如表1所示。结合该表及本文、文献［20，23］中各公式的图示法验证结果看出，本文公式具有较高精度，可以满足实际需要。

大量资料验证结果表明，尽管不去考虑概念难以判明且难以定量描述的黏结力及附加下压力的复杂影响，也不考虑流速沿水深的分布计算公式如何表达，但得到的公式同样可以同时适用于黏性细颗粒泥沙及散粒体泥沙的起动流速计算。

表1 各家公式计算值与实测资料相比的平均相对误差

鉴于黄河下游冲刷只有把床面细沙冲起进入悬移运动才显得有效，黄河下游取D=0.03mm、H=1.5 m，根据汛期糙率n=0.011、非汛期糙率n=0.015，利用式（29）计算得到黄河下游汛期起动流速Vc=0.49 m/s、非汛期起动流速Vc=0.413 m/s。再以D=0.05 mm、H=2.5 m、汛期糙率n=0.013、非汛期糙率n=0.017 计算宁蒙河段，可得汛期起动流速Vc=0.414 m/s、非汛期起动流速Vc=0.37 m/s。由此说明，引入糙率n 的起动流速公式可说明黄河沙质河段汛期起动流速都比流速（约为1.6～2.5 m/s）小很多，非汛期流速（约为0.8～1.5 m/s）减小较多，但起动流速比汛期小，使起动流速仍比流速小不少，试图揭示出黄河沙质河段往往非汛期冲刷但不出现显著冲刷的河流动力学背景，为了解非汛期小流量输沙潜力提供了基础。

4 结论

本文在分析了前人泥沙起动流速公式基础上，认为李保如泥沙起动流速公式的物理图形既考虑了泥沙颗粒所处边界的水力学差异，又能体现近壁层流层对细颗粒泥沙的综合“隐蔽”影响，只是由于当时缺乏黏性细颗粒泥沙的试验资料，使得公式仅适用于黏结力不占主导的泥沙颗粒。以李保如泥沙起动流速公式研究思路为基础，将摩擦系数f的函数修正为近壁层流层厚度对细颗粒泥沙“隐蔽”影响项及重力作用项两项组成，把前者表示为近壁层流层厚度δ 与泥沙粒径D之比δ/D 的指数形式，后者为常数，并利用实测泥沙起动流速资料确定出待定系数，即得出新的泥沙起动流速公式，克服了原公式分4个区间给出泥沙起动流速表达形式的缺点。分析该式结构发现，泥沙粒径较大时，可忽略前者，泥沙粒径较小时，可忽略后者。利用大量泥沙起动实测资料验证表明，尽管不去考虑概念难以判明且难以定量描述的黏结力及附加下压力的复杂影响，且也不引入目前还不成熟的流速沿水深的分布公式，但只要合理考虑近壁层流层厚度对细颗粒泥沙起动的影响，建立的公式同样能够既适用于黏性细颗粒泥沙又能适用于散粒体泥沙的起动流速计算，同时可以反映出泥沙起动流速随水深增加而增大的变化规律。