王家林，段紅俊

(重慶交通大學(xué) 土木工程學(xué)院，重慶 400074)

基于TetGen的箱梁三維四面體網(wǎng)格自動(dòng)劃分分析

王家林，段紅俊

(重慶交通大學(xué) 土木工程學(xué)院，重慶 400074)

以VC++為平臺(tái)，利用四面體網(wǎng)格生成接口軟件TetGen，研究了單箱單室變截面箱梁的高質(zhì)量四面體網(wǎng)格生成方法。用多連通平面多邊形面片描述橫截面、用字符串函數(shù)描述箱梁截面上各特征點(diǎn)的坐標(biāo)變化規(guī)律，將平面橫截面面片映射為空間橫截面面片，再自動(dòng)生成相鄰空間面片的連接面片，從而獲得箱梁的空間多面體。利用四面體網(wǎng)格劃分接口軟件TetGen對(duì)箱梁空間多面體進(jìn)行四面體網(wǎng)格剖分，自動(dòng)生成高質(zhì)量的四面體網(wǎng)格；通過(guò)將網(wǎng)格數(shù)據(jù)寫(xiě)為inp文件導(dǎo)入Abaqus 軟件進(jìn)行可視化檢查和有限元計(jì)算驗(yàn)證了網(wǎng)格的有效性。

橋梁工程；TetGen；四面體網(wǎng)格；三維實(shí)體網(wǎng)格劃分

隨著我國(guó)經(jīng)濟(jì)的飛速發(fā)展，作為基礎(chǔ)設(shè)施的橋梁工程也在不斷發(fā)展，箱型梁橋也日益受到關(guān)注。箱梁有截面抗扭剛度大，截面效率高，有效地抵抗正負(fù)彎矩等優(yōu)點(diǎn)[1]。基于以上優(yōu)點(diǎn)，箱型梁橋在我國(guó)得到廣泛的使用。

在橋梁計(jì)算中，大多采用有限元軟件(橋梁博士[2]、GQJS[3]、Midas Civil[4]等)對(duì)橋梁進(jìn)行建模。這些軟件的建模與分析一般采用梁?jiǎn)卧M(jìn)行。但梁?jiǎn)卧Ｐ碗y以揭示箱梁的多種復(fù)雜變形效應(yīng)，例如扭轉(zhuǎn)畸變、翹曲、剪力滯等。隨著計(jì)算機(jī)硬、軟件技術(shù)的發(fā)展和分析要求逐步地提高，采用空間實(shí)體單元對(duì)橋梁進(jìn)行高精度分析逐漸成為一種新的發(fā)展趨勢(shì)。趙安等[5]對(duì)橋梁結(jié)構(gòu)空間模型進(jìn)行了有限元分析；交通部公路科學(xué)研究院利用空間實(shí)體單元對(duì)橋梁三維預(yù)應(yīng)力空間進(jìn)行分析，并研發(fā)出針對(duì)這一問(wèn)題的BridgeKF系統(tǒng)[6]；劉艷萍[7]利用空間實(shí)體單元模型對(duì)橋梁預(yù)應(yīng)力鋼筋混凝土結(jié)構(gòu)局部進(jìn)行了ANSYS二次開(kāi)發(fā)；余錢(qián)華等[8]利用實(shí)體單元模擬了連續(xù)剛構(gòu)橋的施工過(guò)程。

基于VC++面向?qū)ο蠹夹g(shù)，筆者通過(guò)研究單箱單室箱型橋梁的結(jié)構(gòu)特點(diǎn)，首先利用典型截面自動(dòng)生成空間多面體，再利用成熟的四面體網(wǎng)格劃分接口軟件TetGen[9]對(duì)箱梁的多面體進(jìn)行四面體剖分，自動(dòng)生成高質(zhì)量的四面體網(wǎng)格，最后導(dǎo)入Abaqus 軟件完成有限元計(jì)算來(lái)驗(yàn)證網(wǎng)格的有效性。

1 TetGen的網(wǎng)格劃分原理

TetGen能對(duì)三維實(shí)體結(jié)構(gòu)四面體剖分，其目的是為了產(chǎn)生合適的四面體網(wǎng)格。TetGen的代碼采用VC++編寫(xiě)，這些代碼可被編譯成直接運(yùn)行的程序，也可被編譯成靜(動(dòng))態(tài)連接庫(kù)。目前的操作系統(tǒng)基本都能運(yùn)行TetGen，如Macos，Windows等。TetGen使用的算法約束delaunay四面體剖分。

1.1 約束delaunay四面體剖分

約束delaunay四面體剖分是用一個(gè)個(gè)四面體網(wǎng)格將三維區(qū)域Ω進(jìn)行劃分。三維區(qū)域Ω可能有著很多形狀上的變化，可能包含內(nèi)邊界，或包含洞。

點(diǎn)集V中的兩個(gè)點(diǎn)p，q可見(jiàn)，但若有一個(gè)約束平面f，且點(diǎn)p，q兩點(diǎn)分別位于平面f的兩邊，p，q間的線(xiàn)段交叉于平面f(圖1)，那么頂點(diǎn)p，q將不可見(jiàn)。那么約束delaunay四面體則就是由這個(gè)點(diǎn)集V的頂點(diǎn)組成。像這樣一個(gè)四面體剖分被稱(chēng)作是約束delaunay四面體剖分。

圖1 約束delaunay四面體剖分示意Fig.1 Constraint delaunay tetrahedral subdivision schemes

1.2 四面體網(wǎng)格質(zhì)量

TetGen處理四面體網(wǎng)格質(zhì)量問(wèn)題的算法是利用shewchuk的delaunay優(yōu)化算法。給出一系列三維的頂點(diǎn)、約束段、不少于90°的輸入角和小平面，則利用delaunay優(yōu)化算法就能夠產(chǎn)生delaunay四面體網(wǎng)格，并且這個(gè)delaunay四面體網(wǎng)格的半徑邊緣比不大于2。該算法產(chǎn)生的網(wǎng)格普遍比較合適，并且對(duì)消除最小或最大二面角的四面體非常有效。

網(wǎng)格劃分中需要有質(zhì)量衡量標(biāo)準(zhǔn)。一般情況下，用邊界縱橫比來(lái)衡量網(wǎng)格劃分的質(zhì)量。邊界縱橫比η(t)是最大邊的長(zhǎng)度lmax與最小高度hmin的比值，即η(t)=lmax/hmin。判斷四面體網(wǎng)格質(zhì)量，邊界縱橫比要愈小愈好。在TetGen衡量網(wǎng)格劃分質(zhì)量的標(biāo)準(zhǔn)中則采用Miller，Teng，Wang Talmor 和Walkington提出的半徑邊緣比的標(biāo)準(zhǔn)。

一個(gè)四面體t有且只有一個(gè)外接球，令r為四面體外接球的半徑，d為四面體最短邊的長(zhǎng)度，則這個(gè)四面體的半徑邊緣比ρ(t)為

(1)

式中：θmin為四面體t的最小面角。

在判斷四面體網(wǎng)格質(zhì)量的應(yīng)用中，半徑邊緣比要愈小愈好。

2 箱梁空間多面體的建立

2.1 平面橫截面面片幾何特征分析

用多連通平面多邊形面片描述橫截面(以單箱單室梁橫截面為例)，可分為外圍平面多邊形與內(nèi)圍平面多邊形，將其特征點(diǎn)依次編號(hào)。

筆者利用VC++編程對(duì)此截面面片節(jié)點(diǎn)進(jìn)行編號(hào)，步驟如下：

1)利用AutoCAD畫(huà)出箱梁截面面片，將其保存成dxf形式。

2)由于dxf文件存取著所畫(huà)截面的線(xiàn)段信息，故編寫(xiě)了函數(shù)，用來(lái)從dxf文件中讀取線(xiàn)段信息到線(xiàn)鏈表中。

3)利用函數(shù)將各線(xiàn)段生成平面面片，再將其節(jié)點(diǎn)進(jìn)行編號(hào)。

(1)定義一個(gè)具有x,y兩個(gè)坐標(biāo)的點(diǎn)列表；

(2)將線(xiàn)鏈表中的線(xiàn)段連入點(diǎn)列表中；

(3)讀取一條線(xiàn)段作為多邊形的首尾，讀取下一段線(xiàn)段時(shí)，若這條線(xiàn)段與多邊形首尾的點(diǎn)有相同的坐標(biāo)，將其添加到多邊形中，連接的線(xiàn)段將從線(xiàn)段列表中刪除；若這條線(xiàn)段與多邊形首尾的點(diǎn)沒(méi)有相同的坐標(biāo)，則接著讀取再下一條線(xiàn)段，此線(xiàn)段將不從線(xiàn)段列表中刪除；

(4)若點(diǎn)列表已閉合，說(shuō)明一個(gè)多邊形已生成。則需生成一個(gè)新的多邊形，如此循環(huán)(3)，(4)兩步，直至線(xiàn)段列表為空；

(5)計(jì)算多個(gè)多邊形的面積，若出現(xiàn)面積為負(fù)值，則表示多邊形為順時(shí)針連接而成，取面積為其相反數(shù)，并將多邊形調(diào)整為逆時(shí)針連接；

(6)將多邊形的面積從大到小排列，并從大到小對(duì)各個(gè)多邊形的各個(gè)點(diǎn)進(jìn)行編號(hào)，并存儲(chǔ)其所對(duì)應(yīng)的坐標(biāo)(圖2)。

圖2 截面面片節(jié)點(diǎn)編號(hào)Fig.2 Section node number

2.2 箱梁截面特征及處理方法

對(duì)等截面而言，直接利用函數(shù)將全局坐標(biāo)從二維變?yōu)槿S，則所有特征點(diǎn)均發(fā)生變化。而對(duì)變截面而言，首先將全局坐標(biāo)從二維變?yōu)槿S，然后將用字符串函數(shù)描述箱梁截面上各關(guān)鍵點(diǎn)的y坐標(biāo)變化規(guī)律。y坐標(biāo)變化規(guī)律包括箱梁根部高度的變化曲線(xiàn)、箱梁底板厚度變化曲線(xiàn)。例如變化的點(diǎn)包括1，2，17，18，19，20；而1，2點(diǎn)的變化規(guī)律為箱梁根部高度的變化曲線(xiàn)；17～20點(diǎn)的變化規(guī)律為箱梁底板厚度變化曲線(xiàn)。

2.3 將平面橫截面面片映射為空間橫截面面片

利用函數(shù)在z坐標(biāo)每隔一定距離得到一個(gè)空間面片(僅包含兩個(gè)空間面片)，如圖3。

圖3 空間面片示意Fig.3 Space diagram

2.4 生成連接面片

以VC++為平臺(tái)利用函數(shù)自動(dòng)生成相鄰空間面片的連接面片，即可生成箱梁的空間多面體，如圖4。

圖4 箱梁空間多面體示意Fig.4 Box girder polyhedron space diagram

3 四面體網(wǎng)格的自動(dòng)生成

在Windows系統(tǒng)中，將TetGen編譯為靜態(tài)庫(kù)“Tetgen.lib”。創(chuàng)建一個(gè)“Tetgenio”的輸入對(duì)象“Win”，將數(shù)據(jù)導(dǎo)入輸入對(duì)象“Win”，調(diào)用函數(shù)“Tetrahedralize()”創(chuàng)建四面體網(wǎng)格，并將結(jié)果輸出到對(duì)象“Wout”中，實(shí)現(xiàn)了三維實(shí)體結(jié)構(gòu)網(wǎng)格劃分。

3.1 Tetgenio結(jié)構(gòu)的輸入與輸出數(shù)組

用多連通平面多邊形面片描述橫截面(以單箱單室梁橫截面為例)，可分為外圍平面多邊形與內(nèi)圍平面多邊形，將其特征點(diǎn)依次編號(hào)。Tetgenio結(jié)構(gòu)的作用是傳入空間實(shí)體的數(shù)據(jù)并輸出四面體及其數(shù)據(jù)。Tetgenio結(jié)構(gòu)的輸入與輸出數(shù)組用來(lái)保存點(diǎn)、多邊形、面片和四面體等內(nèi)容。它為一個(gè)C++類(lèi)，包括一些數(shù)據(jù)域和函數(shù)。Tetgenio結(jié)構(gòu)的輸入輸出數(shù)組如下:

1)“REAL *pointlist”頂點(diǎn)的坐標(biāo)數(shù)組

第一個(gè)頂點(diǎn)的x坐標(biāo)是index[0]，y坐標(biāo)是index [1]，z坐標(biāo)是index[2]；第二個(gè)頂點(diǎn)的x坐標(biāo)是index[3]，y坐標(biāo)是index [4]，z坐標(biāo)是index[5]，依次類(lèi)推，每個(gè)頂點(diǎn)占3個(gè)REALS。

2)“REAL *pointattributelist”頂點(diǎn)的屬性數(shù)組

每一個(gè)頂點(diǎn)的屬性占用一個(gè)“numberofpointat tributes”的REALS。

3)“int *pointmarkerlist”頂點(diǎn)標(biāo)識(shí)符數(shù)組

每一個(gè)點(diǎn)的標(biāo)識(shí)符為一個(gè)整數(shù)。

4)“facet *facetlist”面片數(shù)組

每一個(gè)面都是一個(gè)面類(lèi)型的對(duì)象。

5)“int *facetmarkerlist”面片標(biāo)識(shí)符數(shù)組

每個(gè)面的一個(gè)標(biāo)識(shí)符占一個(gè)整數(shù)。

6)“REAL *holelist”孔洞的數(shù)組

第一個(gè)洞的x坐標(biāo)是index[0]，y坐標(biāo)是index [1]，z坐標(biāo)是index[2]；第二個(gè)洞的x坐標(biāo)是index[3]，y坐標(biāo)是index [4]，z坐標(biāo)是index[5]，依次類(lèi)推，每個(gè)洞占3個(gè)REALS。

7)“int *tetrahedronlist”四面體角數(shù)組

第一個(gè)四面體的第一個(gè)角為index[0]，一個(gè)四面體一共有4個(gè)角。每個(gè)四面體占用numberofcorners整型。

8)“REAL *tetrahedronvolumelist”四面體體積列表

每一個(gè)四面體占一個(gè)REAL。

3.2 TetGen命令行選項(xiàng)

TetGen有多個(gè)命令行選項(xiàng)，它通過(guò)對(duì)命令行選項(xiàng)設(shè)置參數(shù)來(lái)實(shí)現(xiàn)四面體網(wǎng)格的自動(dòng)劃分。以下為常用的命令行選項(xiàng)：

1)-p：對(duì)三維實(shí)體模型進(jìn)行約束delaunay四面體剖分；

2)-q：網(wǎng)格控制參數(shù)，主要設(shè)置最小半徑邊緣比，默認(rèn)為2；

3)-a：控制產(chǎn)生四面體的最大體積；

4)-r：重構(gòu)和優(yōu)化之前產(chǎn)生的網(wǎng)格文件。

4 變截面箱型梁劃分實(shí)例

該連續(xù)剛構(gòu)橋采用單箱單室。箱梁高度是由距離根部高7.5 m處按照1.9次拋物線(xiàn)變化至端部及跨中高3.0 m處。箱梁底板厚度是由距離根部厚100 cm處按照二次拋物線(xiàn)變化至端部及跨中厚30 cm處。

根據(jù)圖紙，利用AutoCAD畫(huà)出單箱單室梁橫截面圖。利用VC++編程對(duì)此單箱單室梁橫截面面片節(jié)點(diǎn)進(jìn)行編號(hào)，如圖5。

圖5 截面面片節(jié)點(diǎn)編號(hào)Fig.5 Section node number

以VC++為平臺(tái)，利用映射點(diǎn)的原理，根據(jù)變截面變形曲線(xiàn)，映射出新的單箱單室截面面片；其中1，2點(diǎn)變形曲線(xiàn)為y—z1.9×3.519 7×10-5；17～20變形曲線(xiàn)為y—z1.9×3.519 7×10-5+z2×2.31 4×10-6。

將網(wǎng)格數(shù)據(jù)寫(xiě)為inp文件導(dǎo)入到Abaqus[8]中，便得到了三維實(shí)體網(wǎng)格效果見(jiàn)圖6。

圖6 箱梁的三維網(wǎng)格Fig.6 Three-dimensional grid section of box girder

為檢查三維網(wǎng)格劃分的準(zhǔn)確性，筆者將導(dǎo)入的空間四面體進(jìn)行有限元計(jì)算。將梁的一端設(shè)為固定端，另一端的中間位置施加垂直向下的荷載1 000 N，計(jì)算結(jié)果如圖7。

圖7 一跨梁的有限元計(jì)算結(jié)果Fig.7 Finite element calculation results for girder of one span

5 結(jié) 語(yǔ)

基于VC++面向?qū)ο蠹夹g(shù)，筆者通過(guò)研究單箱單室箱型橋梁的橫截面特征和截面的變化規(guī)律，建立了一種變截面箱型橋梁空間四面體自動(dòng)生成方法；利用成熟的四面體網(wǎng)格劃分接口軟件TetGen對(duì)箱梁的空間多面體進(jìn)行四面體剖分，自動(dòng)生成高質(zhì)量的四面體網(wǎng)格；通過(guò)將網(wǎng)格數(shù)據(jù)寫(xiě)為inp文件導(dǎo)入Abaqus 軟件，進(jìn)行可視化檢查；并用有限元計(jì)算驗(yàn)證了網(wǎng)格的有效性。

筆者以單箱單室箱梁為例檢驗(yàn)了基于TetGen的箱梁三維四面體網(wǎng)格自動(dòng)劃分的思想，但文中思想并不局限于單箱單室，對(duì)于雙箱等也適用，僅僅是特征點(diǎn)增多。

[1] 郭金瓊,房貞政,鄭振.箱形梁設(shè)計(jì)理論[M].北京:人民交通出版社,2008. GUO Jinqiong, FANG Zhenzheng, ZHENG Zhen.DesignTheoryofBoxGirder[M]. Beijing: China Communications Press,2008.

[2] 周水興,王小松,田維峰,等.橋梁結(jié)構(gòu)電算[M].北京:人民交通出版社,2013. ZHOU Shuixing, WANG Xiaosong, TIAN Weifeng, et al.ComputingforBridgeStructures[M]. Beijing: China Communications Press,2013.

[3] 交通部公路科學(xué)研究所.公路橋梁結(jié)構(gòu)設(shè)計(jì)計(jì)算系統(tǒng)CQJS用戶(hù)手冊(cè)[M].北京:交通部公路科學(xué)研究所,1999. The Ministry of Communications Highway Research Institute.HighwayBridgeStructureDesignandCalculationSystemCQJSUserManual[M]. Beijing: The Ministry of Communications Highway Research Institute,1999.

[4] 邱順冬.橋梁工程軟件Midas Civil應(yīng)用工程實(shí)例[M].北京:人民交通出版社,2011. QIU Shundong.BridgeEngineeringApplicationSoftwareMidasCivilEngineeringPractice[M]. Beijing: China Communications Press,2011.

[5] 趙安,呂建鳴.橋梁試驗(yàn)的結(jié)構(gòu)空間有限元分析—橋梁荷載試驗(yàn)分析系統(tǒng)BLT的研究[J].土木建筑工程信息技術(shù),2013,5(1):33-40. ZHAO An, LV Jianming. Finite element analysis of spatial structure in bridge testing—study of bridge load testing analysis system BLT[J].JournalofInformationTechnologyinCivilEngineeringandArchitecture,2013,5(1):33-40.

[6] 楊昀,向?qū)W建,王磊.連續(xù)剛構(gòu)橋預(yù)應(yīng)力仿真分析系統(tǒng)BridgeKF的研究[C]//中國(guó)公路學(xué)會(huì)橋梁和結(jié)構(gòu)工程分會(huì).中國(guó)公路學(xué)會(huì)橋梁和結(jié)構(gòu)工程分會(huì)2004年全國(guó)橋梁學(xué)術(shù)會(huì)議論文集.北京:人民交通出版社,2004:764-772. YANG Yun, XIANG Xuejian, WANG Lei. Prestressed concrete continuous rigid-frame bridges owing to the simulation analysis system of BridgeKF research [C]// The Chinese Institution of Highway Bridge and Structures Engineering Branch.ChinaRoadBridgeandStructuresEngineeringBranch2004NationalConferenceonBridgeAcademic. Beijing: China Communications Press,2004:764-772.

[7] 劉艷萍.橋梁預(yù)應(yīng)力鋼筋混凝土結(jié)構(gòu)局部有限元分析的ANSYS二次開(kāi)發(fā)[D].武漢:華中科技大學(xué),2005. LIU Yanping.PartialPrestressedConcreteBridgeStructureFiniteElementAnalysisofANSYSSecondaryDevelopment[D]. Wuhan: Huazhong University of Science and Technology,2005.

[8] 余錢(qián)華,朱鵬志,周偉,等.連續(xù)剛構(gòu)橋?qū)嶓w建模方法介紹與結(jié)果分析[J].中外公路,2007,27(2):69-73. YUN Qianhua, ZHU Pengzhi, ZHOU Wei, et al. Entity modeling method is introduced and the result analysis of long span continuous rigid frame[J].JournalofChina&ForeignHighway,2007,27(2):69-73.

[9] SI Hang.TetGen:AQualityTetrahedralMeshGeneratorandA3DDelaunayTriangulator(Version1.5:User’sManual)[R]. Berlin: Weierstrass Institute for Applied Analysis and Stochastics, 2013.

Automatic Meshing of Box Girder Based on TetGen 3D Tetrahedral Element

WANG Jialin，DUAN Hongjun

(School of Civil Engineering，Chongqing Jiaotong University,Chongqing 400074，P.R.China)

Based on VC++ platform, the tetrahedral mesh generation interface software—TetGen was used to study high quality tetrahedral element generation method of box girder with varied cross section of generation single box with single cell. The girder cross section was described with facet of multiply connected plane polygon and letter and symbol string function was used to describe coordinate variation at respective feature points on the box girder cross section. First of all, it was mapped from plane cross-section facet to space cross-section facet. Secondly the connected facet of adjacent space facet was generated automatically, thus spacial polyhedron of box girder was achieved. The tetrahedral mesh generation interface software—TetGen was used to subdivide spacial polyhedron. High quality tetrahedron element was generated automatically. The network data was written in inp document which was subscribed in software Abaqus which allowed for visual inspection and verification of effectiveness of this grid.

bridge engineering;TetGen； tetrahedral mesh； the mesh of three-dimensional entity

10.3969/j.issn.1674-0696.2016.05.01

2015-05-20；

2015-11-03

王家林(1968—)，男，重慶人，教授，博士，主要從事土木工程、結(jié)構(gòu)力學(xué)方面的研究。E-mail:jialinwang@163.com。

段紅俊(1991—)，女，山東煙臺(tái)人，碩士，主要從事土木工程、結(jié)構(gòu)力學(xué)方面的研究。E-mail:1091782733@qq.com。

U448.21+3

A

1674-0696(2016)05-001-04