尹改梅，汪新慶，蔣子陽

（1.中煤科工集團(tuán)武漢設(shè)計研究院，武漢 430064；2.中國地質(zhì)大學(xué)(武漢) 資源學(xué)院，武漢 430074；3.重慶地質(zhì)礦產(chǎn)研究院，重慶 400042）

復(fù)雜邊界的FLAC二維地質(zhì)模型構(gòu)建

尹改梅1，汪新慶2，蔣子陽3

（1.中煤科工集團(tuán)武漢設(shè)計研究院，武漢 430064；2.中國地質(zhì)大學(xué)(武漢) 資源學(xué)院，武漢 430074；3.重慶地質(zhì)礦產(chǎn)研究院，重慶 400042）

FLAC3D是一款優(yōu)秀的巖土力學(xué)分析軟件，但其建模效率較低，尤其面對復(fù)雜地質(zhì)條件時，其應(yīng)用受到一定程度的限制。本文分析了FLAC建模命令的特點，提出一種新的較簡單的建模方法，并開發(fā)了相應(yīng)的計算機(jī)輔助建模程序，該方法通過坐標(biāo)變換、復(fù)雜邊界點坐標(biāo)均勻化處理、逐列定義子網(wǎng)格坐標(biāo)，可生成具有復(fù)雜邊界的地質(zhì)模型，且便于賦材料參數(shù)及應(yīng)力邊界條件參數(shù)，有效提高工作效率。

坐標(biāo)變換；FLAC；復(fù)雜模型

0 引言

唐FLAC是美國ITASCA公司開發(fā)的一款用于二維、三維固體力學(xué)計算的顯式有限差分程序，該程序采用快速拉格朗日求解模式，在計算中網(wǎng)格隨單元的變形而更新。FLAC是一款中高端產(chǎn)品，除了能分析一般的固體力學(xué)問題外，特別適合于巖土介質(zhì)的大變形、破壞行為的研究[1-2]。FLAC具有強(qiáng)大的計算分析能力以及良好的二次開發(fā)性能，目前已成為巖土力學(xué)領(lǐng)域應(yīng)用最為廣泛的數(shù)值模擬軟件之一，F(xiàn)LAC的建模功能相對較弱，在實際應(yīng)用中經(jīng)常需要運(yùn)用FLAC內(nèi)置的建模命令對規(guī)則的模型進(jìn)行網(wǎng)格調(diào)整，造成工作效率低下，甚至很難建立符合需求的復(fù)雜的地質(zhì)模型。在FLAC二維模型構(gòu)建方面，王希寶（2008）分析了DXF格式文件包含的信息，通過編程實現(xiàn)了由DXF文件直接轉(zhuǎn)為FLAC命令文件，提高了建模效率，較好地解決了不規(guī)則邊界邊坡建模問題。但這些方法大多沒有考慮到賦材料參數(shù)、應(yīng)力邊界條件及后處理的便利性。本文以含一條舒緩波狀斷層的平面問題建模為例，分析復(fù)雜邊界FLAC二維模型建模方法，并開發(fā)相應(yīng)的計算機(jī)輔助建模程序，以實現(xiàn)復(fù)雜地質(zhì)模型的快速構(gòu)建。

1 FLAC內(nèi)置的建模命令分析

地質(zhì)工程數(shù)值模擬的基本步驟可包括：首先進(jìn)行地質(zhì)條件分析，建立概化地質(zhì)模型，然后賦予地質(zhì)材料參數(shù)、初始條件和邊界條件，進(jìn)而進(jìn)行數(shù)值計算、模型調(diào)整及后處理分析。其中，由于地質(zhì)體空間形態(tài)的復(fù)雜性，建模工作往往很困難，且費(fèi)時費(fèi)力[10]。

FLAC可接受的合法的模型網(wǎng)格形式只有四邊形，可以接受四邊形退化而來的三角形，但不可接受凹四邊形，而且網(wǎng)格長寬之比不宜大于5。規(guī)則網(wǎng)格的生成命令為“grid i,j”，其中i為列數(shù)目，即x方向上的網(wǎng)格劃分?jǐn)?shù)，j為行數(shù)目，即y方向上的網(wǎng)格劃分?jǐn)?shù)。在FLAC里，“grid”命令只能用一次，即用來確定總的網(wǎng)格數(shù)，然后根據(jù)實際條件調(diào)整網(wǎng)格形態(tài)，可以用“model null”命令修剪掉不需要的網(wǎng)格。調(diào)整網(wǎng)格的命令包括“gen line, arc, circle, table”等，這些命令可以畫出邊界并調(diào)整網(wǎng)格幾何形態(tài)。但有時這些命令也難以得到想要的幾何形態(tài)，而且由于網(wǎng)格單元是按照網(wǎng)格點的列、行編號（i,j）確定的，即使運(yùn)用“gen line”或“gen table”命令得到了不規(guī)則的內(nèi)、外邊界，但給后面賦材料參數(shù)、邊界條件也帶來了困難。

2 復(fù)雜模型構(gòu)建方法

本文以某平面中存在的舒緩波狀斷層建模為例（圖1），通過以下步驟進(jìn)行復(fù)雜平面模型的構(gòu)建。該斷層破碎帶相對較窄，延伸方向NE，由“gen line”或“gen table”等命令很難得到理想的建模效果。

（1）獲得復(fù)雜邊界線坐標(biāo)

對于矢量格式如CAD格式基礎(chǔ)數(shù)據(jù)，可以在命令欄輸入list命令，即可輸出邊界點的坐標(biāo)。對于位圖格式的圖片，可以應(yīng)用surfer等軟件，對邊界線進(jìn)行數(shù)字化，得到邊界點坐標(biāo)。

（2）確定模型網(wǎng)格范圍參數(shù)

根據(jù)研究的具體需求確定模型邊界范圍及分析尺度，即確定模型整體區(qū)域、網(wǎng)格的尺寸。由此，可初步確定需要劃分網(wǎng)格總的數(shù)目，包括列數(shù)i和行數(shù)j。用“grid i, j”命令生成模型總網(wǎng)格。

圖1 不規(guī)則邊界地質(zhì)模型示意圖

（3）確定模型內(nèi)部不同地質(zhì)體的子網(wǎng)格參數(shù)

FLAC按照組（group）劃分具有不同力學(xué)性質(zhì)的地質(zhì)體并賦予材料參數(shù)，每一組可根據(jù)上、下邊界線上的點坐標(biāo)劃分子網(wǎng)格的行數(shù)和列數(shù)。如圖1中，整個區(qū)域范圍大致為長25000m，寬25000m，區(qū)域網(wǎng)格共劃分為84列×84行，初始網(wǎng)格大小約為300m×300m，分A、B、C3個組，3個組以網(wǎng)格行號區(qū)別，A、B、C組對應(yīng)的網(wǎng)格分別為44～84行，42～43行，1～41行?？梢钥闯?，B組的上、下邊界均為不規(guī)則的曲線，且距離較近，依靠FLAC內(nèi)置網(wǎng)格調(diào)整命令來劃分B組子網(wǎng)格較困難。為了使得每行每列盡可能大小均勻，要將每個組的上下邊界線的點進(jìn)行均勻化處理，使之在橫向上均勻分布（圖2）。

圖2 邊界線均勻化處理

a——原始數(shù)據(jù)點 b——均勻化處理后的數(shù)據(jù)點

（4）借助計算機(jī)程序生成建模命令流

根據(jù)均勻化的上、下邊界，結(jié)合該組的寬度大小，生成FLAC命令流，并定義組名。C組均勻化處理后的上、下邊界各有85個數(shù)據(jù)點，因而C組從左至右依次劃分為84個子網(wǎng)格，定義每個子網(wǎng)格的FLAC命令為：“gen x1,y1 x2,y2 x3,y3 x4,y4 i=i1,i2 j=j1,j2”。（x1,y1）、(x2,y2)、(x3,y3)、(x4,y4)為順時針排列點坐標(biāo)，其中（x1,y1）為左下角的點坐標(biāo)。FLAC中x方向或y方向上網(wǎng)格點數(shù)目等于網(wǎng)格數(shù)加1。“j=j1, j2”表示該條語句的子網(wǎng)格點行編號范圍，對于C組而言，均為“j=1,42”。 “i=i1,i2”為子網(wǎng)格點列數(shù)，對于本例中的C組而言， 84個子網(wǎng)格依次為：“i=1,2”、“i=2,3”、 “i=3,4”…… “i=84,85”。所以C組網(wǎng)格命令流生成后，即可以定義組：“group 'granite 1' i =1,85 j=1,41”。圖中A、B組網(wǎng)格劃分、組定義方法與C組方法相同，不同在于對應(yīng)的網(wǎng)格行數(shù)不一樣，即B組的網(wǎng)格點行號范圍：“j=42,44”，A組的網(wǎng)格點行號范圍：“j=44,85”。

本文應(yīng)用Visual C#編程實現(xiàn)了上述處理過程，該程序讀入文本文件格式的邊界線坐標(biāo)，經(jīng)過線性插值處理后輸出新的均勻化的邊界線坐標(biāo)，保存為文本文檔格式，依次得到各組的均勻分布的上、下邊界線數(shù)據(jù)點。再由該組上、下邊界線坐標(biāo)點，結(jié)合用戶在交互界面輸入的網(wǎng)格點的行編號參數(shù)（圖3），生成該組的子網(wǎng)格，將各組子網(wǎng)格命令流合并，輸入FLAC程序，生成的網(wǎng)格如下（圖4，5）。

圖3 FLAC輔助建模程序

圖4平面模型網(wǎng)格

圖5 地質(zhì)材料分組

3 坐標(biāo)變換

本文的案例中，復(fù)雜邊界線是橫向展布的，依據(jù)上、下邊界線建立不同地質(zhì)體子網(wǎng)格，與FLAC現(xiàn)有網(wǎng)格調(diào)整命令相比，其優(yōu)點是可以根據(jù)網(wǎng)格點的行編號快速、方便地定義各組的空間位置。如果復(fù)雜邊界線為縱向展布，設(shè)置模型參數(shù)時，可以將x，y方向加以互換。如果復(fù)雜邊界展布方向與x,y方向有一定夾角，則可以對整個區(qū)域進(jìn)行坐標(biāo)變換。

將原坐標(biāo)系xoy繞原點沿逆時針方向旋轉(zhuǎn)θ度，變成座標(biāo)系 sot。設(shè)平面上某點p，在原坐標(biāo)系中的坐標(biāo)為 (x, y), 旋轉(zhuǎn)后的新坐標(biāo)為(s, t)。經(jīng)坐標(biāo)變換后，線L在sot坐標(biāo)系下就成為沿著S方向展布（圖6）。新舊坐標(biāo)系變換公式為：

圖6 坐標(biāo)旋轉(zhuǎn)示意圖

應(yīng)力邊界條件往往是影響數(shù)值模擬的重要影響因素，一般情況下，斷裂延伸方向與最大主應(yīng)力方向基本一致，但與大地坐標(biāo)系的x軸、y軸方向可能有一定夾角，可以通過坐標(biāo)變換使最大主應(yīng)力方向與新坐標(biāo)系s軸或t軸一致。經(jīng)過坐標(biāo)變換后，不僅便于劃分相對均勻的網(wǎng)格，而且有利于施加應(yīng)力邊界條件、初始條件。

4 結(jié)束語

對具有復(fù)雜內(nèi)、外邊界的平面應(yīng)力、平面應(yīng)變問題，F(xiàn)LAC自帶的模型網(wǎng)格調(diào)整命令往往不能有效滿足需求，筆者提出通過坐標(biāo)變換，邊界均勻化處理，逐列生成子網(wǎng)格，并開發(fā)了相應(yīng)的計算機(jī)程序，可快速有效地進(jìn)行網(wǎng)格劃分以及對具有不同材料參數(shù)的地質(zhì)體進(jìn)行分組賦參數(shù)。

對于含有不規(guī)則斷裂的模型，坐標(biāo)變換后，不僅有利于網(wǎng)格劃分，而且便于對模型施加地應(yīng)力邊界條件。

對于不適宜進(jìn)行坐標(biāo)變換的情形，如邊坡問題，有3種方法建模：①用文獻(xiàn)[1]的方法進(jìn)行處理；②參照本文方法，進(jìn)行左、右邊界的均勻化處理，逐行生成子網(wǎng)格；③采用FLAC內(nèi)置命令“gen line”或“gen table”生成封閉的復(fù)雜邊界，再借助region命令修剪多余網(wǎng)格。

[1]劉 波，韓彥輝. FLAC原理、實例與應(yīng)用指南[M]. 北京：人民交通出版社，2005.

[2]Itasca Consulting Group Inc. FLAC user’s manual [M]. Minnesota∶ State University of Minnesota, USA, 1997.

[3]王希寶，李天斌，許 勇. FLAC二維模型簡化建模方法及應(yīng)用[J]. 中國地質(zhì)災(zāi)害與防治學(xué)報, 2008, 19(1)∶82～85.

[4]鄭文棠，徐衛(wèi)亞，童富果，石安池.復(fù)雜邊坡三維地質(zhì)可視化和數(shù)值模型構(gòu)建[J]. 巖石力學(xué)與工程學(xué)報, 2007, 26(8)∶ 1633～1644.

[5]崔芳鵬，胡瑞林，等.基于surfer平臺的FLAC3D復(fù)雜三維地質(zhì)建模研究[J].工程地質(zhì)學(xué)報，2008，16(5)∶699～70.

[6]汪吉林, 丁陳建, 吳圣林. 基于FLAC3D的復(fù)雜地貌三維地質(zhì)建模[J].地質(zhì)力學(xué)學(xué)報，2008,14(2)∶149～157.

[7]廖秋林, 曾錢幫, 劉 彤, 等. 基于ANSYS 平臺復(fù)雜地質(zhì)體FLAC3D模型的自動生成[ J]. 巖石力學(xué)與工程學(xué)報,2005, 24 (6)∶ 1010～1013.

[8]姬保靜, 劉偉韜, 任 強(qiáng). FLAC3D中復(fù)雜模型的三維建模技巧[J]. 安全與環(huán)境學(xué)報,2006,6(s)∶69～72.

[9]徐文杰, 胡瑞林, 李厚恩, 等. CAD軟件在工程地質(zhì)三維建模中的應(yīng)用[J]. 工程地質(zhì)學(xué)報, 2007, 15( 2)∶279～283.

[10]唐輝明, 晏鄂川, 胡新麗. 工程地質(zhì)數(shù)值模擬的理論與方法[M]. 武漢∶ 中國地質(zhì)大學(xué)出版社, 2001.

FLAC 2-D Modeling Geological Bodies with Complex Boundaries

YIN Gaimei1，WANG Xinqing2，JIANG Ziyang3

（1.Wuhan Design & Research Institute of China Coal Technology & Engineering Group, Wuhan 430064;2．School of Earth Resources, China University of Geosciences, Wuhan 430074; 3．Chongqing Institute of Geology and Mineral Resources, Chongqing 400042）

FLAC is the excellent software for rock and soil mechanics numeric simulating, but it does not perform well in modeling, especially to complex geological model, so its application is limited to some extent. In this paper, the modeling commands in FLAC software were analyzed and a new simple modeling method for FLAC was proposed and corresponding computer-aided modeling code was developed. To transform original coordinates of complex boundaries, get the boundaries uniform distributed transversely, and defi ne the coordinates of each sub-grid column by column, then a model with complex boundaries can be built and it is easy to assign the material parameters and geo-stress boundary conditions. Our method can contribute to the numeric simulation of FLAC.

Coordinate transformation；FLAC；Complex model

TU45; O241

A

1007-1903（2013）02-0049-04

中央高?；究蒲袠I(yè)務(wù)費(fèi)專項基金項目（編號：CUGL120258）

尹改梅(1980－ )，女，碩士，工程師，從事地質(zhì)勘查工程；E mail：ygmman@tom.com