水電站不同布置下的攔污浮箱流固耦合計(jì)算分析

黃振峰 ，徐朋威，李金明，肖震岳

（廣西大學(xué)機(jī)械工程學(xué)院，廣西 南寧 530004）

水電站進(jìn)水口兩岸植被差，覆蓋范圍廣，各種垃圾雜物都隨水流進(jìn)入庫(kù)區(qū)，所以經(jīng)常堵塞上游的攔污柵，輕者水頭降低，影響水電站綜合效益；重則損壞攔污柵的結(jié)構(gòu)，進(jìn)而危及下游的發(fā)電機(jī)組。各水電站為了解決這一問題，根據(jù)自身?xiàng)l件增設(shè)不同結(jié)構(gòu)形式的攔污排，作為進(jìn)水口防護(hù)系統(tǒng)的第一道保障。對(duì)不同地區(qū)正在使用的攔污排研究分析后，發(fā)現(xiàn)兩端錨固、多個(gè)浮箱相互鉸接這種結(jié)構(gòu)攔污性能較好，但仍有污物繞過浮箱流過，并且部分浮箱受到嚴(yán)重?fù)p壞。因此，有必要對(duì)不同布置的攔污浮箱進(jìn)行分析研究，以便于采取更好的優(yōu)化措施。

在實(shí)際水利工程中，固體和流體通常是相互作用的，為了得到更準(zhǔn)確的研究結(jié)果，常常需要用到流固耦合分析。Sunghan 等采用大渦模擬對(duì)Re5500 ～41300 區(qū)間的圓柱近尾流場(chǎng)進(jìn)行了研究分析，發(fā)現(xiàn)隨著Re 數(shù)值的不斷增加，在亞臨界流動(dòng)狀態(tài)下的時(shí)均速度脈動(dòng)峰值更加靠近圓柱背壓點(diǎn)。Kwang-Jun Paik 等基于CFD 理論，分別采用單向和雙向流固耦合方法對(duì)船體剛性和彈性模型進(jìn)行動(dòng)力學(xué)分析，最后，則通過相關(guān)實(shí)驗(yàn)來(lái)驗(yàn)證計(jì)算結(jié)果。朱鳳霞等基于CFD 數(shù)值模擬理論和k?ε湍流模型，分析了攔污柵的不同截面形式對(duì)水流特性的影響，結(jié)果表明，圓形截面的攔污柵水流過渡最為平緩，水頭損失最小。

綜上所述，國(guó)內(nèi)外學(xué)者針對(duì)物體的流固耦合特性開展了諸多研究，但是，關(guān)于水電站不同布置下攔污浮箱的研究還是處于空白區(qū)域。本文首次基于ANSYS FLUENT 多場(chǎng)耦合分析平臺(tái)，對(duì)不同流向角下的水體-浮箱耦合模型進(jìn)行求解，分析不同工況下各攔污浮箱的流固耦合特性變化規(guī)律，為水電站改善攔污效果以及浮箱的結(jié)構(gòu)設(shè)計(jì)優(yōu)化提供參考。

1 流固耦合數(shù)學(xué)模型

1.1 控制方程

水體-浮箱耦合問題屬于流固耦合問題，在機(jī)組正常發(fā)電過程中，水體和浮箱之間有明確的耦合界面，浮箱位移比較有限。流固耦合問題遵循流體力學(xué)的三大定律，對(duì)于不可壓縮的牛頓流體，其基本方程包括流體控制方程、結(jié)構(gòu)控制方程和流固耦合方程。

流體控制方程。表達(dá)式為

式中：ρ為流體的密度；uj為j方向流體的速度；Uj為j方向流體的網(wǎng)格速度；p為壓力；μ為流體的動(dòng)力黏性系數(shù)；fi為i方向流體體積力。

（2）結(jié)構(gòu)控制方程。結(jié)構(gòu)整體的運(yùn)動(dòng)是由單元的運(yùn)動(dòng)方程式按照一定的方式疊加得到，表達(dá)式為

式中：[M]為結(jié)構(gòu)質(zhì)量矩陣；[C]為結(jié)構(gòu)阻尼矩陣；[K]為結(jié)構(gòu)剛度矩陣；{ }αR為流體附加動(dòng)力引起的節(jié)點(diǎn)載荷；{ }βR為結(jié)構(gòu)的其余荷載；{ }δ為結(jié)構(gòu)總體節(jié)點(diǎn)位移列陣。

（3）流固耦合方程。在流固耦合交界位置，應(yīng)該滿足流體與固體應(yīng)力、位移的相等或者守恒，表達(dá)式為：g表示來(lái)自流體的變量，fs 表示流體結(jié)構(gòu)的相互作用，e表示來(lái)自結(jié)構(gòu)的變量。

1.2 流固耦合數(shù)學(xué)模型

攔污排在進(jìn)水口水域運(yùn)動(dòng)比較緩慢，考慮到攔污排使用性能以及自身特點(diǎn)，可以忽略整體結(jié)構(gòu)的慣性效應(yīng)。在穩(wěn)定狀態(tài)下，浮箱與水體流向各不相同，為了便于研究，把坐標(biāo)系固定在所研究的浮箱上。如圖1 所示，以浮箱底面中心為坐標(biāo)原點(diǎn)，X 軸指向左側(cè)銷孔，Y 軸指向浮箱向水面，Z 軸垂直于水平面向上，ωV和ωθ分別表示水體隨浮箱坐標(biāo)系的流速和流向角。

由于水電站進(jìn)水口水流特性比較復(fù)雜，為了簡(jiǎn)化分析，在浮箱定位仿真中，假設(shè)水流的方向和速度是恒定的，因此，可得作用在浮箱上的流力計(jì)算公式為：

式中，Xω是浮箱在橫向所受到的流力，ωY是浮箱在縱向所受到的流力，Mω是浮箱在垂搖方向所受到的流力矩，Cxω、Cyω、Cmω分別為相應(yīng)的流載荷無(wú)因次系數(shù)，L為浮箱的入水長(zhǎng)度，εh為浮箱的入水深度，ρ為水密度。

圖1 浮箱坐標(biāo)系

2 數(shù)值模型建立

2.1 工程實(shí)例

以某水電站進(jìn)水口胸墻式泄洪閘水體-攔污排為研究對(duì)象。在水電站正常運(yùn)行過程中，攔污排的張力和形狀是一定的，但各個(gè)浮箱與水體之間的流向角是不同的。由于攔污排各浮箱距離較大，忽略相互之間對(duì)流場(chǎng)的影響。應(yīng)用三維建模軟件CATIA 建立不同的水體-浮箱耦合有限元模型，并將該模型導(dǎo)入數(shù)值分析軟件ANSYS，研究不同工況下各布置浮箱的流固耦合特性。

該水電站日均入庫(kù)流量753m3/s，日均出庫(kù)流量758m3/s，水頭差8.16m，日均發(fā)電量1241600kWh。浮式攔污排平面布置如圖2 所示，在泄流口前方水域，兩端分別通過安全拉桿錨固在上支墩和下支墩，升降設(shè)備均為隨水位變幅升降的全自動(dòng)控制卷?yè)P(yáng)式升降設(shè)備。按照懸鏈線理論計(jì)算最大弧垂25m，整體設(shè)計(jì)軸線長(zhǎng)度為237.8m，其主要結(jié)構(gòu)是32 個(gè)6.7m 長(zhǎng)的浮箱，外面完全密封，浮箱之間以若干個(gè)圓柱銷和長(zhǎng)1.2m 連桿裝配而成，在壩上進(jìn)水口水域形成整體的柔性結(jié)構(gòu)。

2.2 有限元模型及邊界條件

攔污排各浮箱在水庫(kù)正常運(yùn)行工作過程中，河流、空氣和浮箱是一個(gè)典型的氣-液-固多相流問題。由于模型體積較大，劃分網(wǎng)格、計(jì)算收斂等仿真過程不僅對(duì)電腦的硬件有很高的要求，還需要花費(fèi)很多時(shí)間和精力。為了節(jié)約時(shí)間成本，突出水載荷對(duì)各浮箱影響分析的重點(diǎn)，采用疊模分析將系統(tǒng)簡(jiǎn)化為流固耦合模型，僅考慮浮箱與水下計(jì)算域2 部分。

取水電站的進(jìn)水口為計(jì)算流體，賦予其自由液面以下的相應(yīng)水壓力，整個(gè)流體區(qū)域長(zhǎng)60m，寬30m，高6m，不考慮溫度和能量損失，只需給出水的密度和粘度；浮箱入水深度0.59m，與圓柱銷1 和銷2 進(jìn)行裝配，所用材料均為高強(qiáng)度低合金結(jié)構(gòu)鋼，彈性模量泊松比μ=0.3。經(jīng)過不斷試算，根據(jù)重力相似原則對(duì)模型進(jìn)行縮放，縮放比例尺為L(zhǎng)λ=10。采用非結(jié)構(gòu)四面體網(wǎng)格有限元模型進(jìn)行網(wǎng)格劃分，遠(yuǎn)離流固耦合交界面區(qū)域的網(wǎng)格較為稀疏，在交界面區(qū)域的網(wǎng)格相對(duì)密集。

邊界條件：FLUENT 計(jì)算域設(shè)置為高雷諾數(shù)流場(chǎng)，左右壁面使用對(duì)稱邊界條件，對(duì)于該邊界的垂直方向，流動(dòng)變量梯度的值為零。下表面使用無(wú)滑移壁面，上表面使用固壁邊界并設(shè)置相應(yīng)的滑移條件，確保該面有水流方向的速度。流場(chǎng)入口采用速度進(jìn)口的邊界條件，出口采用壓力邊界出口，相對(duì)壓力為0Pa。湍流模型采用模型，該模型對(duì)自由液面的捕捉精度和計(jì)算可信度較為理想，相比標(biāo)準(zhǔn)k?ε模型，湍流粘度增加了旋轉(zhuǎn)和曲率的內(nèi)容，使流動(dòng)更好地符合湍流的物理定律。

圖2 水電站攔污排平面布置

3 水體-攔污浮箱耦合流場(chǎng)特性分析

3.1 浮箱近壁水體流速特性分析

圖3 所示為90°流向角下入口2m/s 對(duì)應(yīng)的三維流場(chǎng)和對(duì)稱面流速矢量。從圖3 可以看出該工況下，流場(chǎng)中浮箱周圍流體的速度關(guān)于中心線對(duì)稱，且流場(chǎng)分布變化比較劇烈。流體從入水口流向浮箱，在向水面駐點(diǎn)壓力達(dá)到最大值，速度迅速減小為零，動(dòng)能全部轉(zhuǎn)化為壓能。之后，流體繞過浮箱向下游流動(dòng)時(shí)，部分壓能轉(zhuǎn)化為動(dòng)能，最大流速出現(xiàn)在浮箱底部位置。當(dāng)流體流經(jīng)浮箱背水面時(shí)，其中一部分向下游流去，而由于浮箱周圍高負(fù)壓區(qū)域的存在，使得此處的另一部分流體壓力急劇減小，產(chǎn)生回流并在背部區(qū)域形成漩渦。

圖3 90°流向角下入口2m/s流速矢量（單位： m/s ）

通過實(shí)際調(diào)研發(fā)現(xiàn)，污物主要堆積在90°流向角浮箱附近。為進(jìn)一步了解90°流向角下浮箱繞流流場(chǎng)的速度大小變化情況，在浮箱周邊選取4 條近壁水體節(jié)線，生成節(jié)線速度位移變化曲線，浮箱近壁節(jié)線選取情況如圖4 所示。

圖4 浮箱近壁節(jié)線選取

圖5 為不同入口流速下4 條節(jié)線速度位移變化曲線。從圖中可以得出：①針對(duì)不同的入口流速，流體的速度在近浮箱區(qū)域都發(fā)生了突變，變化趨勢(shì)大致相同。②在靠近浮箱流域節(jié)線1、3 速度“陡降急增”，反映了流體繞流造成的浮箱邊緣渦流空化情況。而節(jié)線2、4 速度“陡增急降”，并且陡增的速度較大，甚至大于入口速度。由此說明，漂浮的污物很容易繞過浮箱底部堆積在上游攔污柵，水頭減少、威脅增大、攔污排形同虛設(shè)。因此，有必要為90°流向角浮箱增設(shè)一排入水欄桿，入水深度取決于流速過渡區(qū)域的大小。③對(duì)比不同入口流速下浮箱底部的流場(chǎng)特性，從1 ～3m/s，節(jié)線1 最大降幅依次為86%、84%、88%，而節(jié)線2 最大增幅依次為17%、18%、21%，說明流速過渡區(qū)域與入口流速大小無(wú)關(guān)。

圖6 為入口流速2m/s 近壁流體中軸線（沿Z 軸方向）速度位移分布曲線。結(jié)合模型的縮放比例尺，可以看出該類型結(jié)構(gòu)攔污浮箱實(shí)際增設(shè)入水欄桿，深度h 以1.4m 最優(yōu)。

3.2 浮箱近壁水體壓力特性分析

圖7 為入口流速2m/s 不同流向角下近壁水體壓力分布云圖（X0Y 平面，Z=0.05m）。由圖中可以得出：壓力最大點(diǎn)位于浮箱迎流面?zhèn)?，在該位置浮箱容易受到損壞，流體到達(dá)駐點(diǎn)時(shí)，動(dòng)能逐漸轉(zhuǎn)化為壓能。流向角越小，駐點(diǎn)位置越靠近浮箱迎流面?zhèn)鹊倪吘?，此時(shí)，會(huì)在該位置形成一定范圍的負(fù)壓區(qū)域；流向角越大，浮箱的迎流面積越大，近壁水體高壓范圍也顯著變廣。由于此處是浮箱水體的交界面，當(dāng)流體繞過浮箱時(shí)，流體仍會(huì)在交界面邊緣處發(fā)生擾動(dòng)現(xiàn)象，在浮箱背水面的負(fù)壓區(qū)域出現(xiàn)小范圍的正壓區(qū)。由于自由液面的波動(dòng)和浮箱交界面的相互影響，在浮箱近背水面的尾流場(chǎng)中，可以看到形狀不規(guī)則的渦結(jié)構(gòu)。

圖6 入口2m/s 近壁流體中軸線速度位移變化

圖7 入口流速2m/s 不同流向角下近壁水體壓力分布

將該工況下各流向角的最大、最小壓力值匯集起來(lái)，如表1 所示?？梢钥闯觯?0°流向角下的正壓力最大，而75°流向角下的負(fù)壓力和壓力差最大，說明該流向角下的浮箱在長(zhǎng)期運(yùn)行的過程中極易受到損壞，在攔污浮箱的設(shè)計(jì)優(yōu)化過程要重點(diǎn)考慮75°流向角下的水壓力影響，以避免不必要的損失。

表1 不同流向角下最值壓力

4 攔污浮箱結(jié)構(gòu)振動(dòng)分析

在FLUENT 模塊入口流速2m/s 的求解基礎(chǔ)上，文中在ANSYS-Workbench 協(xié)同仿真平臺(tái)進(jìn)行了水體-攔污浮箱耦合數(shù)據(jù)載荷的傳遞，進(jìn)一步對(duì)不同浮箱結(jié)構(gòu)振動(dòng)固有頻率進(jìn)行分析。

4.1 浮箱結(jié)構(gòu)分析

讀取結(jié)構(gòu)網(wǎng)格以及FLUENT 浮箱表面壓力數(shù)據(jù)，加載到水體-浮箱耦合位置，60°流向角浮箱耦合位置加載壓力載荷如圖8 所示。結(jié)合實(shí)際情況，設(shè)置求解信息并對(duì)兩端銷軸施加固定約束，然后，分別對(duì)不同流向角下的浮箱進(jìn)行求解，提取各自的最大變形量以及最大應(yīng)力值，不同流向角下最大變形量和最大應(yīng)力值曲線如圖9 所示。隨著流向角的增大，最大變形量和應(yīng)力值也逐漸增加，然后，再逐漸減小，在75°流向角下達(dá)到最大值，再次驗(yàn)證了前面得出的結(jié)論。

圖8 90°流向角浮箱耦合位置加載壓力載荷

圖9 不同流向角下浮箱最大變形量及應(yīng)力值

4.2 浮箱預(yù)應(yīng)力模態(tài)分析

將不同流向角下的壓力載荷加載到浮箱的耦合位置作為預(yù)應(yīng)力，施加邊界條件，提取前六階模態(tài)分析結(jié)果。60°流向角下浮箱前三階模態(tài)振型如圖10 所示，振型位移主要集中在浮箱的中部位置，該位置相對(duì)更容易遭受振動(dòng)的沖擊損壞。然后，與無(wú)載荷加載情況下浮箱自由振動(dòng)情況的結(jié)果作對(duì)比，匯集成表如表2 所示。從表中可以對(duì)于不同流向角下的浮箱，在受到外部水載荷作用下的模態(tài)頻率與自由模態(tài)相差無(wú)幾，同時(shí)，可以看到，浮箱在受到外部約束的情況下，其自由振動(dòng)模態(tài)頻率前6 階數(shù)值并不為0，說明剛性模態(tài)被抑制。

5 FLUENT 數(shù)值模擬驗(yàn)證

為驗(yàn)證數(shù)值模擬的準(zhǔn)確性，首先，通過經(jīng)驗(yàn)公式計(jì)算不同傅汝德數(shù)下的浮箱模型，得出不同流速時(shí)自由液面情況下的浮箱黏性阻力系數(shù)，然后，與應(yīng)用FLUENT 計(jì)算得到的黏性阻力系數(shù)進(jìn)行比較。

經(jīng)驗(yàn)公式計(jì)算中，黏性阻力系數(shù)：

摩擦阻力系數(shù)(無(wú)量綱）采用二因次傅汝德法來(lái)進(jìn)行求?。?/p>

雷諾數(shù)：

黏壓阻力系數(shù)采用巴普米爾計(jì)算公式：

圖10 60°流向角下浮箱前三階模態(tài)振型

表2 不同流向角及無(wú)載荷加載下各階模態(tài)頻率/Hz

其中，Dh為水力直徑，mA為浮箱中剖面面積，rL為去流段長(zhǎng)度，Rpv為黏壓阻力，S為浮箱的濕表面積。

圖11為兩種方法計(jì)算所得的黏性阻力系數(shù)曲線對(duì)比圖。從圖中可以看出，隨著傅汝德數(shù)的增大，黏性阻力系數(shù)逐漸減小，而FLUENT 計(jì)算得到的黏性阻力系數(shù)略大于經(jīng)驗(yàn)公式摩擦阻力系數(shù)和黏壓阻力系數(shù)之和，這和文獻(xiàn)[8]所得的結(jié)論相同，因此，F(xiàn)LUENT 數(shù)值模擬結(jié)果是合理可信的。而誤差產(chǎn)生的原因主要是計(jì)算黏壓阻力系數(shù)時(shí)只考慮了去流段的影響，沒有考慮隨著雷諾數(shù)的增大對(duì)黏壓阻力的影響。

圖11 黏性阻力系數(shù)曲線

6 結(jié)語(yǔ)

通過對(duì)不同流向角下的攔污浮箱流固耦合特性進(jìn)行計(jì)算分析以及數(shù)值驗(yàn)證，得到以下結(jié)論：（1）對(duì)于實(shí)際攔污效果不好的90°流向角下攔污浮箱，底部流速過渡區(qū)域與入口速度無(wú)關(guān)，應(yīng)當(dāng)增設(shè)入水深度1.4m 的入水欄桿以攔截污物于第一道防線；（2）隨著流向角的增大，浮箱最大變形量及最大應(yīng)力值先增大后減小,而75°流向角下的壓力差值最大，要重點(diǎn)優(yōu)化該布置下的浮箱結(jié)構(gòu)，以避免在水載荷的作用下造成嚴(yán)重?fù)p壞；（3）不同流向角下的浮箱預(yù)應(yīng)力模態(tài)頻率與無(wú)載荷加載時(shí)近乎一致，不同布置下的浮箱動(dòng)力學(xué)特性只與自身結(jié)構(gòu)特性有關(guān)。