田素凱， 寧 濤， 陳志同

(北京航空航天大學(xué)機(jī)械工程及自動(dòng)化學(xué)院，北京 100191)

基于曲面精確表示的距離極值點(diǎn)的計(jì)算及在刀具干涉檢測(cè)中的應(yīng)用

田素凱， 寧 濤， 陳志同

(北京航空航天大學(xué)機(jī)械工程及自動(dòng)化學(xué)院，北京 100191)

針對(duì)基于曲面精確表示的剛體碰撞檢測(cè)中裁剪曲面距離極值點(diǎn)的求解問題，提出了利用平面向量場(chǎng)估計(jì)初始曲面距離極值點(diǎn)的方法，避免了曲面過度細(xì)分，討論了距離極值點(diǎn)滿足的微分幾何條件，給出了解析曲面/參數(shù)曲面、參數(shù)曲面/參數(shù)曲面、點(diǎn)/參數(shù)曲面和曲線/參數(shù)曲面的距離極值點(diǎn)迭代算法。實(shí)例驗(yàn)證分析了該算法的高效性和可靠性。

碰撞檢測(cè)；平面向量場(chǎng)；距離極值點(diǎn)；迭代

碰撞檢測(cè)(collision detection)是指判斷多個(gè)剛體(零件、機(jī)械設(shè)備等)所占據(jù)的空間是否有重疊區(qū)域，一直是計(jì)算機(jī)仿真、干涉檢查、數(shù)控加工、機(jī)器人控制和虛擬現(xiàn)實(shí)等學(xué)科的關(guān)鍵技術(shù)?，F(xiàn)有文獻(xiàn)中的碰撞檢測(cè)算法分為基于圖形[1]和基于圖像[2]兩大類檢測(cè)算法?；趫D形的碰撞檢測(cè)包括模型分解(modal decomposition)和空間分解(space decomposition)兩種，前者常用方法有 AABB(axis aligned bounding box)層次樹[3]、包圍球(bounding ball)層次樹、OBB(oriented bounding box)層次樹和K-DOPs法(discrete orientation polytopes)[4]等；后者利用GPU硬件的加速功能對(duì)多面體模型進(jìn)行凸分解，進(jìn)而提高碰撞檢測(cè)的效率。這些算法都是基于模型離散或空間分解的方法，在精確性和實(shí)時(shí)性方面有一定的欠缺，無法解決最大過切、欠切位置計(jì)算等問題。近年來，針對(duì)曲面精確表示

的剛體碰撞檢測(cè)算法的研究有一定的進(jìn)展，Choi[5]給出了空間橢球體的碰撞條件；Odehnal[6]研究了兩個(gè)圓環(huán)面的共線法矢；Lopes等[7]提出了二次曲面和超二次曲面的碰撞檢測(cè)數(shù)學(xué)框架；Zeng等[8]給出了支持多種幾何體的碰撞檢測(cè)算法庫。由于CAD系統(tǒng)支持的曲面有十多種之多，而以上方法僅適用于幾種簡(jiǎn)單解析曲面，因此基于曲面精確表示的剛體碰撞檢測(cè)算法尚待研究。

在曲面精確表示的剛體碰撞檢測(cè)算法中，曲面之間的距離極值點(diǎn)對(duì)應(yīng)最大干涉(過切、欠切)位置，相關(guān)計(jì)算是碰撞檢測(cè)算法的關(guān)鍵?，F(xiàn)有曲面間最小距離的求解主要有離散法和解析法，算法的實(shí)現(xiàn)效率及其精確度是兩個(gè)重要的性能指標(biāo)。Gilbert和Foo[9]的算法首先將曲面離散為三角片，然后再計(jì)算兩三角片的最近距離，進(jìn)而獲取兩曲面的近似距離。當(dāng)精度要求比較高時(shí)，該算法的效率較低。Johnson和Cohen[10]提出了一種計(jì)算曲面間最小距離的框架，該框架以包圍盒細(xì)分為基礎(chǔ)，算法效率較低。Sohn等[11]提出了一種線性幾何(line geometry)的方法求解曲面間距離的計(jì)算，該算法減少方程組和未知量的數(shù)量，提高了算法效率，尤其適用于求解二次曲面距離極值點(diǎn)。劉浩等[12]提出了一種計(jì)算雙二次 NURBS曲面的最近距離點(diǎn)。文獻(xiàn)[11-12]都沒有給出一般曲面間距離求解的具體算法。陳麗萍等[13]從極值點(diǎn)的幾何條件出發(fā)，利用賦范空間投影法迭代，求取自由曲線到曲面的極值距離，但該算法并沒有給出初始點(diǎn)的計(jì)算方法。

本文結(jié)合平面向量場(chǎng)提出了一種求解裁剪曲面的初始距離極值點(diǎn)的算法[14-15]，該算法利用裁剪曲面的三角剖分網(wǎng)格，依據(jù)平面向量場(chǎng)的Poincare公式，估算出單張三角片內(nèi)部是否有距離極值點(diǎn)。分別針對(duì)解析曲面/參數(shù)曲面、參數(shù)曲面/參數(shù)曲面等情況，給出了距離極值點(diǎn)迭代求精算法。同時(shí)，還討論了點(diǎn)到曲面距離、曲線到曲面距離的求解算法。

1 初始距離極值位置點(diǎn)的估計(jì)

在數(shù)控加工幾何仿真等領(lǐng)域，需要計(jì)算平面、二次曲面及圓環(huán)面等代數(shù)曲面與參數(shù)曲面距離極值位置。如何利用裁剪曲面的三角網(wǎng)格、基于平面向量場(chǎng)技術(shù)估計(jì)初始距離極值點(diǎn)，需引入幾個(gè)基本概念：

(1) 平面向量場(chǎng)[15]v是連續(xù)函數(shù) v:D →R2，即對(duì)于平面區(qū)域D中的任何點(diǎn)p，v定義平面向量v(p)與之對(duì)應(yīng)。向量場(chǎng)在自然界廣泛存在，例如引力場(chǎng)、重力場(chǎng)、電磁場(chǎng)和流體的速度場(chǎng)等。在微分動(dòng)力系統(tǒng)中，微分方程的穩(wěn)定性用向量場(chǎng)刻劃[16]。對(duì)于平面區(qū)域D中一點(diǎn) p，若v(p)=0，則稱p為平面向量場(chǎng)v的臨界點(diǎn)(critical point)。假設(shè)p為平面向量場(chǎng)v的孤立臨界點(diǎn)，N為p的領(lǐng)域，且在N及其邊界?N上沒有其他臨界點(diǎn)。對(duì)于邊界?N 上的每個(gè)點(diǎn)q，平面向量場(chǎng)v都定義非0向量v(q)與之對(duì)應(yīng)。假設(shè)一個(gè)手持指南針的人，沿逆時(shí)針方向在邊界上行走，指南針總保持與向量場(chǎng)同向。當(dāng)回到出發(fā)點(diǎn)時(shí)，指南針亦回到初始位置。所以在此過程中，指南針相對(duì)于底盤必轉(zhuǎn)動(dòng)整數(shù)圈。逆時(shí)針轉(zhuǎn)一圈記為1，順時(shí)針轉(zhuǎn)一圈記為–1。逆時(shí)針和順時(shí)針?biāo)D(zhuǎn)圈數(shù)的代數(shù)和稱為旋轉(zhuǎn)數(shù)，如圖 1所示。向量場(chǎng)的旋轉(zhuǎn)數(shù)可以用 Poincare公式[16]計(jì)算：，其中C是平面區(qū)域中的封閉曲線，在 C上向量場(chǎng)▽?duì)?= (φu,φv)非退化。有定理：設(shè)v是連續(xù)的平面向量場(chǎng)，C為 v上封閉簡(jiǎn)單光滑曲線，若旋轉(zhuǎn)數(shù)W(v，C)不為0，則閉曲線C內(nèi)至少有一個(gè)臨界點(diǎn)。臨界點(diǎn)和旋轉(zhuǎn)數(shù)反映向量場(chǎng)v的拓?fù)浣Y(jié)構(gòu)。當(dāng)向量場(chǎng)v發(fā)生微小變動(dòng)時(shí)，其拓?fù)浣Y(jié)構(gòu)一般保持不變，所以向量場(chǎng)的拓?fù)浣Y(jié)構(gòu)被認(rèn)為是一種穩(wěn)定的結(jié)構(gòu)，能反映向量場(chǎng)的內(nèi)在性質(zhì)。

圖1 平面向量場(chǎng)的旋轉(zhuǎn)數(shù)

Kriezis等[14]用有向距離函數(shù)定義平面向量場(chǎng)。如圖2所示，假設(shè)點(diǎn)p是曲面 S1(u,v)上的任意一點(diǎn)，則q為曲面 S2(s,t)上距離點(diǎn)p最近的一點(diǎn)。當(dāng)曲面 S1(u,v)上的點(diǎn) p變動(dòng)時(shí)，得到映射M(p )=q。對(duì)于曲面 S1(u,v)上任意一點(diǎn)p，S1(u,v)

有向距離函數(shù)定義為：φ(u,v) =(n2(M (p)), p -M(p) )，其中( , )表示內(nèi)積n2，n2表示曲面S2(s,t)的單位法向量。有向距離函數(shù) φ( u,v)的梯度定義了u、v參數(shù)域上的向量場(chǎng)，在 φ( u,v)的有效定義域內(nèi)，有 ▽?duì)?(u,v) =((n2,S1u),(n2,S1v))。所以對(duì)于 φ( u,v)的任何臨界點(diǎn)(u, v)，有 (n2,S1u) =0, (n2 ,S 1v)=0。

圖2 有向距離函數(shù)的定義

(2) 討論距離極值點(diǎn)的估算方法。在以上有向距離函數(shù)的定義中，假設(shè) S1(u,v)、 S2(s,t)是裁剪參數(shù)曲面， ?1=(V1,T1)、?2=(V2,T2)分別是兩裁剪曲面的三角片集合，其中 Vi表示所有頂點(diǎn)集合，Ti表示所有三角片集合，i=1,2。利用 S2(s,t)的三角片集合，使用點(diǎn)/三角片最近距離算法計(jì)算出頂點(diǎn)集合 V1中每個(gè)頂點(diǎn)對(duì)應(yīng)的 φ(u,v)。對(duì)于三角片集合T1中的任何三角形，假設(shè)v1,v2,v3是其3個(gè)頂點(diǎn)對(duì)應(yīng)的平面向量場(chǎng)φ(u,v)的3個(gè)矢量，如果對(duì)某個(gè)i=1,2,3成立，則該三角形內(nèi)可能存在距離極值點(diǎn)，否則計(jì)算v1和v2，v2和v3，及v3和v1的夾角之和α，若，則該三角形內(nèi)可能存在距離極值點(diǎn)。這里“v1和 v2夾角”是指向量v1按逆時(shí)針轉(zhuǎn)到另向量v2的夾角，如圖3所示。

如果 S2(s,t)退化為一個(gè)三維點(diǎn)，則無須利用?2=(V2,T2)，可直接計(jì)算出有向距離φ( u,v)在三角片網(wǎng)格 ?1=(V1,T1)上的取值，進(jìn)而得到對(duì)應(yīng)的平面向量場(chǎng)；如果 S2(s,t)是二次曲面或圓環(huán)面，則由于 S1(u,v)上的任意點(diǎn)p到 S2(s,t)上最近點(diǎn)q容易計(jì)算，所以可以用解析方法得到▽?duì)?u,v) =((n2,S1u),(n2,S1v))。后續(xù)對(duì)距離極值點(diǎn)的估計(jì)同上給出的參數(shù)曲面到參數(shù)曲面距離極值點(diǎn)的算法。

圖3 向量場(chǎng)中向量夾角的計(jì)算

2 距離極值點(diǎn)的迭代算法

以上提出了一種針對(duì)代數(shù)曲面和裁剪曲面的初始極值位置點(diǎn)的估計(jì)方法。表 1針對(duì)具體的代數(shù)曲面，構(gòu)造對(duì)應(yīng)的迭代公式，可實(shí)現(xiàn)初始極值位置點(diǎn)的迭代求精。

說明：

(1) 圖示中S(u,v)代表參數(shù)曲面S(u,v)，n(u,v)代表曲面上某點(diǎn)處的法向量，G(t)為參數(shù)曲線。

(2) 對(duì)于球面與參數(shù)曲面的情況，給定球心為O(x0,y0,z0)，及參數(shù)曲面S(u,v)，曲面上各點(diǎn)處單位法矢為 n( u,v)。曲面 S(u,v)的法矢方程為L(zhǎng)(s) = S(u,v) +s ·n( u,v)，要找到這樣的u, v參數(shù)，使L(s)通過球心O。在L(s)上找一點(diǎn)距離O最近，由于

所以L(s)與球心O的最近距離為 (V,V )-(n,V )2，其中 V =S-O?；谝陨戏治隽畹繕?biāo)函數(shù)為f(u,v) = (V,V ) -(n,V )2，求導(dǎo)得

表1 迭代公式

(3) 對(duì)于圓柱面與參數(shù)曲面的情況，令給定的圓柱面軸線為 P(t)=P+tD，其中P是圓柱軸線上一點(diǎn) (x0,y0,z0)，D是圓柱軸線單位方向，給定參數(shù)曲面 S(u,v)，曲面上各點(diǎn)處單位法矢為 n(u,v)。曲面 S(u,v)的法矢方程為 L(s) = S(u,v) +s ·n(u,v)，要使L(s)也是圓柱面的法矢，則L(s)必須與圓柱面軸線P(t)交于一點(diǎn)，且滿足(n(u,v),D) =0?？紤]到兩直線L(s)和P(t)上兩點(diǎn)距離之平方，令

其中， V =S-P ，對(duì)g(s,t)求偏導(dǎo)得

所以兩直線的最近距離點(diǎn)對(duì)滿足

解此方程得

其中， Δ=1 -(D,n)2，將上式中的s, t代入g(s, t)得直線L(s)與P(t)的最近距離為

所以迭代目標(biāo)函數(shù)定義為

f(u,v)為L(zhǎng)(s)與P(t)最近距離的平方加上L(s) 與P(t)夾角余弦的平方。實(shí)際上，當(dāng)(n,D )→0時(shí)，Δ→ 1，于是迭代目標(biāo)函數(shù) f(u,v)可簡(jiǎn)化為

對(duì)于 f(u,v)求偏導(dǎo)得

利用上面的導(dǎo)數(shù)構(gòu)造迭代公式：

(4) 對(duì)于圓環(huán)面與參數(shù)曲面情況，令給定的圓環(huán)面中心圓方程為： Q (θ) =O+ rX cosθ + rY sinθ，給定參數(shù)曲面 S(u,v)，曲面上各點(diǎn)處單位法矢為n( u,v)。由于圓環(huán)面與參數(shù)曲面的極值位置點(diǎn)對(duì)連線垂直于曲面一階偏導(dǎo)數(shù)Su和Sv，且與中心圓的切矢垂直，所以定義目標(biāo)函數(shù)為

3 距離極值點(diǎn)算法應(yīng)用

在數(shù)控加工幾何仿真領(lǐng)域，為計(jì)算過切量和確定過切位置，需要計(jì)算平面、二次曲面及圓環(huán)面等代數(shù)曲面與參數(shù)曲面極值位置。本文實(shí)現(xiàn)了平面、圓柱面、球面和圓環(huán)面對(duì)參數(shù)曲面的極值位置計(jì)算算法，并結(jié)合平面向量場(chǎng)技術(shù)和裁剪曲面剖分，實(shí)現(xiàn)了三坐標(biāo)加工過切仿真功能。測(cè)試實(shí)例為：B樣條曲面u, v參數(shù)的次數(shù)都是5次，控制頂點(diǎn)是 30×30個(gè)，x, y, z最小最大范圍分別為[–134,134]×[–133,136]×[–172,77]，單位為mm，如圖4所示。以φ10球頭刀為刀具采用變步長(zhǎng)方法計(jì)算刀位軌跡，其逼近 B樣條曲面等距面的精度為0.001，刀位點(diǎn)總計(jì)64 976個(gè)。裁剪曲面的三角片總計(jì)8 818個(gè)，三角片結(jié)點(diǎn)共4 458個(gè)。

圖4 B樣條曲面及其裁剪部分

在仿真算法中，兩個(gè)相鄰刀位點(diǎn)之間的刀具包絡(luò)體分為球體和圓柱體兩部分處理，如圖 5所示。由于此算法中三角片和刀位點(diǎn)數(shù)據(jù)量很大，為提高搜索與兩個(gè)相鄰刀位點(diǎn)之間的刀具包絡(luò)體鄰近的三角片，在加工坐標(biāo)系xoy坐標(biāo)面上劃分正方形網(wǎng)格單元(邊長(zhǎng)等于刀具半徑)，構(gòu)建三角片相鄰關(guān)系索引表。該方法大幅減少了三角片的查詢計(jì)算量和時(shí)間。

圖5 單條插補(bǔ)段的刀具運(yùn)動(dòng)包絡(luò)體

仿真環(huán)境是Windows 7，處理器是Intel Core Duo P8600，主頻是2.4 GHz，內(nèi)存2 G。仿真結(jié)果如圖6~7所示。最大過切量是0.000 86 mm，耗時(shí)

46.4 s，見表2，過切基本出現(xiàn)在參數(shù)曲面曲率大的位置，且為刀具包絡(luò)的圓柱部分與零件曲面干涉，刀具包絡(luò)的球面部分與零件無干涉。仿真結(jié)果表明本文實(shí)現(xiàn)的代數(shù)曲面/參數(shù)曲面極值位置計(jì)算算法穩(wěn)定可靠，能夠?qū)崿F(xiàn)加工過切仿真等CAD/CAM功能。

圖6 裁剪曲面及其刀位軌跡

圖7 放大10 000倍的過切情況

表2 測(cè)試結(jié)果

4 結(jié) 論

現(xiàn)有成熟的碰撞檢測(cè)算法都是基于模型離散或空間分解的算法，難以解決給出最大干涉量和干涉位置等問題，而基于曲面精確表示的碰撞檢測(cè)則可以解決這樣的問題，但相關(guān)算法更為復(fù)雜。本文研究了點(diǎn)、線、面到裁剪曲面的距離極值點(diǎn)求解問題，給出了一種通用算法思路。該算法基于裁剪曲面的三角片網(wǎng)格、利用平面向量場(chǎng)估算初始距離極值點(diǎn)，再通過Newton迭代精確的距離極值點(diǎn)。對(duì)于二次曲面及圓環(huán)面/參數(shù)曲面的情況，由于二次曲面及圓環(huán)面的幾何特性，相關(guān)平面向量場(chǎng)以及Newton迭代的計(jì)算大為簡(jiǎn)化。后續(xù)研究?jī)?nèi)容主要包括二次曲面及圓環(huán)面之間的距離極值點(diǎn)的求解問題，以及基于曲面精確表示的碰撞檢測(cè)系統(tǒng)的框架設(shè)計(jì)。

[1] Lin M C, Gottschalk S. Collision detection between geometric models: a survey [C]//Proceedings of IMA Conference on Mathematics of Surfaces, Cambridge: Andalucia En La Historia, 1998: 37-56.

[2] 范昭煒, 萬華根, 高曙明. 基于圖像的快速碰撞檢測(cè)算法[J]. 計(jì)算機(jī)輔助設(shè)計(jì)與圖形學(xué)學(xué)報(bào), 2002, 14(9): 805-809.

[3] 陳 華. 確定任意形狀物體最小包圍盒的一種方法[J].圖學(xué)學(xué)報(bào), 2010, 31(2): 49-53.

[4] Klosowski J T, Held M, Mitchell J S B, et al. Efficient collision detection using bounding volume hierarchies of k-DOPs [J]. IEEE Transactions on Visualization & Computer Graphics, 1998, 4(1): 21-36.

[5] Choi Y K. Collision detection for ellipsoids and other quadrics [D]. Hong Kong: The University of Hong Kong, 2008.

[6] Odehnal B. Common normals of two tori [J]. Journal for Geometry & Graphics, 2005, (9): 51-65.

[7] Lopes D S, Silva M T, Ambrósio J A, et al. A mathematical framework for rigid contact detection between quadric and superquadric surfaces [J]. Multibody System Dynamics, 2010, 24(3): 255-280.

[8] Zeng L, Ning T, Xi P. Research and application of general collision detection simulation platform [M]//The 19th International Conference on Industrial Engineering and Engineering Management. Berlin: Springer Berlin Heidelberg, 2013: 1315-1323.

[9] Gilbert E G, Foo C P. Computing the distance between general convex objects in three-dimensional space [J]. IEEE Transactions on Robotics & Automation, 1990, 6(1): 53-61.

[10] Johnson D E, Cohen E. A framework for efficient minimum distance computations [C]//IEEE International Conference on Robotics & Automation. New York: IEEE Press, 1998: 3678-3684.

[11] Sohn K A, Juttler B, Kim M S, et al. Computing distances between surfaces using line geometry [C]// Computer Graphics and Applications, 2002. Proceedings. Pacific Conference on. New York: IEEE Press, 2002: 236-245.

[12] 劉 浩, 唐月紅, 廖文和. 雙二次NURBS曲面間的最短距離[J]. 計(jì)算機(jī)輔助設(shè)計(jì)與圖形學(xué)學(xué)報(bào), 2003, 15(10): 1298-1302.

[13] 陳麗萍, 陳 燕, 胡德金. 一種快速完備的自由曲線和曲面間最短距離求取算法[J]. 上海交通大學(xué)學(xué)報(bào), 2003, 37(S1): 41-44.

[14] Kriezis G A, Patrikalakis N M, Wolter F E. Topological and differential-equation methods for surface intersection [J]. Computer-Aided Design, 1992, 24(1): 41-55.

[15] 寧 濤, 馬德昌, 王亞平, 等. 平面向量場(chǎng)與曲率分析在曲面求交中的應(yīng)用[J]. 計(jì)算機(jī)學(xué)報(bào), 1997, 20(12): 1074-1080.

[16] 張錦炎, 錢 敏. 微分動(dòng)力系統(tǒng)導(dǎo)引[M]. 北京: 北京大學(xué)出版社, 1991: 41-47.

On Computation of Distance Extremum Points Based on Exact Surface Representation and Its Application in Tool Interferance Detection

Tian Sukai, Ning Tao, Chen Zhitong

(School of Mechanical Engineering and Automation, Beihang University, Beijing 100191, China)

Aiming at the problem of solving distance extreme points on trimmed surfaces in collision detection of rigid bodies which are represented by exact surface, a method is presented to obtain the initial distance extreme points using the plane vector field, which can avoid the excessive subdivision of the surface. The differential geometry conditions of the distance extreme points are discussed, and iterative formulas for obtaining distance extreme points of the analytic surface/parametric surface and parametric surface/parametric surface are given in this paper. The efficiency and reliability of the algorithm are verified by an example.

collision detection; plane vector field; distance extreme points; iteration

V 260.5；TP 391

10.11996/JG.j.2095-302X.2016050620

A

2095-302X(2016)05-0620-06

2016-03-24；定稿日期：2016-05-13

國(guó)家科技重大專項(xiàng)–高檔數(shù)控機(jī)床與基礎(chǔ)制造裝備科技重大專項(xiàng)課題(2013ZX04011031)

田素凱(1989–)，男，山東棗莊人，碩士研究生。主要研究方向?yàn)橛?jì)算機(jī)輔助設(shè)計(jì)、CAD技術(shù)。E-mail：tskdqq@163.com

寧 濤(1966–)，男，遼寧丹東人，副教授，博士。主要研究方向?yàn)橛?jì)算機(jī)輔助設(shè)計(jì)、CAD/CAM技術(shù)。E-mail：ningtao@buaa.edu.cn