厚覆冰導線升力突變機理及數值模擬研究

樓文娟，潘 晨，孫建平

(浙江大學 結構工程研究所，杭州 310058)

0 引 言

架空輸電線路導線在表面覆冰后會形成非圓形的截面形狀，其風致氣動力可誘發導線舞動，嚴重威脅輸電線路的安全[1-2]。

李萬平[3-4]、樓文娟[5]、閻東[6]、李新民[7]等對新月形、D形覆冰導線進行風洞試驗研究全風攻角下氣動力系數。Hisato[8]等對不同厚度的三角覆冰和圓角三角覆冰導線進行表面測壓試驗和天平測力試驗，驗證兩種試驗結果的一致性，并獲得了覆冰導線氣動力數據。

基于計算流體力學（CFD）的數值模擬方法能夠展現完備的流場信息，可與風洞試驗互補，便于更深層次地研究流場結構的演化。李新民[9]用SA湍流模型對23 mm厚新月形覆冰導線進行數值模擬，發現在風攻角15°附近存在回流渦現象。Ishihara[10]利用大渦模擬（LES）對圓角三角形覆冰導線進行模擬，模擬結果與風洞試驗數據相吻合。林巍[11]、呂翼[12]分別用k-ωSST湍流模型和SA湍流模型對覆冰6.7 mm（導線直徑26.82 mm）和覆冰18 mm （導線直徑32.76 mm）的新月形覆冰導線進行數值模擬，并得到了較好的模擬結果。符玉珊[13]運用格子玻爾茲曼方法對18 mm覆冰的新月形覆冰導線（導線直徑32.76 mm）進行大渦模擬，在較低耗時情況下得到了較高精度的模擬結果。

現有對覆冰導線氣動力特性的數值模擬研究中，覆冰厚度多數小于1D（D為導線直徑），而根據歷史觀測數據，覆冰厚度達到1D及以上的情況在覆冰嚴重地區經常出現[14-15]。圖1根據文獻[5-6]的風洞試驗數據繪制，為不同厚度覆冰導線升力系數隨風攻角變化情況的示意圖。當覆冰厚度較大時，覆冰導線的升力系數曲線在風攻角15°附近區間出現“尖峰”，并且隨著冰厚的增大更加明顯。升力系數的突變使得在該風攻角區間易出現覆冰導線的風致舞動現象[6,16-17]，準確模擬和理解這一升力系數的“尖峰”，對于覆冰導線風致舞動的驗算和防治具有特別意義。

圖1 不同厚度新月形覆冰升力系數（I為湍流度）Fig. 1 Lift coefficients of the crescent-shape ice with different thickness (I indicates the turbulence intensity)

目前覆冰導線氣動力系數的CFD模擬大多局限于無升力“尖峰”的較薄覆冰，通常選用雷諾時均法進行氣動力模擬。但對于厚覆冰導線的模擬結果并不理想，尤其是風攻角15°附近升力系數“尖峰”模擬較為困難；另外，對于厚覆冰導線升力系數在風攻角15°附近出現極大值的原因也尚不明確，鮮有討論。本文分別采用ANSYS Fluent軟件內基于k-ωSST湍流模型的雷諾時均法和大渦模擬的數值方法，通過網格和時間步長的優化提高模擬精度，得到風攻角10°～20°范圍內1D新月形厚覆冰導線氣動力參數和繞流風場結果，并利用風洞試驗結果驗證數值方法的適用性。通過對覆冰導線表面壓力分布和繞流風場結構分析，進一步探究厚覆冰導線升力系數在風攻角15°附近出現“尖峰”的原因。

1 數值方法和計算模型

1.1 大渦模擬

大渦模擬通過空間濾波將湍流分成大尺度渦和小尺度渦。其中，大尺度渦對湍流輸運的動量、質量、能量等物理量起到決定性作用，因此對大尺度的渦進行了不可壓縮Navier-Stokes方程直接模擬。小尺度的渦趨于各向同性，其對大尺度渦的影響通過亞格子模型模擬。

本文的大渦模擬選用局部渦黏度壁面自適應（WALE）亞格子模型。WALE亞格子模型是在Smagorinsky-Lilly亞格子模型的基礎上發展過來的，它修正了Smagorinsky-Lilly亞格子模型在近壁面渦黏系數不為0的缺陷[18]。其亞格子應力張量具體形式如下：

其中，μt為亞格子黏度系數，其定義式如下：

其中，

式中，Cw為WALE常數，取0.325；Δ為過濾之后的尺寸，為網格尺寸。

1.2 基于k-ω SST模型的雷諾時均法

雷諾時均法（RANS）將流場分為系綜平均量和湍流脈動量，并通過湍流模型來封閉N-S方程。其中，基于湍流動能k和比耗散率ω的k-ω模型[19]是一種較為常用的湍流模型。k-ω模型在湍流黏度的定義中考慮了湍流剪切應力的傳遞，并通過設定混合函數將標準k-ε模型和k-ω模型進行結合，在近壁面區域內使用k-ω湍流模型求解，在邊界層外則使用k-ε模型。既能充分發揮k-ε模型在模擬遠離壁面的充分發展湍流的優勢，又能適應于各種壓力梯度下的邊界層模擬，在模擬壓力梯度和分離流上具有較大優勢，常被用于模擬鈍體繞流。

1.3 計算模型

本文研究對象為1D新月形覆冰導線，裸導線直徑D為26.82 mm。計算域入口邊界至導線中心距離為30D，出口邊界距導線中心50D；計算域的橫風向寬度為60D，展向高度為2D，最大阻塞率為2.0%。計算域入口設為速度入口（入口風速U= 10 m/s），出口為自由出流，覆冰導線表面為光滑無滑移壁面，展向邊界為周期性邊界，其余為對稱邊界。幾何模型和計算域如圖2所示。計算網格為使用Pointwise生成的結構網格，近壁面第一層網格高度設為0.001D，以滿足LES計算時y+≈1要求，網格增長率為1.05。圖3為計算網格局部示意圖。

圖2 幾何模型和計算域Fig. 2 Schematic of the model and the computational domain

圖3 近壁面網格Fig. 3 Near-wall grid

大渦模擬采用SIMPLEC算法求解壓力-速度耦合方程，壓強項為二階格式，動量項為有界中心差分格式，時間離散方式為二階隱式格式。在使用基于k-ωSST模型的雷諾平均法在求解黏性子層時同樣需要滿足y+≈1的網格要求，選用與LES相同的計算網格。k-ωSST模型計算時湍動能和耗散率項均設置為二階迎風格式，其余和LES設置相同。

1.4 網格與時間步長無關性驗證

為分析網格數和時間步長對大渦模擬計算結果的影響，分別計算了風攻角15°下的4個算例，結果列于表1。

表1 不同計算工況設置Table 1 Parameters for different computation conditions

由表1可得，周向節點數對氣動力模擬結果有較大影響。隨覆冰厚度增大，需要更多節點來捕捉覆冰導線壁面弧度細節。針對本文工況，600個周向節點能較好表現覆冰導線截面形狀特性。綜合考慮計算精度與耗時，本文LES選用Case2的網格與時間步長設置；k-ωSST湍流模型選取同樣的網格，但由于本文算例中k-ωSST湍流模型對時間步設置不敏感，因此參考文獻[10]中的取值，無量綱時間步設置為0.186。

2 結果對比及分析

2.1 氣動力系數

覆冰導線升力系數、阻力系數以及扭矩系數定義為：

式中，空氣密度ρ為1.225 kg/m3，導線外徑D為26.82 mm，展向長度H為2D，U為試驗風速。升力FL、阻力FD、扭矩M的方向如圖4所示。

圖4 氣動力方向及風攻角定義Fig. 4 Definition of the aerodynamic force directions and wind attack angles

圖5為風攻角10°～20°范圍內1D新月形覆冰導線的氣動力系數。對比k-ωSST和LES的計算氣動力系數可以發現，LES模擬結果與試驗結果非常接近，在風攻角15°下，升力系數、扭轉系數均出現“尖峰”，并且對角度一階求導后得到的Den Hartog系數和Nigol系數（圖6）在負斜率部分也和試驗值非常接近；而k-ωSST模擬結果在風攻角15°～20°范圍內升力系數沒有降低，沒有出現明顯的“尖峰”。因此，LES模型比k-ωSST模型更能準確模擬覆冰導線的氣動力，能滿足覆冰導線舞動分析的需求。

圖5 各風攻角下1D新月形覆冰導線氣動力系數Fig. 5 Aerodynamic coefficients of a 1D crescent-shape iced conductor at various wind attack angles

圖6 1D新月形覆冰導線Den Hartog系數和Nigol系數Fig. 6 Den Hartog and Nigol coefficients of a 1D crescent-shape iced conductor

2.2 平均流場特性

根據1D新月形覆冰導線表面的壓力分布特征，將其分為四個區域：覆冰端上側（0 ＜x/D＜ 1.5）、導線端上側（1.5 ＜x/D＜ 2）、導線端下側（2 ＜x/D＜ 2.5）和覆冰端下側（2.5 ＜x/D＜ 4），并分別將其記為Ⅰ區、Ⅱ區、Ⅲ區、Ⅳ區，如圖7所示。

圖7 覆冰導線分區及周向坐標定義Fig. 7 Partitions of the iced conductor with the circumferential coordinates

圖8為兩種數值方法計算得到的平均流線圖。整體而言，兩種算法的平均流場在Ⅱ區和Ⅲ區均有較大尺度的回流，Ⅱ區逆時針大尺度渦緊貼著導線壁面，Ⅲ區渦為上側流體繞過Ⅱ區渦后再次附著在導線表面的二次渦，且Ⅱ區回流區域面積大于Ⅲ區回流區域面積。

k-ωSST模型計算的平均回流區長度在風攻角10°與12.5°較為接近，在15°附近最大，并在15°至20°范圍內隨風攻角增大而減小。LES計算得到的平均回流區長度則在10°～15°下沒有明顯變化，在15°～20°下隨風攻角增大而明顯增大。

相比較k-ωSST模型而言，LES能夠模擬出Ⅰ、Ⅱ區緊貼壁面逆時針的分離泡。隨著風攻角的增大，這些分離泡的位置逐漸前移靠近導線覆冰端頂點，覆冰導線上側分離點也向覆冰頂端靠近，且在Ⅲ區存在小尺度渦結構。Ⅲ區壁面附近渦結構僅在15°風攻角下為一對方向相反的小尺度渦，其余風攻角下僅存在一個順時針渦。

2.3 表面風壓系數

從圖5可看到，k-ωSST模型和LES模擬的升力系數在15°～20°有明顯區別，圖9為該范圍覆冰導線壁面平均風壓系數分布圖。由圖9可得，1D新月形覆冰導線表面平均風壓系數在Ⅱ、Ⅲ區均保持在0附近，靠近導線覆冰端頂點的平均壓力系數存在大幅度突變。平均壓力系數均在覆冰頂端的上表面和下表面附近（約x/D= 0.1與x/D= 3.9處）達到最小值和最大值。k-ωSST模型和LES計算結果差異主要體現在Ⅰ區。

圖8 不同風攻角下平均流線圖Fig. 8 Averaged streamline patterns under various wind attack angles

圖9 風攻角15° ～ 20°下的平均風壓系數Fig. 9 Averaged wind pressure coefficient at wind attack angles ranging from 15° to 20°

圖10為各風攻角的平均速度和壓力云圖?？梢钥吹?，LES結果中最大平均速度都在15 m/s左右，k-ωSST模型計算得到的最大平均速度和LES基本一致，僅在20°風攻角下突然增大至17 m/s。在各風攻角下均存在兩個高流速區域（x/D= 0.1和x/D= 2.5附近），此處流場對覆冰導線的吸力最大，對應壓力曲線上的兩個極小值。

在覆冰端高速區域內，k-ωSST模型的高速區域面積更大，呈半圓形；LES模型的高速區域呈長條狀，對應于圖9所示的LES模型在覆冰端上側的平均壓力分布明顯小于k-ωSST模型。LES方法捕捉到的小尺度渦脫引起覆冰導線上側出現風速介于最大速度和0之間的中速區，高速區向上抬升，減小了導線在此范圍內的吸力。

k-ωSST和LES模擬結果存在區別，其主要原因是LES能夠模擬出緊貼覆冰導線上側壁面的多個小尺度渦，而這些小尺度渦對覆冰導線壁面壓力有著不可忽視的影響；k-ωSST模型僅模擬出了Ⅰ區和Ⅱ區的大尺度回流渦，無法模擬出Ⅰ區壁面上的小尺度渦，夸大了Ⅰ區壁面的風吸力。并且隨著覆冰厚度的增大，覆冰端上側的小尺度渦結構會更復雜，對覆冰導線的影響更大。因此對于較厚覆冰導線，LES是數值模擬的首選。

圖10 各風攻角下的平均速度云圖和平均壓力系數Fig. 10 Averaged velocity contours and averaged wind pressure coefficients at various wind attack angles

3 升力系數突變原因

由圖8可知，風攻角10°至20°是覆冰導線尾流處的渦不斷擴大并附著于覆冰導線上壁面的過程；并隨著覆冰厚度增大，Ⅰ、Ⅱ區的小尺度渦結構更加復雜。在該過程中，覆冰導線周圍流場變化致使導線氣動力出現較大變化。由圖5可知，LES方法能更有效地模擬1D新月形覆冰導線氣動力系數隨風攻角的變化情況，其表面壓力分布結果更準確。下面根據LES方法模擬所得平均流場結果，分析升力系數在風攻角15°出現“尖峰”的流體動力學機制。

表2為各風攻角下壁面四個區域的升力系數分量，從中可以看到覆冰導線升力主要來源于Ⅳ區正壓以及Ⅰ區負壓。對Ⅳ區，在10°至20°范圍內隨著風攻角的增大，Ⅳ區前沿駐點的位置遠離覆冰頂點，駐點處的壁面曲率減小，來流與壁面夾角增大，根據動量定律，壁面壓力增大；同時，隨著風攻角的增大，Ⅳ區正壓沿升力方向的分量減小，沿阻力方向的分量增大。如表2所示，Ⅳ區正壓提供的升力在15°風攻角下達到最大值。對Ⅰ區，在10°～20°范圍，隨著風攻角增大，分離點向來流方向前移，并出現緊貼壁面的小尺度渦，其與尾流區域渦旋結構的相互作用致使Ⅰ區壁面吸力先減后增。同時，在風攻角15°下，覆冰導線Ⅲ區附近的一對方向相反的小尺度渦增大了沿壁面法向的壓力梯度，使Ⅲ區出現較大的正壓。綜合覆冰導線上側吸力略微減小、下側（Ⅲ、Ⅳ區）的壁面壓力因駐點所在處的壁面曲率、來流夾角和壁面切線方向等因素影響達到最大值，致使導線上下表面壓力差明顯大于其余風攻角，1D厚覆冰導線的升力系數曲線在15°風攻角附近出現“尖峰”。

表2 各風向角下各壁面分區的升力系數分量Table 2 Components of the lift coefficient in each surface partition under various wind attack angles

4 結 論

本文分別采用基于k-ωSST湍流模型的雷諾時均法和大渦模擬的數值方法模擬1D新月形覆冰導線在風攻角10°～20°范圍內的氣動力特性和繞流風場，并結合風洞試驗數據對結果進行對比分析。主要結論如下：

1）相比于覆冰導線氣動力模擬中廣泛應用的雷諾時均法，大渦模擬能夠準確捕捉壁面附近的小尺度渦結構，有較準確的模擬結果。對于較厚覆冰導線的氣動力數值模擬，應首選LES。

2）1D新月形覆冰導線在上側壁面處的渦結構影響整體流場，并在下側壁面曲率、來流夾角和壁面切線方向共同作用下導致風攻角15°時的升力系數突變。

根據本文的模擬結果以及出現升力“尖峰”的原因闡釋，得到了較厚新月形覆冰的流場特征，但具體的防舞措施還需進一步探究嘗試。