康蘇海 ，金 生 ，劉 洋

（1.大連理工大學海岸和近海工程國家重點實驗室，大連116024；2.交通運輸部天津水運工程科學研究所工程泥沙交通行業(yè)重點實驗室，天津300456；3.天津市南港工業(yè)區(qū)開發(fā)有限公司，天津300456）

在二維顯示設備上表現(xiàn)三維流場，實際上是對三維流場場景進行平面投影。如果將三維流場以示蹤粒子方式表示，經過投影后得到的平面顯示流場中的示蹤粒子必然是前后重疊、上下遮蓋，無法清晰地表現(xiàn)三維流場，達不到模擬三維流場的效果。以往，為解決該問題，通常采用截面流場形式對三維流場進行模擬［1-3］，并分為水平截面流場和垂直截面流場兩種方式。無論哪種方式，都是首先根據(jù)給定的條件在三維場景中獲取截面，然后通過正投影方式將處于截面上的示蹤粒子的三維流場數(shù)據(jù)向截面進行投影，并在截平面內以二維流場的方式繪制，且流場繪制是在單獨的二維畫布上進行的。Mike21、Delft3D等水動力數(shù)值模擬軟件都是采用這種方式實現(xiàn)對三維流場模擬的。該方式在獨立窗口進行截面投影流場繪制，雖然通過多次獲得多個截面達到認識三維流場的目的，但是由于脫離所在三維場景而相互獨立模擬，不能實現(xiàn)多截面流場的同時對比比較，而且直觀性不強、表現(xiàn)不生動。

截面是空間三維曲面，在三維場景中具有明顯立體感，設想將獲得的截面及流場直接在其所在的三維場景中以三維方式進行繪制，可以解決上述通常截面流場顯示方式的弊端，不僅使截面流場有空間位置感，截面流場表現(xiàn)也更為真實，而且通過三維場景中多個截面流場的同時模擬，能夠更為方便地比較和認識三維流場特性，表現(xiàn)方式也更生動。下面首先論述復雜邊界及地形條件下的截面生成算法，然后討論截面流場的投影算法，最后以WPF為三維圖形平臺，實現(xiàn)三維截面流場的仿真模擬。

1 立體截面獲取算法

在三維場景中繪制截面，需根據(jù)給定的截面參數(shù)（水平截面的高程值或垂直截面的平面走向拐點坐標集合）求出截面和三維地形曲面的交點坐標集合，進而構造并生成截面構造點及單元鏈接關系，并以三維曲面方式在三維場景中繪制截面。同時需要記錄截面單元節(jié)點與原有節(jié)點的對應關系，以實現(xiàn)新舊節(jié)點對應關系，并由此得到構成截面的節(jié)點處各物理量值。

1.1 水平截面生成

三角形是最簡單的平面多邊形，由它生成的三棱柱體與某平面相交時，出現(xiàn)的相交情況最少，便于交點集合的判斷與截面生成?？紤]計算域平面上采用三角單元網(wǎng)格（四邊形或其他多邊形可以進一步細分為若干三角型），垂向分層方式（目前三維水動力計算基本采用垂向分層方式生成三維網(wǎng)格）構造三維計算網(wǎng)格，于是每個三維單元實際上是上下底不一定平行的三棱柱體。

在進行水平截面與三維計算域相交集合判斷時，因為水平面上采用了三角形網(wǎng)格，所以對每個水平單元所具有的3個節(jié)點分別進行相交情況的判斷（每個點的高程都有可能大于、等于、小于截面高程值），根據(jù)排列組合原理會出現(xiàn)27種組合方式，將27種組合方式進行歸類，最終得到10種本質上不同的組合方式（表1）。表1中L、E、G代表地形節(jié)點高程和指定截面高程的大小關系為小于、等于、大于，“代表標志”為標志字符串，列舉了節(jié)點描述情況中的某一種情形，例如第9種截面交點描述中的“LLE”包含了“LEL”、“ELL”2 種組合方式，代表了3種組合形式。

分析表1，發(fā)現(xiàn)當標志字符串中不含L時，該單元不必計入有效單元，否則當標志字符串中不含有G時，該單元直接計入有效單元，其他情況需要通過插值新節(jié)點構造新的有效單元。插值節(jié)點重新構造有效單元時，根據(jù)大、小、等出現(xiàn)的順序，分不同情況構造新單元，保證單元節(jié)點按逆時針排序。

根據(jù)上述原理，首先循環(huán)所有地形節(jié)點，比較節(jié)點高程與Zm的大小關系，得到標志字符（G、E、L），然后對模型域內全部地形單元進行循環(huán)判斷，根據(jù)單元3個節(jié)點的標志字符得到3個字母構成的單元標志字符串，只對其中含有L的單元進行新單元構建運算，因為沒有L的單元表示截面截到了地形以下，是無效的，而部分含有E的單元也是截到了原有節(jié)點上，該節(jié)點會在相鄰的單元循環(huán)時考慮到，因此仍然可忽略。針對單元標志字符串中既含有L又含有G或E的，需要進行交接點插值，構建新的節(jié)點，根據(jù)指定水平截面與各單元節(jié)點的有效相交點重新構造水平截面三角單元的節(jié)點及單元連接關系，從而完成水平截面的幾何構建。

表1 地形曲面節(jié)點高程與水平截面高程關系Tab.1 Relationship of terrain surface node elevation vs.horizontal cross-section elevation

1.2 垂直截面生成

垂直截面生成可以通過指定x值截取平行于y軸的垂直平面、指定y值截取平行于x軸的垂直平面或者在水平面內用戶通過鼠標操作指定任意折線構成的垂直截面。

在垂直截面生成時，由于折線點可能位于域外、域內、島內的任何位置，加上折線會穿過島嶼，甚至穿過域外部分，所以構建生成垂直截面的函需要考慮上述多種因素。

通過分析，垂直截面與域內的地形曲面網(wǎng)格三角形單元相交，會出現(xiàn)4種情況：（1）交線只通過單元某一頂點；（2）通過某一頂點和該頂點對面的邊；（3）和某邊重合；（4）穿過單元的2個邊。情況（1）必然和臨接三角形有（2）、（3）類型的相交，考察臨接單元時再進行統(tǒng)計，所以參與統(tǒng)計的情形是（2）、（3）、（4）。這 3 種情況下，對三棱柱的垂直切割必然會產生大小不一的矩形，這些矩形都因為分層而被劃分為NL個（層數(shù)）小4邊形（垂向兩邊平行，左右高程不一致導致上下不平行），圖1為垂直截面的某一部分示意圖。

將每一個小四邊形的左下角頂點和右上角頂點連接，便產生一系列三角形，它們構成垂直截面的三維網(wǎng)格曲面幾何要素。圖1-a表示指定垂直截面在平面走向與平面網(wǎng)格局部相交的示意圖，圖1-b為截面垂向與三維分層單元的相交情況。垂向被分為兩層，按截面定義時截線折線走向，確定垂直截面的法向方向，于是按右手法則確定三角單元的構成，如圖1中ABC、ACD、DCE、DEF為其中部分構成單元。

生成垂直截面的關鍵是如何根據(jù)指定的截面走向折線順序查找相交三角單元，求解交點序列，構建截面節(jié)點集合，從而生成最終的垂向截面。從建立通用的垂直截面生成出發(fā)，考慮到截面走向折線的隨意性，對折線中每個線段單獨處理，最后進行合成，于是問題轉化為如何求解任意線段截取垂面的問題。

求解線段與平面顯示域三角單元邊的交點序列，可以采用文獻［4］中關于方向追蹤法的描述方法，確定折線經過的所有單元，并用線段相交算法求出交點。問題的關鍵是采用該方法必須考慮到如何處理各類特殊情況，包括線段起始點不在域內、線段穿越了內部島嶼、線段中途穿越域外區(qū)域等。當線段起始點不在域內時，如果在域外，則求出線段與外邊界線的交點并記錄交點所在的單元號準備下一步方向追蹤法處理之用；如果起始點在內部島嶼內，則求出線段與該島嶼邊界線的交點并記錄交點所在的單元號。當線段穿越內部島嶼時，標記該段垂直截面的結束并查找線段與島嶼邊界的另一個交點，然后繼續(xù)下一段垂直截面的交點計算。當線段中途穿越域外區(qū)域時，標記該段垂直截面的結束并查找線段與外邊界的下一個臨接交點，然后繼續(xù)下一段垂直截面的交點計算。

自由表面三維流場的水面高程在整個模擬時間段往往會隨時間的推進而上下變化，因此生成垂直截面時，需要指定整個過程中最高水面高程作為垂直截面的上邊界，垂直截面流場模擬過程中根據(jù)水面變化值動態(tài)修改垂直截面的上邊界，達到自由水面追蹤之目的。

下面給出垂直截面生成步驟：

（1）初始化先判斷指定平面折線各拐點P1，P2，…，PM所在位置并進行記錄：落在計算域外賦值為-1、落在計算域賦值為0、落在內部島嶼內賦值為島嶼索引號（島嶼索引號必然為大于0的數(shù)值）；

（2）若完成所有段計算則轉到步（4），否則取線段PnPn+1進行計算；

（3）首先確定Pn點所在單元索引號，若Pn不在有效區(qū)域內，查找PnPn+1與邊界交點坐標及交點所在單元索引號。若Pn落在域內，根據(jù)坐標值利用分域法找到所在單元索引號；若Pn落在域外，分塊法定位最近交點并利用分域法［4］定位所在單元索引號；若Pn落在島嶼內部，對所在島嶼邊界線循環(huán)查找交點位置并利用分域法定位所在單元索引號；

（4）用遞歸法獲取線段PnPn+1在平面域內網(wǎng)格線的交點序列，同時記錄被島嶼或陸地分割的位置；

（5）分段構造垂直截面網(wǎng)格單元，落在島嶼內部或陸地上的截面部分不被記錄和顯示；

（6）返回（2）繼續(xù)下一線段的計算；

（7）合并所有分割的垂直截面，形成一個整體三維曲面集合在場景中顯示。

給定坐標值截取平行于x或y軸的垂直截面的情況是上述任意走向垂直截面生成的特殊情況，不再贅述。

2 截面流速投影

獲取截面后，采用對截面上每節(jié)點方式或在截面內隨機布點方式布設流場示蹤粒子，參見文獻［4］。確定示蹤粒子位置后，需要在三維空間內對粒子位置處的流速進行差值，參見文獻［5］。水平截面的流速投影較為簡單，實際上就是x、y方向的流速分量的合成流速，在此不做討論。

垂直截面上的流速投影需要經過若干次坐標變換。三維流速矢量向垂直截面投影的原理圖如圖2所示，其中xyz坐標系是模型現(xiàn)場坐標系將原點平移至粒子坐標位置的局部坐標系，ABCD為截取的垂直截面，截面與x軸夾角為α（逆時針方向為正）。

繪制流速矢量箭頭時，先將三維流速u、v、w 3個分量向ABCD截面投影，則u、v在ABCD的投影合成就是 Uc=u·cosα+v·sinα，w 是截面投影流速矢量c的另一個垂直分量，于是可以在截面局部坐標系（以示蹤粒子所在位置為坐標原點，截面走向為x′軸，平行于 z 軸方向為 y′軸的平面坐標系 x′o′y′）內繪制出截面流速矢量箭頭。通過坐標變換，可以較為容易地將局部坐標系中繪制的流速矢量箭頭換到現(xiàn)場坐標，最終完成截面流場質點流速矢量的繪制。

3 WPF實現(xiàn)截面流場繪制

截面流場要繪制在三維場景中，首先要選擇三維圖形平臺建立三維場景，微軟WPF是輕量級三維圖形平臺，可以實現(xiàn)三維場景的搭建，另外WPF集成與MS.NET框架下，隨微軟操作系統(tǒng)一起發(fā)布，具有較廣的應用環(huán)境［6-9］。

3.1 截面流場繪制

截面流場模擬時，只能描述截面二維流場特性，垂直于截面的另一速度分量不會得到體現(xiàn)，如果此時仍采用拉格朗日法模擬截面流場，已沒有實際意義，因為某時刻在截面上的粒子經過一時間段運動，完全有可能已“逃離”截面，因此已不可能得到其位于截面上的流速。既然截面流場只能代表位于本截面上的粒子的速度投影，采用歐拉法描述截面流場更為科學。

對于水平截面流場，通過前面的水平截面生成算法，得到了截面上所有節(jié)點對應的原有節(jié)點的索引號，因此不難在平面上找到對應的物理量，關鍵是如何在垂向插值得到截面節(jié)點位于哪一個分層或單元。考慮垂向平均分層的三維網(wǎng)格構建情形，首先讀取二維水面高程，減去相應點底床高程，得到凈水深，于是通過截面高程所在靜水深的位置，計算得到該點所在垂向層的位置，然后通過垂向插值計算，得到截面節(jié)點物理量值。對于通過位置插值得到的新的一類截面節(jié)點（即在構建水平截面時，截面與底部曲面相交而插值得到的點），在平面上需要進行二次插值以保證節(jié)點物理量的精度。

對于垂直截面流場，由于構成截面的大部分節(jié)點不和模型節(jié)點重合，需要進行空間差值?？臻g差值分平面差值和垂向差值兩步進行，首先根據(jù)垂向截面與平面網(wǎng)格的交點關系，在水平方向上差值得到通過粒子所在位置的垂線與上下兩個分層交點的物理量，然后根據(jù)粒子所在垂向層的位置比例，插值得到最終的物理量值。垂向截面流場模擬時，需要考慮模擬過程中由于水位變化造成動垂向動網(wǎng)格而引起的垂向插值變化，準確計算物理量的空間插值，還要特別注意到水位漲落引起的流速質點的顯示與隱藏。

3.2 截面疊加標量場

在截面流場模擬顯示的同時，有時需要對截面上某標量物理量如壓強、水位、鹽度、溫度、某污染物濃度值等同時進行顯示，這可以通過截面紋理進行表現(xiàn)。模擬開始前建立標量場對應的色譜表，根據(jù)各截面節(jié)點位置處標量值大小，通過插值得到對應的顏色值，將顏色值以填充畫刷的顏色體現(xiàn)在截面紋理上可以實現(xiàn)截面標量場的顯示。

4 多截面流場模擬及工程應用

為便于分析比較，在同一場景中同時繪制若干個水平和垂直截面，每一個截面可以有不同的截面填充設置和不同顏色的示蹤粒子，如此可對局部三維流場進行不同方向、不同深度的比較，更容易認識三維流場特性。文中選擇“長江、沱江匯合口三維水流特性研究”項目作為驗證模擬效果的工程實例。

江河主流和支流匯合區(qū)域水流稱為匯合口，由于兩股水流的匯流比不同，相互頂托作用有強有弱，攜帶泥沙含量不同，使得匯合口區(qū)域水流三維特性十分明顯，更使得該區(qū)域的泥沙沖淤變化復雜。長期以來對其三維水流結構研究較少、認識不足，因此研究人士匯合口區(qū)域的水流三維特性十分必要。三維計算關注的重點是沱江和長江匯合口區(qū)域三維水流結構，為泥沙沖淤理論分析提供依據(jù)。從水平、垂直截面流場分析區(qū)域三維水流特征。

圖3為不同高程值為215 m、222 m、230 m的3個水平截面流場，圖中白色代表最上一層截面流場，黑色代表最下一層截面流場，灰色代表中間層截面流場?？梢钥闯?，受水流流動慣性力及河道地形共同影響，底層截面流場比表層截面流場偏轉早、偏轉角度大，可見彎道水流表層和底層水平流速大小不等、流向存在明顯差異。

圖4以多截面顯示方式給出了沿程若干個垂直橫截面的投影流場。截面流場由兩側指向匯流交界面，即所謂“剪切層”，匯流交界面處有明顯的垂向流速，流速方向向下，越到深槽，垂直向下流速分量值越大；在匯合口上游，垂向存在兩個漩渦流，即靠近主流側的逆時針漩渦流和靠近支流一側順時針漩渦流，到水流恢復區(qū)截面流場只出現(xiàn)垂向逆時針漩渦流。

通過2種截面流場對比分析，匯合區(qū)流場表現(xiàn)出明顯的三維特性：匯合區(qū)存在橫向流速分量，橫向流方向由兩岸指向匯流交界線；匯合區(qū)深槽和左側淺灘區(qū)域水流垂向運動特點十分明顯，垂直流速方向在回流交界線向下，在左岸附近向上，形成左岸附近水體的翻滾現(xiàn)象；深槽內水體在沿江主流和垂向流動共同作用下呈螺旋前進態(tài)勢；匯合區(qū)右側淺灘水流總體變化平緩，垂向流動不明顯。

5 結論

通過對適應于復雜地形及邊界條件下的水平及垂直截面獲取算法的研究，找到了高效的、適應性強的截面生成算法并將截面繪制于三維場景中；在此基礎上，通過合理的方式在截面上布設示蹤粒子，并利用垂直投影方法在截面上繪制投影流場；同時對位于場景中的多個截面進行流場模擬，可方便地比較局部流場、分析流場的三維特性，比起原有截面流場模擬方法更為直觀、實用。標量場可作為材質顏色繪制在截面上，實現(xiàn)截面標量場的模擬。

［1］孔令勝，南敬實，荀顯超.平面三維顯示技術的研究現(xiàn)狀［J］.中國光學與應用光學，2009（2）：112-118.KONG L S，NAN J S，XUN X C.Research status quo of flat 3-D display technology［J］.Chinese Journal of Optics and Applied Optics，2009（2）：112-118.

［2］辛文杰，陳志昌，羅小峰.河口海岸數(shù)值模擬可視化系統(tǒng)［C］//中國海洋工程學會.第十二屆中國海岸工程學術討論會論文集.北京：海洋出版社，2005.

［3］陳立華，梅亞東，王現(xiàn)勛，等.基于 OpenGL 三維河網(wǎng)地形與數(shù)據(jù)場的可視化［J］.武漢大學學報：工學版，2007，40（3）：34-37.CHEN L H，MEI Y D，WANG X X，et al.Visualization of 3D river network terrain and data field based on OpenGL［J］.Engineering Journal of Wuhan University，2007，40（3）：34-37.

［4］康蘇海.利用 GDI+實現(xiàn)二維動態(tài)流場實時模擬［J］.水道港口，2009，30（6）：390-393.KANG S H.Technology of real time animation on 2D flow with GDI+［J］.Journal of Waterway and Harbor，2009，30（6）：390-393.

［5］KANG S H，JIN S.Application of Video Anaglyph Maker for 3-D Flow Simulation［J］.Journal of Hydrodynamics，2010（2）：289-294.

［6］楊珺，王繼成，劉然.立體圖像對的生成［J］.計算機應用，2007，27（9）：2 106-2 109.YANG J，WANG J C，LIU R.Stereo pairs creation［J］.Journal of Computer Applications，2007，27（9）：2 106-2 109.

［7］Adam Nathan.WPF揭秘［M］.北京：人民郵電出版社，2008.

［8］黎華，肖偉，黃海峰，等.三維真實感地形生成的關鍵技術研究［J］.測繪科學，2006，31（4）：57-60.LI H，XIAO W，HUANG H F，et al.Study on the key technology of generation of 3D realistic terrain［J］.Science of Surveying and Mapping，2006，31（4）：57-60.

［9］賴儀靈.WPF 全景體驗［J］.程序員，2007（3）：98-101.LAI Y L.WPF panoramic experience［J］.Programmer，2007（3）：98-101.