艦船綜合水動(dòng)力分析虛擬試驗(yàn)系統(tǒng)中試驗(yàn)結(jié)果的可視化關(guān)鍵技術(shù)研究

倪 昱，金建海，單威俊

(中國(guó)船舶科學(xué)研究中心，江蘇 無錫214082)

目前，對(duì)于虛擬實(shí)驗(yàn)結(jié)果的可視化，一般通過商業(yè)CAE仿真后處理軟件來實(shí)現(xiàn)，包括MSC.Patran，EnSight，Tecplot，F(xiàn)emap等。這些商業(yè)軟件數(shù)據(jù)可視化功能已經(jīng)比較全面，但是其可視化的具體方法和代碼并不開放，不能很好地應(yīng)用于艦船水動(dòng)力虛擬試驗(yàn)中個(gè)性化的功能要求，同時(shí)針對(duì)水動(dòng)力虛擬試驗(yàn)而言，其可視化顯示效率還有改進(jìn)的余地。

因此，在中國(guó)船舶科學(xué)研究中心“艦船綜合水動(dòng)力虛擬水池試驗(yàn)技術(shù)系統(tǒng)”的開發(fā)過程中，運(yùn)用自主設(shè)計(jì)的相關(guān)算法或改進(jìn)已有算法，實(shí)現(xiàn)了水動(dòng)力虛擬試驗(yàn)中大規(guī)模有限元模型及其結(jié)果信息的快速顯示。其中試驗(yàn)結(jié)果可視化的關(guān)鍵技術(shù)主要包括虛擬試驗(yàn)?zāi)Ｐ涂梢暬夹g(shù)、標(biāo)量試驗(yàn)結(jié)果可視化技術(shù)及矢量試驗(yàn)結(jié)果可視化技術(shù)。本文分別就這三方面，討論虛擬試驗(yàn)結(jié)果的可視化關(guān)鍵技術(shù)及其實(shí)現(xiàn)算法。

1 虛擬試驗(yàn)?zāi)Ｐ涂梢暬?/h2>

CFD的結(jié)果數(shù)據(jù)是在建立網(wǎng)格的基礎(chǔ)上得出來的，因此，進(jìn)行數(shù)據(jù)可視化就需要己知流場(chǎng)中的網(wǎng)格數(shù)據(jù)。因此虛擬模型可視化是試驗(yàn)結(jié)果可視化的基礎(chǔ)。

在計(jì)算機(jī)圖形學(xué)中，對(duì)于三維幾何模型的可視化，一般是通過OpenGL函數(shù)繪制點(diǎn)、線和多邊形等基本的幾何圖形來實(shí)現(xiàn)的。但在艦船虛擬試驗(yàn)中幾何模型一般比較大［1］，本文使用“三維模型外表面提取實(shí)現(xiàn)算法”和“網(wǎng)格消隱算法”，對(duì)于三維模型的內(nèi)部單元和被遮擋的外部單元進(jìn)行消隱，以提高模型的顯示速度和效果。

1.1 三維模型外表面提取實(shí)現(xiàn)算法

三維網(wǎng)格模型是由三維網(wǎng)格單元組成的，由于在三維模型可視化過程中需要顯示的只是所選擇的計(jì)算區(qū)域外表面的信息，因此三維網(wǎng)格數(shù)值模擬計(jì)算結(jié)果可視化快速實(shí)現(xiàn)問題，就轉(zhuǎn)化為快速的判斷在三維網(wǎng)格中哪些單元面是計(jì)算域外表面，以及哪些外表面可見的問題。

在三維模型網(wǎng)格中，內(nèi)部單元面屬于并只屬于兩個(gè)單元;外部單元面只屬于某一個(gè)單元。形成結(jié)構(gòu)外部單元面的一般方法為先形成結(jié)構(gòu)中所有單元的單元面，然后對(duì)各個(gè)單元面進(jìn)行比較，節(jié)點(diǎn)組成相同的面即為內(nèi)部單元面，剩余的即為外部單元面。隨著網(wǎng)格單元數(shù)目的增加，這種方法的運(yùn)算量將呈幾何級(jí)數(shù)增加。

本文使用了一種最小節(jié)點(diǎn)相關(guān)面表的方法，可快速得到位于結(jié)構(gòu)表面的單元面的集合(且稱之為結(jié)構(gòu)邊界面集合)。

步驟1。建立一個(gè)節(jié)點(diǎn)相關(guān)面表H(node)=(H1，H2，…，Hn)。

步驟2。按照后處理模型數(shù)據(jù)結(jié)構(gòu)的設(shè)計(jì)要求，形成單元的各個(gè)單元面;將各個(gè)單元面分別加入節(jié)點(diǎn)相關(guān)面表;各面的第一個(gè)節(jié)點(diǎn)的整體編號(hào)即為其在節(jié)點(diǎn)相關(guān)面表中所屬的節(jié)點(diǎn)。對(duì)所有單元進(jìn)行上述處理后，就形成了結(jié)構(gòu)的最小節(jié)點(diǎn)相關(guān)面表。

步驟3。剔除內(nèi)部單元面，形成外部單元面即結(jié)構(gòu)邊界面。事實(shí)上，由于單元面的第一個(gè)節(jié)點(diǎn)是該單元面中節(jié)點(diǎn)號(hào)最小的節(jié)點(diǎn)，因此用上述方法形成的節(jié)點(diǎn)相關(guān)面表中，各表項(xiàng)是相互獨(dú)立的，即在不同的節(jié)點(diǎn)的相關(guān)面中沒有重合的單元面，重合的單元面只可能出現(xiàn)在相同節(jié)點(diǎn)的相關(guān)面中。對(duì)于某一節(jié)點(diǎn)而言，其相關(guān)面的數(shù)目是比較少的，這樣就大大減少了單元面重合判斷的計(jì)算量。對(duì)非空的Hi中的各相關(guān)面進(jìn)行比較，剔除其中包含的節(jié)點(diǎn)的整體編號(hào)完全相同的單元面，剩余的即為結(jié)構(gòu)的邊界面的集合。

1.2 網(wǎng)格消隱算法

目前主要的消隱算法有:Z緩沖區(qū)算法、畫家算法、掃描線算法和區(qū)域采樣算法。在這些經(jīng)典的消隱算法中，Z緩沖區(qū)(Z.buffer)算法比較適合于處理結(jié)構(gòu)化網(wǎng)格的消隱問題。

Z緩沖器中每個(gè)單元的值是對(duì)應(yīng)像素點(diǎn)所反映對(duì)象的z坐標(biāo)值。Z緩沖器中每個(gè)單元的初值取成z的極小值，幀緩沖器每個(gè)單元的初值可放對(duì)應(yīng)背景顏色的值。圖形消隱的過程就是給幀緩沖器和Z緩沖器中相應(yīng)單元填值的過程。在把顯示對(duì)象的每個(gè)面上每一點(diǎn)的屬性(顏色或灰度)值填入幀緩沖器相應(yīng)單元前，要把這點(diǎn)的z坐標(biāo)值和Z緩沖器中相應(yīng)單元的值進(jìn)行比較。只有前者大于后者時(shí)才改變幀緩沖器的那一單元的值，同時(shí)Z緩沖器中相應(yīng)單元的值也要改成這點(diǎn)的z坐標(biāo)值。如果這點(diǎn)的z坐標(biāo)值小于Z緩沖器中的值，則說明對(duì)應(yīng)像素已經(jīng)顯示了對(duì)象上一個(gè)點(diǎn)的屬性，該點(diǎn)要比考慮的點(diǎn)更接近觀察點(diǎn)。對(duì)顯示對(duì)象的每個(gè)面上的每個(gè)點(diǎn)都做上述處理后，可得到消除了隱藏面的圖。

虛擬試驗(yàn)?zāi)Ｐ途W(wǎng)格顯示及消隱的算法可以描述為

1)將三維模型外表面數(shù)據(jù)讀入動(dòng)態(tài)分配的數(shù)組，即上文建立的結(jié)構(gòu)的邊界面的集合;

2)建立集合中各個(gè)單元面的法向量;

3)利用OpenGL基本庫(kù)提供的圖元繪制函數(shù)構(gòu)造幾何圖元;

4)利用OpenGL庫(kù)函數(shù)實(shí)現(xiàn)實(shí)體的各種變化操作，包括幾何變換、投影變換、裁剪視區(qū)變換等，利用Z緩沖區(qū)算法來實(shí)現(xiàn)消隱;

5)建立光照模型。OpenGL把光照細(xì)分為4個(gè)獨(dú)立組成部分:環(huán)境光、發(fā)射光、散射光、反射光，在這4部分的共同作用下形成所要求的光照效果。

2 標(biāo)量結(jié)果可視化技術(shù)實(shí)現(xiàn)

2.1 概述

大部分的流場(chǎng)數(shù)據(jù)，都是通過數(shù)值求解流體模型方程N(yùn)avier-Strokes方程、Euler方程等得到的，其中絕大部分是標(biāo)量，如溫度、壓強(qiáng)等。對(duì)于流場(chǎng)標(biāo)量數(shù)據(jù)的可視化處理，最重要的是云圖的顯示方法［2］。

艦船虛擬試驗(yàn)生成的結(jié)果文件，一般數(shù)據(jù)量都比較大，甚至可能是大規(guī)模模型，而且可能會(huì)在某些區(qū)域的數(shù)據(jù)場(chǎng)變化劇烈。因此對(duì)于云圖的顯示，本文提供了兩種方法:①顏色是連續(xù)變化的，采用基于三角形的單元網(wǎng)格填充法;②顏色是分段離散的，采用改進(jìn)的基于三角形的區(qū)域填充法。

2.2 云圖顯示算法

2.2.1 顏色是連續(xù)變化的

1)查找表繪制。采用如圖1所示的線性的物理量-顏色映射方法。其中F min表示場(chǎng)量最小值，F(xiàn) max表示最大值，F(xiàn) min用藍(lán)色表示，F(xiàn) max用紅色表示，其中間值F mid一般用綠色表示。(F min，F(xiàn) max)整個(gè)區(qū)間內(nèi)的值使用的顏色都是采用線性插值。

圖1 顏色連續(xù)線性變化

2)模型云圖繪制。在繪制模型網(wǎng)格時(shí)，采用基于三角形的單元網(wǎng)格填充法，把每個(gè)網(wǎng)格節(jié)點(diǎn)上的物理量轉(zhuǎn)化成顏色表中對(duì)應(yīng)的某一種顏色后，單元內(nèi)部每個(gè)空間點(diǎn)上的值，運(yùn)用OpenGL相關(guān)函數(shù)，自動(dòng)使用向內(nèi)Lagrange線性插值的方式完成，見圖2。最終實(shí)現(xiàn)的效果圖如圖3。

2.2.2 顏色是離散分段的

1)查找表繪制。采用如圖4所示的線性的物理量-顏色映射方法，其中F min與F max分別表示場(chǎng)量的最小值與最大值，將F min用藍(lán)色color(0，0，255)表示，將F max用紅色color(255，0，0)來表示，其中間值F mid一般用綠色表示。(F min，F(xiàn) max)整個(gè)區(qū)間內(nèi)的按照劃分段數(shù)，按照一定的階梯分成有限個(gè)單一顏色。

圖4 顏色按梯隊(duì)逐漸變化

2)模型云圖繪制。在繪制模型網(wǎng)格時(shí)，采用改進(jìn)的基于三角形的區(qū)域填充法，依次獲得三角形網(wǎng)格(四邊形可以當(dāng)做兩個(gè)三角形處理)。每個(gè)節(jié)點(diǎn)上的物理量轉(zhuǎn)化成顏色表中對(duì)應(yīng)的某一種顏色后，如果它們?cè)谕粋€(gè)顏色范圍內(nèi)，則直接填充，此時(shí)三角形網(wǎng)格就一種單一的顏色;如果它們處于不同的顏色區(qū)間，求出交點(diǎn)，并依次用不同顏色填充三角形被分割之后的兩個(gè)部分，這兩種不同顏色的分界線，其實(shí)也就是等值線。最后實(shí)現(xiàn)的效果圖見圖5。

3 矢量結(jié)果可視化技術(shù)研究

在艦船虛擬試驗(yàn)的結(jié)果數(shù)據(jù)中，除了標(biāo)量數(shù)據(jù)，還有大量的矢量結(jié)果數(shù)據(jù)。矢量的可視化方法與標(biāo)量不盡相同，它既要顯示物理量的大小，又要顯示物理量的方向［3］。

圖5 顏色離散變化效果

矢量可視化方法一般可以分為點(diǎn)圖標(biāo)法、矢量線法、特征法和紋理法。按照艦船虛擬試驗(yàn)系統(tǒng)的要求，本文提供點(diǎn)圖標(biāo)法、矢量線法兩種最常用方法的矢量可視化技術(shù)。

3.1 點(diǎn)圖標(biāo)法

點(diǎn)圖標(biāo)是最直觀的顯示矢量數(shù)據(jù)的方法。對(duì)于每一個(gè)采樣點(diǎn)，采樣點(diǎn)處矢量的大小和方向可以通過具有大小和方向的圖標(biāo)(一般是箭頭)來表示。對(duì)矢量場(chǎng)進(jìn)行映射時(shí)，用箭頭的長(zhǎng)短和方向來表示矢量的大小和方向。

在二維矢量場(chǎng)中，箭頭可以很好地反映出矢量的方向和大小，在計(jì)算流體力學(xué)中箭頭表示法應(yīng)用較多。在三維矢量場(chǎng)中，對(duì)于數(shù)據(jù)量大、采樣密集的矢量場(chǎng)，過多的箭頭會(huì)導(dǎo)致像雜亂無章，顯示結(jié)果不理想，但如果顯示的箭頭太少，又會(huì)丟失一些數(shù)據(jù)信息。

有兩種方法可以解決這個(gè)問題，一種方法是使用低分辨率對(duì)矢量場(chǎng)重新采樣，減少需要顯示的箭頭數(shù)量;另一種方法是只對(duì)關(guān)注的局部矢量場(chǎng)數(shù)據(jù)進(jìn)行可視化處理。另外，流場(chǎng)的矢量場(chǎng)，通常在某些區(qū)域變化很大，導(dǎo)致這些區(qū)域的部分箭頭重疊，這樣的矢量圖在細(xì)節(jié)上不準(zhǔn)確，且難以辨識(shí)?？梢酝ㄟ^將矢量的大小非線性映射到箭頭長(zhǎng)短上的方法來解決這個(gè)問題。

綜上所述，矢量場(chǎng)點(diǎn)(箭頭)映射算法可描述如下。

1)計(jì)算出流場(chǎng)控件網(wǎng)格數(shù)據(jù)點(diǎn)上的矢量值;

2)搜索所有網(wǎng)格數(shù)據(jù)點(diǎn)矢量的最大值和最小值;

3)將矢量大小非線性映射為箭頭長(zhǎng)短;

4)計(jì)算箭頭坐標(biāo)的旋轉(zhuǎn)方向，旋轉(zhuǎn)箭頭。

其實(shí)現(xiàn)的效果見圖6。

圖6 點(diǎn)圖標(biāo)法效果

3.2 矢量線法

矢量線法源于試驗(yàn)型流場(chǎng)的可視化技術(shù)［4］。與試驗(yàn)方法類似，矢量線包括流線、時(shí)線、脈線和跡線。艦船虛擬試驗(yàn)數(shù)據(jù)，一般都需要進(jìn)行流線的顯示，因此本文將給出流線的顯示技術(shù)［5-6］。

流線指的是某一瞬時(shí)的這樣一條曲線，曲線上的每一個(gè)點(diǎn)的速度矢量，總在該點(diǎn)與曲線相切。數(shù)值積分法是生成流線最常用的方法，整個(gè)過程可以看成是選定初始位置，然后通過數(shù)值積分依次獲取所有的質(zhì)點(diǎn)，連接所有的質(zhì)點(diǎn)就能夠獲取隨時(shí)間變化的曲線。

流線的計(jì)算包括3個(gè)主要過程。

1)速度插值。在求解流線時(shí)，由于速度是定義在離散的網(wǎng)格點(diǎn)上，當(dāng)需要得到空間中任意一點(diǎn)的速度時(shí)，需要對(duì)速度插值。常用的插值算法主要有線性基函數(shù)插值方法、三次樣條法、面積加權(quán)法等。由于面積加權(quán)法一般用于三角網(wǎng)格單元中，精度可以控制的比較高，同時(shí)速度又快，所以本文采用面積加權(quán)方法。

面積加權(quán)法的原理如圖7所示。

圖7 面積加權(quán)法

設(shè)三角形內(nèi)一點(diǎn)P，要得到該點(diǎn)處的數(shù)值。通過將三角形的各頂點(diǎn)與該點(diǎn)相連形成三個(gè)三角形，然后利用它所對(duì)應(yīng)的小三角形的面積占總面積比例決定這個(gè)三角形的每個(gè)節(jié)點(diǎn)的權(quán)值，即

式中:wA，wB，wC——P點(diǎn)在△ABC中的面積坐標(biāo)。

要使得面積坐標(biāo)的3個(gè)分量(wA，wB，wC)在0～1之間，只有讓P點(diǎn)位于三角形△ABC內(nèi)。

2)數(shù)值積分。常用的流線積分方法有Euler法、二階Runge-Kutta法和四階Runge-Kutta法等?？梢砸罁?jù)計(jì)算精度和計(jì)算效率來進(jìn)行選擇流線積分的階數(shù)和步長(zhǎng)，為了解決流線跟蹤過程中的效率與精度之間的問題，本文采用變步長(zhǎng)的自適應(yīng)積分方法，見圖8，可以得到公式

圖8 流線積分過程中的自適應(yīng)步長(zhǎng)

積分(步長(zhǎng)h)有

假設(shè)A、|V0|分別為流線起點(diǎn)所在三角形單元以及流線起始點(diǎn)的速度大小，則可以推導(dǎo)出初始積分步長(zhǎng)為

3)三角網(wǎng)格中的點(diǎn)定位。三角網(wǎng)格中的點(diǎn)定位十分重要，其算法效率的高低直接關(guān)系到整個(gè)流線生成的算法的效率。假設(shè)一點(diǎn)P的坐標(biāo)，該點(diǎn)坐標(biāo)在物理空間內(nèi)，為了獲取P點(diǎn)所在的三角形單元以便插值計(jì)算得到該點(diǎn)的速度值，需要進(jìn)行一次點(diǎn)定位操作。這樣導(dǎo)致在整個(gè)流線構(gòu)造過程中每一個(gè)點(diǎn)都要執(zhí)行一次點(diǎn)定位算法。因此有必要對(duì)該過程進(jìn)行簡(jiǎn)化。本文充分利用Delaunay網(wǎng)格三角化的優(yōu)勢(shì)，Delaunay三角網(wǎng)格共一條邊最多只有兩個(gè)三角形，在設(shè)計(jì)數(shù)據(jù)結(jié)構(gòu)時(shí)，將三個(gè)頂點(diǎn)的指針和該三角形公共邊所在的另一個(gè)三角形的指針記錄在三角形數(shù)據(jù)結(jié)構(gòu)中。當(dāng)該三角形的邊位于邊界上時(shí)，則將相對(duì)三角形指針設(shè)為NULL。本文采用基于公共邊點(diǎn)定位快速算法，來判斷流場(chǎng)空間中的一點(diǎn)位于哪個(gè)三角形單元中。

基于公共邊點(diǎn)定位快速算法的描述如下。

步驟1。在三角形單元幾何的鏈表T上任意取一個(gè)三角形A，作為搜索出發(fā)點(diǎn)。

步驟2。取流線上的任意一點(diǎn)P。

步驟3。計(jì)算點(diǎn)P在A中的面積坐標(biāo):根據(jù)坐標(biāo)的正負(fù)可以判斷出點(diǎn)P與三角形對(duì)應(yīng)的位置關(guān)系;如果為正，則說明點(diǎn)在三角形內(nèi)，停止搜索，輸出A，結(jié)束。如果為負(fù)，說明點(diǎn)在此三角形外，并且能夠獲知公共邊信息。根據(jù)鄰接拓?fù)潢P(guān)系找到此三角形共此邊的鄰接三角形，并以鄰接三角形為出發(fā)點(diǎn)A繼續(xù)搜索，轉(zhuǎn)到步驟3。

步驟4。算法結(jié)束。

假設(shè)三角剖分后共由n個(gè)三角形單元組成，若利用傳統(tǒng)算法判斷點(diǎn)位于哪個(gè)三角形中的算法復(fù)雜度為O(n)。若領(lǐng)用本文所采用基于公共邊點(diǎn)定位快速算法，算法效率能夠達(dá)到通過比較可以得知，若三角形單元的個(gè)數(shù)越大，則本文所提出的算法的效率就越高，使用該算法可以大大加快流線的追蹤速度。

整個(gè)流線生成的具體算法如下。

步驟1。給定流線生成的種子點(diǎn)。

步驟2。對(duì)于每個(gè)種子點(diǎn)P，調(diào)用基于公共邊點(diǎn)定位快速算法，找出包含P點(diǎn)的三角形網(wǎng)格單元。

步驟3。當(dāng)P點(diǎn)不是臨界點(diǎn)時(shí)并且P點(diǎn)在流場(chǎng)邊界內(nèi)，則①P點(diǎn)處的速度矢量值通過插值求得，同時(shí)根據(jù)P點(diǎn)處的速度以及P點(diǎn)所在的網(wǎng)格大小確定需要積分的步長(zhǎng);②采用自適應(yīng)步長(zhǎng)數(shù)值積分方法計(jì)算下一步位置點(diǎn)P'，同時(shí)存儲(chǔ)P'點(diǎn)坐標(biāo)，令P=P';③采用基于公共邊點(diǎn)定位快速算法，找出包含P點(diǎn)的三角形網(wǎng)格單元。

步驟4。算法結(jié)束。

4 結(jié)論

在中國(guó)船舶科學(xué)研究中心開發(fā)的“艦船綜合水動(dòng)力虛擬水池試驗(yàn)技術(shù)系統(tǒng)”基礎(chǔ)上，研究了艦船虛擬實(shí)驗(yàn)結(jié)果的可視化關(guān)鍵技術(shù)，重點(diǎn)介紹了其中的關(guān)鍵算法，包括三維模型顯示、云圖顯示、矢量顯示等。通過這些關(guān)鍵技術(shù)的突破和有機(jī)整合，不僅成功實(shí)現(xiàn)了虛擬水池試驗(yàn)結(jié)果可視化的基本功能要求，而且能獲得更快的顯示速度。

上述可視化研究主要針對(duì)虛擬試驗(yàn)結(jié)果數(shù)據(jù)的后處理，今后將在試驗(yàn)結(jié)果的交互式方面做更多研究工作，使用戶可以通過交互方式修改原始數(shù)據(jù)、邊界條件或其它參數(shù)，實(shí)現(xiàn)用戶對(duì)虛擬實(shí)驗(yàn)過程的交互控制和引導(dǎo)。

［1］LEE SY，LEE S，NISHIDO S.Development of a new simulation method of mold filling based on a body-fitted coordinate system［J］.Cast Metals，2002，15(4):389-393.

［2］李建波，陳健云，林 皋.針對(duì)三維有限元數(shù)據(jù)場(chǎng)的精確后處理算法［J］.計(jì)算機(jī)輔助設(shè)計(jì)與圖形學(xué)學(xué)報(bào)，2004，16(8):1169-1175.

［3］WARD M O.A taxonomy of glyph placement strategies for multidimensional data visualization［J］.Information Visualization.2002(1):194-210.

［4］ROBERT SL.Interactive 3D flow visualization based on textures and geometric primitives［D］.Vienna，Austria:Vienna University of Technology，2004.

［5］吳曉莉.面向空間遙科學(xué)實(shí)驗(yàn)的流場(chǎng)可視化技術(shù)研究［D］.長(zhǎng)沙:國(guó)防科學(xué)技術(shù)大學(xué)，2007.

［6］崔加勝.流線可視化技術(shù)的研究與應(yīng)用［D］.鎮(zhèn)江:江蘇大學(xué)，2008.