蔡昊鵬，馬騁，陳科，錢正芳，2，張赫

確定螺旋槳側(cè)斜分布的一種數(shù)值優(yōu)化方法

蔡昊鵬1，馬騁1，陳科1，錢正芳1，2，張赫3

（1海軍裝備研究院，北京100161；2海軍工程大學(xué)，武漢430033；3中國艦船研究院，北京100192）

合理的側(cè)斜分布可以使螺旋槳各半徑剖面不同時(shí)進(jìn)入船艉非均勻流場的高伴流區(qū)域，可有效降低螺旋槳在船艉非均勻流場運(yùn)轉(zhuǎn)所引發(fā)的振動(dòng)噪聲。文章利用基于螺旋槳面元法的非定常水動(dòng)力性能預(yù)報(bào)工具配合智能尋優(yōu)算法，對螺旋槳側(cè)斜分布問題進(jìn)行迭代優(yōu)化設(shè)計(jì)。側(cè)斜分布曲線采用B樣條曲線表達(dá)，可保證曲線在優(yōu)化過程中的光順性；應(yīng)用改進(jìn)粒子群算法，結(jié)合螺旋槳非定常水動(dòng)力性能預(yù)報(bào)工具，以降低軸向非定常力幅值為目標(biāo)對側(cè)斜分布進(jìn)行了優(yōu)化及相關(guān)的分析探討；編制了基于OpenMP語言的并行優(yōu)化算法以在多核高性能計(jì)算機(jī)上提高運(yùn)算效率。數(shù)值優(yōu)化結(jié)果表明，采用文中方法優(yōu)化側(cè)斜分布，可降低螺旋槳的非定常力幅值，該文的側(cè)斜優(yōu)化方法是可行、有效的，并行算法的引入可明顯提高優(yōu)化效率。

側(cè)斜分布優(yōu)化；改進(jìn)粒子群算法；非定常面元法；并行計(jì)算

1 引言

船舶螺旋槳工作在船后非均勻的三向伴流場中，螺旋槳旋轉(zhuǎn)一周內(nèi)，槳葉各半徑剖面的來流攻角隨著當(dāng)?shù)匕榱鞯牟煌鴷r(shí)刻變化，導(dǎo)致槳葉上承受周期性變化的力，這種螺旋槳葉上的非定常載荷會(huì)使螺旋槳產(chǎn)生激振力；該激振力可分為兩類：非定常載荷通過軸系傳遞到船體的所謂軸承力和螺旋槳誘導(dǎo)的壓力場通過水傳遞至船體表面的所謂表面力。大量理論與試驗(yàn)研究發(fā)現(xiàn)，良好的側(cè)斜分布可以使螺旋槳葉片各半徑剖面不同時(shí)進(jìn)入船后不均勻流場的高伴流區(qū)域，從而有效降低螺旋槳工作在非均勻流場產(chǎn)生的激振力，不僅有利于降低螺旋槳本身的直發(fā)聲，而且可減小螺旋槳激振力所引發(fā)的船艉振動(dòng)及輻射噪聲。

螺旋槳側(cè)斜分布的選擇要結(jié)合船后伴流場的各階諧調(diào)分量一起考慮，使側(cè)斜分布與伴流場相互配合，得到較好的減振、降噪效果。目前，確定螺旋槳側(cè)斜分布的主要步驟如下：在螺旋槳設(shè)計(jì)階段，設(shè)計(jì)者通過對船后伴流進(jìn)行協(xié)調(diào)分析，初步確定螺旋槳的側(cè)斜分布，然后通過理論或試驗(yàn)研究分析，多次調(diào)整側(cè)斜分布至最終側(cè)斜分布方案。該“修改—驗(yàn)證—修改”的人工嘗試方法得到的側(cè)斜分布設(shè)計(jì)方案可能在目標(biāo)流場中達(dá)到了較佳的減振、降噪效果，但不一定是最佳的效果；而且該方法較為繁瑣、需要大量經(jīng)驗(yàn)，會(huì)給螺旋槳設(shè)計(jì)工作帶來較大不便。

本文采用一種基于部分變異策略的粒子群優(yōu)化（PSO algorithm with partial mutation strategy，PMPSO）算法，結(jié)合基于面元法理論的螺旋槳非定常水動(dòng)力性能預(yù)報(bào)方法，以降低軸向非定常力幅值為目標(biāo)，在給定伴流場中對螺旋槳側(cè)斜分布進(jìn)行了數(shù)值優(yōu)化計(jì)算及相關(guān)的分析探討；并進(jìn)一步編制了基于OpenMP語言的并行算法，以在多核高性能計(jì)算機(jī)上提高本文側(cè)斜分布優(yōu)化方法的效率。

2 基本理論

2.1 基于變異策略的粒子群優(yōu)化算法

粒子群優(yōu)化算法（particle swarm optimization，簡稱PSO）是Knenedy和Eberhart發(fā)展的源于群智能和人類認(rèn)知學(xué)習(xí)過程的一種智能優(yōu)化算法[1]。該算法源于鳥類捕食行為的模擬。同遺傳算法類似，PSO算法首先初始化一群隨機(jī)粒子，每個(gè)粒子都代表著優(yōu)化問題的一個(gè)可能解。粒子有自己的位置和速度（決定飛行的方向和距離），位置對應(yīng)的目標(biāo)函數(shù)值作為粒子的適應(yīng)度。在每次迭代中，各個(gè)粒子記憶、追隨當(dāng)前最優(yōu)粒子在解空間中搜索，通過跟蹤兩個(gè)“極值”來更新自己，一個(gè)是粒子本身找到的最優(yōu)解，即個(gè)體極值；另一個(gè)是整體種群的最優(yōu)解，稱為全局極值。找到這兩個(gè)最優(yōu)值后，粒子按照某種算法結(jié)合這兩個(gè)最優(yōu)值迭代更新自己的速度和位置。

算法基于群體的優(yōu)化技術(shù)的，亦即可擁有多條搜索軌道，具有良好的并行性；PSO算法也無需梯度信息，只利用目標(biāo)的取值信息，具有很強(qiáng)的通用性；此外，PSO算法原理簡單、操作方便和易用性極強(qiáng)。PSO算法已在諸多應(yīng)用領(lǐng)域取得了成功。但是，PSO算法在進(jìn)化后期也存在著搜索精度、效率較低，易陷入局部最優(yōu)等缺點(diǎn)。尤其在高維度的優(yōu)化問題上，對于PSO算法，粒子在進(jìn)化過程中一旦早熟，其狀態(tài)近似停滯，很難逃出局部最優(yōu)。本文針對該問題提出了一種基于部分變異策略的PSO（PSO algorithm with partial mutation strategy，簡稱PMPSO）算法，用以提高粒子陷入局部最優(yōu)后的逃逸能力。具體思想為，在進(jìn)化的早、中期就提前篩選早熟的粒子并對其進(jìn)行位置信息變異操作，而不是等到粒子群整體收斂后對所有粒子變異；該方法有望以較高的優(yōu)化效率獲得更多的搜索到全局最優(yōu)的機(jī)會(huì)，達(dá)到避免算法早熟、幫助停滯粒子逃出局部最優(yōu)的目的。PMPSO的具體算法及相關(guān)討論詳見文獻(xiàn)[2]，不再冗述。

2.2 基于面元法理論的螺旋槳非定常水動(dòng)力性能預(yù)報(bào)方法

本文利用基于速度勢的低階面元法對螺旋槳的非定常水動(dòng)力性能進(jìn)行預(yù)報(bào)，采用計(jì)算較為簡捷的關(guān)于擾動(dòng)速度勢的基本積分微分方程；應(yīng)用雙曲面形狀的面元以消除面元間的縫隙，并選用線性尾渦模型模擬尾渦形狀；通過Newton-Raphson迭代來滿足槳葉隨邊等壓Kutta條件；選取某一葉片作為主葉片進(jìn)行面元上奇點(diǎn)強(qiáng)度的求解，其它葉片上的速度勢分布可由主葉片推導(dǎo)求得。該螺旋槳非定常水動(dòng)力性能預(yù)報(bào)方法的詳細(xì)描述見文獻(xiàn)[3]，在此省略。

2.3 螺旋槳側(cè)斜分布曲線的表達(dá)

通常，有兩大類方法來表達(dá)一條曲線，一類是用點(diǎn)集合表達(dá)；一類是用一個(gè)或多個(gè)形狀函數(shù)表達(dá)。本文采用B樣條曲線表達(dá)螺旋槳的側(cè)斜（角）分布曲線，是一種點(diǎn)集合表達(dá)方法。

B樣條曲線具有表示和設(shè)計(jì)自由曲線曲面的強(qiáng)大功能，用較少的控制點(diǎn)就能表達(dá)整個(gè)曲線或曲面，是目前應(yīng)用最廣泛的形狀數(shù)學(xué)描述的主流方法之一[4]，采用B樣條曲線表達(dá)螺旋槳的側(cè)斜角分布，在任意調(diào)整、修改分布曲線形狀的同時(shí)能夠最大程度上保證曲線的光順性，而且B樣條曲線本身較少的控制點(diǎn)也意味著優(yōu)化變量維度的減少，將更易于優(yōu)化方法得到最優(yōu)解。相比于傳統(tǒng)的采用若干個(gè)需要經(jīng)驗(yàn)確定的形狀函數(shù)組來表達(dá)曲線形狀的方法[5]，B樣條曲線方法具備更強(qiáng)的便利性和易用性。

B樣條曲線在具體應(yīng)用時(shí)，分為兩種算法：一種是由控制多邊形（多面體）頂點(diǎn)計(jì)算B樣條曲線（面），即正算過程；另一種是由原始數(shù)據(jù)點(diǎn)反推控制多邊形（多面體）頂點(diǎn)，即反算過程。B樣條曲線基本原理及具體算法詳見文獻(xiàn)[4]。

下圖為分別采用不同控制點(diǎn)個(gè)數(shù)表達(dá)的同一條側(cè)斜角分布曲線，從中可以看出，取6個(gè)控制點(diǎn)生成的B樣條曲線來表達(dá)螺旋槳側(cè)斜角分布，其精度是可以接受的。

圖1 不同B樣條控制點(diǎn)個(gè)數(shù)表達(dá)的側(cè)斜角分布曲線Fig.1 Skewed angles distribution curves with different numbers control points of B spline

2.4 螺旋槳側(cè)斜優(yōu)化原理及方法

本文的螺旋槳側(cè)斜分布的優(yōu)化方法，是利用2.1節(jié)的優(yōu)化算法配合2.2節(jié)的非定常水動(dòng)力性能預(yù)報(bào)方法，對側(cè)斜分布問題進(jìn)行迭代優(yōu)化設(shè)計(jì)的。該優(yōu)化問題具體提法如下：優(yōu)化目標(biāo)：

限制條件：

優(yōu)化變量：

目標(biāo)函數(shù)：

式中：Wi為函數(shù)權(quán)重，有為螺旋槳葉頻、倍葉頻軸向非定常力和力矩，在具體計(jì)算中，為簡化計(jì)算、取得更好的優(yōu)化結(jié)果，只取一倍、二倍葉頻的軸向非定常力和力矩；i為葉頻倍數(shù)。

限制條件中，σT為優(yōu)化后螺旋槳推力與原螺旋槳推力的誤差，εT為誤差限，即為可接受的推力損失上限；此限制條件可保證，經(jīng)過側(cè)斜優(yōu)化后的螺旋槳，其帶來的推力損失在設(shè)計(jì)者可接受的范圍內(nèi)，否則為減小非定常脈動(dòng)力而犧牲太多螺旋槳推力，也不是合理的設(shè)計(jì)理念。

優(yōu)化變量Sj為螺旋槳葉各剖面?zhèn)刃被騻?cè)斜角度，其變化范圍由上限SU和下限SL控制。本文采用B樣條曲線來表達(dá)螺旋槳葉各剖面?zhèn)刃苯堑姆植?，所以，此時(shí)優(yōu)化變量Sj代表了B樣條曲線的控制點(diǎn)的位置。

側(cè)斜優(yōu)化的具體流程可簡單概括為：通過面元法程序預(yù)報(bào)螺旋槳的葉頻、倍葉頻軸向非定常力和力矩，構(gòu)成優(yōu)化目標(biāo)σf，由改進(jìn)的PSO算法進(jìn)行優(yōu)化計(jì)算，根據(jù)優(yōu)化目標(biāo)及限制條件的滿足情況，修改優(yōu)化變量Si的值形成新的設(shè)計(jì)方案，并令面元法程序?qū)υ撔路桨高M(jìn)行預(yù)報(bào)計(jì)算，再構(gòu)成優(yōu)化目標(biāo)σf，由優(yōu)化算法再進(jìn)行評價(jià)、修改優(yōu)化變量Si，如此循環(huán)直至獲得最優(yōu)解。

2.5 基于OpenMP的螺旋槳側(cè)斜優(yōu)化的并行方法

因?yàn)榱Ｗ尤核惴ㄖ械牧Ｗ釉谶M(jìn)行適應(yīng)度函數(shù)的計(jì)算時(shí)是彼此獨(dú)立進(jìn)行的，存在天然的并行性，所以可在此階段引入并行算法，在單位時(shí)間內(nèi)提高粒子的進(jìn)化效率即提高整個(gè)算法的優(yōu)化效率。

OpenMP（Open Multi-Processing）是一套已被廣泛接受的、支持跨平臺(tái)共享內(nèi)存并行系統(tǒng)的多線程程序設(shè)計(jì)的指導(dǎo)性編譯處理方案。OpenMP使用C、C++和Fortran語言，可以在大多數(shù)的處理器體系和操作系統(tǒng)中運(yùn)行。在具體應(yīng)用中，OpenMP可以通過直接的編譯指令實(shí)現(xiàn)并行執(zhí)行。OpenMP提供了一個(gè)“fork-and-join”的執(zhí)行模型，在這個(gè)模型中，一個(gè)程序以線程的形式被執(zhí)行。線程依次執(zhí)行直到發(fā)現(xiàn)一個(gè)并行區(qū)域中的并行指令。此時(shí)，該線程構(gòu)造出一系列線程并且成為主線程，所有的線程同步執(zhí)行聲明，直到并行區(qū)域結(jié)束。

本文應(yīng)用OpenMP并行語言庫[5]開發(fā)了針對側(cè)斜分布優(yōu)化問題的并行PSO優(yōu)化算法，在多核CPU計(jì)算機(jī)上調(diào)用優(yōu)化程序，可明顯減少優(yōu)化設(shè)計(jì)時(shí)間。

3 螺旋槳側(cè)斜分布優(yōu)化結(jié)果及分析

3.1 螺旋槳非定常力預(yù)報(bào)精度校核

本文的螺旋槳側(cè)斜分布優(yōu)化方法，是建立在對螺旋槳在非定常流場中的軸向非定常力準(zhǔn)確預(yù)報(bào)的基礎(chǔ)上的，所以本節(jié)首先檢驗(yàn)螺旋槳非定常水動(dòng)力的預(yù)報(bào)精度。

第14屆ITTC曾組織了螺旋槳非定常力的比較計(jì)算，計(jì)算對象是一只4葉B系列槳，該槳的詳細(xì)參數(shù)及給定的伴流場情況可參見文獻(xiàn)[6]。文獻(xiàn)[6]給出了14家的計(jì)算結(jié)果與實(shí)驗(yàn)值的比較情況。本文按2.2節(jié)的方法應(yīng)用面元法理論計(jì)算了該槳的軸向非定常力，計(jì)算時(shí)采用了同時(shí)計(jì)入軸向和切向伴流的組合伴流。如圖2-3所示為計(jì)算結(jié)果的比較情況，縱坐標(biāo)為相應(yīng)的力或力矩系數(shù)，橫坐標(biāo)標(biāo)號(hào)代表不同的計(jì)算結(jié)果，具體為：（1）Hanaoka計(jì)算值；（2）Tsakonas計(jì)算值；（3）Jacobs計(jì)算值；（4）Isay計(jì)算值；（5）Rank計(jì)算值；（6）Verbrugh計(jì)算值；（7）實(shí)驗(yàn)測試值；（8）本文計(jì)算值。

圖2 螺旋槳一倍葉頻軸向力對比Fig.2 Comparison of 1st order blade rate axial force

圖3 螺旋槳二倍葉頻軸向力對比Fig.3 Comparison of 2nd order blade rate axial force

圖2-3的對比表明，本文的基于面元法理論的螺旋槳非定常脈動(dòng)力預(yù)報(bào)方法精度良好，可作為側(cè)斜分布優(yōu)化的性能計(jì)算基礎(chǔ)。

3.2 螺旋槳側(cè)斜分布優(yōu)化結(jié)果及分析

選擇某五葉大側(cè)斜槳Seiun-Maru HSP作為側(cè)斜優(yōu)化的目標(biāo)槳[7]。HSP螺旋槳的主要幾何參數(shù)、伴流條件及其非定常軸向脈動(dòng)力計(jì)算結(jié)果分別如表1、圖4和圖5所示，由計(jì)算結(jié)果可再次驗(yàn)證本文螺旋槳非定常脈動(dòng)力預(yù)報(bào)方法的可靠性。

表1 HSP螺旋槳的幾何參數(shù)Tab.1 Main parameter of HSP

圖4 HSP槳的伴流場Fig.4 Flowfield of HSP

取HSP槳為母型槳、HSP槳的伴流場為給定伴流，以其現(xiàn)有的側(cè)斜分布的最大值、最小值為優(yōu)化變量的上、下限，并固定葉根和葉梢處側(cè)斜角的值；設(shè)優(yōu)化算法的最大迭代次數(shù)iter=200，粒子數(shù)m=60，推力損失上限εT=2.5%，以螺旋槳一倍、二倍葉頻的軸向脈動(dòng)力和力矩作為優(yōu)化目標(biāo)，優(yōu)化其側(cè)斜分布，優(yōu)化結(jié)果如圖6所示。經(jīng)計(jì)算，優(yōu)化后的螺旋槳一倍葉頻軸向脈動(dòng)力和力矩較HSP槳分別下降了2.1%和1.5%，二倍葉頻軸向脈動(dòng)力和力矩分別下降了1.1%和0.8%，而推力值也較HSP槳下降1%；圖6還給出了未給定SL即側(cè)斜角變化下限的優(yōu)化結(jié)果，對應(yīng)此優(yōu)化結(jié)果螺旋槳一倍葉頻軸向脈動(dòng)力和力矩較HSP槳分別下降了約10%和5%，二倍葉頻軸向脈動(dòng)力和力矩分別下降了約6%和4%，但從圖中可以看出，側(cè)斜角有變大的趨勢，已使所優(yōu)化結(jié)果在實(shí)際應(yīng)用和制造上不成立，且側(cè)斜過大，螺旋槳的強(qiáng)度還需額外考慮，所以，規(guī)定合理的側(cè)斜角變化范圍是得到合理的優(yōu)化側(cè)斜分布的前提。

圖5 HSP螺旋槳主葉片非定常力比較Fig.5 Comparison of unsteady force of HSP main blade

圖6 優(yōu)化后的側(cè)斜角分布Fig.6 Skew angle distribution after optimizaiton

側(cè)斜優(yōu)化計(jì)算工作在聯(lián)想“Think Station”工作站平臺(tái)完成，其CPU處理器為雙6核“Intel Xeon E5-2620”，表2所示為采用并行算法后優(yōu)化算法完成優(yōu)化工作所用的時(shí)間及加速比和計(jì)算機(jī)使用效率，以12個(gè)CPU并行工作所用時(shí)間來看，可比單CPU串行工作節(jié)省時(shí)間約88%?？梢?，采用多核CPU計(jì)算機(jī)運(yùn)行并行優(yōu)化算法，可明顯減少優(yōu)化時(shí)間，從而提高計(jì)算、優(yōu)化效率。

表2 并行計(jì)算的時(shí)間、加速比及效率Tab.2 Time,speed-up ratio and efficiency of parallel computing

4 結(jié)論

本文給出了一種基于改進(jìn)PSO算法和非定常面元法的螺旋槳側(cè)斜分布數(shù)值優(yōu)化方法，并將該算法并行化。針對HSP螺旋槳，以降低軸向非定常力幅值為目標(biāo)對側(cè)斜分布進(jìn)行了優(yōu)化及相關(guān)的分析探討。數(shù)值優(yōu)化結(jié)果表明，采用本文方法優(yōu)化側(cè)斜分布，可降低螺旋槳的非定常力幅值；側(cè)斜優(yōu)化過程中，要確定好側(cè)斜角變化的上下限，以避免強(qiáng)度校核和應(yīng)用加工方面的問題；并行算法的引入可明顯提高優(yōu)化效率，12個(gè)CPU并行工作可比單CPU串行工作節(jié)省時(shí)間約88%；本文的側(cè)斜優(yōu)化方法是可行、有效的。

參考文獻(xiàn)：

[1]Kennedy J,Ebehtart R,Shi Y H.Swarm intelligence[M].San Francisco:Morgan Kuafann Publishers,2010.

[2]Cai Haopeng,Su Yumin,Zhang He.A particle swarm optimization algorithm with partial mutation strategy and its application[C]//2009 IEEE International Conference on Intelligent Computing and Intelligent Systems.Shanghai,China,2009: 522-526.

[3]蘇玉民,黃勝.船舶螺旋槳理論[M].哈爾濱:哈爾濱工程大學(xué)出版社,2003.

[4]施法中.計(jì)算機(jī)輔助幾何設(shè)計(jì)與非均勻有理B樣條[M].北京:高等教育出版社,2001.

[5]王獻(xiàn)孚.船用翼理論[M].北京:國防工業(yè)出版社,1998.

[6]Hermanns M.Parallel Programming in Fortran 95 using OpenMP[M].Universidad Politecnica de Madid,Spain,2002.

[7]Jessup S.Measurement of multiple blade rate unsteady propeller forces[R].David Taylor Research Center,Ship Hydromechanics Department Research and Development Report,DTRC-90/015 May 1990.

[8]Experimental data for unsteady panel calculations and comparisons(SEIUN-MARU HSP)[C]//Propeller RANS/Panel Method Workshop.22nd ITTC Propulsion Committee.Grenoble,France,1998.

[9]王國強(qiáng),董世湯.船舶螺旋槳理論與應(yīng)用[M].哈爾濱:哈爾濱工程大學(xué)出版社,2007.

A numerical optimization method to determine propeller skew distribution

CAI Hao-peng1,MA Cheng1,CHEN Ke1,QIAN Zheng-fang1,2,ZHANG He3
(1 Institute of Naval Academy of Armament,Beijing 100161,China;2 Naval Univ.of Engineering,Wuhan 430033, China;3 China Ship Research and Development Academy,Beijing 100192,China)

Better skew distribution is benefit for propeller blade sections to access the nonuniform ship wake nonsynchronously,so the vibration and noise induced by propeller operating in the ship wake can be decreased.In this paper,an artificial intelligence optimization algorithm and a prediction method of unsteady hydrodynamic characteristics of propeller based on the surface panel method are provided to optimize the propeller skew distribution.The propeller skew distribution curve is described by B-spline,and the curve fairing is ensured;Aiming at reducing unsteady axial forces,an improved particle swarm optimization combined with the prediction method of propeller unsteady hydrodynamic performance is utilized to optimize and analyze propeller skew distribution.The parallel optimization algorithm based on OpenMP is also developed for multi-cup computers in order to improve the computation efficiency.Numerical optimization results show that the propeller skew optimization method provided in this paper can reduce propeller unsteady thrust amplitude,the method is effective,reliable and practical,and the parallel computing can improve the optimizing efficiency greatly.

skew distribution optimization;improved PSO(particle swarm optimization); unsteady surface panel method;parallel computing

U661.313

A

10.3969/j.issn.1007-7294.2014.07.005

1007-7294（2014）07-0771-07

2014-04-15

國家自然科學(xué)基金資助項(xiàng)目（NO.51079158）；中國博士后科學(xué)基金面上資助（2012M512140）

蔡昊鵬（1983-），男，海軍裝備研究院在站博士后，E-mail：cdcaisiqi@gmail.com；

馬騁（1963-），男，研究員，博士生導(dǎo)師。