基于數值計算方法的水下合成射流推力特性

賈連超 耿令波 胡志強 楊 翊 王 超

1.中國科學院沈陽自動化研究所機器人學國家重點實驗室，沈陽,1100162.東北大學機械工程與自動化學院，沈陽,110819 3.中國科學院大學,北京，100049

0 引言

傳統的水下機器人推進方式主要依靠螺旋槳與舵的配合，當航速較低時很難有效地調節機器人的姿態，并且長時間低速航行會產生低速顫振現象進而產生噪聲，機器人隱蔽性較差。為了實現機器人在低速航行時的姿態控制，一種較為常見的解決方案是設置多部噴水推進機構，這又給機器人帶來了阻力大、質量增加、系統復雜等不利因素。未來的水下機器人由于工作任務復雜，需要具備低噪聲、高隱蔽性、高機動性、結構簡單、體積小、水阻小等特點。傳統的操縱系統已經不能滿足下一代機器人的上述需求，因此，需要一種新型的推進及操縱系統，以滿足水下機器人全航速范圍內良好操縱性的要求。

合成射流本質上屬于非穩態射流。國外研究人員已經通過實驗證實，合成射流的推進效率明顯高于穩態射流，并且發現射流中渦的出現可以明顯提高推進效率[1-8]。WHITTLESEY等[5]的實驗結果顯示，相比穩態射流，含有渦的非穩態射流使得推進效率提升40%，且與阻力相關的水動力性能提升70%以上。水下的合成射流受到的關注越來越多[9-16]。KRUEGER[1]通過實驗發現，渦的產生使得射流推力增大20%以上；KRIEG等[12-13]仿照水母的運動原理設計了一種活塞型合成射流推進器，并對該合成射流推進器進行了一系列推力實驗，實驗結果顯示，雖然合成射流在一個周期內的質量變化為0，但其動量變化不為0，可以產生凈推力[12-13]，證明了合成射流用于水下推進、姿態調節的可行性。

合成射流驅動器具有體積小、質量小、結構簡單、效率高等優點，且可以完全集成在機器人本體內部，對機器人流體外型影響很小，是一種非常有潛力的新型水下姿態調節方式。目前關于水下合成射流激勵器的研究大多集中在實驗方面，數值計算方法方面的研究很少，且系統研究網格尺度、邊界層以及湍流模型對推力計算精度的影響未見報道。另外，目前對合成射流推力增強性質的機理認識不夠清晰。為此，本文研究了水下合成射流的數值計算方法，系統研究了網格尺度、邊界層尺度以及湍流模型對推力及流場的影響。建立了合成射流推力的數學模型，并設計了一種出口可偏轉的合成射流激勵器，通過對該激勵器進行數值實驗，并結合本文建立的數學模型，闡釋合成射流推力增強性質的內在機理。

1 數值計算方法

THOMAS等[17]通過CFX軟件對一種采用非對稱型激勵函數的射流激勵器流場進行了數值仿真，并研究了“伴流”對射流誘導動量的作用。KUMAR等[18]在FLUENT軟件中采用大渦模擬，對三維射流激勵器進行了數值仿真。

合成射流本質上屬于一種湍流現象，對湍流的數值計算一般包括直接模擬、大渦模擬及雷諾時均方法模擬，直接模擬與大渦模擬計算量很大，所以在工程應用中普遍采用雷諾時均方法對湍流進行模擬，可以取得理想的精度，因此，本文也采用雷諾時均方法對湍流進行數值模擬。

按照不同的封閉形式，雷諾時均方程又可分為單方程模型(Spalart-Allmaras)、兩方程模型(k-ε、k-ω)以及雷諾應力模型，本文主要對比k-ε模型、k-ω模型以及大渦模擬湍流模型對合成射流數值計算的影響。

本文所研究的合成射流激勵器的作動部件為活塞，活塞的運動指定為正弦型：

y(t)=asin(2πft)

式中，y為活塞位置；a為活塞振幅，取5 mm；f為活塞振動頻率，取10 Hz。

為了盡可能接近物理實際，本文采用動網格對活塞運動進行模擬，網格更新規律如下：

式中，u為流體的速度向量；ug為動網格的速度向量；Γ為擴散系數；SΦ為Φ的源相；A為控制體的面積向量；ρ為密度。

本文所研究的激勵器外形如圖1所示，激勵器的幾何參數見表1。本文采用基于動網格的有限體積法進行數值計算[19]，流場邊界條件設置如圖2所示，激勵器底部設置為動邊界，激勵器的其他邊界均設置為無滑移壁面條件。與激勵器出口相對的外流場邊界設置為壓力出口，其余邊界設置為壓力入口，入口壓力與出口壓力為環境參考壓力，其值為零。

圖1 激勵器的結構和參數Fig.1 Structure and parameters of actuator

mm

圖2 邊界條件Fig.2 Boundary conditions

流場設置為黏性、不可壓、非定常，流體設置為水，空間、壓力、能量均采用二階精度離散。為保證求解精度，每個周期劃分為200個時間步，每個時間步長迭代50次，收斂準則為殘差10-5。

為了保證數值方法的可靠性，本文分別研究了不同網格尺度、不同邊界層尺度以及不同湍流模型對激勵器推力以及流場的影響。除了湍流模型對比之外，本文其他計算均默認采用SSTk-w湍流模型。

1.1 網格尺度的影響

本文研究了不同網格尺度下合成射流推力及流場的變化，為了在保證精度的同時減少整體網格數量，在進行數值計算時，整個流域采用4種不同的網格尺度進行劃分：激勵器出口網格尺度s1、出口邊界層網格尺度s3、激勵器腔體網格尺度s4、外流域網格尺度s2，其定義如圖3所示。

圖3 不同流域內的網格尺寸Fig.3 Mesh sizes in different flow fields

對3種不同尺度的網格分別進行了數值計算，不同網格尺度下的推力曲線對比如圖4所示，流場對比如圖5所示，其中，網格尺度定義為s1-s3-s4-s2。推力通過對激勵器所有表面的壓力應力進行積分得到，由圖4及圖5可知，不同網格尺度下推力及流場基本不變，因此，為減小計算量，本文后續計算均采用0.4-0.04-0.8-10的網格尺度。

圖4 不同網格尺寸下的推力Fig.4 Thrust at different mesh sizes

圖5 不同網格尺寸下的渦量云圖Fig.5 Vorticity cloud chart in different mesh sizes

1.2 邊界層尺度的影響

圖6 不同邊界層網格尺寸下的Yplus值Fig.6 Yplus values in different boundary layer mesh sizes

圖7 不同邊界層網格尺寸下的推力Fig.7 Thrust in different boundary layer mesh sizes

圖8 不同邊界層網格尺寸下的渦量云圖Fig.8 Vorticity cloud chart at different boundary layer grid sizes

為了研究邊界層尺度對合成射流激勵器推力及流場的影響，本文對比了4種不同邊界層尺度網格的推力及流場計算結果，4種不同的邊界層所對應的第一層網格的高度h1分別為0.01 mm、0.04 mm、0.08 mm、0.20 mm。4種不同邊界層網格所對應的Yplus值在一個周期內的變化曲線如圖6所示。推力及流場的對比結果如圖7、圖8所示，可以看出，推力及流場對邊界層網格尺度的變化并不敏感。

1.3 不同湍流模型的影響

本文研究了不同湍流模型對合成射流推力及流場的影響，所研究的湍流模型分別是Realizablek-ε模型、標準k-ε模型、SSTk-ω模型、大渦模擬(LES)模型。4種不同湍流模型所對應的推力及流場計算結果如圖9、圖10所示，可以看出，不同湍流模型所得到的推力結果基本一致，但流場細節有差異，大渦模擬與SSTk-ω模型得到的流場非常接近，而Realizablek-ε及標準k-ε計算得到的流場在渦的結構上與前兩者不同。Realizablek-ε及k-ε模型計算得到的渦強度明顯低于SSTk-ω及LES模型的渦強度。

圖9 不同湍流模型下的推力Fig.9 Thrust in different turbulence models

圖10 不同湍流模型下的渦量云圖Fig.10 Vorticity cloud chart in different turbulence models

如果只關注激勵器的推力性質，可以采用標準k-ε模型，該模型對邊界層要求較低，可以減小計算量。如果要研究激勵器的流場性質，就需要采用低雷諾數模型，如SSTk-ω或大渦模擬，由于大渦模擬計算量較大，故本文默認采用SSTk-ω模型。

1.4 與實驗結果對比

為了驗證本文數值計算方法的正確性，將數值計算結果與KRIEG等[13-14]給出的實驗結果進行對比。KRIEG等并沒有給出激勵器的具體幾何參數，只是提到活塞行程為6.35 mm時排出41.5 mL體積的水，根據這一信息推算出激勵器的活塞直徑L=91 mm，再根據激勵器的L/D=4，推算出此時激勵器的出口直徑D=23.6 mm。采用前文的湍流模型以及網格尺度，對該激勵器在活塞振動頻率為11 Hz條件下的推力進行了數值計算，并對得到的推力進行正則化處理，并將結果與KRIEG等的實驗結果進行對比，結果顯示兩者具有很好的一致性，如圖11所示。

圖11 數值計算和實驗數據的對比Fig.11 Comparison between numerical and experimental data

2 推力模型

2.1 合成射流推力

由關于水環境下合成射流的相關研究可知，與穩態射流相比，合成射流具有較高的推進效率[1-3]。為了對合成射流的推力性質尤其是推力增強性質給出合理的解釋，本文建立了合成射流的推力模型。在建立推力模型時，本文將激勵器內部流體作為控制體，如圖12所示。具體建模過程詳見文獻[20]。

圖12 控制體Fig.12 Control volume

2.2 穩態射流推力

本文將具有相同質量流量的穩態射流的推力作為基準值，通過將合成射流推力與穩態射流推力進行對比，來比較不同條件下下合成射流的推進效率。為此，本文建立了穩態射流作用力的數學模型。具體建模過程詳見文獻[20]。

假設激勵器排出的流體可以用一個以激勵器出口為底面的圓柱表示，則圓柱高度

(1)

式中,v為流經激勵器出口的流體的平均速度；t為時間;a為活塞的位移;k為出口直徑與活塞面積之比；T為激勵器振動周期。

簡化得到，凈動量變化

ΔH=ρπ3D2fl2/16

式中，l為激勵器一個行程排出的流體柱高度；D為激勵器出口直徑；ρ為密度；f為激勵器活塞振動頻率。

則相同質量流量穩態射流在一個周期內的平均作用力

(2)

KRIEG等[12-13]將式(6)作為合成射流的平均推力計算公式，由本文的分析可知，該推力為相同質量流量下穩態射流的推力，合成射流的平均推力比此值偏高。

3 L/D對推進效率的影響

合成射流推進效率易受各種因素影響，其中一個重要因素就是L/D。量綱一參數L/D可用作合成射流效率的一個量度，當L/D值為4左右時，合成射流的效率最高；當L/D值偏離4時，合成射流的效率降低[4,7-8]。

研究圖1所示的合成射流激勵器在不同L/D時的推力變化，為了便于比較，將推力量綱一化：

不同L/D下的量綱一推力變化如圖13所示，圖13中的基準線為f2，通過比較數據點與基準線的高低，即可知合成射流推力與穩態射流推力的相對大小。由圖13可以看出，L/D為4.3或5.8時，合成射流的推力位于基準線之上，其他L/D條件下推力均位于基準線之下。這說明，L/D為4.3或5.8時合成射流的效率較高，推力大于穩態射流的推力；而L/D<4或L/D>6時合成射流效率較低，推力小于穩態射流的推力。

圖13 不同行程比下的平均推力Fig.13 Average thrust at different stroke ratios

對合成射流的研究均表明，合成射流的推進效率受L/D影響較大，普遍認為，當L/D在4附近時，合成射流推進效率最高[4,7-8]，此時合成射流的相對推力大于穩態射流，這與本文得到的結果一致。隨著L/D遠離最優值，合成射流的效率降低，其推力小于相同條件下穩態射流的推力，這也與本文得到的結果一致。

KRIEG等[13]也給出了不同L/D條件下的推力變化，但給出的結果中，隨著L/D的增大，推力呈現單調遞減的趨勢，這與本文得到的結果不一致。且當L/D>9時，增大頻率對推力沒有影響，KRIEG等給出的解釋是空化的影響，本文數值計算時沒有考慮空化，這可能是本文結果與KRIEG等的結果存在差異的原因。

4 機理闡釋

圖14 噴嘴可偏轉的激勵器Fig.14 The actuator with nozzle deflection

本文對圖14所示的激勵器在偏角為0°及90°條件下的推力進行了數值計算。激勵器主腔體的回轉軸與y軸方向平行，因而，當出口不偏轉時，射流的噴射方向為y軸方向；當出口偏轉90°時，射流的噴射方向為x軸方向。兩種條件下，激勵器x、y軸方向的推力曲線如圖15所示。

(a)偏角為0°

(b)偏角為90°圖15 噴嘴在偏角為0°和90°時x和y方向上的推力Fig.15 Nozzles at 0 and 90 degrees, thrust in the x and y directions

由圖15可以看出，當出口不偏轉時，y向推力不對稱度較大，而x向推力可以忽略；當出口偏轉90°時，y向推力接近正負對稱，此時x向推力不對稱度較高。曲線不對稱程度越高，說明一個周期的平均推力越大，對激勵器的推力增強性質貢獻越大。

兩種條件下x、y方向的平均推力見表2，可以看出，當激勵器不偏轉時，推力主要集中在y向，此時，推力增強也集中在y向。而當激勵器偏轉90°時，推力主要集中在x向，因而，推力增強也主要集中在x向。

表2 噴嘴偏角在0°和90°時x和y方向上的平均推力Tab.2 Average thrust under nozzle deflection angle of 0 and 90 degree

|f|引起合成射流推力增強的機理解釋如下：當激勵器噴出流體時，激勵器外部流體處于相對靜止狀態，此時激勵器噴出的射流相對于外部靜止流體具有較大的相對速度，因此，噴射流體會受到外部流體的阻滯作用，將此噴射流體柱視為剛性，則其受力如圖16所示，可以看出，外部流體對控制體的作用力方向與激勵器對控制體的力的方向相反，則有

圖16 激勵器噴出時的受力圖Fig.16 Force of the actuator ejected

當流體流入激勵器時，吸入流體速度相對于激勵器內部流體具有一定的相對速度，因此，吸入流體會受到激勵器內部流體的阻滯作用，將吸入流體柱視為剛性，則其受力如圖17所示，可以看出，外部流體對于控制體流體的作用力的方向與激勵器對控制體的作用力同向，則有

圖17 激勵器吸入時的受力圖Fig.17 Force of the actuator inhalation

由此可以看出，當激勵器噴出流體時，激勵器受到的作用力比基準值偏大，而當激勵器吸入流體時，激勵器受到的作用力在數值上比基準值偏小，因而激勵器作用力曲線整體上有一個向上平移的趨勢。

5 結論

(1)網格尺度及邊界層尺度對推力及流場的計算結果影響較小。

(2)不同湍流模型所得到的推力結果基本一致，研究激勵器的推力特性可以采用k-ε模型，研究激勵器流場特性可以采用SSTk-ω模型。

(3)通過對比分析激勵器不同驅動條件下的推力變化可發現，當L/D接近4時，合成射流的推進效率最高；當L/D遠離4時，推進效率開始降低。

(4)本文根據推力的數學模型，并通過分析揭示了合成射流的推力增強性質來源于外部流體的阻滯作用。

[1] KRUEGER P S. The Significance of Vortex Ring Formation and Nozzle Exit Over-pressure to Pulsatile Jet Propulsion[D]. Pasadena：California Institute of Technology,2001.

[2] MOSLEMI A A, KRUEGER P S. Propulsive Efficiency of a Biomorphic Pulsed-jet Underwater Vehicle[J]. Bioinspiration&Biomimetics,2010,5(3):036003-14.

[3] KRUEGER P S, GHARIB M.Thrust Augmentation and Vortex Ring Evolution in a Fully Pulsed Jet[J]. AIAA Journal,2005，43(4):792-801.

[4] WHITTLESEY R W, DABIRI J O. Optimal Vortex Formation in a Self-propelled Vehicle[J] J. Fluid Mech.,2013,737(4)：78-104.

[5] WHITTLESEY R W, DABIRI J O. Vortex-enhanced Propulsion[J]. Journal of Fluid Mechanics,2011,668(4):5-32.

[6] KRUEGER P S, GHARIB M. The Significance of Vortex Ring Formation to the Impulse and Thrust of a Starting Jet[J]. Physics of Fluids,2003,15(5):1271-1281.

[7] GHARIB M, RAMBOD E, SHARIFF K. A Universal Time Scale for Vortex Ring Formation[J]. Journal of Fluid Mechanics,1998,360:121-140.

[8] KRUEGER P S, DABIRI J O, GHARIB M. Vortex Ring Pinchoff in the Presence of Simultaneously Initiated Uniform Background Co-flow[J]. Physics of Fluids,2003,15(7):L49-L52.

[9] MOHSENI K. Pulsatile Vortex Generators for Low-speed Maneuvering of Small Underwater Vehicles[J].Ocean Engineering,2006,33(16):2209-2223.

[10] THOMAS A P, MILANO M, G’SELL M G, et al. Synthetic Jet Propulsion for Small Underwater Vehicles[C]//IEEE International Conference on Robotics and Automation.Barcelona,2005:181-187.

[11] KRIEG M, PITTY A, SALEHI M, et al. Optimal Thrust Characteristics of a Synthetic Jet Actuator for Application in Low Speed Maneuvering of Underwater Vehicles[C]//Proceedings of the 2005 IEEE International Conference on Robotics and Automation. Barcelona,2005:2342-2347.

[12] KRIEG M, MOHSENI K. Dynamic Modeling and Control of Biologically Inspired Vortex Ring Thrusters for Underwater Robot Locomotion[J].IEEE Transactions on Robotics,2010,26(3):542-554.

[13] KRIEG M, MOHSENI K.Thrust Characterization of a Bio-inspired Vortex Ring Thruster for Locomotion of Underwater Robots[J].IEEE Journal of Oceanic Engineering,2008,33,(2):123-132.

[14] KRIEG M, COLEY C, HART C, et al. Synthetic Jet Thrust Optimization for Application in Underwater Vehicles[C]//Proc. 14th Int. Symp. on Unmanned Untethered Submersible Technology (UUST).Durham,2005:1-10.

[15] THOMAS A M P. Exploration into the Feasibility of Underwater Synthetic Jet Propulsion[D]. Pasadena：California Institute of Technology,2007.

[16] WHITTLESEY R W. Dynamics and Scaling of Self-Excited Passive Vortex Generators for Underwater Propulsion[D]. Pasadena：California Institute of Technology,2013.

[17] THOMAS A M, ABRAHAM J P. Numerical Simulation of Circular Synthetic Jets with Asymmetric Forcing Profiles[J].The Open Mechanical Engineering Journal,2010,4(1):1-7.

[18] KUMAR A, PANDA P K, KUMAR V, et al. Combined Experimental and Numerical Study of Synthetic Jet in Quiescent Flow[C]//37th National and 4th International Conference on Fluid Mechanics and Fluid Power. Chennai,2010:1-11.

[19] 張瑞，姜峰，楊晉，等.基于動網格的液壓缸雙向流固耦合分析[J].中國機械工程,2017,28(2):156-162.

ZHANG Rui, JIANG Feng, YANG Jin，et al. Fluid-structure Coupling Analysis of Hydraulic Cylinder Based on Dynamic Mesh Technology[J]. China Mechanical Engineering,2017,28(2):156-162.

[20] 耿令波, 胡志強, 林揚,等.不對稱激勵函數對水環境下合成射流激勵器推力影響及其機理研究[J].機械工程學報,2017,57(3):13-22.

GENG Lingbo, HU Zhiqiang, LIN Yang, et al. Thrust Characteristics of Synthetic Jet Actuator Underwater with Asymmetric Forcing Profile and Its Working Mechanism [J]. Journal of Mechanical Engineering,2017,57(3):13-22.