白玉川，徐國強

（1.天津大學 水利工程仿真與安全國家重點實驗室，天津 300072；2.天津大學 河流海岸工程泥沙研究所，天津 300072）

kinoshita派生彎道水流三維數(shù)值模擬

白玉川1,2，徐國強1,2

（1.天津大學 水利工程仿真與安全國家重點實驗室，天津 300072；2.天津大學 河流海岸工程泥沙研究所，天津 300072）

本文采用雷諾方程和雷諾應力紊流模型相結(jié)合建立起彎道的三維數(shù)學模型，并結(jié)合shukry彎道實驗結(jié)果對模型進行驗證。利用此模型對 kinoshita型非對稱彎道進行模擬研究。結(jié)果表明：非對稱型彎道的高流速區(qū)比對稱型彎道更加輻聚，水流動力軸線更加的不順滑，主流區(qū)在彎頂?shù)捻敍_作用更強；對稱型彎道的斷面環(huán)流強度不斷衰減，而非對稱型彎道則不同；隨著偏斜度的增大，斷面環(huán)流強度在彎道前半段不斷減小，后半段則相反。

非對稱型彎道；三維數(shù)學模型；kinoshita曲線；偏斜度

引 言

天然彎曲河流的平面形態(tài)多呈現(xiàn)出非對稱的平面形態(tài)，這對彎道水流有著重要的影響。現(xiàn)有的研究多是針對對稱型彎道進行研究，而對于非對稱型彎道的研究較少。kinoshita派生曲線與非對稱型河流的平面形態(tài)吻合較好，本文根據(jù)kinoshita派生曲線建立多組非對稱型彎道，并采用雷諾應力模型對其進行數(shù)值模擬，研究非對稱型彎道的水流運動特性。

1 非對稱型彎道的三維數(shù)學模型

1.1 控制方程

描述不可壓縮流體的控制方程是雷諾方程（Reynolds equation），包括連續(xù)性方程和動量方程。為了使雷諾方程封閉，本文采用Reynolds應力紊流模型。

連續(xù)性方程：

動量方程：

Reynolds應力運輸方程：

k方程：

ε方程：

12C1ε=1.44，C2ε=1.92，Cμ=0.09，σk=0.82，σε=1.0。

1.2 模擬自由表面的VOF方法

VOF方法的基本思想：用體積率函數(shù)表示不同流體在計算單元中所占的比例。對于兩相流來說，分別定義aw(x,y,z,t)和aa( x,y,z,t)表示水和空氣在計算單元中的體積比。在任何一個計算單元中水和空氣的體積比之和均為 1。因此任何一個單元會有全是水、充滿空氣或水氣都有三種情況。

通過求解下面的連續(xù)性方程來跟蹤水氣界面：

引入VOF后，自由表面處的密度ρ和動力粘度μ不在是常數(shù)了，而是體積分數(shù)的函數(shù)了。它們可以由下式計算得到：

式中：μw和ρw分別表示水的動力粘度和密度，μa和ρa分別表示空氣的動力粘度和密度。

1.3 邊界條件與數(shù)值計算方法

本文采用有限體積法對方程進行離散，為了提高計算精度，網(wǎng)格采用結(jié)構(gòu)化六面體網(wǎng)格。網(wǎng)格首先在s-n方向剖分，網(wǎng)格垂直于彎道中心線，然后在垂向（z方向）上加密生成結(jié)構(gòu)化六面體網(wǎng)格，網(wǎng)格三維效果如圖1。對于壓力和速度場采用PISO算法進行耦合計算。

圖1 彎道網(wǎng)格剖分三維效果

彎道進口分為水速度進口和空氣壓力進口，水速度進口采用常數(shù)分布，彎道出口采用壓力出口，彎道底部和側(cè)壁均采用固定邊壁，粘性底層采用壁面函數(shù)處理。

1.4 模型的驗證

采用Shukry[1]的彎道實驗對模型進行驗證。圖2為模型彎道縱向垂線平均速度的實測值（圖2a）和計算值（圖2b）的對比圖。從圖2中可以發(fā)現(xiàn)，實測值和計算值的流速分布吻合較好，高流速區(qū)在彎頂凸岸附近，彎頂下游低流速區(qū)也吻合較好。本文采用的模型適合于彎道水流的模擬，本文采用上述模型對非對稱型彎道進行模擬，研究非對稱型彎道的水流特性。

圖2 等值線實測值與計算值對比

2 非對稱彎道的數(shù)值模擬

kinoshita派生曲線[2]考慮了河流的偏斜度和豐盈度，描述彎道的平面參數(shù)更加的完備，可以通過調(diào)整最大偏角、彎道的波長、豐盈度系數(shù)、偏斜度系數(shù)四個參數(shù)來描述彎曲河流的平面形態(tài)，kinoshita派生曲線能夠很好的吻合天然彎曲河流的平面形態(tài)。

式中：Js為偏斜度系數(shù)；Jf為豐盈度系數(shù)；L為彎道的波長；s為彎道水流方向的坐標。

偏斜度系數(shù)對彎道的非對稱特性有著重要的影響，因此本文通過改變偏斜度系數(shù)，建立不同的非對稱型彎道，探究其水流特性。具體的模擬工況如表1，圖3為各個工況下彎道中心線，圖中可以看出，S1為對稱型彎道，S2、S3為非對稱型彎道。

表1 不同工況的參數(shù)設(shè)置

圖3 各個工況下彎道中心線

2.1 縱向流速

縱向流速對彎道水流有著重要作用，為了分析非對稱彎道下的縱向流速，分別繪制各個工況彎頂斷面的縱向流速。

圖4 各個工況下彎頂斷面縱向速度

從圖中可以看出，縱向速度在橫向上存在明顯的坡降，高流速區(qū)（2倍的平均速度）主要集中在靠近凸岸的0.1倍河寬范圍內(nèi)，凹岸的速度明顯小于凸岸速度。S1工況下彎頂處斷面的高流速中心主要集中 0.4～0.6倍的水深處，并且比較對稱。S2、S3彎頂斷面的高流速中心則是不斷的下移，靠近彎道水槽底部。非對稱型彎道的彎頂高流速區(qū)中心比對稱型彎道更加靠近水槽底部?？梢酝茰y，隨著偏斜度的增大，水流會不斷的侵蝕凸岸的底部，可能會在此形成崩岸，促使河灣由曲變直，使河灣加快消亡。

2.2 水流動力軸線

水流的縱向垂線平均速度是彎道水流一個重要特征，可以整體反映水流的特性。

圖5為各個工況下的縱向垂線平均速度等值線圖，圖中實線為水流的動力軸線。從圖中可以看出對稱型彎道和非對稱型彎道均會在彎頂凸岸附近形成高流速區(qū)，而在凹岸形成低流速區(qū)。隨著偏斜度的增大，高流速區(qū)在彎頂凸岸越來越輻聚，流速也在不斷的增大，而緊靠彎道凸岸下游流速會較小，當偏斜度達到一定程度時，此區(qū)域會形成回流。

S1工況下的水流動力軸線相對比較順滑，形態(tài)上比較對稱。隨著偏斜度的增大，水流動力軸線越來越不順滑，與彎頂凸岸的夾角越來越大，此時彎道主流與彎頂?shù)捻敍_作用較強。可以推測，回流區(qū)的形成會使泥沙在凸岸進一步的淤積，有利于凸岸邊灘的發(fā)展，進一步的促進彎道的發(fā)展，當彎道發(fā)展到一定的程度時，高流速區(qū)域直接切割凸岸邊灘，形成裁彎。

圖5 縱向垂線平均速度及水流動力軸線

2.3 斷面環(huán)流強度

將彎道沿彎道中心線均分成31個斷面，S1工況斷面布置如圖6。橫向環(huán)流強度沒有明確的定義，由于彎道中的水流運動是一種螺旋流，引入氣象學中的螺旋度概念表示環(huán)流強度。螺旋度（H）大小反映的是旋轉(zhuǎn)與沿旋轉(zhuǎn)軸方向的強弱程度[3]。螺旋度的計算公式如下：

為了更好的分析環(huán)流強度沿程的變化，將環(huán)流強度的絕對值在斷面上進行積分，得到斷面的環(huán)流強度（DH）。

圖6 斷面布置

圖7 斷面環(huán)流強度沿程分布

從圖中可以看出，S1工況下，斷面環(huán)流強度會出現(xiàn)三個峰值且峰值是不斷衰減，沿程呈現(xiàn)出一種類阻尼現(xiàn)象，而非對稱型彎道則不同。對比S1、S2和 S3可以看出，隨著偏斜度的增加，彎道斷面環(huán)流強度在彎道的前半段是不斷地減小，而在彎頂下游，尤其是出口處，斷面環(huán)流強度逐漸地變大。

2.4 床面切應力

彎道的床面切應力對于泥沙的輸運有著重要的作用，研究彎道的床面切應力分布可以預測泥沙在彎道中的沖淤情況。床面切應力按照下面的公式計算[4]：

式中：Uw為流速；y為到床面的距離；ks為壁面粗糙度厚度。

圖8 床面切應力分布

從圖中可以看出，S1彎道的切應力的范圍為0.1至0.85之間，在彎頂處切應力在橫向方向上存在明顯的坡降，床面切應力在彎頂凸岸較大，而在彎頂?shù)陌及肚袘t相對較小，切應力在彎道的兩側(cè)基本呈現(xiàn)出對稱型分布。S3彎道在進口右岸的切應力較低，切應力在彎道的上半段相對于明顯小于下半段，在彎頂凸岸的下游存在著較小的切應力區(qū)域。非對稱型彎道，切應力在彎頂凸岸處明顯比S1彎頂大。在彎道的前半段，S3切應力在橫向上存在明顯的坡降，而 S1彎道雖然存在著坡降，但沒有S3的坡降大。對比S1、S2、S3可以看出，隨著彎道偏斜度的增大，彎頂處凸岸的切應力不斷的增大，并且低切應力的區(qū)域范圍也在不斷增大。對比圖5和圖8，可以看出，彎道床面切應力分布與縱向垂線平均流速的分布基本一致，彎頂?shù)母吡魉賲^(qū)域切應力相對也較大，而低流速區(qū)域切應力相對也較小。

3 結(jié) 論

本文采用kinoshita派生曲線建立多組非對稱型彎道，并采用數(shù)值模擬的方法對非對稱型彎道進行研究，根據(jù)模擬結(jié)果分析，得出以下結(jié)論：

1）對稱型彎道和非對稱型彎道在彎頂凸岸處均存在高流速區(qū)，隨著偏斜度的增加，高流速區(qū)在彎頂凸岸會更加輻聚，非對稱型彎道的高流速區(qū)比對稱型彎道更加集中。此外水流動力軸線越來越不順滑，主流區(qū)與彎頂?shù)捻敍_作用較強?？梢酝茰y，隨著偏斜度的增大，彎頂凸岸附近的高流速區(qū)，可能切割凸岸邊灘，形成切灘型裁彎。

2）隨著偏斜度的增加，彎道斷面環(huán)流強度在彎道的前半段是不斷地減小，而在彎道的后半段，斷面環(huán)流強度逐漸地變大。

3）彎道床面切應力分布與縱向垂線平均流速的分布基本一致。對稱型彎道的床面切應力分布基本對稱，而非對稱彎道的床面切應力則不是，隨著偏斜度的增大，彎道前半段的切應力小于彎道后半段的切應力。

[1]Shukry A.Flow around bends in an open flume[J].Transactions of the American Society of Civil Engineers,1950,115:751-779．

[2]Kinoshita R．Investigation of channel deformation in Ishikari River[R]．Report to the Bureau of Resources,1961:1-174．

[3]陸慧娟,高守亭.螺旋度及螺旋度方程的討論[J].氣象學報,2003,61(6):684-691.

[4]Olsen N R B.A three-dimensional numerical model for simulation of sediment movements in water intakes with multiblock option[J].User’s Manual,The Norwegian University of Science and Technology,Trondheim,2006.

3D Numerical Simulation of Water Current along Kinoshita-generated Bend

Bai Yuchuan1,2，Xu Guoqiang1,2
(1.State Key Laboratory of Hydraulic Engineering Simulation and Safety,Tianjin University,Tianjin 300072,China；2.Institude of Sediment on River and Coast Engineering,Tianjin University,Tianjin 300072,China)

Reynolds equation and Reynolds stress turbulence model are combined to build a 3D mathematical model of bend,and the model is verified according to shukry bend testing results.Then this model is used to study kinoshita asymmetric bend.The results show that high flow velocity zone of asymmetric bend is more concentrated than that of symmetric bend,water current dynamic axis of the former one is not smooth as the latter one,and the impact force of main flow is greater at the apex of bend.The intensity of circulating flow on section of symmetric bend weakens continuously.With the increase of skewness coefficient,the intensity of circulating flow on the section of asymmetric bend decreases in the first half of bend and increases in the second half of bend..

asymmetric bend; 3D mathematical model; Kinoshita-generated curve; skewness

U611

：A

：1004-9592（2016）06-0007-05

10.16403/j.cnki.ggjs20160602

2015-12-13

自然科學基金資助項目（51279124）

白玉川（1967-），男，教授，主要從事泥沙運動力學及海岸動力學研究。