NURBS流曲線造型新方法及其在船舶設計中應用

張 彥 儒, 林 焰,2, 陸 叢 紅*, 紀 卓 尚

( 1.大連理工大學 船舶工程學院, 遼寧 大連 116024； 2.大連理工大學 工業裝備結構分析國家重點實驗室, 遼寧 大連 116024 )

NURBS流曲線造型新方法及其在船舶設計中應用

張 彥 儒1, 林 焰1,2, 陸 叢 紅*1, 紀 卓 尚1

( 1.大連理工大學 船舶工程學院, 遼寧 大連 116024； 2.大連理工大學 工業裝備結構分析國家重點實驗室, 遼寧 大連 116024 )

針對當前流曲線造型方法僅被應用于已知型值點及其切矢的問題，在其他約束(面積、面積心、流場特征值)問題中的應用，以及與NURBS設計方法的結合仍有待發展研究的現狀，將NURBS與流曲線造型技術相結合，提出了NURBS流曲線單元的概念以及在其他幾何約束下，應用NURBS流曲線造型的新方法．該方法將待設計流曲線在已知切矢的型值點處分段以確定NURBS流曲線單元的個數，將NURBS流曲線單元控制參數的并集作為設計變量，將幾何約束處理為曲線設計的目標函數，建立普適性的NURBS流曲線造型框架．該框架模型可用進化算法進行求解．通過設計實例，驗證了該方法在相關曲線設計方面的可行性和有效性．

船體曲面；流曲線；NURBS；功能曲線；幾何約束

0 引 言

在運動物體(如汽車、飛機、船舶等)外形設計中，運動物體在空氣、水流等流體中作相對運動，流體對其產生阻力，因此其外形設計變得十分重要．運動物體外形的光滑與否將直接影響其運動性能．而具有流線形外形的運動物體不僅外觀漂亮，并且能極大地減少前進過程中流體對物體的阻力．針對運動物體的外形設計，以流體力學為背景的流曲線曲面的造型方法被提出，這種方法稱為流曲線曲面造型技術．縱觀現有文獻，可以將流曲線曲面造型技術分為兩類：第一類是用切矢的積分曲線進行曲線曲面造型的技術[1-6]，其本質是曲線切矢的拼接，是以弧長為參數的三維插值曲線．國內研究者在國外研究的基礎上，在旋轉曲面[7]、非均勻細分[8]、雙三次非均勻B樣條曲面的G1混合[9]等方面進行了研究．但是這一類流曲線的切矢并非流場中的速度矢量，因此本文定義該類方法為廣義的流曲線曲面，它并非是流線形設計中的功能曲線．第二類是基于流體力學原理，將流體力學中流函數的概念引入CAD來構造流曲線曲面．1996年盧小林等給出了流曲線曲面的概念及多項式和B網構造流曲線曲面的基本方法[10]，1997年游世輝等給出了流體曲面設計的一般流體光順條件[11]．2006年劉國偉提出了來流為直勻流的流曲線造型方法，并于2010年應用于高速列車頭部曲線設計[12-13]．2017年徐陽等提出一種面向葉型和翼型的考慮流線場約束的NURBS曲線擬合方法[14]．上述文獻的造型方法，僅局限在流曲線上型值點及其切矢的設計意圖；而流曲線設計，作為一種功能曲線設計，不只是考慮曲線上型值點和切矢的幾何條件，有時還需要考慮曲線所圍面積、面積心等幾何約束，但未見相關報道．另外，1991年國際標準化組織(ISO)頒布的工業產品數據交換標準STEP中，把NURBS作為定義工業產品幾何形狀的唯一數學方法．因此，研究如何將NURBS曲線與流曲線相結合，對流曲線造型的發展將會有非常積極的推動作用．并且，隨著計算機技術、數值計算方法和優化方法的發展，船舶的外形設計由常規的正向設計逐漸走向反向設計[15]，因此研究如何對周圍流場信息進行提取，反向指導功能曲線設計具有現實的工程意義．總體來說，流曲線曲面造型技術的基本理論、數學模型和相關算法仍然不夠完善，還有待進一步研究擴展．

本文在第二類流曲線造型技術的基礎上，將待設計流曲線用NURBS流曲線單元組合設計，并給出在曲線幾何約束(型值點、面積、面積心、流場中特征點特征值)條件下，NURBS流曲線的數值造型方法，以彌補流曲線造型只能處理型值點、切矢等的造型方法的不足，并可應用于反向設計．

1 流曲線數值求解模型

第二類流曲線基本是以理想不可壓縮流體平面定常無旋運動為力學模型構造柔性樣條，將滿足流函數的全部或部分流線作為待求的功能性曲線．它是具有流線特征、可描述物理上可能的流動、可攜帶流場物理信息的功能曲線[12]．

由流體力學可知二維無旋不可壓縮定常流動的流函數ψ(x,y)滿足拉普拉斯方程：

(1)

流函數的等值線ψ(x,y)=C(C為常數)就是流場中的流線，流線上點的切矢與流場中的速度矢量相一致．如圖1所示，同一個物體在不同的流場a和b中運動有不同的流線，并且即使是在同一個流場中，常數C取不同的值也對應著彎曲程度不同的流線．因此基于拉普拉斯方程構造的流曲線具有良好的可塑性．另外，由物面的不可滲透邊界條件及流線的流體不可穿越性可知，流場中的任意流線均可用同位置、同形狀的物面代替．因此流曲線造型的實質是求解滿足給定流動狀態的拉普拉斯方程以得到在特定位置、滿足一定形狀約束的流線．

(a) 流場a

(b) 流場b

求解拉普拉斯方程，一種是給定物體的外形和繞流形式的前提下求解流函數，這一類稱之為正問題，可以應用常規的求解偏微分方程的方法進行求解．而另一種，是在給定部分流函數的前提下，求解物體形狀，這一類稱之為反問題．由前面敘述可知，流曲線造型屬于反問題．

1.1 流曲線造型的數值解法

流曲線造型數值解法的實質是拉普拉斯方程第一邊值問題的逆問題求解．因此先給出其求解公式，然后給出幾何約束條件的處理，以及需要注意的幾個補充方程．本文采用五點差分格式，對拉普拉斯方程進行求解．五點差分格式為

ψi,j=(ψi,j+1+ψi-1,j+ψi+1,j+ψi,j-1)/4

(2)

1.1.1 計算域 如圖2所示，給出來流為直勻流條件下，坐標系和邊界條件的定義．曲線OEF為待求流曲線．AB、BC、CD、AO、FD為固定邊界，OF為控制邊界，流曲線OEF的形狀完全由控制邊界OF的邊界值確定．h為x、y方向上的網格步長．○點稱為網格內點，×點稱為邊界點．x軸定義在計算域邊界AD上，指向固定邊界CD為正．y軸指向固定邊界BC為正．坐標原點由待求流曲線OEF的幾何特點決定．設來流方向平行于x軸，沿x軸正向．

圖2 計算域

文獻[13]未給出固定邊界的數值確定方法．經過反復比較和分析，證實如果在曲線設計階段，給出了流場中網格點值或者流場中特定點切線向量時，將進流段AO長度和去流段FD長度，以及計算域寬度AB按給定的流場區域進行計算即可；否則可以分別按圖2所示，相對OF長度的對應比例選取．

1.1.2 邊界條件 設來流為直勻流，平行于x軸，沿x軸正向，速度大小為vs．因固定邊界處流場信息與待求流曲線形狀無關，所以固定邊界處邊界條件可按直勻流流場特征給出．邊界條件可定義為

(3)

其中yB表示邊界點B的y坐標．從圖2和式(3)可知，邊界條件中只有控制邊界OF段為未知．文獻[13]給出了已知零流曲線(物面)上部分型值點及其切矢的條件下，將反問題轉化為正問題，運用五點差分法求解拉普拉斯方程的數值造型方法．該方法雖然給出了切矢轉化為網格點上流函數值的方法，但并未給出如何在滿足計算精度的網格步長下，使型值點落到網格點上，而且未知控制邊界上的邊界點的個數只能等于已知內點個數，因此其網格不能細分，計算精度不高．本文不考慮控制邊界OF段的邊界值，直接將待求零流曲線OEF用NURBS表示，進行流曲線造型設計．給出補充方程：

ψ(x,y)=0； (x,y)∈OEF

(4)

這樣，只要零流曲線給定，拉普拉斯方程就可以求解．考慮到待求的零流曲線作為控制邊界時，邊界可能不通過網格點，在給定特征點后劃分網格，邊界點也可能不通過網格點，給出上述情況的臨近邊界點網格點的計算方法如下．如圖3所示在網格點(i,j)處的導數，可按下述方法計算．

圖3 邊界不通過網格點

根據Taylor展開公式：

(5)

(6)

消去一階導數?ψi,j/?x，可以得到網格點(i,j)處關于x的二階導數：

(7)

同理，可以得到網格點(i,j)處關于y的二階導數：

(8)

因此拉普拉斯方法在網格點(i,j)處可近似為

(9)

其中ψA和ψB為拉普拉斯方程的邊界值，在本文中因為是零流曲線，所以取值為0．

1.2 零流曲線的NURBS表示

一條p次NURBS曲線定義為

(10)

經過推導，可以得到NURBS曲線在端點處(u=0，u=1)的一階導矢，其表達式如下：

(11)

考慮到設計流曲線一般不會發生折轉，經過反復比較和分析，證實了利用三次NURBS曲線及5個控制頂點來設計一般自由曲線，可以滿足工程精度和靈活修改等要求．將這個自由曲線段稱為NURBS流曲線單元．

如圖4所示P0是流曲線單元的起點，P4是流曲線單元的終點．假設xi、yi和wi分別是控制頂點Pi(i=0，1，…，4)的縱向坐標、橫向坐標及該控制頂點對應的權因子．自由段曲線的首尾控制點權因子w0和w4可設置為1．節點矢量U可以設置為(0 0 0 0 0.5 1 1 1 1)．首末控制頂點P0和P4為已知物面(零流曲線)的型值點，y1和y3可以由式(11)以及w1和w3聯合確定，因此流曲線單元的設計變量為x1、x2、x3、y2、w1、w2、w3．

圖4 流曲線單元控制頂點分布

將已知的零流曲線上的型值點分為兩種：一種為型值點和型值點處曲線切矢已知，另一種為型值點已知而型值點處切矢未知．用第一種型值點將零流曲線進行分段，每一段用一個流曲線單元表示．第二種型值點作為幾何約束進行處理．

1.3 幾何約束的處理

本文將幾何約束轉化為設計目標進行處理，即將設計曲線相關幾何量與幾何約束間的相對誤差的線性加權值最小作為設計目標．即

(12)

式中：n為幾何約束的個數，αi為各個幾何約束的權重系數，Fi為幾何約束的目標函數．下面分別敘述兩種典型的幾何約束處理方式．

1.3.1 流場中特征值的處理 流曲線造型，不能只考慮物面(零流曲線)上的切矢處理，更應該考慮周邊流場信息的提取，并依此來指導物面(零流曲線)的設計．流場中的特征值一般分為特征點的速度矢量和特征點的流函數兩類．由流線的性質可知：流線上任意一點的速度矢量在該點必與流線相切．如圖5所示，流場中網格點(i,j)處的切矢就是實際流場中該點的速度矢量v．因此可以將特征點上的速度矢量按如下方式轉化為特征點流函數值進行處理．

vx、vy與流函數之間的關系及差分方程為

(13)

圖5 特征點上速度矢量處理

網格點(i,j)處的流函數值ψi,j與網格步長h已知，當vx、vy已知時，可以由式(13)推出相鄰網格點的流函數值．這樣，可以將速度矢量幾何約束的設計目標轉化為最小化對應網格點處流函數值的相對誤差的最大值，設ψdi,j+1和ψdi+1,j為通過本文方法設計的流曲線流場中網格點(i,j+1)和(i+1,j)的流函數值，則設計目標可以表示為

minF1=max{|ψdi,j+1-ψi,j+1|/|ψi,j+1|,

|ψdi+1,j-ψi+1,j|/|ψi+1,j|,

|ψdi,j-ψi,j|/|ψi,j|}

(14)

1.3.2 未知切矢型值點的處理 假設欲使設計曲線通過型值點列Pod=(Pod1Pod2…Podn)，流曲線上對應型值點x坐標的y坐標為(Pd1.yPd2.y…Pdn.y)．因此設計目標為最小化各型值點y坐標值與對應的流曲線上點的y坐標值的相對誤差中的最大值，即

minF2=max{|Pd1.y-Pod1.y|/|Pod1.y|,

|Pd2.y-Pod2.y|/|Pod2.y|,…,

|Pdn.y-Podn.y|/|Podn.y|}

(15)

1.3.3 流曲線面積和面積心的處理 設Awo為所設計的曲線和控制邊界圍成的區域要求的面積，(Laox,Laoy)為該區域面積的面積心坐標，Awd和(Ladx,Lady)分別為通過本文設計方法得到的流曲線的相應參數，目標函數設為最小化設計曲線對應區域的面積相對誤差和面積心坐標的相對誤差之和，即

minF3=min{|Awd-Awo|/|Awo|+

|Ladx-Laox|/|Laox|+

|Lady-Laoy|/|Laoy|}

(16)

2 流曲線造型步驟

本文采用以下方法做NURBS流曲線造型：

Step1研究待設計流曲線的幾何性質特征，給出最少控制頂點的NURBS設計模型即流曲線單元．

Step2給出已知切矢的型值點．然后用已知切矢的型值點將待設計流曲線進行分段，分別對每一段用流曲線單元進行設計，最后匯總給出總的設計變量．

Step3給出其他幾何約束條件，根據1.3約束處理方法將幾何約束轉化為設計目標．

Step4選用進化算法，將Step2中確定的設計變量和Step3中確定的目標函數代入進化算法進行求解．

Step5輸出設計得到的流曲線及其周邊流場信息．

3 數值算例

根據本文的造型方法，采用人機交互的遺傳算法[16]對7 000 t 散貨船(設計航速18 kn)設計水線進行方案設計的數值算例．為了配合圖2計算域的坐標，本船坐標系設置為坐標系原點在舯橫剖面與中縱剖線的交點，x軸指向船尾為正，y軸指向左舷為正．按水線特征將水線分為進流段和去流段．表1列出了已知條件和設計目標．其中Awo表示水線進流段或去流段與船中縱剖線和所設計曲線自由段的端點處橫剖面所圍面積，(Laox,Laoy)表示面積心坐標，Pst表示流曲線單元起點坐標，vst表示起點處切矢，Ped表示流曲線單元終點坐標，ved表示終點處切矢，ψC為計算流場中特征點處的流函數值．

圖6給出了本文方法得到的設計水線結果以及控制頂點分布，圖7給出了周圍流場信息．

表1 設計水線的已知條件和設計目標

Tab.1 The known conditions and design objective of the design waterline

已知條件Pst/mPed/mvstved進流段(-50．764，0．405)(-23．257，8．600)(79．4847．12)(397．5 0)去流段(25．5，8．6)(53．799，3．348)(322．8 0)(1013．4-262)設計目標Awo/m2(Laox，Laoy)/m特征點/mψC進流段168．128(-33．958，3．557)(-50．995，1)4．407去流段187．359(37．656，3．515)(53．055，4)10．255

圖6 設計水線的結果

圖7 周圍流場信息

表2列出了本文方法得到的水線控制頂點坐標．其中，x和y分別是控制頂點的x、y坐標，w是對應控制點的權值．表3列出了算法得到的水線的特征參數以及與設計目標的相對誤差|Er|．由于進流段端部存在拐點，相比去流段的誤差較大．但對于各項幾何設計目標，已基本達到了工程上的初始設計精度要求(工程精度要求一般為0.5%)．

表2 設計水線控制頂點坐標

表3 設計水線設計結果與設計目標間的差異

4 結 語

本文將流曲線造型方法進行擴展，結合NURBS和進化計算方法，給出了在給定約束條件下，利用拉普拉斯方程第一邊值問題的逆問題數值解進行流曲線造型的新方法．研究結果表明：本文可以運用周圍流場的信息來指導流曲線的設計，并且可以在給定來流的條件下初步對其進行評估．

本文對于流曲線造型方法的擴展僅做了初步探討，在以下方面仍需進一步研究：如何將二維設計轉化為三維設計；來流為非直勻流時，流曲線的造型方法；如何提取周圍三維流場信息以指導二維特征線設計；采用并行計算方法以加快流曲線設計速度．總之，流曲線曲面造型技術有待進一步研究，應用領域仍有待擴展．

另外，算例中，進流段誤差相比去流段誤差較大，是因為設計目標水線的進流段端部存在拐點，逼近較困難，造成設計結果與設計目標間面積、面積心及流函數值的誤差較大．在五點差分方法中，對計算域進行網格離散后，網格點值對邊界值比較敏感，微小的邊界值波動，也會引起較大的流函數值的誤差．因此流曲線設計結果取決于相關幾何特性的誤差以及設計結果形狀是否滿足設計預期．而流函數值的作用是作為指導約束，與相關幾何特性一起指導流曲線的設計，使其設計形狀逼近設計預期．其允許誤差與網格點所在位置和設計目標形狀有關，它的具體規律將是后續的研究內容．

[1] MIURA K T, WANG L, CHENG F. Streamline modeling with subdivision surfaces on the Gaussian sphere [J].ComputerAidedDesign, 2001,33(13): 975-987.

[2] SUZUKI Y, TANAKA I, MIURA K T,etal. Streamline modeling based on potential flow [C] //Proceedings-9thPacificConferenceonComputerGraphicsandApplications,PacificGraphics2001. Tokyo: IEEE Computer Society, 2001:305-312.

[3] MIURA K T, CHENG Fuhua, WANG Lazhu. Fine tuning: curve and surface deformation by scaling derivatives [C] //Proceedings-9thPacificConferenceonComputerGraphicsandApplications,PacificGraphics2001. Tokyo: IEEE Computer Society, 2001: 150-159.

[4] MIURA K T, SAKIYAMA N, KANEKO T. Streamline modeling-Designing free-form surfaces by blending tangent vectors [J].IFIPAdvancesinInformationandCommunicationTechnology, 2001,75:235-248.

[5] SAKIYAMA N, MIURA K T, TAKAHASHI T,etal. Fair curve and surface design system using tangent control [C] //PG′98Proceedingsofthe6thPacificConferenceonComputerGraphicsandApplications. Washington D C: IEEE Computer Society, 1998: 27-36.

[6] MIURA K T. Unit quaternion integral curve: A new type of fair free-form curves [J].ComputerAidedGeometricDesign, 2000,17(1): 39-58.

[7] 李秀娟,劉 浩,何 鋼,等. 基于流曲線的旋轉曲面的生成和形狀調整[J]. 機械科學與技術, 2008,27(3):326-329.

LI Xiujuan, LIU Hao, HE Gang,etal. Generation and shape adjustment of revolution surface based on stream curve [J].MechanicalScienceandTechnologyforAerospaceEngineering, 2008,27(3):326-329. (in Chinese)

[8] 李秀娟,劉 浩,何 鋼,等. 基于非均勻細分的流曲線曲面[J]. 機械科學與技術, 2008,27(10):1216-1219, 1224.

LI Xiujuan, LIU Hao, HE Gang,etal. Stream curve and surface based on non-uniform subdivision [J].MechanicalScienceandTechnologyforAerospaceEngineering, 2008,27(10): 1216-1219, 1224. (in Chinese)

[9] 李秀娟,廖文和,劉 浩,等. 基于流曲線曲面的雙三次非均勻B樣條曲面G1混合[J]. 中國科學院大學學報, 2008,25(4):554-559.

LI Xiujuan, LIAO Wenhe, LIU Hao,etal.G1continuity blending between bicubic non-uniform B-spline surfaces with stream curve and surface [J].JournaloftheGraduateSchooloftheChineseAcademyofSciences, 2008,25(4):554-559. (in Chinese)

[10] 盧小林,馬利莊,何志均. 流曲線曲面的概念和構造方法[J]. 計算機工程, 1996,22(3):26-31.

LU Xiaolin, MA Lizhuang, HE Zhijun. The conception and construction method for fluid curve [J].ComputerEngineering, 1996,22(3):26-31. (in Chinese)

[11] 游世輝,卿啟湘,卿 鈞. 計算機輔助流體曲面優化設計的探討[J]. 湖南大學學報, 1997,24(1):71-75.

YOU Shihui, QING Qixiang, QING Jun. Research of computer aided optimal design for fluid surface [J].JournalofHunanUniversity, 1997,24(1):71-75. (in Chinese)

[12] 劉國偉. 流曲線概念及造型方法[J]. 計算機輔助設計與圖形學學報, 2006,18(12):1891-1896.

LIU Guowei. Conception and modeling method of stream curve [J].JournalofComputer-AidedDesign&ComputerGraphics, 2006,18(12):1891-1896. (in Chinese)

[13] LIU Guowei, LIU Weixiao. Stream curve and numerical modeling for high-speed train [C] //InternationalConferenceonOptoelectronicsandImageProcessing. Piscataway:IEEE Computer Society, 2010: 407-410.

[14] 徐 陽,劉 強. 考慮流線場約束的NURBS曲線擬合方法[J]. 計算機輔助設計與圖形學學報, 2017,29(1):137-144.

XU Yang, LIU Qiang. NURBS curve fitting method considering constraint of streamline field [J].JournalofComputer-AidedDesign&ComputerGraphics, 2017,29(1):137-144. (in Chinese)

[15] 趙 峰,李勝忠,楊 磊,等. 基于CFD的船型優化設計研究進展綜述[J]. 船舶力學, 2010,14(7):812-820.

ZHAO Feng, LI Shengzhong, YANG Lei,etal. An overview on the design optimization of ship hull based on CFD techniques [J].JournalofShipMechanics, 2010,14(7):812-820. (in Chinese)

[16] 王運龍,王 晨,彭 飛,等. 基于人機結合遺傳算法的船舶管路三維布局優化設計[J]. 中國造船, 2015,56(1):196-202.

WANG Yunlong, WANG Chen, PENG Fei,etal. Intelligent layout design of ship pipes based on genetic algorithm with human-computer cooperation [J].ShipbuildingofChina, 2015,56(1):196-202. (in Chinese)

NewNURBSstreamcurvemodelingmethodanditsapplicationtoshipdesign

ZHANGYanru1,LINYan1,2,LUConghong*1,JIZhuoshang1

(1.SchoolofNavalArchitecture&OceanEngineering,DalianUniversityofTechnology,Dalian116024,China；2.StateKeyLaboratoryofStructuralAnalysisforIndustrialEquipment,DalianUniversityofTechnology,Dalian116024,China)

The stream curve modeling method has only been applied to the problems with data points and their tangent vectors being given. The application to the problems with other constraints such as area, centroid of area, the eigenvalue of the flow field, etc. and the integration with NURBS need more research and development. NURBS and stream curve modeling method are integrated and the concept of NURBS stream curve unit is proposed. Furthermore, a new NURBS stream curve modeling method based on other geometric constraints is proposed. In this method, the stream curve is segmented at the data points with the given tangent vectors, and the number of the NURBS stream curve cell is determined. The union sets of the controls parameter of the NURBS stream curve cell are designated as the design variables. The geometric constraints are regarded as the design objective function. The universal framework of NURBS stream curve modeling is built. This framework model can be solved by an evolutionary algorithm. The design instances indicate that this method is effective and feasible in the related curve design.

hull surface; stream curve; NURBS; functional curve; geometric constraint

1000-8608(2017)06-0564-07

U662.2

A

10.7511/dllgxb201706003

2017-04-05；

2017-09-25．

國家自然科學基金資助項目(E091002-51109033).

張彥儒(1985-)，男，博士生，E-mail:03401011@163.com；林 焰(1963-)，男，博士，教授，博士生導師，E-mail:linyanly@dlut.edu.cn；陸叢紅*(1972-)，女，博士，副教授，E-mail:lchcad@dlut.edu.cn.