水下航行器阻力參數化分析與優化

劉 峰，梁 旭，苗怡然，屠超華，趙彥凱

(1.哈爾濱工程大學 船舶工程學院，黑龍江 哈爾濱，150001；2.中國船舶工業系統工程研究院，北京，100094；3.中國核動力研究設計院 核反應堆系統設計技術國家重點實驗室，四川 成都，610041)

0 引言

水下航行器是人類用于海洋開發的重要裝備之一，在軍用、民用等領域均有著廣泛的應用與需求[1]。隨著各國開發海洋力度的加大和步伐的加快，人們對于水下航行器的總體性能提出了更高的要求。良好的阻力性能可使水下航行器具備更高的航速、更強的續航能力等，是水下航行器的重要性能指標之一[2]。計算流體力學(computational fluid dynamics，CFD)技術因其精度高[3]、成本低等優點，被廣泛應用于水下航行器阻力計算領域[4-5]。然而，這些阻力計算大多都在一定的范圍即設計空間內選擇參數進行水動力計算，由于比較計算結果需要反復修改參數，因此大量重復性工作導致了設計效率低下[6]。

參數化分析通過對數字參數的修改，實現了模型的重建和自動分析，可提高計算與分析效率。例如，楊卓懿等[6]針對水下航行器進行了格蘭韋爾艇型參數化分析，基于全局探索遺傳算法進行了阻力優化；龐永杰等[7]將二維和三維2 種網格形式計算的精度和效率進行了對比，研究了單獨改變首部或尾部形狀時阻力的變化規律，實現了Myring 型回轉體參數建模與分析。

然而，若直接應用參數化分析進行水下航行器阻力優化，依然需要耗費大量時間。為彌補這一不足，近似模型法采用數學模型對于物理模型近似擬合的方法，在保證模型精度的同時，進一步提高了設計效率，該方法已在阻力分析中得到了廣泛應用[8]。

文中以一型水下航行器為研究對象，基于Java 語言對STAR-CCM+軟件進行了二次開發，實現了水下航行器阻力的參數化建模與分析；采用最優拉丁超立方方法進行了樣本點的選取與分析，并針對設計參數對該水下航行器阻力的影響展開分析。為簡化分析過程，提高阻力分析效率，利用響應面對于樣本點進行了擬合，得到了具有較高擬合精度的阻力近似模型，進一步建立了以直航阻力最小為優化目標、排水體積為約束的水下航行器阻力優化模型，并進行了優化求解，取得了明顯的優化結果。

1 阻力分析模型

1.1 阻力求解方法

利用雷諾平均Navier-Stokes(Reynolds average Navier-Stokes，RANS)方程求解阻力，其在時均化的N-S 方程中，以某種模型的方式將瞬態脈動量體現出來，雷諾時均化后的不可壓連續性方程和RANS 方程分別為

式中：ui、uj為速度分量時均值；p為壓力時均值；ρ為流體密度；μ為動力粘性系數為雷諾應力項。

式中：Gk和Gω為由平均速度梯度引起的湍流動能；Yk和Yω為關于k和ω的湍流耗散項。

1.2 基于STAR-CCM+軟件的阻力分析

利用STAR-CCM+軟件進行阻力計算，通過設計參數的設置可實現幾何模型的參數化。水下航行器外形見圖1。

圖1 水下航行器外形圖Fig.1 Outside view of undersea vehicle

圖中航行器主要參數包括長度x1、半徑x2和尾部長度x3，外形為對稱模型，以中縱剖面為界計算1/2 模型在水中運動的阻力。將原點設置于模型的尾部(尾部斷點)，計算區域前后均為20 m，上、下、側面均為5 m，設置網格基礎尺寸為0.2 m，表面曲率120。在模型首部與尾部設置網格加密，為實現參數化建模與分析，采用切割體網格生成器與棱柱層網格生成器，并啟用表面重構與自動表面修復功能。加密域尺度設置如表1 所示。

表1 加密域尺度設置Table 1 Encryption domain scale settings

將水下航行器置于恒密度靜水中，計算水深為400 m，得到計算壓強為4 MPa，速度為30 m/s，采用湍流模型，湍流強度為0.01，湍流粘度比為10。設置計算域前、側、上及下平面均為速度進口，模型中縱剖面所在平面為對稱平面，計算域后平面為壓力出口，計算域邊界條件設置如圖2 所示。

圖2 計算域網格與邊界條件設置示意圖Fig.2 Calculation domain grid and boundary conditions

利用STAR-CCM+軟件對模型的1/2 進行計算得到阻力為6 873.02 N，換算得到的整體阻力為13 746.04 N。與文獻[10]中同等初始條件下的阻力14 508.20 N 相對比，計算值誤差為5.2%，驗證了計算結果的準確性。

2 阻力參數化分析

利用STAR-CCM+軟件Java 宏錄制功能，對Java 宏文件進行二次編譯。Java 宏文件準備完畢后通過*.bat 批處理文件直接運行Java 宏文件，系統調用cmd.exe 并運行批處理文件中的各條命令。水下航行器阻力參數化分析流程如圖3 所示。

圖3 阻力參數化分析流程圖Fig.3 Parameterized analysis process of resistance

圖中，通過編寫*.bat 文件調用STAR-CCM+，初始設置以*.sim 格式保留在文件中，Java 宏文件進行外形參數設置的保留。分析過程所設計的控制參數和計算結果等借助于Isight 軟件進行后處理，將計算所生成的結果文件存放于相應的工作目錄下。樣本點的生成方法采用最優拉丁超立方方法，得到81 組樣本并進行直航阻力參數化分析，樣本點及所得到的響應值見表2。表中，Rx為水下航行器直航阻力。

表2 阻力樣本點Table 2 Sample points of resistance

圖4 為設計變量靈敏度分析直方圖。由圖可知，對Rx影響最大的是x1，其次是x2，阻力隨這2 個設計變量的增大而增大，x3對結果的影響相對較小。

圖5 為設計變量主效應分析曲線圖。由圖可知，在計算區間內x1對結果的影響基本呈線性關系；x2對結果的影響呈一定非線性關系；x3對結果的影響呈明顯非線性關系。

圖4 設計變量靈敏度分析直方圖Fig.4 Histogram of design variables sensitivity analysis

圖5 設計變量主效應分析曲線Fig.5 Curves of design variables main effect analysis

3 近似模型建立

響應面模型(response surface method，RSM)將試驗設計與數理統計學分析相結合，是一種構建設計變量與響應之間映射關系的方法。其本質是進行多項式求解，具有透明性高、簡單、高效等優點。RSM 按照階數可分為1 階、2 階、3 階和4 階。其中，2 階RSM 的表達式為[11]

式中：為響應面近似值；xi為設計變量；n為設計變量個數；a0，bi，cii和di分別為常數項、一次、二次和交叉項的待定系數。

近似模型擬合精度利用復相關系數法的R2進行分析，即

式中：n為樣本數；yi為真實響應；為響應面擬合值；為真實響應均值。

R2是介于0～1 之間實數，越接近1 說明擬合效果越好，工程上要求其數值在0.9 以上即可初步滿足擬合要求。采用響應面模型進行擬合，近似模型的R2為0.97，滿足工程需要。抽取其中10個預測值(見圖中三角形)與計算值進行對比，如圖6 所示。

圖6 水下航行器直航阻力預測值與計算值對比Fig.6 Comparison of forecast and calculated values of straight running resistance for undersea vehicle

圖6 中，10 個預測值與計算值均非常接近。最終得到近似模型

4 阻力優化

4.1 優化方法

非線性二次規劃法(nonlinear programming quadratic line search，NLPQL)以2 階泰勒級數對目標函數進行展開，同時線性化約束條件，通過求解二次規劃來獲得下一個設計點，再根據2 個可供選擇的優化值進行一次線性搜索，是一種基于順序二次規劃法(sequence quadratic program，SQP)理論的非線性算法。其實質是對式(8)的拉格朗日二次函數的二次規劃子問題近似求解過程。

4.2 優化模型建立及求解

確定直航阻力Rx為水下航行器外形優化的目標函數，即以Rx最小為目標，以排水體積V為約束函數，建立直航阻力Rx優化模型

式中，DV表示設計變量。部分參量單位分別為：Rx為N；V為m3；x1、x2、x3為m。

采用NLPQL，設置迭代步數為200，收斂精度為10-6，相對步長為10-4。迭代16 步時得到優化結果，迭代過程如圖7 所示，圖中直線可以看出求解趨勢。

圖7 直航阻力優化求解過程Fig.7 Optimization solution process of straight running

將優化方案利用STAR-CCM+軟件重新計算，以驗證近似模型計算的準確性，表3 為初始方案、近似模型優化和CFD 仿真優化結果對比。

表3 初始方案與優化結果對比Table 3 Comparison of initial scheme and optimization results

表3 中，相比初始方案，優化方案的直航阻力降低了10.78%，3 個設計變量中x1的變化幅度最大，與設計變量靈敏度分析結果相對應，近似模型計算結果與CFD 仿真計算結果相差1.01%。

5 結論

文中開展了水下航行器阻力參數化分析及優化研究，得出結論如下：

1) 對于STAR-CCM+軟件進行了二次開發，通過Isight 對于STAR-CCM+進行調用，實現了水下航行器阻力的參數化建模與分析，可避免阻力分析模型的頻繁建立與修改，提高了水下航行器阻力分析效率；

2) 針對設計參數對于響應的影響程度進行了分析，采用2 階響應面模型對樣本點進行擬合，得到了精度滿足工程需要的近似模型，為水下航行器阻力分析效率的提高提供參考；

3) 利用非線性二次規劃法進行優化求解，直航阻力降低了10.78%，提高了該型水下航行器的阻力性能。同時設計變量中中段長度的變化幅度最大，說明在該水下航行器的設計過程中，應對該設計變量給予重視。

后續需開展優化后的阻力測試，以驗證文中方法的準確性，還需開展復雜艇型和填加附體的阻力參數化分析工作。