笛卡爾網(wǎng)格下精確高效的壁面距離計(jì)算方法

孟 爽，周 丹，李雪亮，畢 林，3，*

(1.中南大學(xué) 軌道交通安全教育部重點(diǎn)實(shí)驗(yàn)室，長(zhǎng)沙 410075；2.空氣動(dòng)力學(xué)國(guó)家重點(diǎn)實(shí)驗(yàn)室，綿陽(yáng) 621000；3.中國(guó)空氣動(dòng)力研究與發(fā)展中心 計(jì)算空氣動(dòng)力研究所，綿陽(yáng) 621000)

符號(hào)說(shuō)明

0 引言

笛卡爾網(wǎng)格可以不依賴物面網(wǎng)格直接生成，因具有自動(dòng)化程度高、復(fù)雜外形適應(yīng)性好、非定常/多尺度流動(dòng)結(jié)構(gòu)捕捉能力強(qiáng)等優(yōu)勢(shì)而廣受關(guān)注[1]。

壁面距離的精確高效計(jì)算對(duì)實(shí)現(xiàn)笛卡爾網(wǎng)格方法具有重要意義：一方面，由于笛卡爾網(wǎng)格的非貼體特性，通常采用虛擬單元法處理物面邊界[2-6]，需要根據(jù)物面確定虛擬單元的壁面距離，尋找參考點(diǎn)并插值獲得其物理量，壁面距離精確計(jì)算對(duì)虛擬單元物面處理精準(zhǔn)度有較大影響；另一方面，對(duì)于非定常流動(dòng)或運(yùn)動(dòng)物體，笛卡爾網(wǎng)格會(huì)根據(jù)流場(chǎng)的變化或物體的運(yùn)動(dòng)進(jìn)行自適應(yīng)[7-8]，自適應(yīng)后的網(wǎng)格壁面距離需重新計(jì)算，壁面距離計(jì)算效率成為制約流動(dòng)計(jì)算效率的關(guān)鍵因素之一。

目前，壁面距離的計(jì)算方式主要分為求解偏微分方程[9-12]和采用幾何方法直接計(jì)算[13-14]兩類。第一類方法，將最小距離轉(zhuǎn)換為關(guān)于波傳播問(wèn)題的偏微分方程數(shù)值求解[15]，需要額外的計(jì)算花費(fèi)，壁面距離計(jì)算精度受數(shù)值離散精度的影響，且對(duì)于復(fù)雜外形適應(yīng)性較差。第二類方法是根據(jù)空間點(diǎn)與物面離散網(wǎng)格的幾何關(guān)系計(jì)算壁面距離。通常為簡(jiǎn)單起見(jiàn)，在物面網(wǎng)格足夠密的情況下，求解空間網(wǎng)格點(diǎn)的壁面距離時(shí)，一般采用空間點(diǎn)到物面網(wǎng)格中心的距離代替最小距離。趙鐘等[16]發(fā)現(xiàn)，這種近似距離會(huì)引起計(jì)算誤差，不僅影響計(jì)算結(jié)果的精度，而且對(duì)計(jì)算穩(wěn)定性造成影響。為提高計(jì)算精度，有學(xué)者發(fā)展了根據(jù)空間點(diǎn)關(guān)于物面網(wǎng)格所在平面投影點(diǎn)與物面網(wǎng)格的幾何關(guān)系計(jì)算壁面距離的3D 投影法[16]、通過(guò)坐標(biāo)變換轉(zhuǎn)化為二維問(wèn)題的2D 法[17]等，但這些方法涉及大量向量計(jì)算和關(guān)系判斷，計(jì)算量較大。

在幾何方法求解壁面距離時(shí)，由于物面網(wǎng)格數(shù)目較大，計(jì)算效率的核心問(wèn)題是如何快速定位到最短距離相關(guān)的物面網(wǎng)格，物面網(wǎng)格數(shù)據(jù)結(jié)構(gòu)是關(guān)鍵。

最常見(jiàn)的方法是將物面網(wǎng)格順序存入數(shù)組，采用遍歷法搜索物面點(diǎn)得到壁面距離。網(wǎng)格規(guī)模較小時(shí)，計(jì)算量尚在承受范圍之內(nèi)。但是對(duì)于三維復(fù)雜幾何外形的流場(chǎng)計(jì)算，物面網(wǎng)格數(shù)量可達(dá)到O(105)～O(106)量級(jí)，空間網(wǎng)格點(diǎn)可達(dá)到O(108)～O(109)量級(jí)，采用遍歷法計(jì)算量達(dá)到O(1013)～O(1015)量級(jí)。此種規(guī)模的計(jì)算量對(duì)整個(gè)計(jì)算周期的影響是巨大的。Boger[18]提出采用ADT (alternating digital tree)二叉樹(shù)數(shù)據(jù)結(jié)構(gòu)存儲(chǔ)物面網(wǎng)格，計(jì)算效率大幅提高，成為目前最常用的物面網(wǎng)格數(shù)據(jù)存儲(chǔ)方式。但對(duì)于復(fù)雜外形，ADT 平衡性下降，計(jì)算效率仍然差強(qiáng)人意。郭中州等[14]采用KDT (K-dimensional tree)存儲(chǔ)物面網(wǎng)格，解決了ADT 平衡性問(wèn)題，計(jì)算效率進(jìn)一步提升，但對(duì)于離物面較遠(yuǎn)的空間點(diǎn)，在運(yùn)用KDT 最近鄰搜索算法時(shí)，超球面與分割面位置關(guān)系判斷失準(zhǔn)，需要反復(fù)回溯，計(jì)算量大。

本文針對(duì)上述問(wèn)題，引入三角形參數(shù)化方法，將空間點(diǎn)到三角形物面離散網(wǎng)格的最小距離問(wèn)題轉(zhuǎn)換為約束條件下一維極值問(wèn)題；同時(shí)發(fā)展嵌套包圍盒概念的KDT 物面網(wǎng)格數(shù)據(jù)存儲(chǔ)結(jié)構(gòu)，優(yōu)化KDT 最近鄰搜索算法中距物面較遠(yuǎn)數(shù)據(jù)點(diǎn)回溯過(guò)程，實(shí)現(xiàn)最小距離對(duì)應(yīng)的三角形的快速定位。

1 點(diǎn)到三角形最小距離精確算法

點(diǎn)到空間三角形最小距離的問(wèn)題[19]可描述為點(diǎn)P和三角形T(s,t)=B+sE0+tE1之間的最小距離，其中(s,t)∈D={(s,t):s≥0,t≥0,s+t≤1}，B為三角形的一個(gè)頂點(diǎn)，E0和E1分別為此頂點(diǎn)對(duì)應(yīng)的三角形兩條邊。點(diǎn)P和三角形內(nèi)任一點(diǎn)T(s,t)距離的平方為橢圓函數(shù)：

其 中，(s,t)∈D，a=E0·E0，b=E0·E1，c=E1·E1，d=E0·(B-P)，e=E1·(B-P)，f=(B-P)·(B-P)。點(diǎn)P和三角形之間的關(guān)系如圖1 所示。

圖1 點(diǎn)P 和三角形之間的關(guān)系Fig.1 The relationship between P and the triangle

最小距離問(wèn)題轉(zhuǎn)換為在D內(nèi)求連續(xù)可微函數(shù)Q(s,t)的極值問(wèn)題。令：

圖2 用三角域D 劃分s-t 平面以及不同等高線Fig.2 Partition of a s-t plane by triangle domain D and distance contours

圖3 以區(qū)域③為中心的等高線與三角形接觸的兩種情況Fig.3 Two contour levels with centers in region ③touching the triangle

在等值線上的任何一點(diǎn)，-Q的方向均朝向橢圓內(nèi)部。為簡(jiǎn)單起見(jiàn)，我們根據(jù)-Q沿頂點(diǎn)(0,1)的兩條邊的方向?qū)?shù)的符號(hào)進(jìn)行判斷：

1）對(duì)于邊s+t=1，頂點(diǎn)(0,1)處的方向?qū)?shù)為如果值為負(fù)，最小值出現(xiàn)在邊s+t=1上。最小值的情形對(duì)應(yīng)圖3 紅色曲線，處理方式和s∈[0,1]在區(qū)域②類似。

2）對(duì)于邊s=0，頂 點(diǎn)(0,1)處的方向?qū)?shù)為(0,-1)·?Q(0,1)。如果值為負(fù)，最小值出現(xiàn)在邊s=0上。最小值的情形對(duì)應(yīng)圖3 綠色曲線，處理方式和 (,)在區(qū)域②類似。

2 物面單元的存儲(chǔ)

2.1 嵌套包圍盒概念的KDT 參數(shù)定義

與文獻(xiàn)[14]中只存儲(chǔ)物面網(wǎng)格中心點(diǎn)不同，本文發(fā)展了基于嵌套包圍盒概念的KDT，不僅按逆時(shí)針順序存儲(chǔ)物面網(wǎng)格頂點(diǎn)信息（方便笛卡爾網(wǎng)格物面處理），還同時(shí)包含KDT 節(jié)點(diǎn)對(duì)應(yīng)所有物面網(wǎng)格的最小包圍盒，KDT 節(jié)點(diǎn)從上到下形成嵌套包圍盒形式。嵌套包圍盒的引入保證KDT 每個(gè)節(jié)點(diǎn)包含當(dāng)前節(jié)點(diǎn)對(duì)應(yīng)三角形的所有信息，且在最小距離查詢時(shí)，僅需要通過(guò)與剖分面（剖分面的大小由包圍盒確定）坐標(biāo)對(duì)比，就可以快速排除與目標(biāo)點(diǎn)距離大于當(dāng)前最小距離的三角形，效率大幅提升。

2.2 嵌套包圍盒概念的KDT 建立步驟

步驟1：在正式創(chuàng)建KDT 之前，首先需要確定目標(biāo)數(shù)據(jù)集G中三角形元素的數(shù)量N。然后，需要根據(jù)這N個(gè)元素確定G的最小包圍盒，包圍盒的兩個(gè)頂點(diǎn)分為 (Gimin,Gimax)，其中i=1，2，3，其定義如圖4 所示。

圖4 幾何外形包圍盒Fig.4 Geometry bounding box

步驟2：創(chuàng)建一個(gè)空的KDT 備用。

步驟3：若N為1，則將此元素插入到當(dāng)前KDT的節(jié)點(diǎn)中去，同時(shí)，將此節(jié)點(diǎn)的左右子樹(shù)均置空。

步驟4：若N為2，則將第一個(gè)元素插入到當(dāng)前KDT 的節(jié)點(diǎn)中去，同時(shí)規(guī)定第一維為劃分維。此時(shí)，只有第二個(gè)元素尚未被插入到KDT 中。將第二個(gè)元素中心點(diǎn)與第一個(gè)元素的中心點(diǎn)在第一維的值進(jìn)行比較。若第二個(gè)元素大于第一個(gè)元素，則對(duì)當(dāng)前節(jié)點(diǎn)的右子樹(shù)進(jìn)行步驟3 的操作，同時(shí)當(dāng)前節(jié)點(diǎn)的左子樹(shù)置空；否則，對(duì)當(dāng)前節(jié)點(diǎn)的左子樹(shù)進(jìn)行步驟3 的操作，同時(shí)當(dāng)前節(jié)點(diǎn)的右子樹(shù)置空。此步驟中d加1。

步驟5：若N大于2，則首先確定劃分維。為了保證KDT 的平衡，劃分維的確定是計(jì)算N個(gè)三角形的中心點(diǎn)Tcenter在各個(gè)維度的方差，然后找出方差最大的維度即為劃分維。然后將此N個(gè)元素中心點(diǎn)在方差最大的維度上按從小到大的順序進(jìn)行排序。將位于中間的元素插入到當(dāng)前KDT 的節(jié)點(diǎn)中去。小于此元素的將全部位于當(dāng)前節(jié)點(diǎn)的左子樹(shù)，大于此元素的將全部位于當(dāng)前節(jié)點(diǎn)的右子樹(shù)。對(duì)當(dāng)前節(jié)點(diǎn)左子樹(shù)和右子樹(shù)的元素遞歸執(zhí)行步驟3 或步驟4 或步驟5，直至所有元素全部插入到KDT 中。

值得注意的是，在向KDT 中插入節(jié)點(diǎn)的同時(shí)，節(jié)點(diǎn)對(duì)應(yīng)的包圍盒也在進(jìn)行同步劃分。例如KDT 根節(jié)點(diǎn)對(duì)應(yīng)的包圍盒即為數(shù)據(jù)集G的包圍盒 (Gimin,Gimax)。假設(shè)節(jié)點(diǎn)K對(duì)應(yīng)的包圍盒為(Kimin,Kimax)，則K的左子節(jié)點(diǎn)KL的包圍盒(KLimin,KLimax)為：

由于節(jié)點(diǎn)包圍盒表示此空間區(qū)域內(nèi)所有三角形元素集合的最小包圍盒，而在進(jìn)行空間區(qū)域剖分時(shí)的依據(jù)是三角形中心點(diǎn)確定的平面，所以節(jié)點(diǎn)包圍盒與包圍盒之間，剖分面與剖分面之間會(huì)存在相交的部分。這種相交只會(huì)出現(xiàn)在三維及以上的情況。

2.3 創(chuàng)建KDT 時(shí)間復(fù)雜度分析

由于KDT 是一個(gè)二叉樹(shù)，所以對(duì)于給定的N個(gè)元素建成一個(gè)深度為d的最優(yōu)KDT 有：

3 最小距離查詢

計(jì)算空間某點(diǎn)到物面的最小距離需要兩個(gè)步驟。第一，得到此點(diǎn)到物面最小距離的點(diǎn)所在的物面三角形；第二，采用點(diǎn)到三角形最小距離精確算法得到此點(diǎn)到物面的最小距離。本節(jié)基于KDT 最近鄰算法對(duì)第一步進(jìn)行設(shè)計(jì)。

3.1 算法設(shè)計(jì)

1）首先從KDT 的根節(jié)點(diǎn)開(kāi)始查詢，根據(jù)目標(biāo)點(diǎn)與根節(jié)點(diǎn)剖分面的位置關(guān)系，確定與目標(biāo)點(diǎn)T潛在的最近點(diǎn)是在根節(jié)點(diǎn)的左子樹(shù)還是右子樹(shù)。

2）根據(jù)1）確定的子樹(shù)，對(duì)此子樹(shù)的節(jié)點(diǎn)進(jìn)行步驟1）的操作。

3）當(dāng)遞歸訪問(wèn)至葉子節(jié)點(diǎn)時(shí)，記錄T與當(dāng)前節(jié)點(diǎn)之間的距離作為當(dāng)前最近距離lmin，當(dāng)前節(jié)點(diǎn)作為與目標(biāo)點(diǎn)T最近的節(jié)點(diǎn)記為Knearest。

4）進(jìn)行回溯操作。回溯操作從當(dāng)前最近的點(diǎn)Knearest開(kāi)始向上進(jìn)一步查找。構(gòu)造以T為球心，lmin為半徑的球。若球面與Knearest父節(jié)點(diǎn)的剖分面相交，則查找Knearest的兄弟節(jié)點(diǎn)即Knearest父節(jié)點(diǎn)的另外子樹(shù)的節(jié)點(diǎn)，同時(shí)計(jì)算T與Knearest父節(jié)點(diǎn)之間的距離，若小于當(dāng)前最小距離，則Knearest和lmin均進(jìn)行更新；若球面與剖分面不相交，則說(shuō)明Knearest的兄弟節(jié)點(diǎn)中不存在比Knearest與T更近的節(jié)點(diǎn)。

5）執(zhí)行步驟4）直至根節(jié)點(diǎn)。得到的Knearest即為與T最近的點(diǎn)，lmin為Knearest與T的最小距離。至此最小距離查詢模塊結(jié)束。

3.2 算法優(yōu)化

在實(shí)際流場(chǎng)計(jì)算中，存在空間網(wǎng)格點(diǎn)與物面距離較遠(yuǎn)的情況，此時(shí)采用上述最小距離查詢方法將面臨回溯操作效率低下的問(wèn)題。為了直觀說(shuō)明問(wèn)題，本小節(jié)以二維數(shù)據(jù)點(diǎn)為例進(jìn)行說(shuō)明，給定隨機(jī)分布的17 個(gè)點(diǎn)，對(duì)這些點(diǎn)創(chuàng)建KDT，如圖5 所示。給定目標(biāo)單元T，正常查找得到的最鄰近單元是G。在回溯操作中，以T為圓心，T和G之間的距離lTG為半徑的圓與根節(jié)點(diǎn)J的剖分軸m相交。根據(jù)上述回溯步驟，此時(shí)需在J的左子樹(shù)進(jìn)行查詢。然而，實(shí)際情況是圓與J的左子樹(shù)并不相交。因此，后續(xù)的回溯操作均為無(wú)效回溯，這將大大浪費(fèi)計(jì)算資源。為此，本文在回溯操作時(shí)，首先判斷以T為圓心，T和G之間的距離lTG為半徑的圓是否與線段ab相交（ab由剖分軸和當(dāng)前包圍盒確定）。若相交，則進(jìn)行正常的回溯操作；若不相交，則直接排除對(duì)應(yīng)節(jié)點(diǎn)另外子樹(shù)存在最小距離的可能性。采用改進(jìn)前后的回溯法進(jìn)行的回溯操作如圖5（b～c）中綠色箭頭所示。從圖中可以看出，改進(jìn)后的回溯方法可極大減少不必要的回溯過(guò)程，進(jìn)而提高最小距離查詢效率。

圖5 圓與剖分軸相交且回溯無(wú)效的情況Fig.5 Sketch of an optimized backtracking process when the circle intersects the split axis

圖6 為球、導(dǎo)彈和DPW6 三種幾何外形及空間笛卡爾網(wǎng)格。為了定量對(duì)比改進(jìn)前后回溯算法的效率，本文對(duì)這三種不同幾何外形進(jìn)行了壁面距離查詢測(cè)試。對(duì)幾何外形進(jìn)行壁面距離查找時(shí)，存儲(chǔ)所有空間網(wǎng)格點(diǎn)的壁面距離，同時(shí)記錄得到壁面距離所需的查找次數(shù)，然后將查找次數(shù)的總和除以空間網(wǎng)格點(diǎn)數(shù)得到平均查找次數(shù)，結(jié)果如圖7 所示。從圖中可以看出，針對(duì)不同外形，本文改進(jìn)的回溯方法均可有效減小查找次數(shù)。針對(duì)同一幾何外形，物面三角形數(shù)量越多，查找次數(shù)減少的越明顯。

圖6 三種幾何外形Fig.6 Three geometries

圖7 求最小距離所需查找次數(shù)Fig.7 Number of queries for finding the minimum distance

4 算例測(cè)試

本小節(jié)對(duì)本文發(fā)展的壁面距離算法的準(zhǔn)確性和效率進(jìn)行了測(cè)試。測(cè)試的幾何外形仍為球、導(dǎo)彈和DPW6。為了保證測(cè)試結(jié)果的合理性，本文所有工況均在AMD EPYC 7702 的處理器上測(cè)試，測(cè)試平臺(tái)為Window 10 系統(tǒng)及Visual Studio 2012 版本。每個(gè)工況測(cè)試10 次，統(tǒng)計(jì)平均獲得程序運(yùn)行的時(shí)間。

首先對(duì)壁面距離結(jié)果的準(zhǔn)確性進(jìn)行對(duì)比測(cè)試。對(duì)壁面距離的兩種計(jì)算方式進(jìn)行對(duì)比，第一種是近似壁面距離，即空間點(diǎn)到物面三角形中心的最小距離，第二種是精確壁面距離，即采用本文發(fā)展的算法計(jì)算得到的最小距離。計(jì)算結(jié)果如圖8 所示，圖8（a）為近似壁面距離計(jì)算的結(jié)果，圖8（b）為精確壁面距離計(jì)算所得結(jié)果。從圖8（a）中可以明顯看出，球的近似壁面距離云圖及等值線呈現(xiàn)波浪狀，與實(shí)際物理特征明顯不符。采用本文方法計(jì)算得出的精確壁面距離等值線光順平滑，符合物理特征。對(duì)于幾何外形較為復(fù)雜的DPW6，在遠(yuǎn)離壁面區(qū)域，兩種計(jì)算方法得出的壁面距離結(jié)果接近；而在近壁區(qū)域，采用近似壁面距離計(jì)算方法得到的壁面距離為250 的等值線與實(shí)際幾何特征誤差較大，而精確計(jì)算結(jié)果的近壁面等值線可較好反映近壁物面特征。

圖8 不同幾何外形近似壁面距離與精確壁面計(jì)算結(jié)果云圖及等值線Fig.8 Contours of approximate and accurate calculation results for wall distances of different geometric shapes

此外，以球?yàn)檠芯繉?duì)象，我們將本文方法計(jì)算得到的壁面距離與解析解進(jìn)行對(duì)比。結(jié)果表明本文方法的計(jì)算結(jié)果與解析解差異在百萬(wàn)分之一以內(nèi)。因此，本文方法得到的壁面距離精確性較高。

在計(jì)算效率測(cè)試方面，對(duì)遍歷法+精確距離、ADT+精確距離、KDT+精確距離，KDT+近似距離4 種方法的計(jì)算效率進(jìn)行對(duì)比。物面三角形數(shù)量、空間網(wǎng)格數(shù)量以及CPU 運(yùn)行時(shí)間如表1 所示。從表中可以看出，采用KDT 方法的計(jì)算效率相比遍歷法和ADT 方法均有量級(jí)的提升。尤其對(duì)于復(fù)雜的DPW6 來(lái)說(shuō)，采用遍歷法的時(shí)間為3 616.87 s，KDT 只需2.92 s，效率提升超過(guò)1 000 倍。而對(duì)于簡(jiǎn)單的球來(lái)說(shuō)，效率提升僅不到300 倍。因此，壁面距離查詢效率的提升程度與幾何外形密切相關(guān)。

表1 壁面距離查詢所需時(shí)間Table 1 Wall-distance querying time

精確距離的計(jì)算要比近似距離的計(jì)算多調(diào)用一個(gè)點(diǎn)到三角形距離的函數(shù)，此函數(shù)在回溯過(guò)程中頻繁調(diào)用會(huì)增加計(jì)算耗時(shí)。從定量結(jié)果來(lái)看，精確距離計(jì)算耗時(shí)為近似距離計(jì)算耗時(shí)的2 倍左右，但這對(duì)整個(gè)計(jì)算周期的耗時(shí)影響并不明顯。因此，通過(guò)增加少量計(jì)算時(shí)間來(lái)獲取計(jì)算精度收益是值得的。

為了進(jìn)一步測(cè)試本文發(fā)展的方法在大規(guī)模網(wǎng)格計(jì)算時(shí)的效率，我們選用DPW6 外形進(jìn)行物面網(wǎng)格離散并生成空間笛卡爾網(wǎng)格。物面三角形網(wǎng)格數(shù)量跨越萬(wàn)、十萬(wàn)以及百萬(wàn)3 個(gè)量級(jí)；空間笛卡爾網(wǎng)格數(shù)量跨越千萬(wàn)、億以及十億3 個(gè)量級(jí)。效率測(cè)試結(jié)果如表2 所示，結(jié)果表明，對(duì)于同一種物面離散方式，CPU運(yùn)行時(shí)間與空間網(wǎng)格數(shù)量大致呈正比；對(duì)于同一數(shù)量空間網(wǎng)格（由于物面網(wǎng)格數(shù)量不同，加密相同層級(jí)生成的空間笛卡爾網(wǎng)格數(shù)量有略微差異），物面網(wǎng)格數(shù)量從30 021 增加到1 417 181 時(shí)，CPU 運(yùn)行時(shí)間僅增加為原來(lái)的1～2 倍，這再次說(shuō)明了本文基于KDT 的回溯方法的高效性。

表2 DPW6 壁面距離計(jì)算效率Table 2 Wall-distance calculation time for DPW6

為了與文獻(xiàn)[16]中計(jì)算效率進(jìn)行對(duì)比，本文選用JSM（JAXA standard model）標(biāo)模進(jìn)行壁面距離查詢測(cè)試，網(wǎng)格規(guī)模為33 億。本文方法（無(wú)并行）進(jìn)行壁面距離查詢所需時(shí)間大致為0.51 h，文獻(xiàn)[16]中采用并行方法對(duì)進(jìn)行壁面距離查詢所需時(shí)間為0.42 h。因此，除去網(wǎng)格類型、測(cè)試環(huán)境帶來(lái)的誤差，總體來(lái)說(shuō)，本文方法在大規(guī)模網(wǎng)格下仍表現(xiàn)出較高效率。

5 結(jié)論

在本文工作中，我們基于笛卡爾網(wǎng)格發(fā)展了一種壁面距離計(jì)算方法，具備精確性、高效性以及通用性?；诒疚膬?nèi)容，主要有以下3 個(gè)方面結(jié)論：

1）在提高壁面距離計(jì)算的準(zhǔn)確性方面，我們將物面三角形進(jìn)行參數(shù)化（s和t）表示，根據(jù)s和t的取值范圍將s-t平面分為7 個(gè)區(qū)域，對(duì)不同區(qū)域進(jìn)行分情況討論，從而將求最小距離問(wèn)題巧妙地轉(zhuǎn)化為求約束條件下一維極值問(wèn)題，計(jì)算量大幅減少。

2）在提高計(jì)算效率方面，采用平衡的KDT 存儲(chǔ)物面三角形所有頂點(diǎn)信息并建立相應(yīng)的嵌套包圍盒，使得在最小距離查詢時(shí)盡可能快速排除在目標(biāo)范圍之外的單元；同時(shí)優(yōu)化了當(dāng)前超球面與剖分平面的位置判斷法則，解決了在距離物面較遠(yuǎn)處因無(wú)效回溯次數(shù)過(guò)多導(dǎo)致查詢效率下降的問(wèn)題。

3）此外，本文的壁面距離查詢模塊是一個(gè)單獨(dú)的模塊，在計(jì)算壁面距離時(shí)，僅需知道空間網(wǎng)格點(diǎn)的位置信息以及物面三角形信息。因此，本文發(fā)展的壁面距離計(jì)算方法不僅對(duì)笛卡爾網(wǎng)格，對(duì)結(jié)構(gòu)網(wǎng)格、非結(jié)構(gòu)網(wǎng)格以及重疊網(wǎng)格等也可較好兼容，可方便程序的移植。

本文采用3 個(gè)三維幾何外形對(duì)本文方法進(jìn)行測(cè)試，結(jié)果表明計(jì)算得到的壁面距離與解析值的誤差在百萬(wàn)分之一以內(nèi)，有較高的精確性；十億量級(jí)網(wǎng)格規(guī)模下的單核計(jì)算效率與已有文獻(xiàn)[16]的并行計(jì)算效率相當(dāng)。雖然在大規(guī)模網(wǎng)格測(cè)試中表現(xiàn)出較高的計(jì)算效率，但仍有提升空間。下一步我們計(jì)劃添加OpenMP/MPI 并行方法進(jìn)一步提升計(jì)算效率，以滿足下一代超大規(guī)模CFD 湍流模擬的需求。