機動航態下水中航行器邊界層特征參數計算建模研究

謝 華，沈泓萃，張 楠，田于逵

（中國船舶科學研究中心，江蘇 無錫 214082）

0 引 言

流動噪聲是水中航行器表面繞流形成的脈動壓力激勵產生的直接輻射噪聲與結構振動二次輻射噪聲，其噪聲強度隨航速增加而迅速增加，在高航速情況下它將成為水下航行體的主要噪聲源。尤其是在機動航態下，水中航行器周圍流場在時間上和空間上都呈現出強烈脈動特性的湍流，存在不同時間尺度和不同空間尺度的相互干擾，伴隨有強烈的渦發放，充滿了形式各異的大小旋渦，它們的擴散與耗散必然會導致能量的損失，增加了船體阻力，而且與船體、推進器固壁相互作用，誘發船體振動及遠場噪聲輻射。

邊界層參數是水中航行器流激二次聲輻射計算的輸入參數，它的變化直接影響湍流邊界層脈動壓力譜的分布，進而影響流動輻射噪聲的強度，邊界層參數的計算是流動噪聲預報的前期工作，可為水中航行器水動力性能評估提供技術基礎[1]。

1 邊界層參數理論簡介

1904年普朗特在德國海德爾堡第三屆國際數學家學會上宣讀題為“關于摩擦極小的流體運動”的論文，建立了邊界層理論。他根據對水槽中水流實驗的觀測分析，提出了邊界層的概念：粘性極小的流體繞物體流動時，在緊靠物體附近存在著一層極薄的邊界層，其中粘性起著很大的影響。而在邊界層外，流體中的粘性可以忽略不計，可認為是理想流體。由于邊界層極薄，經簡化N－S方程，得出普朗特邊界層方程[2]。

伴隨著電子計算機的廣泛應用，1970年以來邊界層理論發展異常迅猛。近年來邊界層基本方程向縱深發展的方向有兩個方面：一方面由平面、軸對稱流動等二元問題向三元問題發展、并從采用直角坐標系向采用任意正交曲線坐標系、或甚至于用非正交曲線坐標系發展；另一方面逐步深入的考慮了邊界層的高階效應，從古典的“薄”邊界層理論向“厚”邊界層理論、“部分拋物型”邊界層理論等發展。

薄邊界層方程中一個重要假設是在薄邊界層中壓力穿過邊界層不變，所以使方程能夠求解。然而由于船尾流動的復雜特性，致使許多薄邊界層假定失效，必須對其作改進。因此，Bradshaw、Kline[3]運用了復雜剪切流動模型；Patel、Huang等人[4]對厚邊界層進行了詳細測量，并利用軸對稱厚邊界層方程作了計算，這些方程考慮了物體表面縱向和橫向曲率，穿過邊界層的壓力變化和湍流模型，且顯示了計算結果和測量數據吻合良好。

雖然出現了具有復雜形式的各種邊界層方程，例如厚邊界層方程、部分拋物型邊界層方程等，但是薄邊界層方程仍然是有用的。一般都采用混合的方式，流線型物體前端75%仍采用薄邊界層方程來計算，后端25%再改用厚邊界層方程計算。

目前，國內關于水中航行器機動航態下的邊界層特征參數計算尚屬空白，本研究建立的水中航行器機動航態下厚邊界層計算方法將為機動航態下水中航行器流激二次聲輻射計算提供數據輸入，為水中航行器水動力性能評估提供技術基礎。

2 邊界層特征參數計算建模

基于船體的幾何型值和給定的航速，可以對船體進行勢流流場的速度分布和壓力分布預報計算。勢流計算理論在船舶水動力學中已經有廣泛的應用，基于船體表面源分布的Hess-Smith邊界元方法是一種主要的計算方法和勢流流場分析手段。通過對船體主體和附體的表面進行網格劃分的Hess-Smith方法，可以確定船體周圍勢流流場的速度分布和壓力分布信息。船體周圍勢流流場的計算是船體周圍邊界層計算的基礎，它可以為邊界層計算建立船體表面的流線坐標系和確定度量系數。

通常，在進行水中航行器周圍粘性流場的理論計算過程中把航行器周圍的繞流流場分為三個流域：勢流流域、薄邊界層流域和厚邊界層流域。勢流流域是邊界層以外的流域；除了航行體的首尾之外,其75%的范圍可以考慮為薄邊界層流域；航行器尾部和尾流之中可以考慮為厚邊界層流域。

2.1 三維薄邊界層計算建模

薄邊界層計算可以適用于水中航行器繞流的大部分范圍，可用于分析這些范圍上的邊界層發展。在薄邊界層理論中，假定法向壓力梯度為零。邊界層發展可以用微分法和積分法進行計算，通過計算對比來看，積分法具有精度高、速度快和穩定性好的特點。在薄邊界層的積分法中，在水中航行器的表面建立流線坐標系，通過該坐標系中的縱向動量積分方程、橫向動量積分方程、卷吸積分方程和給定的表面摩擦關系式，以及給定的縱向、橫向速度剖面模型，獲得一組封閉的方程組，聯立求解后可以獲得有關的薄邊界層積分參數，薄邊界層發展的計算也為尾部厚邊界層計算提供基礎。以下我們給出了三維薄邊界層計算的詳細的理論計算公式和算法。

在三維薄邊界層假定的基礎上，通過量級分析后可得到流線坐標系下的三維薄邊界層動量積分方程：

式中：h1，h3為度量系數，K13，K31分別表示為 x和 z坐標曲線的曲率，τ120和 τ320分別表示在 x和 z方向的壁面摩擦應力，式中的六個積分厚度分別為：

為建立各積分厚度之間的關系，需要對速度剖面作出適當的假設。試驗表明，縱向速度剖面可用二維的表達式描述，橫向與縱向的速度剖面之間可用經驗關系式建立聯系。這些經驗關系式又稱為速度剖面模型，本文采用Mager速度剖面模型：

把（5）式分別代入（1）式和（2）式，可以看出，原來的七個未知數變成了四個未知數，他們分別是θ11：積分厚度；β2：橫流參數；H：形狀因子；Cf1：壁面摩擦系數。

因此，需要補充兩個方程才能使求解方程組封閉。這里，我們分別選擇了卷吸積分方程和縱向壁面摩擦關系式作為兩個補充方程。

通過求解（6）～（9）式聯立的方程組，可得三維薄邊界層參數分布。

2.2 三維厚邊界層計算建模

水中航行器的尾部流動是很復雜的，但是仍然可以把這個范圍的流動視為湍流邊界層流動，可以采用復雜剪切流處理的FTSL方法和厚邊界層方法來處理。我們采用Patel[5-6]和Larsson[7]的厚邊界層計算方法，在厚邊界層的計算過程中，需要考慮邊界層內的速度、壓力沿法向的變化。在水中航行器表面的流線坐標系中，厚邊界層理論微分方程包括有縱向、橫向和法向的動量方程以及連續方程。通過引入積分厚度參數的假定之后，可以建立考慮了邊界層中速度和壓力沿法向變化的縱向、橫向的動量積分方程、卷吸積分方程和表面摩擦關系式。通過假定縱向、橫向和法向的速度剖面模型可以建立一組封閉的方程組，聯立求解后可以獲得有關的厚邊界層積分參數。以下給出了三維厚邊界層的詳細的理論計算公式和算法。

流線坐標系下的三維厚邊界層動量積分方程如以下所示：

將（5）式代入（10）、（11）式，原方程可變為：

式中有四個未知數，θ11、β2、H和Cf1，因此補充卷吸積分方程和縱向壁面摩擦關系式使方程組封閉。 通過求解（12）、（13）式和（8）、（9）式組成的聯立方程組，可得三維厚邊界層參數分布。

根據上述理論，編寫了相關邊界層特性計算程序。其計算流程圖如圖1所示。

圖1 邊界層特性計算程序流程圖Fig.1 Processing steps of the boundary layer calculation

3 算例分析

3.1 計算模型

計算模型選用Suboff潛艇，模型長4.356 m，最大寬度0.508 m。圖2為潛艇表面計算網格分布，圖3為潛艇表面勢流速度分布。

3.2 計算方法驗證

圖2 Suboff潛艇表面計算網格分布Fig.2 Field-point distribution on surface

圖3 Suboff潛艇表面勢流速度分布Fig.3 Potential velocity distribution on surface

圖4和圖5分別為潛艇表面壓力系數分布以及摩擦系數分布計算結果與試驗結果[8]比較，自由來流風速35 m/s，對應艇長雷諾數1×107。從圖中可以看出，計算結果與試驗值吻合較好，僅在首部和尾部之間略有差異，這是由于計算程序的邊緣效應造成的。從量級與分布趨勢來看，計算得到的潛艇表面時均壓力分布以及摩擦系數分布能夠準確反映物理真實，研究表明本文建立的厚邊界層計算方法對于水中航行器模型邊界層參數計算是可信的。

圖4 Suboff潛艇表面壓力系數分布Fig.4 The pressure coefficient distribution along the hull

圖5 Suboff潛艇表面摩擦系數分布Fig.5 The friction coefficient distribution along the hull

3.3 機動航態下艇體邊界層特征參數計算分析

本文針對Suboff潛艇光體，進行了不同偏航角度和不同雷諾數情況下的三維厚邊界層理論計算，得到了其邊界層特征參數沿艇體分布曲線。計算內容包括：勢流速度、壓力系數、摩擦系數、動量損失厚度、邊界層厚度、形狀因子等。

圖6為機動航態下，艇體表面邊界層特征參數分布。從圖中可以看出，機動航態下潛艇吸力面的動量損失厚度、排擠厚度和邊界層厚度大于直航狀態，形狀因子基本相當，摩擦系數小于直航狀態。潛艇壓力面的動量損失厚度、排擠厚度和邊界層厚度小于直航狀態，形狀因子基本相當，摩擦系數大于直航狀態。

圖7-8為不同偏航角度下，艇體表面邊界層特征參數分布。從圖中可以看出，隨著偏航角度增大，艇體吸力面的動量損失厚度、排擠厚度和邊界層厚度逐漸增大，形狀因子基本相當，摩擦系數逐漸減小；艇體壓力面的動量損失厚度、排擠厚度和邊界層厚度逐漸減小，形狀因子基本相當，摩擦系數逐漸增大。

圖7 不同偏航角度艇體吸力面邊界層特征參數分布（水速6 m/s）Fig.7 The boundary layer characteristic parameters distribution along the suction surface in different yaw angle(flow velocity=6 m/s)

圖8 不同偏航角度艇體壓力面邊界層特征參數分布（水速6 m/s）Fig.8 The boundary layer characteristic parameters distribution along the press surface in different yaw angle(flow velocity=6 m/s)

圖9為不同雷諾數下，艇體表面邊界層特征參數分布。自由來流水速分別為6 m/s、9 m/s和12 m/s，對應雷諾數分別為：2.6×107，3.9×107和5.2×107。從圖中可以看出，隨著雷諾數增大，邊界層特征參數沿艇體表面分布將逐漸減小。

圖9 不同雷諾數艇體表面邊界層特征參數分布Fig.9 The boundary layer characteristic parameters distribution along the hull in different Reynolds number

4 結 論

本文基于邊界層積分求解理論，進行了機動航態下水中航行器厚邊界層計算建模，發展了求解動量積分方程組的數值計算方法。編寫了相關計算程序，并以SUBOFF潛艇零攻角下的表面壓力分布及摩擦系數分布試驗曲線進行了驗證。利用該自主程序，在不同Re和偏航角下，計算獲得了邊界層諸特征參數沿SUBOFF艇體的分布曲線，表明了該軟件具備了評估水下航行體水動力性能的能力。

參 考 文 獻：

[1]Lv S J,Yu M S.Prediction of hydrodynamic radiation noise of underwater vehicle[J].Journal of Hydrodynamics,2007,22：475-482.

[2]Xue Z Y.Boundray layer thoery[M].Beijing：Science Publishing Company,2005：1-5.

[3]Gadd G E.A simple calculation method to assess quality of viscous flow near a ship system[C]//Proc.ISSVR.Goteborg,1978,2：23-33.

[4]Patel V C,Lee Y T.Thick axisymmetric boundary layers and wakes：Experiment and theory[J].International Symposium on Ship Viscous Resistance,1978,5：156-167.

[5]Patel V C,Nakayama A,Damian R.An experimental study of thick turbulent boundary layer near the tail of a body of revolution[R].Iowa Institute of Hydraulic Research,Report 142,1973.

[6]Patel V C.On the equations of a thick axisymmetric turbulent boundary layer[R].Iowa Institute of Hydraulic Research,Report 143,1973.

[7]Larsson L,Chang M Sh.Numerical viscous and wave resistance calculations including interaction[C].Proc．13th Symposium on Naval Hydrodynamics,1980：707-728.

[8]Huang T,Liu H L,Groves N,Forlini T,Blanton J,Gowing S.Measurements of flows over an axisymmetric body with various appendages in a wind tunnel：The DARPA suboff experimental program[C].Proceedings of 19th Symposium on Naval Hydrodynamics,1992.