基于三向Dexel的GPU加速數控仿真方法

李一鵬 賈世宇 潘振寬 王吉強 李 冰

(青島大學計算機科學技術學院 山東 青島 266000)

0 引 言

數控仿真對數控加工有著重要的意義。數控仿真通過在虛擬環境中使用數控加工程序模擬加工過程來驗證其正確性，避免因數控加工程序錯誤而導致的工件質量下降和資源浪費。

現代數控仿真一般使用基于空間分割的方法來表示工件模型實體，其中由Hook[1]提出的Dexel方法廣泛地應用于國內外商業軟件和仿真系統中。基于Dexel的仿真方法的主要問題在于表示幾乎平行于Dexel方向的部分實體時精度較差，存在難以避免的階梯化現象與信息丟失的問題。其改進版本三向Dexel方法則可以在保留其簡單布爾操作以及較小空間占用優點的同時增加建模的精度。三向Dexel方法在數控仿真中應用的主要問題在于復雜的表面重構過程。目前主流的做法是利用三向Dexel模型與體素(Voxel)模型的相似性，將三向Dexel模型轉換成體素模型，再使用成熟的體素模型表面重構算法進行表面重構[2]。于珊[3]通過使用GPU加速的移動立方體(Marching Cube)[4-5]方法進行表面重構提高了數控仿真的速度。Peng等[6]提出的基于輪廓的表面重構方法及Ren等[7]提出的基于網格的表面重構方法相比于將三向Dexel模型轉換為體素模型，其算法復雜度較高，時間消耗較大，且難以使用并行技術進行加速。Jachym等[8]通過使用OptiX光線追蹤引擎來獲取更具真實感的重構表面，但難以實現流暢的實時仿真。Inui等[9]提出的四面柱(Quad Pillars)算法擁有較高的并行度和執行效率，但是僅適用于單向Dexel模型的表面重構。

相比于移動立方體方法，雙輪廓(Dual Contouring，DC)[10]方法通過使用厄米特數據(Hermite data,物體表面點及其法線數據)計算特征點，可以還原由于采樣問題丟失的特征(如尖銳棱角特征)。雙輪廓法在數控仿真中應用的主要限制在于其高昂的計算開銷，在CPU端進行的雙輪廓法往往難以滿足數控仿真對于實時性的要求，因此目前該方法還未在數控仿真中有所應用。

針對上述問題，本文提出一種并行化的基于三向Dexel的數控仿真方法，通過使用三向Dexel建模方法提高仿真精度，采用雙輪廓方法進行表面重構，并通過運用并行加速技術提高計算速度。相較于現存的數控仿真方法，在仿真的精度相同時，本文方法使用的GPU加速的雙輪廓算法能夠更準確地還原出數控加工工件的尖銳棱、角特征，為實際加工提供更加準確的仿真結果。同時通過應用多線程加速和GPU加速技術，使得中等配置(參考配置見仿真結果部分)的計算機也能實現流暢、準確的仿真，降低了對計算機硬件的要求。

1 仿真方法

數控仿真方法的結構如圖1所示。數控仿真程序在CPU端、GPU端有不同的分工。仿真程序從存儲工程文件中讀取信息后，在CPU端根據胚料信息進行實體建模，根據刀具參數計算每行NC代碼控制下刀具所掃過的空間體(下簡稱為掃掠體)；將工件實體與掃掠體進行相減操作；將所得三向Dexel模型轉換成適合進行表面重構的結構后發送到GPU端；由GPU對表面進行重構，并對重構所得多邊形表面進行圖形渲染。如此重復上述過程，直至NC代碼讀取完畢。

圖1 數控仿真方法示意圖

2 仿真算法實現

2.1 三向Dexel實體建模方法

數控仿真使用實體間的布爾減運算來描述數控加工中刀具切削工件的過程，為此需要對工件和刀具進行實體建模。單向與三向Dexel模型示意圖如圖2所示。單向Dexel模型通過記錄一組平行且等距的光線在實體內部的部分來表示實體，這種方式會造成如圖2(a)所示的信息丟失。三向Dexel模型則通過如圖2(b)所示的方式，在三個相互垂直的方向進行取樣來減少信息丟失，從而提高仿真精度。三向Dexel模型的結構由三組單向Dexel模型組成。

圖2 單向Dexel與三向Dexel模型示意圖

2.1.1單向Dexel模型的建模

單向Dexel模型在其XY平面上，沿X軸、Y軸正方向以仿真精度P為間隔分布著向Z軸正向發射的光線。每條光線與實體相交的交點都被保存下來，每兩個交點之間的線段代表著光線在實體內部的部分，稱為一個Dexel。構成Dexel的兩個交點中距離XY平面距離較遠的為Dexel的上端點，較近的為下端點。同一條光線上的Dexel端點依據距離發射點的距離按升序排列形成的鏈狀結構稱為Dexel鏈。這些Dexel鏈可以構成一個二維的網格(Grid)結構。假定網格結構沿X軸、Y軸方向分為nX、nY份，則該Dexel模型的分辨率為nX×nY，空間復雜度為O(nX×nY)。

為了壓縮Dexel模型所需的儲存空間，在記錄端點位置信息時僅記錄端點到XY平面的距離。需要獲取端點的確切坐標時，依據當前Dexel鏈在二維網格中的位置進行計算即可。坐標為(i,j)的Dexel鏈上端點的坐標(x,y,z)為：

x=OX+P×iy=OY+P×jz=OZ+dZ

(1)

式中：{OX,OY,OZ}為單向Dexel坐標系中的原點;P為仿真精度;dZ為端點到XY平面的距離。

2.1.2三向Dexel模型的建模

三向Dexel模型的結構由三個相互垂直方向的單向Dexel模型疊加得到。由于數控加工的仿真過程中僅包含布爾減運算，仿真中可能發生布爾運算的范圍可以限定在工件胚料的包圍盒中。建模時選擇在胚料包圍盒的X、Y、Z方向分別建立編號為0、1、2的三組單向Dexel模型。第i組單向Dexel模型的局部坐標系{O，eX，eY，eZ}為：

O=BmineX=e(i+1)mod3eY=e(i+2)mod3eZ=ei

(2)

式中：O為單向Dexel模型局部坐標系的原點;eX、eY、eZ為單向Dexel模型局部坐標系X、Y、Z軸的單位向量;Bmin為胚料包圍盒的最小點;e0、e1、e2為胚料包圍盒坐標系的X、Y、Z軸單位向量。

三向Dexel模型的分辨率NX×NY×NZ為：

NX=ceil[(Bmax.x-Bmin.x)/P]NY=ceil[(Bmax.y-Bmin.y)/P]NZ=ceil[(Bmax.z-Bmin.z)/P]

(3)

式中：Bmax為胚料包圍盒的最大點;P為仿真精度。

2.1.3三向Dexel模型的布爾運算

文獻[11]詳細描述了三向Dexel的并、交、差布爾運算，數控仿真中主要用到的為三向Dexel的減運算。三向Dexel模型之間進行布爾運算時，對雙方X、Y、Z方向的單向Dexel模型分別進行布爾運算，即為相對應光線上Dexel的減運算。數控仿真中運算根據刀具Dexel和工件Dexel的相對位置關系，共有圖3所示的6種情況。

圖3 Dexel間的相減運算

(1) 刀具Dexel與工件Dexel沒有任何相交，不需要修改工件Dexel。

(2) 與情況1相同，不需要修改工件Dexel。

(3) 工件Dexel的下端點到刀具Dexel上端點之間的部分被切除。

(4) 刀具Dexel的下端點到工件Dexel的上端點之間的部分被切除。

(5) 工件Dexel被完全刪除。

(6) 刀具Dexel的下端點到上端點之間的部分被切除，工件Dexel分裂成兩個Dexel。

可以看出，相減產生的新Dexel端點的位置為相應刀具Dexel端點的位置，法線與該刀具Dexel端點法線方向相反。

2.2 刀具掃掠體生成

刀具掃掠體的生成有兩種方法：第一種為近似的方法[12]，根據操作的軌跡位置、刀軸角度信息以一定的間距進行插值計算，將掃掠體離散為多個刀具實體；第二種方法為精確的方法[13-14]，通過刀具輪廓線和操作起止點信息精確地計算刀具的掃掠體。精確計算的主要優勢在于精度高于前者，劣勢在于計算需要更多時間。因此本文選擇近似的構造方法，以仿真精度P的一半為間隔插值生成刀具實例，將這些實例的集合作為當前操作所生成掃掠體的近似。

工件與刀具的減運算首先需要生成刀具實例的三向Dexel模型。為了減少建模的運算量，在三向Dexel模型的局部坐標系中生成該實例的包圍盒，僅在該范圍內進行模型生成。之后從工件三向Dexel模型中減去所生成的模型。其多線程加速運作機制如圖4所示。假定共有i個輔助線程，將其分別編號為1,2,…,n，創建輔助線程的主線程也作為0號線程參與運算，則第i號線程負責每個刀具實例包圍盒內第k(n+1)+i條Dexel鏈的構建與對應工件Dexel鏈的相減運算。

圖4 布爾減運算的并行運作機制

2.3 工件的表面重構

三向Dexel模型不適用于可視化，為了獲取適合圖形設備渲染的多邊形表面，需要對實體進行表面重構。本文方法在表面重構時，在CPU端將三向Dexel模型轉化為厄米特數據，并在GPU上進行八叉樹的構建、特征點的計算、八叉樹的化簡，以及重構表面的提取。

如圖5所示，由于三向Dexel模型的采樣精度有限，直接使用表面重構算法(如移動立方體算法，以下簡稱MC算法)進行表面重構會產生如圖5(c)所示的特征丟失。雙輪廓(簡稱DC)算法則可以通過使用由三向Dexel模型提取而來的厄米特數據，通過計算每個體素內部的特征點(圖5(d)中的白色點)來還原其丟失特征。表面重構算法的步驟如下：

(1) 在CPU中將三向Dexel數據轉換為適合構建八叉樹的數據，并將數據送入GPU；

(2) 在GPU上生成能夠完全容納網格數據的最小八叉樹；

(3) 并行處理每個葉子節點，對其特征點進行計算；

(4) 對八叉樹進行化簡，并提取重構表面。

圖5 移動立方體與雙輪廓表面重構

2.3.1網格信息與赫米特數據的提取

八叉樹的建立需要包含完整三向Dexel模型的網格信息，因此需要建立與三向Dexel模型同等分辨率的網格，網格格點及厄米特點的數據結構如圖6所示。為了節省顯存，格點數據中不保存完整的厄米特點數據，只保存索引值。索引數組的0、1、2號元素分別代表格點X、Y、Z方向的厄米特點。厄米特點的數據單獨封裝，以一維數組的形式傳遞給GPU。格點信息的提取通過遍歷三向Dexel模型中的Dexel完成，如圖7所示。首先為每個Dexel的端點分配索引值，并將其所經過的格點均設置為內部格點，在最靠近Dexel上下端點的內部格點處記錄對應Dexel端點的索引值。

圖6 格點與厄米特點的數據結構示意圖

圖7 從三向Dexel模型中提取網格信息

2.3.2八叉樹的建立

雙輪廓法結合自適應精度八叉樹可以簡化生成的表面，有效減少表面包含三角形的數量，提高渲染速度。本文使用的在GPU端運行的并行八叉樹生成算法[15]能夠在節約顯存的同時提升算法的運行效率，可以極大縮短八叉樹構建時間。所構建八叉樹的最大深度由傳入網格數據的尺寸決定，大小為能夠完全容納網格的最小八叉樹。由于GPU并不支持動態內存分配，每層節點所需要的存儲空間需要提前分配。該方法首先構造葉子節點，并從底層向上依次進行八叉樹的構造。若一個葉子節點為非空葉節點(8個角點中既包含內部格點又包含外部格點的葉節點)，則構建一個新的節點，使用原子計數器對其進行編號并使用查詢表來存儲子節點與父節點之間的關系。在建立八叉樹的過程中需要統計構建葉節點數量，用于為后續步驟會用到的二次誤差方程數據分配存儲空間。八叉樹構建完成后，使用構建的八叉樹來進行特征點的計算和八叉樹的簡化。

2.3.3特征點的計算

雙輪廓法中特征點的計算是通過求解葉節點的二次誤差方程(Quadratic Error Function，QEF) 得到的。QEF的定義如下：

(4)

式中：pi、ni代表著該葉子節點所對應體素邊上厄米特點的位置與法線，可以通過查詢該體素8個角點對應格點的索引數組得到。將式(4)轉化為矩陣形式并展開為：

E[x]=xTATAx-2xTATb+bTb

(5)

式中：矩陣A代表著體素邊上厄米特點的法線ni;向量b的第i個分量為nipi，求解方程最小值所得的x即是該體素的特征點。方程的具體解法參照文獻[7]。為了方便后續八叉樹的簡化及多邊形表面的提取，由體素構造出的QEF及其計算結果被單獨封裝在QEF類中，存儲于事先分配好的存儲空間中并由原子計數器為其分配索引值。構建完成后，在八叉樹節點中記錄其QEF的索引值以建立對應關系。所有葉節點的特征點計算完成后，對八叉樹進行簡化。

2.3.4八叉樹的化簡與重構表面的提取

對八叉樹進行化簡的目標是使用更少的三角形來表示重構表面的平緩區域，在保持圖形特征不改變的前提下減少渲染三角形的總量，以降低圖形設備的負擔。首先將同屬于一父節點的8個葉節點的QEF中各項相加，若計算相加所得QEF解得的誤差值小于一給定閾值，則可以對這些葉節點進行化簡。化簡后將葉節點的特征點修改為父節點的特征點，并將父節點標記為偽葉子節點，如此往復由底向上逐層進行化簡。

化簡完成后遍歷八叉樹的每個葉節點，對其進行表面提取。將每個葉節點的特征點提取到數組中，使用原子計數器為其分配索引值。最后根據葉節點之間的拓撲關系設置索引，完成多邊形表面的提取。

3 仿真結果

數控仿真系統軟件使用C++編寫CPU端程序，使用OpenCL編寫GPU端程序，使用OpenGL進行圖形渲染，運行硬件環境為Intel Core i7-3770 四核CPU，12 GB內存，2 GB顯存的NVIDIA GeForce GTX 660 Ti顯卡。

仿真實驗1測試不同建模方法的仿真效果，如圖8所示，其中：(a)、(c)為使用單向Dexel建模方法，(b)、(d)為使用三向Dexel建模方法。可以看出，應用三向Dexel建模方法后消除了使用單向Dexel建模方法時的階梯化問題。

圖8 不同建模方法的仿真效果對比

仿真實驗2為雙輪廓法對特征修復的驗證，圖9、圖10分別為該方法對雕像模型和葉輪模型丟失特征的修復效果。圖9(a)、圖10(a)為使用三向四面柱(Tri-Quadpillars，簡稱TQ，將文獻[9]中的算法擴展到三維得到)算法進行表面重構得到的結果圖；圖9(b)、圖10(b)為使用雙輪廓法進行表面重構得到的結果圖。可以看出，使用雙輪廓法后，由特征丟失導致的毛刺、鋸齒問題得到了修復，同時能夠更好地保存尖銳棱角等信息。

圖9 雙輪廓法對雕像模型丟失特征的修復效果

圖10 雙輪廓法對葉輪模型丟失特征的修復效果

仿真實驗3測試仿真方法加速前后的運行效率。表1為仿真實驗3中所使用的實驗對象，實驗共使用兩組數控加工程序：第1組為雕像的加工，第2組為葉輪的加工。

表1 仿真所用實驗對象

表2所示為多線程加速前后表面重構的平均時間。可以看出，多線程加速比單線程雙輪廓法表面重構快232%～289%。由于雙輪廓算法中八叉樹構建所耗費的時間在算法總耗時中占比最大，且該部分并行化程度較高，多線程加速能夠獲得較高的加速比，但由于雙輪廓算法本身時間開銷較大，僅在CPU端進行加速仍然難以滿足實時仿真的要求。

表2 多線程加速前后表面重構平均時間

表3展示了GPU加速前后表面重構平均時間。可以看出GPU加速比CPU下多線程加速雙輪廓法表面重構快5～15倍。完全在GPU端進行的雙輪廓算法在重構對象分辨率升高時，因網格數據與厄米特數據的增大而增加了GPU的讀寫壓力，造成了加速比的下降，但仍然能滿足實時仿真的要求。

表3 GPU加速前后表面重構平均時間

4 結 語

本文實現了一種并行化的基于三向Dexel的數控仿真方法。采用三向Dexel建模方法提升仿真的精度，使用雙輪廓法來獲取更加真實的重構表面。該方法的創新之處在于其通過雙輪廓算法修復因表面特征丟失造成的毛刺、鋸齒問題，并通過將多線程加速、GPU加速與雙輪廓算法結合起來，使得表面重構的速度相對于多線程的實現提升了5～15倍，能夠在滿足仿真實時性要求的前提下，修復由于三向Dexel模型采樣問題造成的特征丟失，從而提高仿真的真實度。