閆博文 黃 俊 劉志勤 王耀彬

西南科技大學(xué)計(jì)算機(jī)科學(xué)與技術(shù)學(xué)院，四川綿陽(yáng) 621000

隨著航空航天技術(shù)的發(fā)展和國(guó)防科技的需求，飛行器的外形設(shè)計(jì)越來(lái)越受關(guān)注，飛行器設(shè)計(jì)的第一步便是外形參數(shù)化，它的優(yōu)劣程度對(duì)于飛機(jī)設(shè)計(jì)的成功與否有著至關(guān)重要的影響，直接影響著總體項(xiàng)目的進(jìn)度。一個(gè)好的參數(shù)化方法應(yīng)該具備以下特性：1)較大的尋優(yōu)空間；2)良好的局部控制能力，以及較少的設(shè)計(jì)變量；3)保證在優(yōu)化過(guò)程中飛機(jī)的幾何外形是光滑的[1-2]。

Kulfan和Bussoletti提出了一種比較優(yōu)秀的參數(shù)化方法：CST技術(shù)。CST技術(shù)同時(shí)使用類別函數(shù)與形狀函數(shù)來(lái)對(duì)外形進(jìn)行控制，類別函數(shù)(Class Function)生成幾何圖形的基本外形，再通過(guò)形狀函數(shù)(Shape Function)對(duì)這個(gè)基本的幾何圖形進(jìn)行修正，從而得到所需要的幾何外形。因該技術(shù)具有設(shè)計(jì)變量少，可調(diào)節(jié)，設(shè)計(jì)空間廣，可產(chǎn)生光滑的幾何外形，以及較好的魯棒性等優(yōu)點(diǎn)而被廣泛應(yīng)用在飛行器外形參數(shù)化的設(shè)計(jì)中[3-5]。由于Bezier多項(xiàng)式對(duì)局部幾何特性描述能力不足，王迅等人提出了使用B樣條函數(shù)替代Bezier多項(xiàng)式構(gòu)造形函數(shù)的改進(jìn)方法，以提高對(duì)外形的局部描述能力[6]，因?yàn)锽樣條方法在數(shù)據(jù)的擬合、平滑和插值等方面有良好的效果，所以在階數(shù)較小的情況下，該方法確實(shí)起到了優(yōu)化作用。但B樣條基函數(shù)為分段函數(shù)，無(wú)法使用統(tǒng)一的解析式表達(dá)，大多使用迭代方法求值，會(huì)出現(xiàn)多次重復(fù)計(jì)算，隨著階數(shù)的增加，計(jì)算量也急劇增加，文獻(xiàn)[6]也只采用了3階來(lái)進(jìn)行證明。因此，只有解決B樣條在高階數(shù)計(jì)算量大的問(wèn)題，找出一種快速求解的算法，提高運(yùn)算效率，才能讓B樣條改進(jìn)CST方法發(fā)揮出更大的價(jià)值。因此，本文主要研究B樣條基函數(shù)快速求值改進(jìn)CST算法。

1 B樣條替代Bezier多項(xiàng)式建模

1.1 B樣條基函數(shù)定義

B樣條的基函數(shù)有許多定義方法，它們中最被大家廣泛認(rèn)可的是de Boor-Cox遞推定義法，因?yàn)閐e Boor-Cox遞推定義充分展示了B樣條基函數(shù)的性質(zhì)，并給出了明確的B樣條基函數(shù)遞推求解算法，故本文將其作為B樣條基函數(shù)的標(biāo)準(zhǔn)定義[7]。

基函數(shù)定義：

(1)

1.2 替代Bezier多項(xiàng)式的CST方法建模

一般機(jī)翼的CST方法數(shù)學(xué)表達(dá)式：

(2)

(3)

修改后的B樣條基函數(shù)CST公式對(duì)求解時(shí)出現(xiàn)的病態(tài)現(xiàn)象有很大的緩解，求解更加穩(wěn)定，具有很好的局部描寫(xiě)能力。

2 B樣條基函數(shù)快速求值算法

2.1 de Boor-Cox算法弊端

在實(shí)際計(jì)算B樣條時(shí),往往不需要將整個(gè)B樣條數(shù)值全部計(jì)算出來(lái)，有時(shí)只需要計(jì)算某一參數(shù)值處的B樣條基函數(shù)。由B樣條局部支撐性質(zhì)，可知當(dāng)參數(shù)u∈[ui,ui+1),只需計(jì)算出Nm,k(u),m=i-k,i-k+1,…,i；而其他k次B樣條在該處均為0。

圖1 B樣條基函數(shù)運(yùn)算路徑

圖1中，B樣條基函數(shù)計(jì)算過(guò)程已用虛線三角形標(biāo)出，可以直觀的看到，每一個(gè)基函數(shù)不僅依賴于它下一層相鄰的2個(gè)基函數(shù)，而且還影響它上一層相鄰的2個(gè)基函數(shù)。顯然，式(1)的計(jì)算方法并不能充分利用這種網(wǎng)格結(jié)構(gòu)，存在大量的重復(fù)計(jì)算，例如2個(gè)虛線三角形重疊的部分就是重復(fù)計(jì)算的部分，故而要加快B樣條基函數(shù)求解效率就要從這里著手。

2.2 準(zhǔn)備工作

為后續(xù)闡述B樣條基函數(shù)快速求值首先定義一種向量擴(kuò)展算法。

Zj=Nj+Qj-1-Qj

(4)

這種遞推式的運(yùn)算有利于使用計(jì)算機(jī)運(yùn)算，稱定義3給出的運(yùn)算為向量遞推擴(kuò)展運(yùn)算。

2.3 快速求解算法證明[8]

首先，對(duì)將要使用的符號(hào)進(jìn)行說(shuō)明：

Nj,m(u):表示m次B樣條第j個(gè)樣條基函數(shù)，其中，m=0,1,2,…,k；j=1,2,3；

(5)

證明B樣條基函數(shù)快速求值前，還需證明引理1。

根據(jù)定義2的向量擴(kuò)展運(yùn)算法，則有：

(6)

其中：0≤j≤m，規(guī)定向量下標(biāo)溢出時(shí)對(duì)應(yīng)元素置0。

(7)

再由定義(1)可知：

定理1提供了實(shí)現(xiàn)B樣條基函數(shù)求值的向量遞推擴(kuò)展方法。其計(jì)算路徑如圖2。由圖2可知，B樣條基函數(shù)向量擴(kuò)展快速求值方法正是一種網(wǎng)狀計(jì)算結(jié)構(gòu)，它是以1維的單位向量{[1]}為主向量，通過(guò)遞推擴(kuò)展，實(shí)現(xiàn)了同次數(shù)B樣條基函數(shù)值的同步計(jì)算(圖2中虛線方框部分)。向量遞推擴(kuò)展方法避免了重復(fù)計(jì)算，提高了計(jì)算效率，這由圖2算法計(jì)算路徑?jīng)]有重疊部分可以充分說(shuō)明。綜上所述，B樣條的向量擴(kuò)展求值算法在計(jì)算路徑上比式(1)的de Boor-Cox算法效率更高。

圖2 B樣條向量擴(kuò)展法基函數(shù)運(yùn)算路徑

2.4 算法效率

同de Boor-Cox遞推定義法相比，式(4)的B樣條基函數(shù)向量擴(kuò)展求值算法效率更高，計(jì)算k次B樣條基函數(shù)非零值，統(tǒng)計(jì)對(duì)比式(1)與式(4)的運(yùn)算次數(shù)，可得：使用向量同步遞推運(yùn)算獲得k+1個(gè)非零基函數(shù)值的計(jì)算效率是使用de Boor-Cox遞推方法的2k+1倍。

表1 k次B樣條的(k+1)個(gè)非零基函數(shù)運(yùn)算量對(duì)比

3 應(yīng)用展示

3.1 B樣條快速求值CST二維建模

圖3為使用B樣條快速求值CST方法生成的二維圖形。

圖3 B樣條向量擴(kuò)展法模擬Clark-Y機(jī)翼

3.2 B樣條快速求值CST三維建模

3.2.1 局部控制展示

圖4中，(a)圖展示使用B樣條快速求值方法控制橫截面從圓形漸變到正方形，(b)圖展示截面函數(shù)指數(shù)NC隨著長(zhǎng)度的增加而進(jìn)行變化，若將該立體圖形的總長(zhǎng)度視為100%，在總長(zhǎng)度的0%～20%區(qū)間內(nèi)，立體圖形的橫截面是標(biāo)準(zhǔn)圓形，NC=0.5；在總長(zhǎng)度的95%～100%區(qū)間內(nèi)，立體圖形的橫截面是標(biāo)準(zhǔn)的正方形，NC=0.005，在總長(zhǎng)度的20%～95%區(qū)間內(nèi)，函數(shù)指數(shù)隨長(zhǎng)度的增加而減小。

圖4 B樣條快速求值局部控制展示

3.2.2 B樣條快速求值CST方法三維機(jī)翼

三維機(jī)翼其實(shí)本質(zhì)上可以看成是二維機(jī)翼的擴(kuò)展，擴(kuò)展時(shí)需要注意控制外形的變化。

圖5 CST算法三維機(jī)翼展示

4 結(jié)論

通過(guò)B樣條改進(jìn)CST參數(shù)化方法，提高了CST算法的局部控制能力，同時(shí)使用向量擴(kuò)展方法實(shí)現(xiàn)B樣條快速求值，提高了B樣條CST算法的計(jì)算效率，增大了設(shè)計(jì)空間，提高了可控性，為提高后續(xù)翼型氣動(dòng)優(yōu)化的精度，最終得到良好的氣動(dòng)外形奠定了基礎(chǔ)。

[1] Isakova N P, Kraiko A A, P’Yankov K S. Direct Method for the Design of Optimal Three-dimensional Aerodynamic Shapes[J]. Computational Mathematics & Mathematical Physics, 2012, 52(11):1520-1525.

[2] 廖炎平, 劉莉, 龍騰. 幾種翼型參數(shù)化方法研究[J]. 彈箭與制導(dǎo)學(xué)報(bào), 2011, 31(3):160-164.(Liao Yanping, Liu Li, Long Teng. The Research on Some Param-eterized Methods for Airfoil[J]. Journal of Projectiles Rockets Missiles & Guidance, 2011, 31(3):160-164.)

[3] Kulfan B. A Universal Parametric Geometry Representa-tion Method - "CST"[C]. AIAA--2007-0062, Jan 2007

[4] Kulfan, B. M, Bussoletti, J. E., Fundamental Parametric Geometry Representations for Aircraft Component Shapes[C]. AIAA-2006-6948, 11th AIAA/ISSMO Multidisciplinary Analysis and Optimization Conference: The Modeling and Simulation Frontier for Multidisciplinary Design Optimization, 6 - 8 September, 2006

[5] Ciampa P D, Zill T, Nagel B. CST Parametrization for Unconventional Aircraft Design Optimization[C].Congress of the International Council of the Aerospace Sciences. DLR, 2010.

[6] 王迅, 蔡晉生, 屈崑,等. 基于改進(jìn)CST參數(shù)化方法和轉(zhuǎn)捩模型的翼型優(yōu)化設(shè)計(jì)[J]. 航空學(xué)報(bào), 2015, 36(2):449-461.(Wang Xun，Cai Jinsheng，Qu K，et al. Airfoil Optimization Based on Improved CST Parametric Method and Transition Model[J]. Acta Aeronautica et Astronautica Sinica，2015,36(2)：449-461.)

[7] Piegl L, Tiller W. The NURBS, Book[M]. Springer Ber-lin Heidelberg, 1997.

[8] 孫海洋, 范大鵬. B樣條快速求值算法及其在數(shù)控插補(bǔ)中的應(yīng)用[J]. 計(jì)算機(jī)工程與應(yīng)用, 2007, 43(34):81-84.(Sun Haiyang, Fan Dapeng.Fast Algorithm to Compute B-spline Functions and Its Application in CNC Interpolation[J]. Computer Engineering and Applications,2007,43(34)：81-84.)