非均勻湍流中螺旋槳激振力寬頻譜研究

蒲汲君, 周其斗, 孟慶昌

(海軍工程大學艦船工程系, 430033, 武漢)

水下螺旋槳激振力產生的噪聲一直以來都受到國內外許多學者的關注和研究[1],該非定常力一般分為兩部分:周期力分量和寬頻力分量,這些非定常力的產生源于螺旋槳與周圍非均勻、非定常流場的相互作用。寬頻激振力一般在一階和二階葉頻處產生較大的波峰,這些波峰的產生主要是由于螺旋槳周圍非定常湍流場的影響。為降低螺旋槳激振力產生的噪音,深入研究非定常湍流場中螺旋槳激振力的寬頻公式是很有必要的。

已有很多學者都對湍流中螺旋槳激振力進行了研究。Blake使用頻譜法研究了螺旋槳激振力的寬頻譜公式,但由于在頻譜法中對螺旋槳旋轉方向的速度相關性進行了過度假設,因此與試驗值相比誤差較大[2]。Sevik則首次使用了相關性法研究螺旋槳的寬頻激振力,但他在研究中只考慮了來流速度對激振力的影響,忽略了螺旋槳旋轉帶來的影響[3],因此他的研究結果只能顯示激振力譜的大致趨勢,而不能預報葉頻處的激振力波峰現象。隨后,Martinez等將螺旋槳的旋轉影響考慮在內,成功地在計算中得到了葉頻處的激振力波峰現象,但波峰值與試驗值相比仍有較大差距[4-8]。Sevik在水洞中研究了10槳葉固定弦長螺旋槳在不同湍流尺寸下的非定常激振力,在試驗中顯示了激振力寬頻譜在葉頻處出現波峰的結果,從試驗中還發現波峰值隨著湍流積分尺度增大而增大的現象[3]。Thompson隨后研究了不同槳葉數的螺旋槳激振力寬頻譜,他的研究結果顯示,隨著葉間距離的增大,激振力寬頻譜的波峰值隨之增大[9]。

在湍流方面,Dryden在試驗中對湍流場的基本特性進行了研究,并引入湍流積分尺度Λ來描述湍流速度的相關性規律,他認為推導得到的相關性公式能表達湍流場的基本特性[10]。Taylor等引入湍流特征尺度λ同樣對湍流速度的相關性公式進行了研究[11-12],結果表明,在湍流中兩點距離相較于湍流積分尺度Λ較小時,使用Taylor等推導的公式是合適的,在距離較大時,Dryden推導的公式則能較好地描述速度的相關性規律[13]。

總的來說,現今對螺旋槳激振力寬頻譜的研究還停留在湍流均勻且各向同性的條件下。本文考慮了在螺旋槳前安置導流片的情況(導流片的存在會改變螺旋槳前的湍流特性),并研究了不同葉數的導流片對激振力譜的影響規律。

1　螺旋槳激振力公式推導

1.1　各向同性湍流

相關性法最早由Sevik提出,并與試驗值進行了對比,但只能預測推力譜分布的大致趨勢[3]。后來,Jiang等在該方法基礎上考慮了螺旋槳槳葉的轉動,成功地預報了葉頻處的波峰現象[5-7]。圖1給出了相關性法的原理示意圖。下面簡要介紹相關性法的原理。單片槳葉在半徑r截面處受到的推力為

(1)

圖1　湍流導致的激振力示意圖

假定湍流場是均勻且各向同性的,于是湍流速度場可由任意兩點的速度相關性表示。通過將不同半徑處的升力相關性做空間上的積分運算,再取傅里葉變換后可以得到螺旋槳推力譜的表達式如下

(2)

將式(1)代入式(2),得到推力譜的表達式如下

φF(ω)=π2ρ2?drdr′cosγcosγ′C(r)C(r′)·

(3)

現在問題集中在確定任意兩點垂直于槳葉表面湍流速度分量的相關性函數QNN上。在實際工程中,實時測量各點的相關性函數幾乎是不可能的,因為這需要在螺旋槳周圍安置大量的探測儀,同時還不能對流場產生影響。于是,在計算中一般使用無螺旋槳下的湍流速度場,該湍流速度場已做過大量實驗[10],有較好的理論和數據基礎。首先引入固定坐標系,將wN的相關性函數QNN用任意方向湍流速度wi的相關性函數Qij來表示

eiN(r,2πS/N)ejN′(r′,2πS′/N+Ωτ)

wj(r′,2πS′/N+Ωτ)=

eiN(R,2πS/N)ejN′(r′,2πS′/N+Ωτ)·

Qij(r,2πS/N,r′,2πS′/N+Ωτ,τ)

(4)

式中:eiN(r,2πS/N)和ejN′(r′,2πS′/N+Ωτ)為槳葉表面垂直方向的方向余弦。

如前所述,湍流場假定為均勻且各向同性的。在該假設下,任意兩點縱向和橫向速度的相關函數可用g(η)和f(η)表示,其中η為兩點間距離。此時,Qij的表達式如下

(5)

式中:δij為Kronecker算子。假設流體是不可壓縮的,于是f(η)與g(η)的關系如下

(6)

(7)

其中兩點間距離分量的表達式如下

ηx=U∞τ

(8)

1.2　各向異性湍流

與各向同性湍流下公式相比,這里主要考慮了槳葉旋轉到不同位置時湍流度和來流速度的變化,于是引入時間t

L(r,t,S)=

(9)

與前面均勻湍流下激振力的推導相同,這里假定湍流速度場可由任意兩點的速度相關性表示。通過將不同半徑處的升力相關性做空間上的積分運算,再取傅里葉變換后可以得到螺旋槳推力譜的表達式如下

?drdr′cosγcosγ′·

(10)

將式(9)代入式(10),得到推力譜的表達式如下

φF(ω)=π2ρ2?drdr′cosγcosγ′U(r,t,S)U(r′,t+
τ,S′)wi(r,2πS/N+Ωt,t)wj(r′,2πS′/N+

Ω(t+τ),t+τ)eiN(r,2πS/N+Ωt)ejN′(r′,2πS′/N+

(11)

式中:t為絕對時間;τ為兩點間的相對時間。首先設定

wj(r′,2πS′/N+Ω(t+τ),t+τ)·

ejN′(r′,2πS′/N+Ω(t+τ))·

U(r,t,S)U(r′,t+τ,S′)·

eiN(r,2πS/N+Ωt)dt

(12)

這里f(r,r′,S,S′,τ)為包含時間t和τ的所有相關項,詳見附錄A。

圖2　4葉導流片產生湍流示意圖

本文計算的是導流片后的螺旋槳激振力,以4葉導流片為例,4葉導流片安裝在螺旋槳正前方,距離螺旋槳盤面一倍螺旋槳直徑。圖2給出了4葉導流片后湍流的示意圖。由圖可見,在導流片后會產生大量的旋渦,使得相關區域湍流度明顯激增,這些不同尺度的渦則是產生寬頻激振力的主要原因,而在其他位置來流速度的湍流度較小,螺旋槳盤面處的湍流實際上是非均勻的。這里可以將螺旋槳盤面分為8個區域,從圖2中觀察可知只有ACEG區域的湍流度較大,所以只有流經這些區域的槳葉葉片才能產生激振力。另外,值得注意的是,ACEG區域中任一區域與其他區域之間都相隔一個湍流度幾乎為0的區域(BDFH),所以若槳葉元位于ACEG中不同的兩個區域,則它們之間湍流速度的相關性依然為0;只有兩槳葉元同時位于ACEG中的任一相同區域時,它們的湍流速度才存在相關性。由于這里主要探究非均勻流場對螺旋槳寬頻激振力的影響規律,所以劃分的湍流區域只要保證理論推導與實際情況近似相同即可。為簡便計算,將湍流區域劃分為角度為20°的扇形,以后在試驗中可以精確測量各種條件下的湍流場。關于湍流場的數學模型見附錄B。

文獻[14]指出,在來流靠近螺旋槳盤面時,湍流渦會拉伸變形。圖3顯示了由于螺旋槳抽吸作用導致的湍流渦拉伸變形的原理。從圖中可以看到,在靠近螺旋槳的過程中,湍流渦沿縱向被拉伸變形,導致了縱向湍流積分長度Λx變大,橫向湍流積分長度Λy、Λz變小,從而減小了x方向湍流速度的相關性,增大了y、z方向湍流速度的相關性。為考慮湍流渦變形的影響,引入拉伸系數ax、ay、az[6],得到如下關系式

(13)

值得注意的是,湍流渦變形程度與湍流渦的尺寸有很大關系。一般來說,湍流渦的尺寸越大,湍流渦的變形程度越大,因此ax越小,ay、az越大。

圖3　湍流渦拉伸變形圖

2　計算和對比

2.1　各向同性湍流下計算結果

Sevik在直徑為1.22 m的水洞中對螺旋槳在湍流中的非定常激振力進行了研究[3],試驗螺旋槳有10片槳葉,無側斜,槳葉弦長恒定且為25.42 mm,槳葉半徑為10.16 cm,槳轂半徑為2.5 cm,具體幾何形狀如圖4所示。用于產生湍流的網格有兩種尺寸,分別為10.16 cm和15.25 cm,螺旋槳安裝的位置分別為網格后24.3 m(20倍網格尺寸)和36.6 m(20倍網格尺寸)。由Dryden的研究結果可知,該種工況下產生的湍流度相同且都為3%,湍流積分尺度為0.2倍網格尺寸。來流速度都保持為4.57 m/s,螺旋槳進速為1.22。

圖5和圖6分別給出了兩種湍流尺度下不同方法計算得到的螺旋槳激振力寬頻譜與試驗值的比較結果,其中計算參數的取值如表1所示。

表1　計算參數的取值

觀察圖5a可以發現,相關性法的計算結果在一階葉頻處出現較小波峰,該現象與試驗結果相同,但波峰值與試驗值相比差距較大,在一階葉頻后,激振力寬頻譜隨著頻率的增大而減小,總的來說,相關性法結果的趨勢與試驗值相比大致一致,但計算值小于試驗值,波峰處的誤差尤其大。圖6a的結果顯示,隨著湍流積分尺度增大,相關性法的計算結果幾乎無明顯改變,而試驗值在二階葉頻和三階葉頻處都出現了較為明顯的波峰,這在計算結果中幾乎沒有得到體現。由此可以發現,相關性法的計算結果在大致趨勢上能與試驗值保持一致,但無法準確地預報葉頻處的波峰現象。

觀察圖5b和圖6b可以發現,改進后相關性法的計算結果與試驗結果的趨勢相同,且在葉頻波峰處能較好地符合試驗結果。特別是在15.24 cm網格下,改進后的相關性法能較好反映二階和三階葉頻處的波峰現象,雖然二階葉頻的波峰值與試驗值相比誤差較大,但取得波峰的頻率點位置大致相同。由此可以得出結論:將新的湍流速度相關性公式和拉伸系數ax、ay、az引入相關性法后,計算得到的螺旋槳激振力寬頻譜能更好地符合試驗結果。

(a)初始公式計算結果

(b)改進后公式計算結果圖5　　　　10.16 cm網格下激振力寬頻譜計算結果和試驗結果[3]

(a)初始公式計算結果

(b)改進后公式計算結果圖6　　　　15.24 cm網格下激振力寬頻譜計算結果和試驗結果[3]

2.2　各向異性湍流下計算結果

在螺旋槳前安置導流片后,螺旋槳盤面的湍流是不均勻且各向異性的,這里分別研究了單葉導流片、兩葉導流片和4葉導流片后10葉、5葉和3葉螺旋槳的激振力寬頻譜,其中其他參數的設置都與前面各向同性的計算相同,湍流積分尺度Λ=0.02。

圖7為10葉螺旋槳在各向異性和同性湍流下的計算結果,可以發現,與各向同性湍流下的螺旋槳激振力寬頻譜相比,導流片的存在會大幅度地提高葉頻處的波峰,特別是在多階葉頻下激振力寬頻譜的波峰現象仍比較明顯,與各向同性湍流下的激振力寬頻譜有明顯區別。這主要是由于在該流場下只有特定區域才能產生激振力,使得激振力的周期性明顯增大。另外,較各向同性湍流下的計算結果,加入導流片后螺旋槳激振力寬頻譜整體較小。這是由于各向同性湍流下各個區域都有湍流作用在螺旋槳槳葉上,而各向異性湍流下只有導流片后才會產生湍流,隨著導流片葉數增加,能產生激振力的區域也隨之增大,因此隨著導流片葉數增加,螺旋槳激振力寬頻譜隨之增大。除此之外,還可以發現,各向異性和各向同性湍流下螺旋槳激振力寬頻譜取得波峰的位置一致。由此得出,導流片并不會影響激振力寬頻譜波峰的位置,該位置依然取決于葉頻的大小。

圖7　　　　不同葉數導流片下10葉螺旋槳激振力寬頻譜計算結果

圖8為不同葉數導流片和螺旋槳下激振力寬頻譜計算結果,其中10葉、5葉和3葉螺旋槳葉頻分別為190、95、57 Hz,觀察發現螺旋槳激振力寬頻譜波峰的位置都近似位于螺旋槳葉頻位置,這也進一步驗證了前面的研究。另外,從圖中還可以發現,隨著螺旋槳葉數增大,激振力寬頻譜波峰的幅度也隨之增大。由此得出,螺旋槳葉數越多,產生的激振力寬頻譜波峰的幅度越大。除此之外還可以發現,在導流片葉數不多(各個導流片產生的湍流都局限在自己的區域,沒有發生交匯)的情況下,各個數量導流片下產生的激振力寬頻譜在各個頻率下的大小相差都相同,激振力寬頻譜隨頻率的變化趨勢都基本相同。

(a)10葉螺旋槳

(b)5葉螺旋槳

(c)3葉螺旋槳圖8　　　　不同葉數導流片和螺旋槳下激振力寬頻譜的計算結果

從以上研究可以得出結論,流場的均勻性、各向同性或各向異性并不會對激振力寬頻譜的波峰位置產生影響,它們只能影響激振力寬頻譜葉頻處的波峰程度,一般來說,流場的不均勻程度越大,激振力寬頻譜葉頻處的波峰程度越大。

3　結　論

本文介紹并推導了計算螺旋槳激振力寬頻譜的相關性法,在相關性法的基礎上,考慮了螺旋槳抽吸作用導致的湍流渦拉伸變形的影響,引入了拉伸系數ax、ay、az,分別推導了在各向同性和各向異性下螺旋槳激振力的寬頻公式,雖然本文并沒有對導流片產生的各向異性湍流場進行準確測量,只是進行了一定的假設和估計,但該計算對研究各向異性湍流場下螺旋槳寬頻激振力的變化規律有一定的指導意義,本文得到了以下結論。

(1)相關性法只能預報螺旋槳激振力寬頻譜的大致趨勢,無法得到二階以上的葉頻處的波峰現象。

(2)改進后的相關性法不僅能得到螺旋槳激振力寬頻譜的大致趨勢,還能較好地預報葉頻處的波峰現象。

(3)流場均勻性、各向同性或各向異性并不會對激振力寬頻譜的波峰位置產生影響,它只能影響激振力寬頻譜葉頻處的波峰程度。

參考文獻:

[1]屈煥成, 張荻, 謝永慧, 等. 汽輪機調節級非定常流動的數值模擬及汽流激振力研究 [J]. 西安交通大學學報, 2011, 45(11): 39-44.

QU Huancheng, ZHANG Di, XIE Yonghui, et al. Numerical simulation of unsteady flow phenomena and excitation force in control stage of steam turbine [J]Journal of Xi’an Jiaotong University, 2011, 45(11): 39-44.

[2]BLAKE W B. Mechanics of flow-induced sound and vibration: IIComplex flow-structure interactions [M]. New York, USA: Academic Press, 1986; 927-931.

[3]SEVIK M. Sound radiation from a subsonic rotor subjected to turbulence: NASA SP-304 Part 11 [R]. Washington, USA: NASA, 1974.

[4]MARTINEZ R, WEISSMAN K. Spatial-domain analysis of the thrust on a propeller cutting through isotropic turbulence: CAA Rept U-1894-358.47 [R][S.l.]: CAA, 1990.

[5]JIANG C W, CHANG M S, LIU Y N. The effect of turbulence ingestion on propeller broadband thrust [R][S.l.]: Taylor Naval Ship Research and Development Center, 1991: 51.

[6]CATLETT M, ANDERSON J, STEWART D. Aeroacoustic response of propellers to sheared turbulent inflow[C/OL]∥18th AIAA/CEAS Aeroacoustics Conference, 4-6 June 2012, Colorado Springs, CO. [2017-06-15]. DOI: 10.2514/6.2012-2137.

[7]ANDERSON J M, CATLETT M, STEWART D O. Modeling rotor unsteady forces and sound due to homogeneous turbulence ingestion [J]AIAA Journal, 2015, 53(1): 81-92.

[8]CHEN Y, WANG L, HUA H X. Longitudinal vibration of marine propeller-shafting system induced by inflow turbulence [J]. Journal of Fluids & Structures, 2017, 68: 264-278.

[9]THOMPSON D E. Propeller time-dependent forces due to nonuniform flow [D]Pennsylvania, USA: The Pennsylvania State University, 1976.

[10] DRYDEN H L, SCHUBAUER G B, MOCK W C J, et al. Measurements of intensity and scale of wind-tunnel turbulence and their relation to the critical Reynolds number of spheres: NACA Rep [R]. [S.l.]: NACA, 1937: 581.

[11] TAYLOR G I. Statistical theory of turbulence [J]. Proc Roy Soc of London: Series A, 1935, 4151(873): 21-478.

[12] VON K T. On the statistical theory of turbulence [J]. Proceedings of the National Academy of Sciences of the United States of America, 1937, 23(2): 98-105.

[13] DRYDEN H L. A review of the statistical theory of turbulence [J]. Quart Appl Math, 1942: 7-42.

[14] MAJUMDAR S J, PEAK N. Noise generation by the interaction between ingested turbulence and a rotating fan [J]. Journal of Fluid Mechanics, 2000, 359: 181-216.

附錄A

這里f(r,r′,S,S′,τ)為包含時間t和τ的所有相關項,繼續簡化得到

f(r,r′,S,S′,τ)=

E(g(r,t,S)g(r′,t+τ,S′)wN(r,2πS/N+

Ωt,t)wN(r′,2πS′/N+Ω(t+τ),t+τ))

式中:g(r,t,S)為關于來流速度U(r,t,S)的函數,其表達式如下

(A1)

由于湍流函數w(r,2πS/N,t)與來流速度的相關函數g(r,t,S)之間的相關性幾乎為0,于是可近似認為

f(r,r′,S,S′,τ)≈E(g(r,t,S)g(r′,t+τ,S′))·
E(wN(r,2πS/N+Ωt,t)·
wN(r′,2πS′/N+Ω(t+τ),t+τ))

(A2)

式中:E為期望函數。將湍流速度和來流速度分別求期望,得到

(A3)

(A4)

E(g(r,t,S)g(r′,t+τ,S′))≈

(A5)

附錄B

按照上文對湍流速度的分析,相關函數有如下定義。

Qij(r,r′,S,S′,τ,t)=

Qij(r,2πS/N,r′,2πS′/N+Ωτ,τ)

此外,Qij(r,r′,S,S′,τ,t)=0,其中b=1,2,…,M,M為導流片總葉數。

[本刊相關文獻鏈接]

陳堯興,李志剛,李軍.非均勻進汽時迷宮密封汽流激振動力特性的研究.2017,51(9):145-152.[doi:10.7652/xjtuxb2017 09021]

李震,琚亞平,張楚華.高精度中心-WENO混合格式在可壓縮各向同性湍流大渦模擬中的應用.2017,51(1):25-30.[doi:10.7652/xjtuxb201701005]

李典,劉小民,楊羅娜.仿鸮翼的三維仿生翼型葉片氣動特性研究.2016,50(9):111-118.[doi:10.7652/xjtuxb201609018]

任蕓,吳登昊,牟介剛,等.考慮旋轉和曲率的湍流模型修正及應用.2016,50(9):25-30.[doi:10.7652/xjtuxb201609004]

蘭艾青,王志恒,于洪石,等.合成射流激勵位置對控制翼型大攻角分離流動影響的數值研究.2016,50(7):110-117.[doi:10.7652/xjtuxb201607017]

李志剛,李軍,豐鎮平.進口防旋板對孔型密封非定常氣流激振特性的影響.2016,50(1):16-21.[doi:10.7652/xjtuxb 201601003]