基于SPH方法的流體粒子與軟體碰撞檢測①

施 鵬,陳 飛,廖晉民

1(中南財經(jīng)政法大學(xué),武漢 430073)

2(中南大學(xué) 湘雅二醫(yī)院,長沙 410008)

3(陸軍勤務(wù)學(xué)院 戰(zhàn)爭經(jīng)濟(jì)實驗室,重慶 404000)

虛擬現(xiàn)實(Virtual Reality,VR)技術(shù)在近年來結(jié)合各個行業(yè)有了跨越式的發(fā)展[1],在醫(yī)學(xué)培訓(xùn)領(lǐng)域中應(yīng)用也十分廣泛,特別是利用VR 技術(shù)開發(fā)的虛擬手術(shù)系統(tǒng)(Virtual Surgery,VS).虛擬手術(shù)系統(tǒng)以其操作可重復(fù)性、逼真的沉浸感而廣受醫(yī)療培訓(xùn)人員的歡迎[2].在虛擬手術(shù)系統(tǒng)中,流血模擬是其重要組成部分,其伴隨著手術(shù)過程頻繁出現(xiàn),是虛擬手術(shù)系統(tǒng)逼真性和完整性的重要體現(xiàn).為了確保流血模擬的真實性,需要解決流血過程中與軟組織的碰撞檢測.

碰撞檢測的對象大致可以分為剛體、軟體和流體.層次包圍盒[3]作為一類經(jīng)典的碰撞檢測算法,大量應(yīng)用于處理剛體間的碰撞檢測,如游戲設(shè)計、計算機(jī)輔助設(shè)計等領(lǐng)域.構(gòu)成層次包圍盒樹的包圍盒類型有很多種,其中較為著名的有軸向包圍盒[4],方向包圍盒[5],包圍球,離散有向多面體包圍盒[6],凸包包圍盒等[7].

近年來,國內(nèi)外研究者對剛體與軟體以及軟體間的碰撞檢測做了大量的研究,Wang TT 等[8]實現(xiàn)了一種可形變四面體之間的連續(xù)碰撞檢測,該算法能處理軟體間精確連續(xù)碰撞檢測的同時,魯棒性和實時性都很好.但由于流體不同于軟體,其拓?fù)浣Y(jié)構(gòu)變化較大,更新結(jié)構(gòu)的代價也較大,不能簡單的使用上述算法.在流體與剛體或軟體碰撞檢測方面,Vanhoey 等[9]提出一種歐拉網(wǎng)格流體與彈簧振子剛體模型之間的碰撞檢測,但其剛體模型不能太精細(xì),否則實時性受到很大的影響.Chentanez 等[10]提出一種SPH 流體粒子與有限元四面體網(wǎng)格之間的碰撞檢測算法,通過在血管壁上附著SPH 粒子并對其進(jìn)行約束從而實現(xiàn)血管與血管壁碰撞檢測的效果,該算法能處理SPH 流體與有限元剛體或軟體之間的碰撞檢測,但對軟體的動力學(xué)模型有相應(yīng)的要求,適應(yīng)性不強(qiáng),特別是對虛擬手術(shù)中流血與軟組織器官的碰撞檢測,限制了軟組織器官動力學(xué)模型的選擇范圍.宗智等[11]將SPH 模擬的流體粒子與質(zhì)點彈簧模型(Mass-Spring Model,MSM)之間進(jìn)行碰撞檢測,實現(xiàn)了血液流動與血管壁形變之間的相互作用,從而模擬出了血液在血管里面流動的現(xiàn)象.但必須使用MSM 模型作為軟件組織的動力學(xué)模型的限制,使其無法在流體與軟體間的碰撞檢測中廣泛推廣.

本文針對流體與軟體碰撞檢測的問題,提出一種基于空間劃分的SPH 流體粒子與軟體碰撞檢測算法,該算法首先對虛擬手術(shù)系統(tǒng)中流血粒子與軟組織器官碰撞檢測的特點進(jìn)行分析,在此基礎(chǔ)上進(jìn)行SPH 流體粒子與軟體組織表面的碰撞檢測研究.該算法克服了以往研究中對軟體組織模型選擇的局限性,而且實現(xiàn)了精確性和實時性上的優(yōu)化,以滿足虛擬手術(shù)流血模型對碰撞檢測要求.

1 虛擬手術(shù)流血碰撞檢測的特點

在虛擬手術(shù)系統(tǒng)中,流血模擬是關(guān)系著整個系統(tǒng)真實性和沉浸感非常重要的一環(huán),而處理流血與器官或者皮膚表面的碰撞檢測是流血模型中必不可少的部分.本文采用SPH 算法作為流血模擬的流體動力學(xué)模型.SPH 算法是一種對流體力學(xué)方程—Navier-Stokes方程進(jìn)行離散化求解的拉格朗日粒子法.在本文中,主要關(guān)注SPH 流血粒子與虛擬軟組織器官表面三角面片之間的碰撞檢測.這類碰撞有以下幾個特點:

1)流血粒子的相鄰關(guān)系和空間位置在每一個時間步內(nèi)變化都很大,要求進(jìn)行碰撞檢測的結(jié)構(gòu)更新算法能快速適應(yīng)和處理流體的拓?fù)浣Y(jié)構(gòu)變化,這就決定了很難采用層次包圍盒樹對拓?fù)浣Y(jié)構(gòu)進(jìn)行更新.

2)虛擬器官軟組織不同于剛體,其形態(tài)和相對位置會隨著形變和切割而改變,這就要求碰撞檢測算法能夠處理和感知到這種變化.

3)流血粒子與虛擬器官軟組織之間的碰撞檢測對精確性也需要相當(dāng)?shù)谋ＷC,否則流血粒子很有可能穿透到虛擬器官軟組織內(nèi)部,模擬效果也大打折扣.

2 算法概述

在處理大規(guī)模三維場景下的碰撞檢測時,不同碰撞檢測算法的共同特點是進(jìn)行逐層次、由粗到細(xì)的碰撞檢測,即先將明顯不可能相交的空間進(jìn)行快速排除,再對已縮小范圍中的物體進(jìn)行更為精確的判斷.本文提出的碰撞檢測算法能結(jié)合SPH 算法中最近相鄰粒子搜索(Nearest Neighboring Particle Search,NNPS)過程進(jìn)行碰撞檢測空間的均勻劃分,并通過均勻劃分的空間快速確定碰撞檢測范圍,在保證空間劃分時間和三角面片定位時間適當(dāng)?shù)耐瑫r,盡可能的減少精確的基元間碰撞檢測的次數(shù).

從上述數(shù)據(jù)來看，配電網(wǎng)重構(gòu)后，該網(wǎng)絡(luò)的功耗從重構(gòu)前的511.4 kW下降到了466.1 kW，從圖5中來看，網(wǎng)絡(luò)中的第9個節(jié)點為電壓最低點，重構(gòu)前為0.968 p.u，重構(gòu)后為0.971 p.u?？梢娊?jīng)過算法重構(gòu)將系統(tǒng)的功耗降低又改善了電壓質(zhì)量。驗證混合GA-PSO對配網(wǎng)重構(gòu)的優(yōu)越性。

2.1 基于空間劃分的碰撞檢測算法

基于空間劃分的碰撞檢測算法是一種確定粗略碰撞檢測范圍的算法:將空間劃分為不同區(qū)域并測試基元是否在同一空間區(qū)域中相交.這種方法大大降低了基元測試的數(shù)量[12].對于三維空間中可形變物體的碰撞檢測,一般采用均勻網(wǎng)格對空間實施覆蓋,其優(yōu)勢在于確定任意點所在的均勻網(wǎng)格單元序號十分簡單:只需利用世界坐標(biāo)系的值除以網(wǎng)格單元長度即可.一般可將基于空間劃分的碰撞檢測算法分為3 個階段,首先確定碰撞檢測空間范圍并在該空間上劃分空間網(wǎng)格單元,再根據(jù)劃分好的網(wǎng)格單元確定需測試基元所在的網(wǎng)格序號,最后遍歷所有網(wǎng)格單元,對網(wǎng)格單元中包含兩種以上基元的網(wǎng)格進(jìn)行基元間的精確碰撞檢測.

2.2 基于空間劃分的SPH 流體粒子與虛擬器官軟組織碰撞檢測算法

本文提出的碰撞檢測算法對傳統(tǒng)的空間劃分碰撞檢測算法進(jìn)行改進(jìn)的地方是:利用SPH 算法NNPS 過程中劃分的均勻網(wǎng)格空間對流體粒子和軟組織三角面片進(jìn)行定位.相對于傳統(tǒng)的基于空間劃分的碰撞檢測算法,本算法節(jié)省下了進(jìn)行空間劃分以及SPH 流血粒子定位的資源消耗,很大程度上的降低了算法的時間復(fù)雜度[13].

在本節(jié)中,首先介紹SPH 算法進(jìn)行流血模擬的過程,并針對SPH 算法進(jìn)行鄰域搜索時間復(fù)雜度過高的問題,提出基于均勻空間劃分的NNPS 過程,再結(jié)合NNPS 過程進(jìn)行SPH 粒子和器官表面三角面片的定位,然后進(jìn)行基元間的碰撞檢測計算.

2.2.1 基于SPH 算法的模擬流體

對于模擬流體效果的流體模擬來說,主要是求解不可壓縮流的Navier-Stokes 方程,其形式如下:

其中,j表示的是r的鄰域中粒子的序號,mj表示的第j個粒子的質(zhì)量,rj是粒子j的位置坐標(biāo),ρ表示其密度,Aj表示粒子j的對應(yīng)值.函數(shù)W(r,h)稱為光滑函數(shù),其中h被稱為光滑長度.

2.2.2 基于空間劃分的NNPS

SPH 算法中,NNPS 過程對算法的實時性影響較大,若采用遍歷搜索法,則其時間復(fù)雜度為O(n2),n為粒子數(shù)目,顯然會成為SPH 算法的瓶頸.基于空間劃分的方法進(jìn)行NNPS 過程能將NNPS 過程的時間復(fù)雜度降低到O(nm),其中,m表示單元三維網(wǎng)格內(nèi)流血粒子的平均數(shù)目,顯然,采用基于空間劃分的方法進(jìn)行NNPS 過程能大大降低其計算資源的消耗.

在基于空間劃分的NNPS 過程中,首先根據(jù)流體的范圍確定三維坐標(biāo)的最大值(xmax,ymax,zmax)和最小值(xmin,ymin, zmin),根據(jù)最大值和最小值劃定的范圍,以光滑長度h對空間進(jìn)行垂直坐標(biāo)系3 個方向上的劃分,劃分后,流體運動的空間就被包含在三維網(wǎng)格中,并伴隨著流體的運動不斷改變,如圖1所示.

圖1 SPH 算法中NNPS 過程劃分的網(wǎng)格空間

基于空間劃分的NNPS 過程建立的三維均勻空間網(wǎng)格可以復(fù)用于本文算法中碰撞檢測范圍確定以及SPH 流血粒子的定位.

2.2.3 空間網(wǎng)格長度劃分依據(jù)以及對碰撞檢測算法精度和性能的影響

光滑長度的確定決定了粒子的影響范圍,光滑長度過長,則影響范圍中的粒子過多,計算復(fù)雜度增加;光滑長度過短,則影響范圍內(nèi)粒子太少,不能表現(xiàn)出粒子的性質(zhì);所以選擇一個合適的光滑長度對SPH 算法以及碰撞檢測算法都非常重要.

由于本算法著重于表達(dá)流體效果,而非精確的物理仿真,對流體的計算不需要非常精準(zhǔn),所以我們認(rèn)為處理的流體密度在任何時候都是均勻的.基于以上設(shè)定,我們可以選擇固定的光滑長度h,其數(shù)值與流體的密度成反比[14].

用于碰撞檢測的空間網(wǎng)格使用的是NNPS 過程中建立起來的空間網(wǎng)格.在光滑長度h 是常量的情況下,我們采用h作為NNPS 網(wǎng)格劃分的單元長度,這樣,對于給定空間中的粒子,其相鄰粒子只能在本身所在的網(wǎng)格或者相鄰的網(wǎng)格中.在一維、二維、三維空間中,尋找相鄰粒子只需要搜索對應(yīng)的3、9、27 個單元即可.

在本文提出的碰撞檢測算法中,空間網(wǎng)格的大小對算法性能的影響較大,單元網(wǎng)格尺寸太大會導(dǎo)致過多流血粒子和三角形包圍到同一網(wǎng)格中,從而增加了基元測試的負(fù)擔(dān),太小又會導(dǎo)致一個三角形占用過多網(wǎng)格,并影響碰撞檢測的精確性.Ericson 等的研究表明,當(dāng)三維網(wǎng)格大小與三角形的AABB 包圍盒尺寸相當(dāng)時為最優(yōu)[15],因此我們?nèi)∷腥切蔚腁ABB 包圍盒尺寸的平均值.

鑒于上述兩種方案都需要影響空間網(wǎng)格的尺寸和數(shù)量,所以我們將兩者計算出來的網(wǎng)格長度取平均值,使得網(wǎng)格長度兼顧流體模擬和碰撞檢測.

空間網(wǎng)格的數(shù)量由流體所分布的空間所決定的,我們使用計算出來的網(wǎng)格長度對空間進(jìn)行分割.比較典型的例子如圖1(b)所示,網(wǎng)格數(shù)量為 9×9×9個,圖1(d)所示的網(wǎng)格數(shù)量為1 7×17×12個,在此數(shù)量的網(wǎng)格數(shù)量下,在保持精度的情況下,時間復(fù)雜度影響也很小.

2.2.4 結(jié)合空間劃分NNPS 過程進(jìn)行SPH 粒子與軟組織表面的碰撞檢測

本文提出的算法能將SPH 算法中基于空間劃分的NNPS 過程用于空間均勻劃分并將流體粒子定位.要進(jìn)行基于空間劃分的碰撞檢測,還需利用劃分出來的均勻空間網(wǎng)格進(jìn)行軟組織表面三角面片的定位.對于軟組織表面的三角面片i的頂點(xi,yi,zi),用其減去SPH 算法中NNPS 過程確定的最小坐標(biāo)值(xmin,ymin,zmin),再除以單元網(wǎng)格的長度h,得到網(wǎng)格序號此序號表示該頂點在劃分的空間中所在的相對位置,若網(wǎng)格序號超過流血確定的序號,說明該頂點超出流體所在的空間,可以排除出碰撞檢測的范圍.反之,則將該三角面片加入到對應(yīng)網(wǎng)格的碰撞檢測序列中.當(dāng)所有三角面片上的點均判斷完畢后,再對每個網(wǎng)格單元進(jìn)行遍歷,若網(wǎng)格單元內(nèi)既存在SPH 流血粒子,又存在軟組織三角面片,則表明該網(wǎng)格中存在碰撞檢測的可能性,這時,就對相應(yīng)的序列進(jìn)行基元間的精確碰撞檢測.這樣,就將SPH 算法中的NNPS 過程應(yīng)用于流體粒子的定位以及三角面片的定位,再進(jìn)行進(jìn)一步的精確碰撞檢測.將該碰撞檢測算法的流程圖如圖2所示.

3 實驗結(jié)果及分析

對本文提出的基于空間劃分的碰撞檢測算法進(jìn)行實驗,試驗環(huán)境為Interi7 9500 K CPU,16 GB 內(nèi)存,NVIDIA GeForce GTX1080 顯卡.使用VS 2017,OpenGL 編程實現(xiàn)算法.實驗的對比算法是傳統(tǒng)的基于空間劃分的碰撞檢測算法.本文的算法中,由于空間劃分過程和流血粒子定位過程均在基于SPH 算法的流體模擬過程實現(xiàn),其相對于傳統(tǒng)的空間劃分算法優(yōu)勢明顯,結(jié)果如圖3、圖4所示,圖3是傳統(tǒng)的空間劃分算法的時間性能表現(xiàn),圖4是本文所提出的算法的時間性能表現(xiàn).時間性能對比如表1和表2所示.

對肝臟軟組織與流血粒子的碰撞檢測進(jìn)行實驗,肝臟表面三角面片數(shù)量為13 070,粒子數(shù)目為4800,模擬的幀率為151 fps,模擬效果如圖3所示.

對流血粒子使用Marching Cubes 算法進(jìn)行表面渲染,模擬肝臟表面流血的效果,結(jié)果如圖4所示.

圖4 模擬肝臟表面流血效果

表1 傳統(tǒng)算法時間性能(單位:ms)

表2 本文算法時間性能(單位:ms)

將本文算法運用到虛擬肝臟手術(shù)中,虛擬肝臟實質(zhì)切開后,肝臟產(chǎn)生形變,傷口內(nèi)部開始滲血,流血能夠根據(jù)虛擬肝臟形態(tài)的變化進(jìn)行自適應(yīng)的調(diào)整,具體結(jié)果見圖5.

圖5 使用本文算法模擬虛擬肝臟手術(shù)流血

4 結(jié)束語

虛擬手術(shù)中的流血模擬對碰撞檢測的實時性和精確性均有較高要求,由于手術(shù)仿真的過程中涉及到軟體模型的形變和切割,給流血模擬的碰撞檢測帶來很大的難度.本文論述了虛擬手術(shù)流血模擬中碰撞檢測問題的特點,并給出了一種基于空間劃分的SPH 流體粒子與軟碰撞檢測算法.該算法不僅能有效提高碰撞檢測實時性,而且能實現(xiàn)SPH 模擬的流體粒子與任意動力學(xué)模型的軟體模型之間的精確碰撞檢測.下一步,將此算法結(jié)合GPU 加速,使算法實時性進(jìn)一步提高.實驗證明,該方法能滿足虛擬手術(shù)流血模擬的要求.