展俊德 李錫文 史鐵林

華中科技大學(xué)，武漢，430074

0 引言

立式捏合機作為固體推進(jìn)劑研制和生產(chǎn)的關(guān)鍵設(shè)備，是固體推進(jìn)劑技術(shù)的重要組成部分。兩槳立式捏合機具有物料混合均勻、生產(chǎn)率高等優(yōu)點，在固體推進(jìn)劑的生產(chǎn)中得到了廣泛的應(yīng)用［1］。立式捏合機的混合效率及混合效果在很大程度上取決于攪拌槳的幾何形狀，故攪拌槳曲面設(shè)計是立式捏合機設(shè)計和制造過程中的一個重要環(huán)節(jié)［2］。兩槳立式捏合機攪拌槳由實心槳與空心槳組成，攪拌槳曲面是由其水平截面輪廓線按照一定數(shù)值的螺旋升角順時針或逆時針掃掠生成的。攪拌槳水平截面輪廓線是由自由曲線、圓弧和直線段等以一定的幾何連續(xù)條件依次組合而成的復(fù)雜曲線，攪拌槳水平截面輪廓線和攪拌槳曲面不能用初等解析函數(shù)來表示。

非均勻有理B樣條（NURBS）能夠?qū)崿F(xiàn)二次解析幾何和自由曲線曲面的統(tǒng)一，通過局部改變控制點或權(quán)因子即可調(diào)整局部的曲線曲面形狀，而不影響其他部分，具有良好的幾何性質(zhì)［3］。因此，本文根據(jù)攪拌槳不同高度的截面輪廓線生成了若干個3次NURBS平面曲線，利用蒙皮法生成了立式捏合機攪拌槳的NURBS曲面；在生成的NURBS曲面誤差分析的基礎(chǔ)上，提出了一種新的控制頂點投影法計算曲面理論坐標(biāo)點到NURBS曲面的最小距離的算法，用以評價NURBS曲面的重構(gòu)精度。

1 NURBS曲線曲面表示

一條k次NURBS曲線可以表示為一分段有理多項式矢函數(shù)［4］：

其中，ωi為權(quán)因子，分別與控制頂點Pi相關(guān)。首末權(quán)因子ω0＞0，ωn＞0，其余ωi≥0。Ni，k（u）是節(jié)點矢量U＝ （u0，u1，…，un＋k＋1）決定的k次規(guī)范B樣條基函數(shù)，滿足德布爾－考克斯遞推公式：

在實際應(yīng)用中，節(jié)點矢量首尾的重復(fù)度取為k＋1次，即u0＝u1＝ … ＝uk＝0，un＋1＝un＋2＝…＝un＋k＋1＝1。這樣使得曲線具有同次有理Bezier曲線的端點幾何性質(zhì)，即在權(quán)因子ω1，ωn－1≠0的情況下，曲線在首末端點處分別與控制多邊形首末邊相切。

類似于NURBS曲線，一張p×q次的NURBS曲面可表示為

其中，Pi，j呈矩形陣列，形成一個控制網(wǎng)格；ωi，j是與控制頂點Pi，j有關(guān)的權(quán)因子，規(guī)定四角控制頂點處的權(quán)因子ω0，0，ωm，0，ω0，n，ωm，n都大于0，其余ωi，j≥0，且順序p×q個權(quán)因子不同時為0；Ni，p（u）和Nj，q（v）分別為沿u向的p次和沿v向的q次規(guī)范B樣條基函數(shù)，同樣由德布爾－考克斯遞推公式確定；節(jié)點矢量U、V與NURBS曲線節(jié)點矢量的定義形式相同，這樣可使曲面角點具有Bezier曲面的角點性質(zhì)。

2 立式捏合機攪拌槳曲面的NURBS表示

2.1 立式捏合機攪拌槳的幾何特征

立式捏合機兩個攪拌槳（實心槳和空心槳）均垂直安裝，其支撐軸位于槳葉的上部，兩攪拌槳軸線均偏離混合釜中心軸線一定距離，即具有一定的偏心距。其中，實心槳的偏心距es小于空心槳的偏心距ek。兩攪拌槳中心距為L，實心槳和空心槳直徑為d，基圓直徑為db，槳葉高度為h，混合釜直徑為D。攪拌槳葉刀尖與混合釜內(nèi)壁的最小間隙為c1，兩攪拌槳葉之間的最小間隙為c2，槳葉底端與混合釜內(nèi)底部之間的最小間隙為c3。在物料混合過程中，兩攪拌槳在混合釜內(nèi)做行星捏合運動，實心槳和空心槳的自轉(zhuǎn)角速度分別為ωs和ωk，公轉(zhuǎn)角速度均為ωH。為了使物料在混合釜內(nèi)自上而下、自下而上充分翻動混合，攪拌槳葉表面采用螺旋曲面，根據(jù)捏合原理，空心槳和實心槳螺旋角方向應(yīng)該相反，根據(jù)兩個槳葉的轉(zhuǎn)速比，螺旋角也應(yīng)滿足一定的比例關(guān)系。立式捏合機攪拌槳外形及主要幾何尺寸如圖1所示。

2.2 立式捏合機攪拌槳的NURBS表示

圖1 立式捏合機工作原理

立式捏合機攪拌槳曲面的NURBS表示過程是：首先根據(jù)攪拌槳水平截面輪廓線上的理論型值點，利用參數(shù)曲線方法反算出控制頂點，同理可以求出槳葉不同高度水平截面輪廓線的控制頂點，進(jìn)而利用蒙皮法反求出槳葉曲面的控制頂點，再用參數(shù)曲面方法構(gòu)造出曲面模型，最后在Visual C＋＋6.0編程環(huán)境下，調(diào)用 OpenGL 接口函數(shù)，生成攪拌槳水平截面輪廓的NURBS曲線和攪拌槳的NURBS曲面。具體構(gòu)造過程如下［5－7］：

（1）根據(jù)攪拌槳水平截面輪廓線上的型值點，依次生成各曲線段（自由曲線、圓弧段、直線段等）的NURBS表示；把階次小于3次的NURBS曲線段升階為3次；依次拾取各曲線段的控制頂點和權(quán)因子，相鄰兩段的控制頂點和權(quán)因子要滿足連接點處的幾何連續(xù)性要求，確定復(fù)合后整條曲線的控制頂點和權(quán)因子；計算各曲線段的控制多邊形長度在復(fù)合曲線的控制多邊形總長中所占有的比例，得出各曲線段在整體節(jié)點域u∈［0，1］中的定義區(qū)間，依次對各曲線段的節(jié)點做整體參數(shù)變換，得到復(fù)合曲線在整體參數(shù)下的節(jié)點矢量。得出的立式捏合機攪拌槳水平截面輪廓線的3次NURBS曲線如圖2所示。圖2中左半部分是實心槳的水平截面曲線，右半部分是空心槳中部的水平截面曲線。

圖2 立式捏合機攪拌槳截面輪廓的NURBS表示

（2）由于攪拌槳曲面是由其水平截面輪廓線以一定的螺旋升角順時針或逆時針掃掠生成的，所以在攪拌槳不同軸向位置處的水平截面曲線形狀是相同的，只是在其徑向平面內(nèi)繞軸向做了不同角度的旋轉(zhuǎn)。所以，在求出空心槳和實心槳某一個水平位置截面輪廓線的NURBS表示的基礎(chǔ)上，依據(jù)NURBS曲線的幾何不變性，可以求出槳葉若干個不同高度水平截面輪廓線的NURBS表示。

（3）在復(fù)雜曲面的外形設(shè)計中，曲面多是由一組給定的截面曲線蒙皮生成的。由于不同高度攪拌槳截面輪廓的NURBS曲線的節(jié)點和階次相同，不需要進(jìn)行截面曲線間的u向相容處理，只需計算相鄰截面NURBS曲線的對應(yīng)控制頂點之間的距離，故采用累積弦長參數(shù)法計算v向節(jié)點矢量。以各截面曲線的NURBS控制頂點為型值點反算出NURBS曲面的控制頂點。對于空心槳，由于上下實心部分曲面與中間空心曲面的交線滿足文獻(xiàn)［8］提出的條件，故在曲面拼接處至少是G1連續(xù)的。利用蒙皮法生成的立式捏合機攪拌槳的NURBS曲面如圖3所示。

圖3 立式捏合機攪拌槳曲面的NURBS表示

3 曲面誤差分析

為了驗證所生成的NURBS曲面的精確性，本文以計算得到的曲面理論坐標(biāo)點到NURBS曲面的最小距離作為評價NURBS曲面的誤差值［9］。最小距離越小，說明生成的NURBS曲面精度越高。關(guān)于點到NURBS曲面的投影和距離的算法，國外相關(guān)學(xué)者已經(jīng)做了很多研究［10－14］。通常地，計算一點P到曲面S（u，v）的最小距離通常可以轉(zhuǎn)化為計算點在曲面上的投影［15］，即要滿足：

由于點P在曲面S（u，v）上可能不存在或存在若干個投影點，即式（3）可能無解或存在若干個解，為此，本文提出了一種新的利用控制頂點投影法計算點到NURBS曲面最小距離的算法。該算法首先將控制頂點投影到NURBS曲面上，得到相應(yīng)的投影點；其次在NURBS曲面上依次生成相鄰?fù)队包c之間的曲線，這些曲線將NURBS曲面分割成若干個至少是G1連續(xù)的曲面片，使得點P在每個曲面片上至多存在一個投影點；最后通過判斷點P與每個曲面片邊界曲線之間的關(guān)系，可以得出點P在曲面片上是否存在投影點，計算點P與投影點（如果存在）之間的距離，即為點P到該曲面片的最小距離；同理可以求出點P到其他曲面片的最小距離，這些最小距離中的最小值即為點P到整個NURBS曲面的最小距離。下面給出該算法的具體實現(xiàn)步驟。

3.1 將控制頂點投影到NURBS曲面上

依次將每行控制頂點向NURBS曲面進(jìn)行投影，得出對應(yīng)的投影點。以控制頂點Pi，j為例，Pi，j在曲面上的投影點Pti，j應(yīng)該滿足式（3），即控制頂點Pi，j與其在曲面上的投影點Pti，j的線段分別和投影點Pti，j的u向、v向切矢之間的點積為0。特別地，對于第一行和最后一行的u向控制頂點，投影到對應(yīng)的NURBS曲面u向邊界曲線上，僅要求控制頂點與投影點的線段和投影點的u向切矢之間的點積為0；同樣地，對于第一列與最后一列的v向控制頂點，也投影到對應(yīng)的NURBS曲面v向邊界曲線上，僅要求控制頂點與投影點的線段和投影點的v向切矢之間的點積為0；如果控制頂點在NURBS曲面上，則該控制頂點與投影點重合；如果曲面上存在u向或v向G0連續(xù)的折痕曲線時，在折痕曲線兩側(cè)，找出一個距離折痕曲線最近的控制頂點，將該控制頂點投影到折痕曲線上，僅要求該控制頂點與投影點的線段與投影點的v向或u向切矢之間的點積為0。

在求投影點的計算過程中，目前一般利用Newton－Raphson迭代方法對式（3）進(jìn)行求解，此方法計算速度快，可以達(dá)到較高的計算精度，迭代公式為

假定上一個控制頂點Pi，j－1在曲面上的投影點Pti，j－1對 應(yīng) 的 節(jié) 點 值 為 （uti，j－1，vti，j－1），以（uti，j－1，vti，j－1）作為迭代初始值，其收斂準(zhǔn)則為

ε為計算精度，這里ε取1×10－10，求出控制頂點Pi，j在NURBS曲面上的投影點對應(yīng)的節(jié)點值（uti，j，vti，j），再利用式（2）求出Pi，j在曲面上的投影點Pti，j。類似地可以求出其他控制頂點在NURBS曲面上的投影點以及投影點對應(yīng)的節(jié)點值。

3.2 NURBS曲面片的生成

由于在曲面尖點、角點等只滿足G0連續(xù)的位置處有一個或多個控制頂點投影點，故在u向或v向折痕曲線上也有多個控制頂點的投影點。這樣，u向或v向相鄰?fù)队包c之間至少滿足G1連續(xù)。以控制頂點Pi，j、Pi＋1，j、Pi，j＋1和Pi＋1，j＋1在 NURBS曲面對應(yīng)的投影點Pti，j、Pti＋1，j、Pti，j＋1和Pti＋1，j＋1為例，投 影 點 對 應(yīng) 的 節(jié) 點 值 分 別 為 （uti，j，vti，j）、（uti＋1，j，vti＋1，j）、 （uti，j＋1，vti，j＋1） 和 （uti＋1，j＋1，vti＋1，j＋1），當(dāng) 節(jié) 點 變 量 （u，v）在 （uti，j，vti，j）與（uti＋1，j，vti＋1，j）之間以微小量（Δu，Δv）變化時，可以得到投影點Pti，j與Pti＋1，j之間在 NURBS曲面上的連續(xù)曲線，由多元函數(shù)的微分公式知，這些曲線至少是G1連續(xù)的。類似地可以得到其他相鄰?fù)队包c之間在NURBS曲面上的曲線。這些曲線把NURBS曲面分割成若干個至少是G1連續(xù)的曲面片，而曲面上尖點、角點或折痕曲線等G0連續(xù)的位置被分割在曲面片的邊界曲線上。

圖4形象化地說明了控制頂點投影點間的曲線將NURBS曲面分割成若干個曲面片的實現(xiàn)過程。

圖4 控制頂點投影法分割NURBS曲面

3.3 計算點P到每個曲面片的最小距離

對 于 曲 面 片Pti，jPti＋1，jPti，j＋1Pti＋1，j＋1， 由NURBS曲面的凸包性質(zhì)及曲面片邊界曲線的變差縮減性質(zhì)可知，曲面片內(nèi)在u向或v向上至多存在一條極值點曲線。由于曲面片內(nèi)的Su和Sv是連續(xù)的，所以曲面片的Su和Sv變化邊界值在曲面片的四條邊線上，又由于邊線也是連續(xù)的，所以曲面片的Su和Sv邊界值在曲面片的四個角點Pti，j、Pti＋1，j、Pti，j＋1和Pti＋1，j＋1處。將曲面片四個角點的Su和Sv分別投影到XY平面、XZ平面和YZ平面，得出四個角點的Su和Sv二維投影向量，依據(jù)連續(xù)函數(shù)極值存在的判別條件，可以點積四個角點Pti，j與Pti＋1，j、Pti，j＋1與Pti＋1，j＋1的Su，Pti，j與Pti，j＋1、Pti＋1，j與Pti＋1，j＋1的Sv在三個投影面上的向量，如果點積值均是正，說明曲面片的Su或Sv是單調(diào)連續(xù)的；如果存在負(fù)值，說明曲面片內(nèi)存在u向或v向的極值曲線，則利用二分法可以在相應(yīng)邊線上找出極值點，在曲面片上連接極值點，極值點間的曲線將曲面片分為兩片u向或v向單調(diào)連續(xù)的小曲面片。類似地，可以使得其他曲面片在u向和v向上是單調(diào)連續(xù)的。

由于曲 面 片Pti，jPti＋1，jPti，j＋1Pti＋1，j＋1或 其 子曲面片在u向和v向上是單調(diào)連續(xù)的，則在曲面片上f（u，v）與g（u，v）是連續(xù)函數(shù)，可以通過計算f（u，v）與g（u，v）在曲面片四個角點的值，判斷點P在曲面片內(nèi)的投影點Pt是否存在。具體過程如下：

f（u，v）如果在Pti，j與Pti＋1，j的值異號，或在Pti，j＋1與Pti＋1，j＋1的值異號，說明點P在曲面片上的u向投影存在；如果都是同號，說明點P在曲面片上的u向投影不存在。同樣地，g（u，v）如果在Pti，j與Pti，j＋1的值異號，或在Pti＋1，j與Pti＋1，j＋1的值異號，說明點P在曲面片上的v向投影存在；如果都是同號，說明點P在曲面片上的v向投影不存在。

如果點P在曲面片上的u向投影和v向投影同時存在，說明點P在曲面片上的投影點存在；如果點P在曲面片上的u向投影和v向投影同時不存在，說明點P在曲面片上的投影點不存在。

如果投影點存在，則點P到曲面片的最小距離就是點P與其投影點Pt之間的距離。如果投影點不存在，則點P到曲面片的最小距離點存在于曲面的四個邊線上。如果僅u向投影存在，則在相應(yīng)的曲線段Pti，jPti＋1，j或Pti，j＋1Pti＋1，j＋1上利用二分法求出u向投影點；同樣地，如果僅v向投影存在，則在相應(yīng)的曲線段Pti，jPti，j＋1或Pti＋1，jPti＋1，j＋1上利用二分法求出v向投影點；點P到曲面片的最小距離就是點P與u向投影點或v向投影點之間的距離。如果u向投影和v向投影都不存在，計算點P與曲面片四個角點之間的距離，其中最小值就是點P到曲面片的最小距離。

在投影點的計算上，以曲面片一個角點作為初始值，利用Newton－Raphson迭代方法求解式（3）。由于曲面片至少是G1連續(xù)的，在u向和v向是單調(diào)連續(xù)的，u向邊界曲線間Su變化不大，v向邊界曲線間Sv變化也不大，所以，以曲面片的角點為初始值，一般都可以使迭代過程收斂。

3.4 計算點P到NURBS曲面的最小距離

比較點P到各曲面片的最小值，找出全局最小值，即為點P到整個NURBS曲面的最小距離。取多組攪拌槳水平截面輪廓線上的理論坐標(biāo)點，利用控制頂點投影法計算理論點到攪拌槳NURBS曲面的最小距離，得出的最大距離為0.1mm，說明立式捏合機攪拌槳的NURBS曲面精度是較高的。

4 結(jié)語

本文在立式捏合機攪拌槳水平截面曲線的NURBS表示的基礎(chǔ)上，利用蒙皮曲面方法，生成了立式捏合機攪拌槳的NURBS曲面。在此基礎(chǔ)上，提出了一種新的控制頂點投影法計算點到NURBS曲面最小距離的方法，驗證了立式捏合機攪拌槳曲面的NURBS表示的精確性和可靠性，為攪拌槳的優(yōu)化設(shè)計打下了較好的基礎(chǔ)。同時，控制頂點投影法計算點到NURBS曲面最小距離的算法具有重要的實用價值，可以在CAD/CAM等領(lǐng)域內(nèi)進(jìn)行進(jìn)一步推廣和應(yīng)用。

［1］葉定友，張德雄.中國航天固體火箭推進(jìn)技術(shù)的發(fā)展［J］.中國航天，2002（12）：24－27.

Ye Dingyou，Zhang Dexiong.Evolution of Chinese Aerospace Solid Rocket Propulsion Technology［J］.Aerospace China，2002（12）：24－27.

［2］王正方，翟瑞清.立式捏合機攪拌槳的設(shè)計［J］.固體火箭技術(shù)，1993（1）：65－69.

Wang Zhengfang，Zhai Ruiqing.The design of Stirring Blades of Vertical Kneading Machine and Parametric Calculations［J］.Journal of Solid Rocket Technology，1993（1）：65－69.

［3］朱心雄.自由曲線曲面造型技術(shù)［M］.北京：科學(xué)出版社，2000.

［4］施法中.計算機輔助幾何設(shè)計與非均勻有理B樣條（CAGD ＆ NURBS）［M］.北京：高等教育出版社，2001.

［5］張樂年，周來水，周儒榮.基于復(fù)雜曲線的NURBS蒙皮曲面［J］.工程圖學(xué)學(xué)報，1996，17（2）：58－63.

Zhang Lenian，Zhou Laishui，Zhou Rurong.Research in Nurbs Skinning Surface Based on Complex Curve［J］.Journal of Engineering Graphics，1996，17（2）：58－63.

［6］Piegl L A，Tiller W.Algorithm for Approximate NURBS Skinning［J］.Computer－aided Design，1996，28（9）：699－706.

［7］Yoshimasa T.Skinning－Surface Generation Based on Spine－Curve Control［J］.The Visual Computer，2000，16（2）：134－140.

［8］于丕強，施錫泉.雙三次NURBS曲面的G1連續(xù)性條件［J］.大連理工大學(xué)學(xué)報，2004，44（3）：330－333.

Yu Piqiang，Shi Xiquan.G1Continuous Conditions for Bicubic NURBS Surfaces［J］.Journal of Dalian University of Technology，2004，44（3）：330－333.

［9］柯映林，李齊敏.基于截面輪廓的NURBS曲面重建［J］.中國機械工程，2005，16（12）：1083－1087.

Ke Yinglin，Li Qimin.Reconstruction of NURBS Surface from Spatial Cross－sections［J］.China Mechanical Engineering，2005，16（12）：1083－1087.

［10］DyllongE，Luther W.Distance Calculation Between a Point and a NURBS Surface［C］//4th International Conference on Curves and Surfaces.Saint－Malo，1999：55－62.

［11］Piegl L A，Tiller W.Parameterization for Surface Fitting in Reverse Engineering［J］.Computer－aided Design，2001，33（8）：593－603.

［12］Ma Y L，Hewitt W T.Point Inversion and Projection for NURBS Curve and Surface：Control Polygon Approach［J］.Computer－aided Geometric Design，2003，20（2）：79－99.

［13］Selimovic I.Improved Algorithms for the Projection of Points on NURBS Curves and Surfaces［J］.Computer Aided Geometric Design，2006，23（2）：439－445.

［14］PieglL A，Rajab K，Smarodzinava V，et al.Point－Distance Computations：A Knowledge－Guided Approach［J］.Computer－aided Design and Applications，2008，5（6）：855－866.

［15］Zhou J，Sherbrooke E C，Patrikalakis N M.Computation of Stationary Points of Distance Functions［J］.Engineering with Computers，1993，9（4）：231－246.