貼壁二維方柱繞流對壁面摩擦應力的影響

張之豪，傅奇星，王慶洋，徐勝金,*

1.清華大學 航天航空學院，北京 100084 2.中國航天空氣動力技術研究院，北京 100074 3.中國汽車工程研究院股份有限公司，重慶 401122

0 引 言

二維方柱繞流問題是典型的分離流問題，在工程中應用廣泛。自由平行來流繞零迎角方柱運動時，根據雷諾數的不同，可能在第一個迎流的棱角處產生分離，出現(xiàn)再附后還可能在尾流附近的棱角處再次分離。在全局不穩(wěn)定性和局部不穩(wěn)定性機制作用下，方柱上下分離的自由剪切層會繞曲形成交替出現(xiàn)的大尺度渦，并在方柱下游產生一定寬度的尾流區(qū)[1-2]。在方柱后緣布置分離板[3-4]，可以抑制上下自由剪切層的繞曲，在近尾流區(qū)不會形成交替出現(xiàn)的大尺度渦。

當二維方柱置于平板邊界層內時（如廣泛存在于建筑、橋梁、燃氣輪機葉片冷卻、海底管道輸運等工程應用中的貼壁方柱繞流），來流速度呈梯度分布，因平板的存在，流動只會在方柱上表面出現(xiàn)分離，形成含有渦量的大尺度流動結構[5-8]。平板的存在，相當于將無限長的分離板置于方柱下方，不僅使方柱下表面無法產生分離流，也會對上表面的分離流產生影響。Panigrahi 和Acharya[9-10]采用熱線與激光多普勒測速技術（LDV）研究了不同雷諾數下的方柱下游流動速度功率譜及特征頻率，發(fā)現(xiàn)大尺度流動結構的脫落頻率與雷諾數呈線性關系，并指出剪切層的Kelvin-Helmholtz 不穩(wěn)定性模式與后臺階流動類似。通過模式識別提取了渦旋運動的隨機部分與相干部分，發(fā)現(xiàn)剪切層外緣的流動由渦旋誘導的大尺度上拋運動所主導，這種上拋運動與流動速度的非高斯分布形式存在重要關聯(lián)。流動繞過方柱后產生2 個差異明顯的區(qū)域[11]：一為分離/回流區(qū)，一為恢復/重建區(qū)。分離/回流區(qū)是指從方柱前緣流動分離點至下游平板流動再附點之間的區(qū)域，該區(qū)域內存在狹長的回流泡，是流動結構生成、發(fā)展的主要區(qū)域，流動具有較強的間歇性；恢復/重建區(qū)位于流動再附點下游區(qū)域，該區(qū)域內流動結構逐漸耗散，流場逐漸恢復為湍流邊界層。進一步研究發(fā)現(xiàn)，在Kelvin-Helmholtz 不穩(wěn)定性作用下，方柱剪切層產生低頻振蕩，在方柱上表面形成大尺度流動結構，并類似渦脫一樣產生流動分離[12]。Shi 等[13]采用時間分辨粒子圖像測速技術（TR-PIV）對方柱下游流場進行了研究，分析了大尺度流動結構對流場非定常特性的影響，發(fā)現(xiàn)大尺度再附渦與低頻振蕩分離泡之間存在相互作用，導致再附點附近流動間歇性因子發(fā)生急劇變化。

本文就貼壁二維方柱繞流對下游壁面摩擦應力的影響機制開展實驗研究，以期為深入理解表面沖蝕、污染物聚集、近壁面湍流耗散等問題的機理提供參考。

1 實驗方案

實驗在低速直流風洞（實驗段長、寬、高尺寸為2 m × 0.5 m × 0.5 m）中進行。如圖1 所示，在距風洞底壁0.16 m 處水平放置光滑平板。平板前緣為楔形，在距前緣20 mm 處放置直徑5 mm 的絆線，以促進邊界層轉捩，獲得充分發(fā)展的湍流邊界層。平板后緣安裝角度可調節(jié)的尾板，實驗時可通過調節(jié)尾板角度使平板邊界層沿流向滿足零壓梯度條件。貼壁二維方柱置于距平板前緣1 m 處，方柱寬度D=10 mm，展長0.5 m，與風洞實驗段截面寬度一致。來流風速U0=15 m/s，基于方柱寬度D 和來流風速U0定義的雷諾數ReD=1.1 × 104。經熱線測量，來流湍流度為0.3%。

圖1 平板及貼壁二維方柱示意圖Fig.1 Diagram of flat plate and wall mounted 2D square cylinder

無方柱情況下，經邊界層熱線測量，以方柱所在位置邊界層動量厚度θ（根據流向平均速度剖面積分得到）定義的雷諾數為Reθ=U0θ/ν=5 270（空氣的運動黏性系數ν=1.51 × 10-5m2/s）。邊界層熱線單點采樣頻率20 kHz，采樣時間10 s。

圖2 平板湍流邊界層流向平均速度剖面Fig.2 Profile of time-averaged streamwise velocity of flat plate TBL

1.1 流場測量

使用2D2C TR-PIV 測速系統(tǒng)對貼壁二維方柱下游流場進行測量。測速系統(tǒng)包括高速相機（分辨率1 280 像素×800 像素，配備最大光圈3.5、180 mm微距鏡頭）、Nd-YAG 激光器（波長532 nm，最大能量40 mJ，頻率10～10 000 Hz）、BNC 575 同步控制器。選用經Laskin 噴嘴霧化的癸二酸二辛脂（DEHS）煙霧顆粒（直徑2 μm）作為示蹤粒子，在直流風洞入口處釋放，經充分混合后進入風洞實驗段。圖3 為流場測量示意圖。

圖3 貼壁二維方柱下游流場PIV 測量示意圖Fig.3 Diagram of PIV measurement of flow field downstream of the 2D wall-mounted square cylinder

如圖3 所示，為在同等像素分辨率下獲得更為精細的流場結構，縮小拍攝視場，沿流向（x 向）分2 個區(qū)域進行測量。區(qū)域1 和2 的法向（y 向）視場范圍均為0 < y/D < 4，流向視場范圍分別為0 < x/D <8、6 < x/D < 14。2 個區(qū)域有寬為2D 的重疊區(qū)域，以便在進行數據統(tǒng)計平均時拼接測量區(qū)域。PIV 使用雙幀雙曝光工作模式，采樣頻率為1 kHz，2 個區(qū)域均采集8 000 對圖像，2 幀圖像時間間隔為50 μs。

1.2 壁面動態(tài)摩擦應力測量

熱線能夠獲得精準的動態(tài)速度信息，但常規(guī)單絲熱線在近壁面應用時會面臨壁面干擾、無法靠近壁面、無法判斷動態(tài)速度方向等困難。本文研發(fā)了平行雙絲熱線傳感器，結合基于機器視覺的高精度定位系統(tǒng)，實現(xiàn)壁面動態(tài)摩擦應力測量。

摩擦應力測量原理如圖4 所示。熱線探頭尖端裝有2 根平行熱絲（直徑5 μm、法向間距0.025 mm、流向間距0.6 mm）。測量時，將平行雙絲熱線（熱線1、熱線2）伸入至黏性底層，熱線探桿垂直于壁面，2 根熱線與壁面平行、與來流方向垂直。根據2 根熱線測點的流速U1(t)、U2(t)及熱線與壁面的距離h1、h2，由牛頓內摩擦定律可得到2 根熱線所對應的壁面動態(tài)摩擦應力τ1(t)、τ2(t)。

圖4 摩擦應力測量原理示意圖Fig.4 Principle of wall shear stress measurement

熱線1 和2 的流向位置不同，上游熱線和連接熱線的支架對流動有阻礙作用，使得下游熱線處的流速略低于上游熱線。當近壁區(qū)域流動方向改變時，2 根熱線測點處的流速差異也會發(fā)生改變。流速變化體現(xiàn)為熱線輸出電壓的變化，因此，可根據2 根熱線輸出電壓差的大小，判斷壁面動態(tài)摩擦應力的方向。圖5 給出了動態(tài)摩擦應力方向判別的具體方法。設定熱線1 位于上游，摩擦應力與主流來流方向相同時為正，即τ(t) > 0；與主流來流方向相反時（回流方向）為負，即τ(t) < 0。測量時，同步獲得2 根熱線的電壓信號E1(t)和E2(t)。在判斷ti時刻的摩擦應力方向時，為增強魯棒性，計算ti時刻附近時間段[ti- ts,ti+ ts]內電壓差的平均值=與判斷閾值ΔEth進行比較：若>?Eth，則τ(ti) > 0；若 0，τ(ti)=τ1(ti)；若τ(ti) < 0，則τ(ti)=-τ2(ti)。

圖5 動態(tài)摩擦應力方向判斷示意圖Fig.5 Diagram of dynamic shear stress direction judgment

根據牛頓內摩擦定律，壁面摩擦應力與熱線測點處流速存在對應關系，可直接標定熱線電壓與摩擦應力之間的函數關系。標定時的熱線位置與測量位置相同，盡量確保壁面影響一致，將壁面影響包含在標定曲線中。利用基于圖像識別的精密定位系統(tǒng)，將熱線定位至湍流邊界層黏性底層內。本文平板鏡像效果良好，可根據熱線圖像及熱線在平板上的鏡像確定熱線與壁面之間的距離，定位精度約4.2 μm。本文將2 根熱線定位至距壁面h1=0.1 mm和h2=0.075 mm 處，在無方柱情況下Reθ=5 270的湍流邊界層中分別對應y1+=4.5、y2+=3.3，確保2 根熱線處于黏性底層內。

為減小熱線和支架的流動干擾導致的標定誤差，分別將熱線1 和2 置于上游進行標定。標定時，將平行雙絲熱線傳感器置于距平板（如圖1 所示，但無絆線）前緣水平距離55 mm 處，利用平板層流邊界層已知的壁面摩擦應力對平行雙絲熱線進行標定。層流邊界層的壁面摩擦應力根據Blasius 解計算得到，標定時的來流風速U0=0 m/s 及U0=3～18 m/s（間隔1 m/s），對應平均摩擦應力標定范圍為0～0.51 Pa。標定結果如圖6 所示。

圖6 平行雙絲熱線標定結果Fig.6 Calibration curves of parallel double hot wire

基于四次多項式擬合建立壁面摩擦應力與熱線電壓之間的關系，根據Coles-Fernholz 經驗模型數據（如圖7 所示），可估計出Reθ< 5 270 范圍內的湍流邊界層平均摩擦應力小于0.4 Pa，說明標定數據滿足湍流邊界層的測量需求。完成標定后，選取=0 標定點對應的2 根熱線電壓差均值作為摩擦應力方向判斷閾值，其物理意義在于：當摩擦應力由正向變?yōu)樨撓蚧蛴韶撓蜃優(yōu)檎驎r，必會經過τ=0，因此τ=0 可視為摩擦應力方向改變的臨界點。

圖7 不同雷諾數下平板湍流邊界層平均摩擦應力測量結果Fig.7 Results of time-averaged shear stress in flat plate TBL at different Reynolds numbers

為驗證平行雙絲熱線測量摩擦應力的可靠性，對距平板（有絆線）前緣1 m 處的零壓梯度平板湍流邊界層壁面摩擦應力進行測量。測量時，將熱線1 和2 分別置于上游位置測得摩擦應力。測量結果如圖7 所示（無量綱平均摩擦應力=/0.5ρU02）。與Coles-Fernholz 經驗模型[14]相比，熱線1 和2 測得的壁面摩擦應力的平均相對誤差分別為3.6%和3.3%，平均相對不確定度分別為2.7%和3.0%（95%置信度），驗證了摩擦應力測量的準確性。

湍流會增大流體熱傳導系數[15]，因此，利用考慮壁面效應影響的層流邊界層標定的熱線測量湍流邊界層時，會略低估流體帶走的熱量，導致測量的速度略偏高，這是圖7 中的實測數據比Coles-Fernholz 經驗模型數據略高的原因。但平均偏差小于3.6%，仍然滿足測試要求。

使用平行雙絲熱線對貼壁二維方柱下游壁面摩擦應力進行測量。如圖8 所示，平行雙絲熱線探桿搭載于二維高精度滑臺上，滑臺可控制熱線探針沿壁面法向（y 向）和流向（x 向）運動，快速掃描測量各流向位置的摩擦應力。熱線1 置于上游，熱線2 置于下游。對x/D=1～15 范圍內（測量間距Δx/D=0.5）的壁面摩擦應力進行測量，熱線采樣頻率為25 kHz，每個測點的采樣時間為15 s，共重復測量4 次，獲得摩擦應力統(tǒng)計特征。然后，結合平行雙絲熱線與TR-PIV，同步測量貼壁二維方柱下游流場和壁面摩擦應力。為避免PIV 激光直接照射熱線，熱線與PIV 測試面沿展向錯開2.5 mm。流場測試區(qū)域為：1.5 < x/D < 8.5，0 < y/D < 4.3。壁面摩擦應力測點位置選取近壁面流向平均速度為0 處（x/D=7）。PIV 采樣頻率為1 kHz，平行雙絲熱線采樣頻率為25 kHz，同步采樣時間為8 s。

圖8 貼壁二維方柱下游平板摩擦應力測量示意圖Fig.8 Diagram of WSS measurement downstream of the wallmounted 2D square cylinder

2 結果與討論

2.1 流動結構特征

在貼壁二維方柱的阻塞作用下，流動會在方柱上游壁面某處開始發(fā)生分離并“爬升”，之后在方柱前緣拐角附近發(fā)生二次流動分離[16]，在方柱下游形成回流區(qū)。圖9（a）為貼壁二維方柱下游流場平均流向速度分布。以近壁面流向平均速度為0 處作為判據（圖9 中白色等值線代表=0），可以得到流動再附點[6]位于壁面x/D=7 處，在x/D=7、y/D=4 處，=0.95U0，說明邊界層厚度大于4D，即40 mm。圖9（b）為貼壁二維方柱下游流場雷諾剪應力分布。可以發(fā)現(xiàn)，高雷諾剪應力區(qū)域主要分布于回流區(qū)與外部流動交界處的剪切層附近，在高雷諾剪應力區(qū)域下方，隨著與壁面距離的減小，雷諾剪應力逐漸降低。理論上，壁面處流體的流向與法向脈動速度為0，因此壁面處流動雷諾剪應力為0。

圖9 貼壁二維方柱下游流場統(tǒng)計特征量分布云圖Fig.9 Contour of statistical values of flow field downstream of the wall-mounted 2D square cylinder

圖10 給出了x/D=7 處雷諾剪應力沿法向的分布。受壁面光污染影響，距離壁面最近的PIV 測量位置為y/D=0.1，對應無方柱湍流邊界層y+=45，與熱線測量位置（y+=4.5）相比，距離壁面較遠。但由圖10 可以推斷：隨著y/D 減小，雷諾剪應力逐漸趨近于0，與流向速度梯度相比，雷諾剪應力為小量。因此，當法向位置無限趨近于壁面時，壁面附近流體受到的總剪應力壁面摩擦應力可直接由黏性底層速度梯度確定：值得注意的是，受PIV 采樣頻率限制，計算出的雷諾剪應力值相當于進行了低通濾波，圖9（b）與圖10 中的數值僅有參考價值，但仍能反映雷諾剪應力變化趨勢。

圖10 x/D=7 處流動雷諾剪應力沿法向分布Fig.10 Distribution of Reynolds shear stress of flow along normal direction at x/D=7

由于本文貼壁二維方柱完全浸沒于湍流邊界層中（δ/D=5），方柱前緣流動分離所產生的大尺度流動結構與眾多小尺度非相干結構摻混在一起。為提取大尺度流動結構，采用降階POD 重構方法[17-19]對方柱下游瞬時流場進行重構。本文使用前20 階模態(tài)對流場進行POD 重構，保證重構后的流場所含有的脈動能量占原流場60%以上，并使用λci渦識別準則[20]對大尺度流動結構進行辨識。

本文重點關注流動分離區(qū)內（0 < x/D < 7）的流動結構特征。該區(qū)域流動結構豐富，且具有較強的渦旋與非定常特性。通過統(tǒng)計大量測試結果，將貼壁二維方柱下游流動結構進行分類。從方柱前緣產生的大尺度流動結構以2 種形式向下游運動，分別如圖11 和12 所示（圖中λci*=λciD/U0，λci基于瞬時速度場計算得到，其正負號與當地渦量一致；箭頭表示瞬時脈動速度矢量）。在圖11 中，黑框中從方柱前緣脫落的流動結構在向下游運動過程中發(fā)生分裂：紅框中的子結構沿水平方向向下游運動，與壁面保持一定距離；綠框中的子結構向壁面運動，最終與壁面接觸。在圖12 中，黑框中從方柱前緣脫落的流動結構沿水平方向向下游運動，在運動過程中保持較為完整的形態(tài)。本文將流動結構分裂產生的向壁面運動的子結構（圖11 綠框中的流動結構）稱為“Ⅰ渦”，未分裂的流動結構（圖12 黑框中的流動結構）及分裂產生的沿水平方向向下游運動的子結構（圖11 紅框中的流動結構）統(tǒng)稱為“Ⅱ渦”。結合圖11 和12 中的瞬時脈動速度矢量圖可知：當Ⅰ渦與壁面接觸時，接觸點下游一段距離內出現(xiàn)強烈的下掃流動；當Ⅱ渦從壁面上方經過時，壁面附近出現(xiàn)局部回流。

圖11 貼壁二維方柱下游流動結構運動過程，4 圖為連續(xù)時間序列，時間差1 msFig.11 Motion process of flow structure downstream of the wall mounted 2D square cylinder,where the four figures are continuous time series with interval of 1 ms

圖12 貼壁二維方柱下游流動結構運動過程，4 圖為連續(xù)時間序列，時間差2 msFig.12 Motion process of flow structure downstream of the wall mounted 2D square cylinder,where the four figures are continuous time series with interval of 2 ms

圖13（a）和（b）分別為圖11（a）和（b）流動結構分裂過程中的渦旋強度λci*等值線圖。從圖中可以發(fā)現(xiàn)，即將分裂的流動結構中存在2 個逐漸遠離的高渦量區(qū)域，流動結構中心區(qū)域逐漸被低渦量流動占據，最終2 個渦量比較強的流體微團發(fā)生分離，該過程可能與流動的黏性耗散、湍流擴散及邊界層外部流動的擾動作用有關。

圖13 流動結構分裂過程中的渦旋強度等值線圖Fig.13 Vorticity intensity contour of flow structure during splitting

2.2 流動結構對測點壁面摩擦應力的影響

圖14（a）為貼壁二維方柱下游壁面平均摩擦應力沿流向的分布，與壁面附近流場的平均流向速度分布（圖9）高度一致。圖14（a）中的平均摩擦應力0 值位置（x/D=7 附近）與圖9 中的流動再附點位置（x/D=7）基本一致。圖14（b）為摩擦應力回流間歇因子γτ（γτ為回流方向動態(tài)摩擦應力出現(xiàn)的時間占比）沿流向的分布。在流動再附點處（x/D=7），γτ接近0.5，說明該處來流方向、回流方向的摩擦應力時間占比相當。由上述結果可知，流動再附點處（x/D=7）的摩擦應力具有很強的非定常特性，動態(tài)摩擦應力在0 值附近頻繁變化（即摩擦應力的方向在來流方向、回流方向之間頻繁變化）。本文選取x/D=7 位置作為壁面摩擦應力測點，有利于分析特征流動結構對壁面摩擦應力的影響機理。

圖14 貼壁二維方柱下游壁面摩擦應力統(tǒng)計值沿流向分布Fig.14 Distribution of statistical value of WSS along the flow direction downstream of the wall-mounted 2D square cylinder

圖15 和16 分別為貼壁二維方柱下游流場、x/D=7 處壁面摩擦應力的同步測量結果。圖15 為某些特征時刻對應的瞬時流場連續(xù)時間序列（時間間隔Δt=1 ms，箭頭表示瞬時脈動速度矢量，使用λci準則進行渦識別）。圖15 中的ta～tf時刻對應圖16中標注的ta～tf時刻，tb、td、te對應Ⅰ渦出現(xiàn)的時刻，ta、tc、tf對應Ⅱ渦出現(xiàn)的時刻。在tb、td、te時刻附近時間段內，Ⅰ渦出現(xiàn)于摩擦應力測點上游并與壁面接觸，在圖14（b）、（d）、（e）紅框中Ⅰ渦的誘導下，測點附近產生很強的下掃流動，使流向速度梯度向來流方向增大，測點壁面摩擦應力出現(xiàn)以下3 種變化：1）摩擦應力原為來流方向，Ⅰ渦的出現(xiàn)使摩擦應力絕對值陡增（tb時刻）；2）摩擦應力原為回流方向，Ⅰ渦的出現(xiàn)使摩擦應力絕對值銳減（td時刻）；3）摩擦應力原為回流方向，Ⅰ渦的出現(xiàn)使摩擦應力方向改變?yōu)閬砹鞣较颍╰e時刻）。在ta、tc、tf時刻附近時間段內，圖15（a）、（c）、（f）紅框中的Ⅱ渦從測點上方經過，在壁面附近誘導出局部回流，使流向速度梯度向回流方向增大，測點壁面摩擦應力出現(xiàn)以下3 種變化：1）摩擦應力原為回流方向，Ⅱ渦的出現(xiàn)使摩擦應力絕對值陡增（ta時刻）；2）摩擦應力原為來流方向，Ⅱ渦的出現(xiàn)使摩擦應力絕對值銳減（tf時刻）；3）摩擦應力原為來流方向，Ⅱ渦的出現(xiàn)使摩擦應力方向改變?yōu)榛亓鞣较颍╰c時刻）。值得注意的是，雖然熱線與PIV 測試面已沿展向錯開2.5 mm，但熱線探桿仍會對PIV 圖像造成光污染，影響PIV 圖像質量，使x/D=7 附近區(qū)域出現(xiàn)了流動結構破碎的假象，如圖15（c）所示。

圖15 貼壁二維方柱下游某些時間段內的瞬時流場Fig.15 Instantaneous flow field downstream of wall-mounted 2D square cylinder during certain time periods

圖16 壁面摩擦應力隨時間變化曲線（x/D=7）Fig.16 Curve between WSS and time (x/D=7)

在壁湍流中，近壁面流向渦[21]也會對壁面摩擦應力造成影響。與上述Ⅰ渦、Ⅱ渦對壁面摩擦應力的影響機制類似，在流向渦的誘導下，渦結構的兩側分別出現(xiàn)上拋流動和下掃流動：上拋流動會使壁面附近出現(xiàn)低速條帶，減小流向速度梯度，使摩擦應力值減小；下掃流動會增大壁面附近的流向速度梯度，使摩擦應力值增大[21]。不同之處在于，Ⅰ渦、Ⅱ渦為展向渦，不僅會改變摩擦應力大小，還會改變摩擦應力方向。

3 結 論

本文利用平行雙絲熱線和TR-PIV 研究了貼壁二維方柱下游流動結構對壁面摩擦應力的影響機制。流動經過貼壁二維方柱，會產生2 種典型的大尺度流動結構：向壁面靠近并接觸壁面的近壁流動結構（Ⅰ渦）；平行壁面沿流向運動的流動結構（Ⅱ渦）。Ⅰ渦、Ⅱ渦的出現(xiàn)改變了壁面附近流動速度的大小和方向、影響了測點壁面摩擦應力：增大了流向速度梯度，導致摩擦應力陡增；減小了流向速度梯度，導致摩擦應力銳減；改變了摩擦應力方向。