曹勝鴻,楊燕初,張?zhí)鞎?張航悅

(1. 中國(guó)科學(xué)院空天信息創(chuàng)新研究院,北京100094;2. 中國(guó)科學(xué)院大學(xué),北京100190)

1 引言

臨近空間(NearSpace)通常是指20～100km的空域,在地球大氣中處于平流層和中層,其在環(huán)境監(jiān)測(cè)、區(qū)域通信、實(shí)時(shí)監(jiān)視、預(yù)警與導(dǎo)彈防御等領(lǐng)域有很大的應(yīng)用價(jià)值[1]。飛艇是一種輕于空氣的航空器,具有垂直起降、長(zhǎng)航時(shí)、載重大、低能耗、噪聲及排放小等特點(diǎn)[2],在臨近空間的探測(cè)方面具有廣泛的應(yīng)用價(jià)值。平流層飛艇是臨近空間飛行器的一種,近幾十年來(lái),隨著航空技術(shù)的發(fā)展,采用新材料和新技術(shù)的大型飛艇得到飛速發(fā)展。

因此,有針對(duì)性地對(duì)平流層飛艇進(jìn)行優(yōu)化設(shè)計(jì),即多方面考慮飛艇的設(shè)計(jì)因素就顯得尤為重要。平流層飛艇總體設(shè)計(jì)與優(yōu)化技術(shù)作為方案設(shè)計(jì)的關(guān)鍵,將會(huì)直接決定飛艇系統(tǒng)的整體性能和操控[3]。國(guó)內(nèi)多家研究機(jī)構(gòu)開(kāi)展了大量的飛艇設(shè)計(jì)優(yōu)化工作,例如,中航通飛研究所的秦何軍針對(duì)飛艇的氣動(dòng)外形進(jìn)行優(yōu)化設(shè)計(jì)[4],中國(guó)特種飛行研究所的張超針對(duì)重在飛艇低阻囊體進(jìn)行外形優(yōu)化[5]。但上述優(yōu)化設(shè)計(jì)僅以阻力系數(shù)作為優(yōu)化對(duì)象,缺少對(duì)飛艇結(jié)構(gòu)因素的考慮。

通過(guò)對(duì)比分析發(fā)現(xiàn),在改變長(zhǎng)細(xì)比時(shí),體面比與長(zhǎng)細(xì)比成正相關(guān),阻力系數(shù)與長(zhǎng)細(xì)比呈現(xiàn)負(fù)相關(guān),因此,本文保留阻力、體面比兩個(gè)飛艇設(shè)計(jì)因素,建立數(shù)學(xué)模型,結(jié)合實(shí)際需求分配兩者權(quán)重,結(jié)合改進(jìn)的復(fù)合形法,構(gòu)建優(yōu)化設(shè)計(jì)框架,通過(guò)Matlab軟件運(yùn)算得出艇體長(zhǎng)度與艇體最大截面直徑,依照實(shí)際需求,提供三種合理的飛艇外形設(shè)計(jì)方案,并結(jié)合實(shí)際算例仿真分析,驗(yàn)證優(yōu)化艇形和優(yōu)化算法的可行性。

2 算法介紹

2.1 飛艇優(yōu)化程序框架

本文在進(jìn)行飛艇優(yōu)化設(shè)計(jì)過(guò)程中,結(jié)合最優(yōu)化設(shè)計(jì)中的復(fù)合形法進(jìn)行改進(jìn),進(jìn)行目標(biāo)函數(shù)的最優(yōu)值求解,首先輸入原始艇形數(shù)據(jù),結(jié)合多學(xué)科優(yōu)化設(shè)計(jì)理念,建立飛艇優(yōu)化目標(biāo)函數(shù),并根據(jù)實(shí)際需求,合理分配各項(xiàng)函數(shù)模型的權(quán)重,利用Matlab軟件結(jié)合改進(jìn)復(fù)合形法探索優(yōu)化值,依照實(shí)際需求,有側(cè)重的針對(duì)算法所得飛艇長(zhǎng)度、飛艇最大截面直徑、飛艇長(zhǎng)細(xì)比進(jìn)行定向飛艇優(yōu)化,并提供后期仿真分析的艇形數(shù)據(jù),其飛艇優(yōu)化流程框架如圖1所示。需要指出的是,本文所提出的改進(jìn)復(fù)合形算法也可與其它最優(yōu)值求解方法相結(jié)合使用,推廣至其它最優(yōu)值求解問(wèn)題。

圖1 飛艇優(yōu)化流程框架

2.2 改進(jìn)復(fù)合形法

在約束問(wèn)題的最優(yōu)化求解過(guò)程中,往往會(huì)出現(xiàn)導(dǎo)數(shù)和梯度較難求得的情況,尤其是多元函數(shù),甚至?xí)](méi)有明確的工程函數(shù)[6]。因此常常采用直接求解的方式,復(fù)合形法的基本求解思想是基于無(wú)約束最優(yōu)化問(wèn)題直接求解方式的單純形法,僅需在每次擴(kuò)張、收縮之后驗(yàn)證新得到的點(diǎn)是否在可行域內(nèi)即可。

單純形法最早由Spendley、Hext和Himsworth提出[7],后由Nelder和Mead改進(jìn)成為變型的單純形法[8],其基本思想是:三點(diǎn)構(gòu)造單純形,沿最差點(diǎn)至剩余點(diǎn)所構(gòu)中心點(diǎn)的方向?qū)ふ易畈铧c(diǎn)的替代點(diǎn),從而構(gòu)造新的單純形,其探索過(guò)程如圖2所示。

圖2 單純形法探索過(guò)程

取三個(gè)點(diǎn)并計(jì)算其函數(shù)值,分析三個(gè)點(diǎn)的函數(shù)值進(jìn)行排列得到最差點(diǎn)Xh,次差點(diǎn)Xg最好點(diǎn)Xl,取次差點(diǎn)Xg與最好點(diǎn)Xl的形心Xc,連接最差點(diǎn)Xh與形心點(diǎn)Xc做第一次等步長(zhǎng)的擴(kuò)張得到反射點(diǎn)Xr,若反射點(diǎn)的函數(shù)值優(yōu)于最差點(diǎn)Xh則再進(jìn)一步擴(kuò)張,否則,將步長(zhǎng)減半,如果新得到的點(diǎn)函數(shù)值依舊差于最差點(diǎn)函數(shù)值,則直接改變探索方向,其中單純形法的探索過(guò)程,還需要將新得到的點(diǎn)與次差點(diǎn)進(jìn)行比較,從而考慮新得到點(diǎn)的刪減。

經(jīng)典單純形法需要計(jì)算剩余點(diǎn)所構(gòu)中心點(diǎn),并將新得到的點(diǎn)與最差點(diǎn)、次差點(diǎn)進(jìn)行比較,且隨著探索迭代次數(shù)的增加,在最優(yōu)值附近的探索速度較慢,探索步驟緩慢、繁瑣,極大地降低了約束問(wèn)題的求解效率。

單純形法求解最優(yōu)化工程函數(shù)問(wèn)題的關(guān)鍵在于單純形的構(gòu)建,而在單純形構(gòu)建的過(guò)程中,選點(diǎn)將會(huì)直接影響探索方向的準(zhǔn)確性,經(jīng)典單純性法求解時(shí)需要計(jì)算次差點(diǎn)與最好點(diǎn)的形心,并與最差點(diǎn)相連,構(gòu)成新的探索方向,楊晉提出可以直接連接最差點(diǎn)與最好點(diǎn)構(gòu)成新的探索方向[9],李慶高在非線(xiàn)性約束問(wèn)題同樣提出直接在最好點(diǎn)與最差點(diǎn)的方向探索得到最差點(diǎn)的替代點(diǎn)[10]。

本文在楊晉單純形法的基礎(chǔ)上進(jìn)一步改進(jìn),提出雙向探索的求解思想,即沿著最差點(diǎn)和最好點(diǎn)的方向探索最差點(diǎn)的替代點(diǎn)的同時(shí),沿著次差點(diǎn)和最好點(diǎn)的方向探索次差點(diǎn)的替代點(diǎn),在兩個(gè)收斂方向同時(shí)求得第二次單純形的構(gòu)造點(diǎn),精簡(jiǎn)單純形構(gòu)造的計(jì)算步驟,提高約束問(wèn)題的求解效率,其探索步驟如圖3所示。

圖3 改進(jìn)算法探索過(guò)程示意圖

改進(jìn)單純形法的計(jì)算過(guò)程如下:

1)取三個(gè)點(diǎn)并計(jì)算其函數(shù)值,分析三個(gè)點(diǎn)的函數(shù)值并進(jìn)行排列,從而得到最差點(diǎn)Xh,次差點(diǎn)Xg,最好點(diǎn)Xl,連接次差點(diǎn)Xg與最好點(diǎn)Xl做第一次等步長(zhǎng)的擴(kuò)張,同時(shí)連接最差點(diǎn)Xh與最好點(diǎn)Xl做第一次等步長(zhǎng)的擴(kuò)張,如式(1)所示

(1)

α：反射系數(shù),單純形法中通常取為1,復(fù)合形法中常取1.3;

Xhr：最差點(diǎn)Xh相對(duì)于最好點(diǎn)Xl的反射點(diǎn);

Xgr：次差點(diǎn)Xg相對(duì)于最好點(diǎn)Xl的反射點(diǎn);

2)分別比較兩個(gè)反射點(diǎn)與最差點(diǎn)、次差點(diǎn)的函數(shù)值,如果反射點(diǎn)函數(shù)值小于最差點(diǎn)、次差點(diǎn)的函數(shù)值,則說(shuō)明探索正確,進(jìn)一步擴(kuò)張,如式(2)所示

(2)

γ：擴(kuò)張系數(shù),一般取為2;

Xh4：反射點(diǎn)Xhr相對(duì)于最好點(diǎn)Xl的擴(kuò)張點(diǎn);

Xg4：反射點(diǎn)Xgr相對(duì)于最好點(diǎn)Xl的擴(kuò)張點(diǎn);

3)如果反射點(diǎn)函數(shù)值比最差點(diǎn)、次差點(diǎn)函數(shù)值要大,則減少擴(kuò)張步長(zhǎng),即收縮。如式(3)所示

(3)

β：收縮系數(shù),一般取為0.5;

Xh2：反射點(diǎn)Xhr相對(duì)于最好點(diǎn)Xl的收縮點(diǎn);

Xg2：反射點(diǎn)Xgr相對(duì)于最好點(diǎn)Xl的收縮點(diǎn);

4)如若收縮點(diǎn)的函數(shù)值仍舊大于最差點(diǎn)、次差點(diǎn)的函數(shù)值,則說(shuō)明探索方向有誤,改變探索方向,即縮小單純形,如式(4)所示

(4)

Xh3：最差點(diǎn)Xh相對(duì)于最好點(diǎn)Xl的壓縮點(diǎn);

Xg3：次差點(diǎn)Xg相對(duì)于最好點(diǎn)Xl的壓縮點(diǎn);

5)探索所得點(diǎn)與第一次單純形的最好點(diǎn)構(gòu)成第二次探索單純形的三個(gè)點(diǎn),圖3所示藍(lán)色線(xiàn)即為第二次單純形的可能情況。單純形構(gòu)造完成后,逐次迭代,直到最好點(diǎn)與形心點(diǎn)幾乎為一個(gè)點(diǎn),終止迭代的判據(jù)如式(5)所示

(5)

針對(duì)收斂性判據(jù),除上述判斷形心點(diǎn)與最好點(diǎn)是否幾近重合之外,也可以判斷最后一次所構(gòu)單純形是否幾近為一點(diǎn),其收斂性判斷的主體思想都是一樣的,即判斷最終構(gòu)造點(diǎn)幾乎為一點(diǎn)。

其改進(jìn)復(fù)合形法的算法流程圖如圖4所示。

圖4 改進(jìn)算法程序框圖

上述單純形法的改進(jìn),省略了經(jīng)典單純形法中對(duì)于形心點(diǎn)的計(jì)算,同時(shí)對(duì)楊晉算法[9]單一方向的探索進(jìn)行改進(jìn),在次差點(diǎn)與最好點(diǎn)、最差點(diǎn)與最好點(diǎn)兩個(gè)方向進(jìn)行計(jì)算和探索,以此求得新的單純形探索構(gòu)建點(diǎn),省略了反射點(diǎn)與次差點(diǎn)的比較過(guò)程,算法簡(jiǎn)單,探索速度快。本文借用一個(gè)案例函數(shù)對(duì)改進(jìn)算法和經(jīng)典算法進(jìn)行比較,驗(yàn)證改進(jìn)算法的可行性。案例目標(biāo)函數(shù)如式(6)所示,收斂性要求0.01。

(6)

兩種方法單次單純形最優(yōu)值與迭代次數(shù)的關(guān)系對(duì)比圖,如圖(5)所示。

通過(guò)圖5的對(duì)比分析,針對(duì)式(6)所示案例目標(biāo)函數(shù),改進(jìn)算法在第3次迭代過(guò)程中,出現(xiàn)大幅度改變,從第6次迭代開(kāi)始進(jìn)入函數(shù)收斂區(qū)域,而經(jīng)典算法在第13次迭代才有大幅度的轉(zhuǎn)變,從第20次迭代開(kāi)始進(jìn)入函數(shù)收斂區(qū)域。改進(jìn)算法的迭代速度相較于經(jīng)典算法提高了58.33%,由此可得,改進(jìn)算法的探索速度高于經(jīng)典單純形法,極大地提高了函數(shù)最優(yōu)值求解的效率。

圖5 經(jīng)典算法與改進(jìn)算法單次最優(yōu)值對(duì)比

3 飛艇優(yōu)化目標(biāo)函數(shù)

以往的飛艇優(yōu)化設(shè)計(jì)著重于阻力系數(shù)的優(yōu)化,忽略其它設(shè)計(jì)因素,本文在進(jìn)行多次設(shè)計(jì)優(yōu)化分析后,針對(duì)飛艇初始優(yōu)化階段,保留阻力系數(shù)和體面比兩個(gè)重要的設(shè)計(jì)因素。其中,阻力系數(shù)直接決定飛艇推進(jìn)功率,進(jìn)而影響能源供給;體面比是艇囊體積與面積的比值,決定著系統(tǒng)浮力的有效利用率[11],其目標(biāo)函數(shù)數(shù)學(xué)模型如式(7)所示

(7)

w1,w2為阻力和體面比的權(quán)重,FCD為阻力,FCD0為初始阻力,λ為體面比,λ0為初始體面比值。

3.1 阻力模型

Hoerner研究了流線(xiàn)型物體的阻力,對(duì)于多數(shù)飛艇,其艇體的氣流是湍流[12],因此阻力系數(shù)關(guān)系如式(8)所示

(8)

其中l(wèi)為艇體長(zhǎng)度,d為艇體最大截面直徑。

Hoerner提出對(duì)于實(shí)際水準(zhǔn)的表面光潔度和Re>5×106,表皮摩擦系數(shù)Cf如式(9)所示

(9)

綜合式(8),式(9),給出阻力系數(shù)公式,如式(10)所示

(10)

其中雷諾數(shù)計(jì)算參見(jiàn)式(11)

(11)

其中μ為粘性系數(shù),ρa(bǔ)為氣體密度,v為風(fēng)速。

故阻力模型函數(shù)如式(12)所示

(12)

其中kD為整艇與裸艇之間的阻力系數(shù)比,由飛艇的實(shí)際構(gòu)型決定,Ve為囊體體積。

3.2 體面比模型

針對(duì)飛艇初始化設(shè)計(jì)階段,王全保[2]僅考慮表面積帶來(lái)的影響,忽略?xún)?yōu)化過(guò)程中囊艇體積的變化,故此,本文引入體積模型,建立體面比影響模型,體面比模型可參照式(13)所示,式(14)-(15)為優(yōu)化前期體積和表面積估算公式,其精確計(jì)算公式參照式(16)-(17)

(13)

其中囊體體積估算公式參照式(14)

(14)

對(duì)于飛艇優(yōu)化初始階段,艇體表面積估算公式,可參照式(15)

Se=(2.14+(l-60)·0.08-(d-15)·0.028×

(15)

其中體積精確值的計(jì)算公式參照式(16)

(16)

對(duì)于飛艇表面積精確值的計(jì)算公式,可參照式(17)

(17)

4 算例優(yōu)化仿真分析

本文選取英國(guó)物理實(shí)驗(yàn)室所提出的一種低阻力飛艇形狀,即后期Skyship飛艇形狀的基礎(chǔ),NPL低阻力系數(shù)飛艇形狀[12],作為飛艇前期優(yōu)化的外形母線(xiàn)函數(shù),其數(shù)學(xué)表達(dá)式如式(18)所示

(18)

l為飛艇長(zhǎng)度,d為飛艇直徑,其N(xiāo)PL低阻力系數(shù)飛艇外形如圖6所示。

圖6 NPL低阻力系數(shù)飛艇外形

本文參照雷諾數(shù)Re>5×106,以及實(shí)際飛艇設(shè)計(jì)案例,取初始飛艇長(zhǎng)度60m,最大截面直徑20m,長(zhǎng)細(xì)比3。取飛艇高度為20km高空,大氣壓力5474.8Pa,風(fēng)速25m/s,空氣密度取0.088kg/m3,根據(jù)實(shí)際設(shè)計(jì)經(jīng)驗(yàn),以一定的最佳長(zhǎng)細(xì)比范圍為約束,此外考慮到飛艇的最大浮力,故以固定的囊體體積為約束,如式(19)所示

(19)

運(yùn)用改進(jìn)復(fù)合形法進(jìn)行目標(biāo)函數(shù)的優(yōu)化探索,借用Matlab軟件生成優(yōu)化后的飛艇外形,如圖6至圖8所示,其中圖6對(duì)應(yīng)阻力權(quán)重為1,體面比權(quán)重為0,以此求得以阻力為優(yōu)化重點(diǎn)時(shí)的優(yōu)化艇形,圖7對(duì)應(yīng)阻力權(quán)重為0,體面比權(quán)重為1,以此求得以體面比為優(yōu)化重點(diǎn)時(shí)的優(yōu)化艇形,圖8對(duì)應(yīng)阻力權(quán)重為0.5,體面比權(quán)重0.5,以此求得阻力與體面比平衡時(shí)的優(yōu)化艇形,需要指出的是,在實(shí)際設(shè)計(jì)優(yōu)化工作過(guò)程中,設(shè)計(jì)者可以根據(jù)實(shí)際需求更改調(diào)配、本文只給出常用的幾組權(quán)重分配比例,此外,也可以針對(duì)目標(biāo)函數(shù)增添其它設(shè)計(jì)影響因素,本文保留兩項(xiàng)隨長(zhǎng)細(xì)比變化呈現(xiàn)相反變化的阻力和體面比因素項(xiàng)。

圖7 權(quán)重為1:0時(shí)飛艇優(yōu)化外形

圖8 權(quán)重為0:1時(shí)飛艇優(yōu)化外形

此外,根據(jù)所得優(yōu)化值點(diǎn),本文給出三個(gè)不同側(cè)重點(diǎn)的優(yōu)化結(jié)果,分別對(duì)應(yīng)每組圖中的除原始艇形外的三條優(yōu)化曲線(xiàn),若設(shè)計(jì)考慮目標(biāo)函數(shù)值最優(yōu),則選擇以?xún)?yōu)化點(diǎn)的最大直徑為恒定值的優(yōu)化結(jié)果,若設(shè)計(jì)考慮長(zhǎng)細(xì)比的范圍,則選擇以?xún)?yōu)化點(diǎn)的長(zhǎng)細(xì)比為恒定值的優(yōu)化結(jié)果,若設(shè)計(jì)考慮飛艇的微小外形優(yōu)化,則以?xún)?yōu)化點(diǎn)的艇長(zhǎng)為恒定值的優(yōu)化結(jié)果。

權(quán)重分配為阻力為1,體面比為0時(shí)的優(yōu)化結(jié)果如圖7所示。

從圖7可以看出,當(dāng)以阻力作為優(yōu)化目標(biāo)時(shí),優(yōu)化后的艇形會(huì)變細(xì)長(zhǎng),即整體艇體長(zhǎng)度變長(zhǎng),最大截面直徑變小,抗阻能力增強(qiáng),其詳細(xì)相關(guān)對(duì)比數(shù)值參見(jiàn)表1。

表1 權(quán)重為1:0時(shí)優(yōu)化艇形與原始艇形對(duì)比

從表1中可以看出,當(dāng)權(quán)重分配為1:0時(shí),從三個(gè)側(cè)重點(diǎn)優(yōu)化所得艇形的阻力得到有效地改善,體面比則會(huì)相應(yīng)的減小。

權(quán)重分配為阻力為0,體面比為1時(shí)的優(yōu)化結(jié)果如圖8所示。

從圖8可以看出,在以體面比為優(yōu)化目標(biāo)時(shí),優(yōu)化后的艇形會(huì)變短胖,即整體艇體長(zhǎng)度變短,最大截面直徑變大,更有利于艇體浮力的利用,其詳細(xì)相關(guān)對(duì)比數(shù)值參見(jiàn)表2。

表2 權(quán)重為0:1時(shí)優(yōu)化艇形與原始艇形對(duì)比

從表2可以看出,當(dāng)權(quán)重因子全部分配給體面比時(shí),所得優(yōu)化艇形的體面比得到顯著提升,相應(yīng)的截面變大,艇長(zhǎng)變短,阻力會(huì)增大。

通過(guò)圖7與圖8的優(yōu)化結(jié)果可以看出,體面比和阻力是優(yōu)化過(guò)程中兩個(gè)呈現(xiàn)相反增量的目標(biāo)量,所以,在實(shí)際的優(yōu)化過(guò)程中,還需要依照實(shí)際需求分配阻力與體面比的權(quán)重。

本文同時(shí)給出權(quán)重分配為1:1時(shí)的優(yōu)化艇形結(jié)果,如圖9,表3所示。

圖9 權(quán)重為1:1時(shí)飛艇優(yōu)化外形

從上述表1、2、3的仿真分析可以看出,優(yōu)化過(guò)程中體面比與阻力呈現(xiàn)相反優(yōu)化增量關(guān)系,故在實(shí)際應(yīng)用中應(yīng)根據(jù)實(shí)際需求分配權(quán)重因子,此外,本文考慮實(shí)際優(yōu)化設(shè)計(jì)的三個(gè)關(guān)鍵變量,最大截面直徑、長(zhǎng)細(xì)比和飛艇長(zhǎng)度,分別有側(cè)重地提出優(yōu)化方案,從優(yōu)化結(jié)果來(lái)看,改進(jìn)算法在飛艇優(yōu)化中是可行的。

5 結(jié)論

本文結(jié)合傳統(tǒng)飛艇優(yōu)化設(shè)計(jì)過(guò)程中單一考慮阻力因素的問(wèn)題,首先引入體面比,構(gòu)建優(yōu)化函數(shù)數(shù)學(xué)模型,構(gòu)建優(yōu)化流程框架,其次針對(duì)傳統(tǒng)最優(yōu)值求解的復(fù)合形法進(jìn)行改進(jìn),依照雙向同步探索思維,提出最優(yōu)值求解的改進(jìn)算法,并將其應(yīng)用至飛艇外形優(yōu)化設(shè)計(jì)中,運(yùn)用Matlab進(jìn)行函數(shù)的計(jì)算和仿真,通過(guò)仿真分析驗(yàn)證了算法的可行性。