李光釗，常 欣，鄧南彥，于鵬垚

(大連海事大學 船舶與海洋工程學院，遼寧 大連 116026)

0 引 言

在自然界中，蜻蜓、蜜蜂等昆蟲的撲翼結構是最常見的串列撲翼系統(tǒng)，學者們從中獲得靈感，采用多翼串列的方式來提升撲翼系統(tǒng)的推進性能。于憲釗等分析了前翼拍動后翼靜止的串列翼系統(tǒng)水動力；Jones等分析了串列水翼的推進性能并進行了相關試驗研究；Broering等采用數值模擬方法研究了串列撲翼的氣動力性能；石成明等采用氣動力試驗方法研究了串列翼型所產生的推力與前后翼距離之間的變化規(guī)律；李永澤等研究了3種串列雙翼的氣動布局，分析了最優(yōu)的前后翼布局方案。許國冬等針對典型推進運動模式與能量采集模式，分析了串列雙翼之間渦的干擾規(guī)律。

盡管國內外針對串列翼的氣動力和水動力性能已開展許多研究，但過往研究中的串列翼運動形式大多為正弦運動。自然界生物的撲翼運動軌跡往往是更為復雜的非正弦運動，目前尚未看到針對非正弦運動串列水翼推進表現的研究。在串列撲翼式水下航行器的研發(fā)驅動下，本文采用計算流體力學方法研究非正弦升沉運動下串列水翼的推進性能，為水下航行器的撲翼運動模式選擇提供指導。

1 計算方法

1.1 撲翼運動的描述

本文中非正弦運動通過橢圓形軌跡實現，其中，撲翼的升沉運動位移方程定義如下：

其中：為撲翼升沉運動的非正弦參數；為升沉運動幅值； ω=2π，為拍動頻率；為時間。=1時，撲翼運動即為正弦運動，不同值對應的撲翼運動軌跡如圖1所示。

圖1 不同S值的撲翼運動軌跡Fig. 1 Flapping trajectories according to different values of S

當兩撲翼前后串列排布時，前翼和后翼的升沉運動可分別表示為：

其中：()為串列排布兩翼中前翼的運動方程，()為后翼的運動方程；為撲翼升沉運動非正弦參數；為前后兩翼升沉運動的相位差。

圖2 串列水翼布置圖Fig. 2 Layout of in-line tandem hydrofoils

1.2 物理量的定義

撲翼在運動時會受到不同方向上的水動力和水動力矩。本文對串列翼系統(tǒng)的前后水翼力學特性分別討論，作用在單獨進行升沉運動水翼上的瞬時推力系數C和瞬時升力系數C可表達為：

其中：F()為作用在單獨水翼上沿軸的瞬時推力；F()為作用在單獨水翼上沿軸的瞬時升力；為流體的密度；為流場來流速度。

則平均推力系數C和平均功率系數C可表達為：

其中：為撲翼升沉運動周期；P為撲翼運動時所需功率。

對于升沉運動水翼，一個運動周期所需的功率可以通過下式計算：

其中：為水翼沿軸軸向速度。

單獨水翼的推進效率 η可表達為：

1.3 數值方法與計算模型

本文采用基于有限體積法的商用CFD求解器STAR-CCM+進行串列水翼的水動力分析。其中，水翼的運動利用重疊網格技術來實現。數值模型的計算域如圖3所示，其中，分析的水翼為NACA0012標準翼型，翼弦長=0.1 m。方形區(qū)域為背景域，尺寸設置為70×70；中間圓形區(qū)域為重疊域，半徑為1。根據Chao等的研究，當撲翼和計算域邊界的距離大于20時，流域尺寸對水翼水動力的影響可以忽略不計，因此，可以判斷本文中的計算域尺寸已足夠大。具體的計算域邊界條件如圖3所示。

圖3 串列雙翼計算域及邊界條件Fig. 3 The domain and boundary conditions of in-line tandem hydrofoils

1.4 數值收斂性分析

性分析，具體計算參數取值為=0.7，=1×10，=0.35 ，=0.25 ，相 位 差 ε=0 。其 中，為雷諾數，其定義為：

針對數值模型中的時間步長和網格尺寸進行收斂

式中：為來流速度，為水的動力粘度。為斯特勞哈爾數，反映流體非定常運動的相似性，其定義為：

式中，為運動頻率。

圖4 不同時間步長下瞬時推力系數時歷曲線Fig. 4 Time history of instantaneous thrust coefficient of in-line tandem flapping foil under different time steps

圖5 不同網格下的瞬時推力系數時歷曲線Fig. 5 Time history of instantaneous thrust coefficient of in-line tandem flapping foil under different grids

圖4為不同時間步長下前后水翼推力系數時歷曲線的比較。計算中選取的時間步長分別為d=/1 000，/2 500，/4 000，對應計算模型的網格數目為24.5萬?？梢钥闯?，不同時間步長下水翼的推力曲線吻合較好，故時間步長d=/2 500用于后續(xù)計算。圖5為不同網格密度下前后水翼推力系數時歷曲線的比較。對應3個不同的網格密度，撲翼表面分別具有200，400，600個節(jié)點，對應計算模型的網格數目分別為6.8萬、24.5萬和54.1萬，計算中采用的時間步長d=/2 500?？梢钥闯?，不同時間步長下水翼的推力曲線總體吻合較好，表明當前數值模型具有較好的網格收斂性，在后續(xù)計算中選取網格數目為24.5萬的計算模型。

1.5 數值算法的驗證

圖6為本文數值結果與Heathcote等試驗結果以及Y.L.Wang等和Young等數值結果的對比。其中，計算翼型為NACA0012水翼，弦長=0.1 m，雷諾數=2×10，運動幅值=0.175，水翼運動形式為單翼的正弦純升沉運動。橫坐標為無量綱所減頻率，定義為=2π/?？梢钥闯?，本文數值結果與文獻[9]試驗結果吻合較好。圖7為本文數值結果與Kinsey等結果的對比。其中，文獻中的計算結果也被許國冬等采用，用于驗證數值模型。其中，計算翼型為NACA0015水翼，弦長=0.1 m ，雷諾數=2×10，運動幅值=1.0，水翼運動形式為串列水翼的俯仰升沉復合運動?？梢钥闯?，本文數值結果與文獻數值結果吻合較好。總體而言，當前數值方法在本文中的應用是合理可行的。

圖6 單翼升沉運動驗證Fig. 6 Plunge motion verification of single hydrofoil

圖7 串列雙翼復合運動驗證Fig. 7 Combined motion verification of in-line tandem hydrofoils

2 結果與分析

主要討論非正弦運動對串列升沉水翼推進性能的影響，其中無特殊說明，計算參數=0.7，=1×10，=0.35，L=7.5，保持不變。

2.1 非正弦升沉系數S對瞬時水動力的影響

采用前后翼相位差 ε=90的結果進行分析。如圖8（a）所示，對于前翼的瞬時推力系數曲線：從峰值來看，非正弦系數從=1增大或減小，都會使前翼瞬時推力系數曲線峰值增加；從周期來看，增大非正弦系數時周期基本不變，而減小非正弦系數時原有周期發(fā)生變化，峰值出現位置由水翼運動的平衡位置（0.25、0.75）變化為水翼運動的極端位置（0、0.5、1）。如圖8（b）所示，對于后翼的瞬時推力系數曲線：從峰值來看，隨著非正弦系數的減小，后翼推力系數曲線峰值增大，隨著非正弦系數的增大，峰值先減后增；從周期來看，由于相位差的存在，后翼周期應比前翼推遲0.25，而當=1時，峰值卻提前出現，這是前翼的尾流作用所致。

圖8 不同非正弦系數下串列撲翼的瞬時推力系數時歷Fig. 8 Time history of instantaneous thrust coefficients of in-line tandem flapping foils under different non-sinusoidal coefficient

圖9 不同非正弦系數下串列撲翼的瞬時升力系數時歷Fig. 9 Time history of instantaneous lift coefficients of in-line tandem flapping foils under different non-sinusoidal coefficient

如圖9（a）所示，對于前翼的升力系數曲線：從峰值和谷值來看，從=1增大或減少非正弦系數均可增大或減小升力系數曲線的峰值或谷值，且從=1減小時相比于增大時峰值增大或谷值減小的程度更大；從周期來看，當＞1時周期基本不變，而當＜1時，峰值或谷值出現位置由水翼運動的平衡位置（0.25、0.75）變化為水翼運動的極端位置（0、0.5、1）。如圖9（b）所示，對于后翼的瞬時升力系數曲線：從峰值和谷值來看，從=1增大或減小非正弦系數同樣會使升力系數曲線峰值更大，谷值更??；從周期來看，由于相位差的存在，后翼的周期性變化比前翼慢0.25。

2.2 非正弦升沉系數S對渦結構的影響

圖10和圖11分別為前后兩翼在=0.5時刻不同非正弦升沉運動時周圍流場的渦量云圖。此時前翼處于最大負向位移位置，后翼處于向下運動的平衡位置。將水翼的上表面稱為吸力面，下表面稱為壓力面，水翼形成升力的原因是吸力面與壓力面的壓強差。

前翼在非正弦系數=1時，一組尾渦由2個正向旋轉、1個負向旋轉的渦組成。在脫落過程中，正向旋轉的渦迅速變弱至消失。水翼首部有一對旋向相反的渦，且?guī)е粋€未與水翼表面接觸的負向旋渦。當這組渦向水翼尾部運動時，渦強度逐漸減弱。由于上表面的渦弱，上下表面壓力差不大，升力無明顯變化。當＜1時，尾渦組成與=1時相似，當減小時，尾渦中一個正向旋轉的渦強度增強，推力明顯增大。當＞1時，一組尾渦由2對旋向相反的渦組成，且脫落后負向旋渦迅速減弱，當增大時，上表面脫落的負向旋渦強度增大。

圖10 前翼運動至最大負向位移位置處(t=0.5T)的渦量云圖Fig. 10 The vorticity contour of the forefoil moving to the position of maximum negative displacement (t=0.5T)

圖11 后翼運動至平衡位置處（t=0.5T）的渦量云圖Fig. 11 The vorticity contour of the hindfoil moving to the equilibrium position (t=0.5T)

對于后翼在非正弦系數=1時，有一個負向旋渦在首端脫落，有一個正向旋渦在尾端脫落。負向旋渦比正向旋渦強度更強，有一組旋向相反的渦在尾流中，觀察到有前翼的尾渦運動到后翼附近。當＜1時，首端負向旋渦分為2個，尾端正向旋渦和尾渦中的負向旋渦隨的減小而減弱。當＞1時，隨著的增大，尾端的正向旋渦強度增強，首端負向旋渦強度減弱。

2.3 非正弦升沉系數S對推進表現的影響

從圖12（a）可知，對于前翼，從非正弦系數=1增大或減少，均能增加前翼的平均推力系數。在=0.5時前翼平均推力系數在相位差為90°處變化明顯，在此非正弦系數下，后翼對前翼的推力影響較大。從圖12（b）可知：在大多數相位差下，從=1增大或減小非正弦系數可以增強后翼平均推力系數，但在某些相位差下，非正弦運動的平均推力系數小于正弦運動的平均推力系數。如相位差為0時，=0.5的平均推力系數小于=1的平均推力系數?？傮w來看，隨著相位差的增加，平均推力系數先增后減。

從圖12（c）可知，對于前翼的推進效率，當＞1時，推進效率隨的增大而減小，當＜1，推進效率隨的減小而增大?？梢杂^察到在=0.5時，前翼推進效率曲線在不同的相位差下，數值發(fā)生較明顯變化，推測是此非正弦系數和相位差時后翼對前翼的影響較大所致。從圖12（d）可知，對于后翼的推進效率同樣在＞1時，推進效率隨的增大而減小，隨相位差的增大先減后增，在＜1時，隨的減小而增大，隨相位差的增大先增后減。

圖12 前后翼在不同相位差下，撲翼推進表現隨非正弦系數的變化Fig. 12 Under different phase differences, the propulsion performance of the fore and hind foils changes with the non-sinusoidal coefficient

3 結 語

本文采用計算流體力學方法和重疊網格技術研究了非正弦運動對串列升沉水翼推進性能的影響，得到如下結論：

1）通過不同時間步長和網格尺寸的數值結果比較，驗證了數值模型的收斂性；通過與公開發(fā)表試驗和數值結果的比較，驗證了數值模型的計算精度。

2）非正弦升沉運動相比于正弦升沉運動會增大撲翼的瞬時推力峰值和升力峰值，而且，非正弦系數的量值與1的偏差越大，瞬時峰值增大程度越明顯。

3）非正弦系數可以改變水翼表面渦的脫落模式，當＜1時，一組尾渦由2個旋向相同、一個旋向相反的渦組成；當＞1時，一組尾渦由2對旋向相反的渦組成。

4）除個別少數相位差對應的運動狀態(tài)，與正弦運動相比，當非正弦系數從=1增大或減小時，非正弦運動下的前翼和后翼平均推力增大。

5）除個別少數相位差對應的運動狀態(tài)，與正弦運動相比，當非正弦系數從=1增大時，前翼和后翼的推進效率減?。划敺钦蚁禂祻?1減小時，前翼和后翼的推進效率增大。除此之外，在此基礎上改變雙翼相位差可以進一步使串列水翼達到最佳推進效率，這將為串列雙翼水下航行器推進系統(tǒng)的設計提供指導。