黃華坤 ，張桂勇*，2，3，孫鐵志 ，宗智 ，2，3

1大連理工大學(xué)船舶工程學(xué)院遼寧省深海浮動(dòng)結(jié)構(gòu)工程實(shí)驗(yàn)室，遼寧大連 116024

2大連理工大學(xué)工業(yè)裝備結(jié)構(gòu)分析國(guó)家重點(diǎn)實(shí)驗(yàn)室，遼寧大連 116024

3高新船舶與深海開(kāi)發(fā)裝備協(xié)同創(chuàng)新中心，上海200240

0 引 言

沖擊射流現(xiàn)象是指具有一定速度的流體通過(guò)撞擊物體表面，使得撞擊區(qū)形成極薄的邊界層，從而得到較高傳熱效率的現(xiàn)象。沖擊射流傳熱原理被廣泛應(yīng)用于工業(yè)生產(chǎn)，如紙張干燥和電子元器件冷卻等。在船舶領(lǐng)域，船舶甲板上的武器發(fā)射、無(wú)人機(jī)發(fā)射和垂直起降飛機(jī)起飛過(guò)程中產(chǎn)生的高溫高速噴射氣流，不僅會(huì)引起甲板的燒蝕，而且還會(huì)帶來(lái)殘余熱應(yīng)力的影響。因此，準(zhǔn)確預(yù)測(cè)沖擊射流的物理現(xiàn)象，對(duì)其工業(yè)應(yīng)用具有重要意義。

許多學(xué)者對(duì)沖擊射流現(xiàn)象開(kāi)展了研究。Gardon等［1］通過(guò)實(shí)驗(yàn)發(fā)現(xiàn)，在較低沖擊距離下，以局部努塞爾數(shù)（Nu）表示的傳熱率表現(xiàn)出第2峰值的特征；當(dāng)沖擊距離較大時(shí)，局部努塞爾數(shù)的第2峰值消失。Colucci等［2］研究了低沖擊距離下（H/D ＜2，H/D為圓形沖擊射流的沖擊距離，H為沖擊高度，D為噴嘴直徑）的沖擊射流現(xiàn)象，指出局部努塞爾數(shù)第1峰值主要受邊界層厚度的影響，第2峰值主要受層流到湍流轉(zhuǎn)捩的影響。隨著計(jì)算流體力學(xué)（CFD）和計(jì)算機(jī)技術(shù)的發(fā)展，CFD方法被廣泛用于模擬沖擊射流現(xiàn)象。在早期的數(shù)值研究中發(fā)現(xiàn)，兩方程湍流模型在滯止點(diǎn)處易產(chǎn)生較高的湍動(dòng)能 k［3-4］，為此，Kato-Launder模型被用于降低滯止點(diǎn)處的湍動(dòng)能［5］。此外，針對(duì)局部努塞爾數(shù)第2峰值的特征，研究表明，SST k-ω模型表現(xiàn)較好［6］。然而，Dutta等［7］指出，基于SST k-ω的低雷諾數(shù)模型雖然在低沖擊距離下（H/B=4，H/B為狹縫沖擊射流的沖擊距離，B為噴嘴寬度）獲得了第2峰值的特征，但在沖擊距離變大時(shí)（H/B=9.2），該模型預(yù)測(cè)出了實(shí)際不存在的第2峰值。為正確反映沖擊射流在不同高度下的特點(diǎn)，Alimohammadi等［8］研究了基于轉(zhuǎn)捩理論的四方程γ-θ模型（也稱Transition SST模型）在圓形沖擊射流中的表現(xiàn)。結(jié)果表明，γ-θ模型可準(zhǔn)確預(yù)測(cè)圓形沖擊射流的傳熱特征；針對(duì)不同的沖擊高度，γ-θ模型表現(xiàn)出較好的預(yù)測(cè)能力。但在圓形沖擊射流中表現(xiàn)良好的湍流模型在狹縫沖擊射流中可能獲得較差的結(jié)果，反之亦然［9］。同時(shí)，γ-θ模型不能預(yù)測(cè)滯止區(qū)附近低湍流強(qiáng)度（Tu＜0.5%）下的橫流失穩(wěn)現(xiàn)象。

因此，本文擬引入基于γ-θ模型的一方程轉(zhuǎn)捩模型來(lái)模擬狹縫沖擊射流現(xiàn)象。該模型耦合了橫流轉(zhuǎn)捩模型，在模擬橫流失穩(wěn)現(xiàn)象和分離流問(wèn)題上具有一定的優(yōu)勢(shì)。同時(shí)，將一方程轉(zhuǎn)捩模型與Kato-Launder模型和SST k-ω模型耦合，并在此基礎(chǔ)上研究一方程轉(zhuǎn)捩模型在預(yù)測(cè)狹縫沖擊射流傳熱規(guī)律上的表現(xiàn)。

1 數(shù)值計(jì)算方法

1.1 控制方程

不可壓縮湍流流動(dòng)的時(shí)均Navier-Stokes方程的質(zhì)量守恒方程、動(dòng)量方程和能量方程分別為：

式中：ρ為密度；μ為動(dòng)力粘性；u和u′分別為平均速度和脈動(dòng)速度；p為壓力；Cp為比熱容；K為傳熱率；T和T′分別為平均溫度和脈動(dòng)溫度；分別為雷諾應(yīng)力張量和湍流熱通量矢量；i和j分別為x，y方向。

1.2 SST k-ω湍流模型

本文采用Menter［10］提出的一種混合標(biāo)準(zhǔn)k-ε模型和k-ω模型的SSTk-ω模型。該模型的湍流粘性系數(shù)μt定義為

式中：k為湍動(dòng)能；ω為比耗散率；α*為阻尼因子；a1=0.31；S為應(yīng)變率；，y為到壁面的距離。

為了封閉控制方程，在模型中增加了湍動(dòng)能k和比耗散率ω方程，表達(dá)式為：

式中，Gk為湍動(dòng)能k的產(chǎn)生項(xiàng)；Gω為比耗散率ω的產(chǎn)生項(xiàng)；Yk，Yω分別為k和ω的耗散；Dω為交叉擴(kuò)散項(xiàng)；Sk，Sω分別為k和ω的源項(xiàng)；σk和σω通過(guò)混合函數(shù)F1進(jìn)行定義，具體為

式中，下標(biāo)1，2分別為k-ω模型和k-ε模型的湍流常數(shù)。

1.3 Kato-Launder模型

Kato-Launder模型修正在Gk中采用渦量代替某個(gè)應(yīng)變率，能避免兩方程模型在滯止點(diǎn)處預(yù)測(cè)的湍動(dòng)能k過(guò)高的問(wèn)題。由于在滯止點(diǎn)附近的流動(dòng)幾乎呈無(wú)旋狀態(tài)，渦量趨于0，因此湍動(dòng)能k變小。該方程的湍動(dòng)能產(chǎn)生項(xiàng)為

Kato和Launder修正后的表達(dá)式為

1.4 一方程轉(zhuǎn)捩模型

在γ-θ模型中引入了轉(zhuǎn)捩臨界動(dòng)量厚度雷諾數(shù)和間歇因子輸運(yùn)方程?；讦?θ模型的一方程轉(zhuǎn)捩模型保留了間歇因子輸運(yùn)方程，同時(shí)使用經(jīng)驗(yàn)關(guān)系式取代了臨界動(dòng)量厚度雷諾數(shù)方程。其中間歇因子γ的輸運(yùn)方程為

式中：σγ=1.0；Pγ和Eγ的定義為

式中：Flength為控制轉(zhuǎn)捩區(qū)的長(zhǎng)度；ca2=0.03，ce2=50；Fonset和Fturb的定義為

1.5 計(jì)算模型及求解設(shè)置

計(jì)算模型如圖1所示，L為沖擊板長(zhǎng)度。由于計(jì)算模型關(guān)于y軸對(duì)稱，為提高計(jì)算效率，建立1/2模型進(jìn)行計(jì)算。表1給出了計(jì)算模型尺度和網(wǎng)格節(jié)點(diǎn)分布情況。

圖1 沖擊射流的計(jì)算模型Fig.1 The computational model for jet impingement flow

本次計(jì)算基于Fluent軟件，使用有限體積法對(duì)計(jì)算域進(jìn)行離散，選擇壓力基穩(wěn)態(tài)求解器和壓力耦合方程半隱式（SIMPLE）算法。梯度離散格式采用最小二乘法，壓力離散格式采用交錯(cuò)壓力格式（Pressure Staggering Option），其他離散格式采用二階迎風(fēng)格式。湍流模型采用SST k-ω模型、Kato-Launder模型和轉(zhuǎn)捩模型（以下簡(jiǎn)稱為“一方程轉(zhuǎn)捩模型”）。網(wǎng)格劃分方案如圖1所示。本文采用結(jié)構(gòu)型網(wǎng)格，根據(jù) Jaramillo［12］的研究，保持第1層網(wǎng)格節(jié)點(diǎn)與壁面間的無(wú)量綱距離y+＜2.5（，其中τw為壁面處的切應(yīng)力），如圖2所示。在壁面和射流出口附近進(jìn)行網(wǎng)格加密，加密方向如圖1箭頭方向所示，以漸疏式網(wǎng)格擴(kuò)展到整個(gè)計(jì)算域，同時(shí)保證網(wǎng)格膨脹率小于2。表1給出了計(jì)算模型尺度以及進(jìn)行網(wǎng)格無(wú)關(guān)性檢驗(yàn)后在x軸和y軸方向的節(jié)點(diǎn)分布情況。

2 結(jié)果分析與討論

2.1 穩(wěn)態(tài)沖擊射流驗(yàn)證

圖2 3種算例下 y+沿沖擊板的分布Fig.2 y+distributions along the impinging plate for 3 cases

表1 計(jì)算模型的尺度和網(wǎng)格節(jié)點(diǎn)分布Table 1 Scale and nodes distribution of computationalmodel

為了驗(yàn)證一方程轉(zhuǎn)捩模型對(duì)穩(wěn)態(tài)沖擊射流換熱的有效性，本文采用了 Ashforth-Frost等［13］的實(shí)驗(yàn)工況，并與其他學(xué)者的數(shù)值計(jì)算結(jié)果進(jìn)行了對(duì)比。實(shí)驗(yàn)研究結(jié)果表明，在H/B=4時(shí)，局部努塞爾數(shù)表現(xiàn)出明顯的第2峰值特征；而在H/B=9.2時(shí)，局部努塞爾數(shù)的第2峰值消失。這2個(gè)算例可以驗(yàn)證一方程轉(zhuǎn)捩模型對(duì)轉(zhuǎn)捩的預(yù)測(cè)能力。Ashforth-Frost等［13］的實(shí)驗(yàn)是在1個(gè)標(biāo)準(zhǔn)大氣壓下進(jìn)行的。射流出口溫度為恒溫300 K，射流出口平均速度雷諾數(shù)基于射流出口寬度，其大小為20 000。沖擊板設(shè)為恒溫310 K，無(wú)滑移固定壁面。上限板設(shè)為恒溫300 K，無(wú)滑移固定壁面。出口采用出流（outflow）邊界條件，出口壓力為1個(gè)大氣壓P=101 325 Pa。計(jì)算的普朗特?cái)?shù)Pr=0.72。

圖3給出了H/B=4時(shí)，得到的局部努塞爾數(shù)與實(shí)驗(yàn)數(shù)據(jù)以及其他學(xué)者的數(shù)值結(jié)果沿沖擊面的分布情況。從圖中可以看出，一方程轉(zhuǎn)捩模型可準(zhǔn)確預(yù)測(cè)局部努塞爾數(shù)第1峰值和第2峰值的大小及位置。受轉(zhuǎn)捩過(guò)程影響，在3≤x/B≤7范圍內(nèi)，傳熱率逐漸升高；當(dāng)轉(zhuǎn)捩完成后，傳熱率逐漸降低。在滯止點(diǎn)附近，流動(dòng)處于低湍流強(qiáng)度狀態(tài)，即準(zhǔn)層流狀態(tài)，這使得滯止點(diǎn)附近相對(duì)于壁面射流區(qū)的傳熱率得到了極的大提高。

圖3 H/B=4時(shí)局部努塞爾數(shù)分布Fig.3 The distribution of local Nusselt number（H/B=4）

受益于Kato-Launder模型修正作用，可用一方程轉(zhuǎn)捩模型準(zhǔn)確預(yù)測(cè)局部努塞爾數(shù)的第1峰值。如前所述，在滯止點(diǎn)處，流動(dòng)幾乎處于無(wú)旋狀態(tài)，因此Kato-Launder模型中渦量的引入使滯止點(diǎn)處的湍動(dòng)能下降；且在第2峰值附近可明顯看到渦的存在（圖4）。雖然Kato-Launder模型無(wú)法對(duì)渦進(jìn)行計(jì)算，但是渦的脫落導(dǎo)致渦量增加，故一方程轉(zhuǎn)捩模型在2次峰值附近預(yù)測(cè)的傳熱率與SST k-ω模型相比有所增加，從而提高了渦附近湍動(dòng)能的預(yù)測(cè)精度［14］。而無(wú)修正模型的RANS k-ω模型［15］和 SST k-ω模型在滯止點(diǎn)處過(guò)高地預(yù)測(cè)了局部努塞爾數(shù)。其中RANS k-ω模型無(wú)法預(yù)測(cè)局部努塞爾數(shù)的第2峰值。SST k-ω模型雖然能夠預(yù)測(cè)第2峰值，但其預(yù)測(cè)結(jié)果低于實(shí)驗(yàn)值，且第2峰值的位置與實(shí)驗(yàn)值相比提前了31%，這也說(shuō)明轉(zhuǎn)捩模型的加入使一方程轉(zhuǎn)捩模型具備了良好的轉(zhuǎn)捩預(yù)測(cè)能力。因粘性底層的流動(dòng)狀態(tài)與湍流轉(zhuǎn)捩相關(guān)，故沒(méi)有粘性底層修正的SST k-ω模型和RANS k-ω模型難以準(zhǔn)確預(yù)測(cè)轉(zhuǎn)捩區(qū)間的位置。RANS/LES模型［15］繼承了大渦模擬（LES）模型在滯止點(diǎn)附近湍流強(qiáng)度幾乎為零的特點(diǎn)，因此正確地捕捉到了局部努塞爾數(shù)的第1峰值，但卻提前預(yù)測(cè)了第2峰值的位置。除一方程轉(zhuǎn)捩模型外，其他數(shù)值模型同樣存在第2峰值預(yù)測(cè)不準(zhǔn)確的現(xiàn)象。整體上，一方程轉(zhuǎn)捩模型獲得了與實(shí)驗(yàn)結(jié)果基本一致的分布趨勢(shì)，這說(shuō)明該模型不僅能準(zhǔn)確模擬滯止點(diǎn)附近的流動(dòng)狀態(tài)，還能準(zhǔn)確預(yù)測(cè)壁面射流區(qū)層流到湍流的轉(zhuǎn)捩過(guò)程。

圖4 H/B=4時(shí)的流線圖Fig.4 The streamlines for H/B=4

圖5給出了H/B=9.2時(shí)的數(shù)值計(jì)算與實(shí)驗(yàn)結(jié)果中局部努塞爾數(shù)沿沖擊面分布的對(duì)比結(jié)果。從圖可看出，局部努塞爾數(shù)的第2峰值消失。圖6顯示了H/B=9.2時(shí)的流線圖，可見(jiàn)回流中心位置較H/B=4時(shí)的遠(yuǎn)，因此在0≤x/B≤10范圍內(nèi)受渦的影響較小。由于轉(zhuǎn)捩受渦的影響，其對(duì)傳熱的影響變?nèi)?，局部努塞爾?shù)的第2峰值消失。RANS/LES模型和一方程轉(zhuǎn)捩模型的預(yù)測(cè)結(jié)果與實(shí)驗(yàn)結(jié)果吻合較好。此外，RANS k-ω模型［15］在滯止點(diǎn)附近預(yù)測(cè)的局部努塞爾數(shù)與實(shí)驗(yàn)值相比更低。這是因?yàn)樵谳^大的沖擊距離下，模型中的應(yīng)力限制器對(duì)流動(dòng)幾乎無(wú)影響，這導(dǎo)致滯止點(diǎn)附近的湍流強(qiáng)度升高，傳熱率下降。而RANS/LES模型［15］則稍好于RANS k-ω模型［14］，但誤差仍然有12%。一方程轉(zhuǎn)捩模型較為準(zhǔn)確地預(yù)測(cè)了滯止點(diǎn)處的局部努塞爾數(shù)。RANS k-ω模型［14］預(yù)測(cè)的局部努塞爾數(shù)第2峰值與實(shí)驗(yàn)值不相符，而RANS/LES模型和一方程轉(zhuǎn)捩模型在趨勢(shì)上則與實(shí)驗(yàn)結(jié)果基本保持一致。與一方程轉(zhuǎn)捩模型相比，SST k-ω模型在滯止點(diǎn)處過(guò)高地預(yù)測(cè)了局部努塞爾數(shù)，這反映出Kato-Launder模型起到了降低滯止點(diǎn)處湍動(dòng)能的作用。但在下游區(qū)域，SST k-ω模型與一方程轉(zhuǎn)捩模型預(yù)測(cè)的局部努塞爾數(shù)基本一致，這說(shuō)明在高沖擊射流情況下，轉(zhuǎn)捩模型對(duì)傳熱預(yù)測(cè)的影響較小，同時(shí)表明一方程轉(zhuǎn)捩模型可正確反映沖擊射流中的轉(zhuǎn)捩強(qiáng)度，從而準(zhǔn)確給出局部努塞爾數(shù)的分布。

圖5 H/B=9.2時(shí)局部努塞爾數(shù)分布Fig.5 The distribution of local Nusselt number（H/B=9.2）

圖6 H/B=9.2時(shí)的流線圖Fig.6 The streamlines for H/B=9.2

2.2 脈動(dòng)沖擊射流驗(yàn)證

脈動(dòng)沖擊射流與穩(wěn)態(tài)沖擊射流的不同之處在于，其提供了一個(gè)擾動(dòng)源。Lytle等［16］的研究表明，平行于沖擊面的速度脈動(dòng)峰值與沖擊面的平均傳熱率峰值相關(guān)，而擾動(dòng)源又與渦的形成相關(guān)。渦影響了射流的發(fā)展以及層流到湍流的轉(zhuǎn)捩過(guò)程，進(jìn)而影響到?jīng)_擊面上的傳熱分布。在此，我們采用了Mladin等［17］的實(shí)驗(yàn)工況。在本算例中，射流出口采用正弦變化的速度分布，其公式為

式中：vavg為入口的平均速度；A為脈動(dòng)幅度；fr為脈動(dòng)頻率；t為時(shí)間。射流進(jìn)口溫度為300 K。沖擊板是厚度為5 mm的鋁鎳材質(zhì)，導(dǎo)熱率為150 W·m-1·K-1，溫度設(shè)為 323 K，無(wú)滑移固定壁面；上限板的溫度設(shè)為恒溫300 K，無(wú)滑移固定壁面。出口設(shè)置為出流邊界條件，出口壓力為1個(gè)大氣壓P=101 325 Pa。計(jì)算得到的沖擊高度為H/B=5，普朗特?cái)?shù)Pr=0.72。

圖7給出了在穩(wěn)定狀態(tài)下局部努塞爾數(shù)沿沖擊面的分布情況。從圖中可以看出，局部努塞爾數(shù)的第2峰值仍然較為明顯。這是因?yàn)榛亓髦行脑趚/B=10.8處，如圖8所示。一方程轉(zhuǎn)捩模型能準(zhǔn)確捕捉到滯止點(diǎn)處局部努塞爾數(shù)的大小，同時(shí)也能準(zhǔn)確預(yù)測(cè)局部努塞爾數(shù)第1個(gè)最低點(diǎn)的值。而k-ω-v2-f（k-ω模型和v2-f模型）混合模型［18］則過(guò)低地預(yù)測(cè)了滯止點(diǎn)和最低點(diǎn)處的局部努塞爾數(shù)；在下游區(qū)，k-ω-v2-f模型［18］的表現(xiàn)較一方程轉(zhuǎn)捩模型差。一方程轉(zhuǎn)捩模型獲得了與實(shí)驗(yàn)基本一致的結(jié)果。這是因?yàn)橐环匠剔D(zhuǎn)捩模型引入了間歇因子輸運(yùn)方程來(lái)模擬轉(zhuǎn)捩，而轉(zhuǎn)捩的發(fā)生與擾動(dòng)相關(guān)，因此對(duì)于脈動(dòng)沖擊射流，一方程轉(zhuǎn)捩模型依然能夠給出較為準(zhǔn)確的結(jié)果。

3 結(jié) 語(yǔ)

本文提出采用一方程轉(zhuǎn)捩模型對(duì)狹縫沖擊射流問(wèn)題進(jìn)行數(shù)值模擬，并研究不同沖擊高度、不同雷諾數(shù)下的穩(wěn)態(tài)沖擊射流和脈動(dòng)沖擊射流傳熱分布規(guī)律。結(jié)果表明，一方程轉(zhuǎn)捩模型不僅能捕捉到局部努塞爾數(shù)的第1峰值，也能準(zhǔn)確預(yù)測(cè)受轉(zhuǎn)捩影響的局部努塞爾數(shù)第2峰值的大小和位置。在4≤H/B≤9.2時(shí)，一方程轉(zhuǎn)捩模型在傳熱預(yù)測(cè)上表現(xiàn)出了良好的魯棒性特征。