基于FLUENT剛性植被莖粗對坡面流特性影響研究

周建森，張升堂，張景洲，徐雪峰

（山東科技大學(xué)地球科學(xué)與工程學(xué)院，山東青島266590）

隨著我國經(jīng)濟(jì)社會的高速發(fā)展，坡面土壤侵蝕已成為我國環(huán)境問題一大挑戰(zhàn)。在減緩坡面土壤侵蝕方面，植被起著至關(guān)重要的作用。由于坡面植被生長狀態(tài)的不同，一定程度上改變了坡面水流內(nèi)部結(jié)構(gòu)，增加動能損失，對坡面流具有阻延、攔蓄作用，使坡面水流具有復(fù)雜多變的特性[1,2]。在河流和生態(tài)過程的影響中，河岸植被在洪水災(zāi)害、河流恢復(fù)和水生態(tài)管理中變得越來越重要[3]。近些年來，國內(nèi)外諸多學(xué)者開始關(guān)注生態(tài)環(huán)境，對含有植被的坡面展開了很多研究，主要是通過植被的行列走向，莖稈粗度，行列間距等方面，從而使人們對這一問題有了逐步深入的了解[4]。CAROPPI G 等用聲學(xué)多普勒測速儀測定了平均流速和湍流結(jié)構(gòu)，并從錄像中研究了植物的動態(tài)運(yùn)動。結(jié)果表明復(fù)雜的灌木狀植被經(jīng)歷季節(jié)性葉片脫落和隨著水動力的增加而導(dǎo)致阻力減小[5]。張升堂等分析在植被淹沒或非淹沒狀態(tài)不同植被覆蓋密度影響下，地表糙率的變化規(guī)律與特征。結(jié)果表明地表糙率取值與植被覆蓋密度成正相關(guān)；同一下墊面情況下，地表糙率取值會隨著不同的水流方向發(fā)生改變[6]。WANG Y T 等通過水槽實(shí)驗(yàn)研究3 種密度（密集、中間和稀疏）和3 個位置（山頂、后坡和腳坡）組合排列的剛性挺水植被的水流動力學(xué)，表明Re不是植被邊坡水力粗糙度的唯一預(yù)測因子。在坡向上，由于植被區(qū)前的時序效應(yīng)和植被區(qū)內(nèi)的快速輸運(yùn)效應(yīng)，植被坡上的所有水動力參數(shù)呈現(xiàn)出上下波動的趨勢[7]。王曉燕等使用PN 激光圖像速度場儀，獲取在相似水流條件下不同剛性植被周圍的流場矢量圖。結(jié)果表明：流場結(jié)構(gòu)會因?yàn)樗壑械闹脖?，出現(xiàn)大量漩渦，紊動強(qiáng)烈，變得越發(fā)復(fù)雜；隨著植被剛度的減小，流場結(jié)構(gòu)慢慢趨向穩(wěn)定[8]。王洪虎等采用聲學(xué)多普勒流速儀分析植物的水流流速沿垂線分布的特征，分析了植被密度對水流阻力的影響，植被密度的增加導(dǎo)致水位升高，水流阻力增大[9]。然而，對于明渠問題，大多是依靠實(shí)驗(yàn)，很少有人通過模擬的方法進(jìn)行研究。

也有部分對數(shù)值模擬的驗(yàn)證及應(yīng)用。魯婧等利用FLUENT與實(shí)驗(yàn)數(shù)據(jù)驗(yàn)證，得出相比于剛蓋假定方法，研究不同寬深比條件下的水流特性，VOF模型試驗(yàn)值更加接近試驗(yàn)值[10]。馮美娟等通過二維建模利用Fluent 中的VOF 模型得出在相同流量、床面粗糙度狀況下，水流平均速度、雷諾數(shù)Re、傅汝德數(shù)Fr以及阻力系數(shù)λ會隨著坡度的增大，有著小幅增長的趨勢，并且與以往的實(shí)驗(yàn)數(shù)據(jù)進(jìn)行對比驗(yàn)證，得出了模型的可靠性[11]。先前的幾位研究人員主要都是驗(yàn)證模型的可行性。ZHAO F 等通過大渦模擬不連續(xù)剛性水下植被斑塊對水流湍流的影響，結(jié)果表明間隙區(qū)域的速度明顯慢于冠層區(qū)域，間隙區(qū)域有利于水生生物的食物供應(yīng)和物理環(huán)境，也有利于沉積物沉積[12]。WANG W 等通過k-e模型對植被多級復(fù)式河道水流進(jìn)行數(shù)值研究，得出k-e模型可以很好捕捉了多級復(fù)合通道中的二次流，主河道底部的剪應(yīng)力大于植被階地。第1階地和第2階地之間的床面剪應(yīng)力沒有顯著差異，植被的存在可以有效地降低河床附近的剪應(yīng)力，避免淤積和減少侵蝕[13]。ANJUM N等通過雷諾應(yīng)力模型分析了明渠中部分分布不連續(xù)剛性植被紊流特性，植被斑塊的遮擋效應(yīng)降低了間隙區(qū)內(nèi)的流速，與單層植被流相比，其對通過垂直雙層植被的水流的影響顯著較高[14]。然而，在一些天然河流的河岸環(huán)境或農(nóng)田灌溉區(qū)，隨著時間的推移，即使坡面上存在相同的植被，莖稈粗度在不同的生長階段也會不斷變化。莖稈直徑對洪水沖刷和灌溉流量影響很大，而且了解植被對坡面流阻力影響規(guī)律，根據(jù)農(nóng)田植被類型控制坡面灌溉水流流速及流量，一方面避免農(nóng)作物出現(xiàn)倒伏情況，使得農(nóng)業(yè)生產(chǎn)力降低。另一方面可以提高水資源的利用效率，對高效灌溉和節(jié)約用水均有重要意義[15-17]。因此，研究植被的根莖大小對水流特性的影響很有必要，而且更貼近自然。

本研究的目的是對矩形明渠中不同根莖大小植被的復(fù)雜三維流動結(jié)構(gòu)進(jìn)行數(shù)值研究。保持在相同密度下，研究3種植被根莖大小，即直徑為3.5、5.0、6.5 mm，闡明由于植被形態(tài)變化而產(chǎn)生的流動結(jié)構(gòu)。采用雷諾應(yīng)力湍流模型，通過計算流體力學(xué)工具FLUENT進(jìn)行模擬。研究了植被形態(tài)的影響，特別是對流速和湍流的影響。在不同位置和斷面進(jìn)行水流特性測量，分析當(dāng)水流流經(jīng)植被時發(fā)生的湍流結(jié)構(gòu)變化。

1 材料和方法

1.1 建模設(shè)置和邊界條件

對160 mm 長和80 mm 寬的計算域進(jìn)行建模。該區(qū)域由縱向連續(xù)的剛性植被組成。x軸、y軸和z軸分別表示縱向、橫向和垂直方向。用圓柱體模擬植被，橫向相鄰植株間距為30 mm，縱向相鄰植株間距為40 mm，圓柱的高度都是60 mm。在坡度為0°的定坡狀態(tài)下進(jìn)行模擬。植被排列的模擬計算域的示意圖(俯視圖)如圖1 所示。采用了2 個縱向截面(Ls1 和Ls2)，即一個通過圓柱之間的間隙，以避免圓柱尾流的顯著影響，另一個直接通過圓柱陣列，以闡明植被的直接影響。此外，還研究了12 個重要位置，即位于x1-x12(位于植被之間和間隙)的垂直水流分布。

圖1 植被排列和指定位置的俯視示意圖（單位：mm）Fig.1 Top view of vegetation arrangement and designated location

水槽的幾何形狀較簡單，因此使用網(wǎng)格質(zhì)量高的結(jié)構(gòu)化三維六面體網(wǎng)格（見圖2）。為了保證計算精度，植被部分使用的是O-Grid 網(wǎng)格形式（見圖3），對明渠底部、邊壁、植被和水面進(jìn)行網(wǎng)格加密，網(wǎng)格數(shù)為3 032 000 個。為了獲得高質(zhì)量的模擬結(jié)果，還進(jìn)行了網(wǎng)格獨(dú)立性的試驗(yàn)，無論怎么增加網(wǎng)格，結(jié)果幾乎沒有發(fā)生改變。

圖2 結(jié)構(gòu)化三維六面體網(wǎng)格Fig.2 Structured 3D hexahedral mesh

圖3 植被O-Grid形式Fig.3 Vegetation O-Grid form

模型中采用的邊界條件是：進(jìn)口采用速度進(jìn)口，出口是自由出流，將模型一分為二，中間面為對稱面。采用基于壓力法的Pressure-Based 求解器，使用的是壓力修正算法，求解的控制方程是標(biāo)量形式的，對于不可壓縮流動有著很好的解。近壁處理采用標(biāo)準(zhǔn)壁函數(shù)，空間離散采用二階迎風(fēng)格式。采用SIMPLE 算法實(shí)現(xiàn)壓力-速度耦合。物理模型選用VOF模型，它合適求解分層流和需要追蹤自由表面的問題。初始化時，將上下分層，下方是水，上方是空氣，避免由于頂部邊界對實(shí)驗(yàn)結(jié)果帶來擾動。當(dāng)所有殘差都低于1×10-6時，認(rèn)為模擬已經(jīng)收斂。隨著進(jìn)一步的迭代，殘差沒有改變，并且還檢查了質(zhì)量流量，在隨后的迭代中不再改變。

1.2 數(shù)值模型的驗(yàn)證

對ZHANG 等[18]的實(shí)驗(yàn)數(shù)據(jù)進(jìn)行模擬驗(yàn)證。在長、寬、高分別為5.0 m、0.4 m、0.3 m 的矩形水槽中，植被株行距為60 mm×60 mm，圓柱直徑為3.0、4.0 和5.0 mm，通過室內(nèi)水槽沖刷實(shí)驗(yàn)記錄實(shí)驗(yàn)數(shù)據(jù)，水力參數(shù)實(shí)驗(yàn)數(shù)據(jù)見表1。

表1 流量參數(shù)Tab.1 Flow parameters

由圖4 可知，隨著水深h的增大，水流阻力系數(shù)f隨之增大，這與張景洲的研究結(jié)論一致。然而，植被莖稈粗度每增加1 mm，f平均增加40.4%。與實(shí)驗(yàn)結(jié)果相比，模擬結(jié)果出現(xiàn)3%的差異，這可能是由于模型的簡化或測速不可避免的測量誤差。實(shí)驗(yàn)結(jié)果與數(shù)值結(jié)果吻合較好，驗(yàn)證了數(shù)值模型的有效性。

圖4 數(shù)值模擬模型的驗(yàn)證圖Fig.4 Validation diagram of numerical simulation model

2 結(jié)果和討論

2.1 平均流動特性

2.1.1 速度的垂直分布

圖5顯示了指定位置的流向速度的垂直分布。平均流向速度參照初始平均速度“U”進(jìn)行無量綱化。植被區(qū)縱向流速呈“S”形分布，通道區(qū)縱向流速符合“J”形分布，坡面上不同位置垂線的無量綱縱向流速的分布與槐文信等的粒子圖像速度儀(PIV)測量結(jié)果相符[19]。在同一水深h下，不同植被莖稈粗度情況下，水流的流速不同，植被莖稈粗度越小對應(yīng)的流速值越大，3種植被莖稈粗度的水流流速關(guān)系為：V3.5mm＞V5.0mm＞V6.5mm。通過線性差值分析，與3.5 mm 粗度相比，粗度為5.0 mm和6.5 mm的平均流速分別降低了6.3%和10.6%。這一規(guī)律說明隨著植被的生長，莖稈粗度的變化對水流流速產(chǎn)生了一定的影響，有助于作物、植物的水肥管理。在水槽底部附近，由于阻力，流速降低到最小值。縱向兩植被之間的區(qū)域，受植被阻擋效應(yīng)慢慢減弱，直到完全擺脫阻擋效應(yīng)以及回流的影響，流速沿程成逐漸增大的趨勢，即Vx3＞Vx2＞Vx1、Vx6＞Vx5＞Vx4。在農(nóng)田灌溉中，更有利于化肥沉積在農(nóng)作物的附近。因?yàn)閤8 處旁邊即是植被，水流沖擊植被，向兩側(cè)擠壓，使得此處流速相較于x7 和x9 最大，x10、x11 和x12 也是如此的規(guī)律。這與ANJUM N 研究結(jié)果一致[20]。由圖5、圖7、圖9 可知，隨著植被粗度的增加，植被后方所選位置平均流速逐漸減小，即Vd=3.5mm,x1＞Vd=5.0mm,x1＞Vd=6.5mm,x1，x1～x6都是如此規(guī)律。植被莖直徑的增加，植被占水流斷面面積不斷增大，減小了水流過流的實(shí)際斷面面積，流經(jīng)植被的水流的濕潤周長增加，導(dǎo)致水流與植被之間的摩擦面積增加，提高了水流阻力，導(dǎo)致植被后方流速降低。由圖6、圖8、圖10 可知，當(dāng)植被莖粗變大時，通道處所選位置的流速逐漸變大，即Vd=6.5mm,x7＞Vd=5.0mm,x7＞Vd=6.5mm,x7，x7～x12 都是如此規(guī)律。這是由于植被的生長，通道處空間相對變小，水流擠壓通過造成。由圖5 可知，在d=3.5 mm 時，植被后方所選位置的流速成劇烈波動的趨勢，可能是因?yàn)樗鞯酿p小，而慣性力發(fā)揮越來越大的作用，流場變得不穩(wěn)定，這與槐文信等PIV結(jié)果相符[19]。

圖5 粗度3.5 mm植被后方指定位置的平均流向速度垂直分布Fig.5 Vertical distribution of average flow velocity at the designated position behind the vegetation with a roughness of 3.5 mm

圖6 粗度3.5 mm通道處指定位置的平均流向速度垂直分布Fig.6 Vertical distribution of average streamwise velocity at the designated position at the channel with a thickness of 3.5 mm

圖7 粗度5.0 mm植被后方指定位置的平均流向速度垂直分布Fig.7 Vertical distribution of average flow velocity at the designated position behind the vegetation with a roughness of 5.0 mm

圖8 粗度5.0 mm通道處指定位置的平均流向速度垂直分布Fig.8 Vertical distribution of average streamwise velocity at the designated position at the channel with a thickness of 5.0 mm

圖9 粗度6.5 mm植被后方指定位置的平均流向速度垂直分布Fig.9 Vertical distribution of average flow velocity at the designated position behind the vegetation with a roughness of 6.5 mm

圖10 粗度6.5 mm通道處指定位置的平均流向速度垂直分布Fig.10 Vertical distribution of average streamwise velocity at the designated position at the channel with a thickness of 6.5 mm

2.1.2 速度等值線分布圖

圖11～圖13顯示了深度為40 mm橫向剖面上速度等值線的空間分布，顯示了當(dāng)水流沖擊植被時流速發(fā)生的顯著變化。水流在經(jīng)過x=20 mm 植被時，由于植被的阻擋，水流會出現(xiàn)回流，使得植被前的流速降低，并且向兩側(cè)擠壓，植被兩側(cè)流速突增。植被的后方出現(xiàn)尾流和渦流，植被間的流速開始遞增，直到完全擺脫植被的阻擋效應(yīng)以及下一個植被的回流的影響，開始變得穩(wěn)定。隨著水流繼續(xù)前進(jìn)，經(jīng)過下一個植被也是呈如此的規(guī)律。通道處的流速成鋸齒狀，植被側(cè)邊的流速達(dá)到峰值，而且隨著植被粗度的增大，通道區(qū)的速度變得越來越大。這是因?yàn)橹脖淮侄仍龃?，通道區(qū)的寬度變小，由于水流擠壓，在剛經(jīng)過植被時流速變得越來越大。流速在經(jīng)過第1個植被時，兩側(cè)的流速達(dá)到最大值。這是由于當(dāng)水流經(jīng)過越多的植被，消耗了水流的部分動能，降低了流速。由于植被阻擋的影響，植被后方的速度比通道區(qū)的速度低。隨著植被的粗度增大，尾流的范圍增大。這種圓柱局部的流動結(jié)構(gòu)不容易通過實(shí)驗(yàn)研究得到，這表明了數(shù)值模擬的重要性。過去的研究表明，植被后面的尾流區(qū)有利于細(xì)顆粒的沉積，這進(jìn)一步刺激了植被的生長和擴(kuò)張[12]。因此，隨著植被的不斷生長，可以提供更多的遮擋，這表明水生生物與沉積物沉積的正反饋。

圖11 粗度3.5 mm在深度40 mm處的橫向剖面上速度等值線Fig.11 Velocity contour on the transverse section with a thickness of 3.5 mm and a depth of 40 mm

圖12 粗度5.0 mm在深度40 mm處的橫向剖面上速度等值線Fig.12 Velocity contour on the transverse section with a thickness of 5.0 mm and a depth of 40 mm

圖13 粗度6.5 mm在深度40 mm處的橫向剖面上速度等值線Fig.13 Velocity contour on the transverse section with a thickness of 6.5 mm and a depth of 40 mm

2.2 湍流特性

2.2.1 湍流動能

圖14 顯示了沿縱截面深度平均湍流動能(TKE)的分布。TKE是與湍流中渦流相關(guān)的每單位質(zhì)量的平均動能。湍流動能的特征是均方根速度的波動?？梢杂^察到，水流在經(jīng)過植被時受到干擾，并表現(xiàn)出不均勻性。3種情況中，對于水流流經(jīng)植被的界面上即Ls1，TKE沿程呈慢慢增大的趨勢，快經(jīng)過第1 個植被的時候即x=20 mm，急劇增大并且達(dá)到了第1 個峰值，隨后急劇減小，這是由于植被產(chǎn)生了較高的阻力，從而產(chǎn)生了較高的湍流動能。在植被之間繼續(xù)呈慢慢增大的趨勢，在水流流經(jīng)第2 個第3 個以及第4 個植被時即x=60 mm、x=100 mm、x=140 mm，TKE達(dá)到了第2到第4個峰值，并且第4個峰值的TKEx=140mm達(dá)到了最大，在流過最后一個植被后，TKE降到了最低并且達(dá)到了穩(wěn)定。這是由于缺少了植被的阻擋，TKE僅僅受底面和壁面的影響而且底面和壁面的影響不發(fā)生變化。流經(jīng)植被，進(jìn)入低流速區(qū)和高流速區(qū)，呈上升和下降趨勢，呈鋸齒狀分布。這一結(jié)果與Gurnell 研究的內(nèi)容相符[21]。植被間TKE整體趨勢也是呈逐漸增大的趨勢即TKEx=100～140mm＞TKEx=60～100mm＞TKEx=20～60mm。對于水流流過通道處即Ls2，在x=0～100 mm 處，TKE緩緩變化，在x=100～160 mm 趨于穩(wěn)定。這是由于通道處的TKE不受植被阻擋直接影響，因此沒有急劇的變化，但是因植被的阻擋產(chǎn)生的尾旋以及水流沖擊植被向通道處擠壓的影響，使得通道處的TKE發(fā)生緩緩的變化。并且隨著植被阻擋和擠壓效應(yīng)慢慢降低，TKE在x=100 mm 之后趨于穩(wěn)定。對于同一粗度的植被，植被區(qū)TKE均明顯高于通道區(qū)，即TKEd=3.5mm,Ls1＞TKEd=3.5mm,Ls2、TKEd=5.0mm,Ls1＞TKEd=5.0mm,Ls2、TKEd=6.5mm,Ls1＞TKEd=6.5mm,Ls2。此外，隨著植被粗度的增大，較粗的植被可以提供更多的阻力導(dǎo)致升降的幅度也變得更大。即TKEd=6.5mm,Ls1＞TKEd=5.0mm,Ls1＞TKEd=3.5mm,Ls1、TKEd=6.5mm,Ls2＞TKEd=5.0mm,Ls2＞TKEd=3.5mm,Ls2。

圖14 沿區(qū)域長度縱向分布的湍流動能Fig.14 Turbulent kinetic energy distributed longitudinally along the length of the region

2.2.2 湍流強(qiáng)度

圖15 給出了沿縱向截面的湍流強(qiáng)度（Urms）分布。Urms是均方根速度波動與平均流速之比，表示坡面流中的湍流，進(jìn)一步揭示了植被對紊流的影響。在所有情況下，由于植被的流動阻力，在植被結(jié)構(gòu)后面觀察到較高的湍流強(qiáng)度，隨后是尾渦，需要足夠的距離才能達(dá)到穩(wěn)定狀態(tài)。這一結(jié)果與以往研究工作中觀察到的結(jié)果一致[22]。與通道處相比，即同一粗度Ls1 與Ls2 處的Urms相比，植被區(qū)域的湍流強(qiáng)度顯著增加，表明植被周圍區(qū)域的湍流強(qiáng)度較高。因此，隨著植被生長，由于水流阻力和植被遮擋而形成的低流速處可能更適合水生生物的營養(yǎng)和物理環(huán)境。此外，由于當(dāng)植被直徑d越大時，水流阻力越大，水流阻力越大，消耗了水流的部分動能，降低了水流流速，從而隨著植被的生長，湍流強(qiáng)度逐漸增加，即Urmsd=6.5mm,Ls1＞Urmsd=5.0mm,Ls1＞Urmsd=3.5mm,Ls1。對于通道處的湍流強(qiáng)度也呈如此規(guī)律。在一年中不同生長階段，植被直徑的增大引起更大的尾流和渦流的產(chǎn)生，有利于農(nóng)田化肥沉積在植被前后。

圖15 沿區(qū)域長度縱向分布的湍流強(qiáng)度Fig.15 Turbulence intensity longitudinally distributed along the length of the region

3 結(jié) 語

隨著我國社會經(jīng)濟(jì)快速發(fā)展，我國越來越重視生態(tài)環(huán)境建設(shè)。植被作為生態(tài)系統(tǒng)的基本要素之一，對坡面徑流具有阻滯作用，具有蓄水保土、改良土壤、改善生態(tài)環(huán)境等諸多功能[23,24]。本研究采用基于RANS 技術(shù)的計算流體力學(xué)工具，利用FLUENT研究了剛性植被不同根莖大小的條件下，明渠中縱向連續(xù)植被的湍流特性。對結(jié)果進(jìn)行了定性和定量分析。實(shí)驗(yàn)結(jié)果與數(shù)值結(jié)果吻合較好，驗(yàn)證了數(shù)值模型的有效性，有助于作物、植物的水肥管理，為灌區(qū)精準(zhǔn)量水的發(fā)展提供了新思路。

研究結(jié)果表明：

（1）植被的遮擋效應(yīng)降低了后方區(qū)域的流速，通道處的流速成鋸齒狀，在植被側(cè)邊的通道處的流速達(dá)到峰值。兩株植被之間的流速逐漸增加。隨著植被粗度的增加，通道處的流速增加，植被后方的流速逐漸減小，以及導(dǎo)致更大的尾流和渦流的產(chǎn)生。

（2）通過線性差值分析，結(jié)果表明隨著植被粗度的不斷增大，平均流速不斷減小，與3.5 mm 粗度相比，粗度為5.0 mm和6.5 mm的平均流速分別降低了6.3%和10.6%。

（3）植被提供的阻力導(dǎo)致湍流動能和強(qiáng)度顯著增加，相反，水流在通道處，水流中的湍流動能和強(qiáng)度減少。湍流動能和強(qiáng)度與植被的粗度呈正相關(guān)關(guān)系。