金屬礦山三維應(yīng)力計(jì)算建模技術(shù)

梅婷婷，張電吉，周春梅，楊 曼,王家國

(1.武漢工程大學(xué)環(huán)境與城市建設(shè)學(xué)院,湖北 武漢 430074；2.湖北宜化集團(tuán)礦業(yè)有限責(zé)任公司，湖北 宜昌 443000)

0 引 言

礦山三維地質(zhì)體建模技術(shù)研究，是“數(shù)字地球”、“數(shù)字礦山”的核心組成部分，是現(xiàn)代礦山信息化研究的熱點(diǎn)和難點(diǎn).近年來，國內(nèi)外學(xué)者吸收科學(xué)可視化、計(jì)算機(jī)圖形圖像學(xué)中的可視化技術(shù)在三維地質(zhì)建模理論與方法方面進(jìn)行了許多有益的探索[1].侯恩科[2]、吳立新、王占剛[3]、齊安文[4]、侯遵譯[5]、張菊明[6]、柴賀軍[7]、毛善軍[8]等學(xué)者對(duì)建模理論和方法的研究和探討，為三維可視化模型的研制和開發(fā)奠定了基礎(chǔ).針對(duì)金屬礦山開挖的建模過程都極為復(fù)雜，尤其是建立三維實(shí)體模型，如果采取傳統(tǒng)的建模方法，不僅需要花費(fèi)大量的時(shí)間，而且還容易出現(xiàn)錯(cuò)誤[9].FLAC3D在解決巖土問題上具有很大的優(yōu)越性，是求解巖土工程問題的最理想工具之一,但是在模型建立以及單元網(wǎng)格劃分等前處理方面卻一直是軟肋.相比之下，ANSYS通用有限元軟件對(duì)于復(fù)雜模型的實(shí)體建模，其前處理功能優(yōu)點(diǎn)尤為突出.筆者以“點(diǎn)→線→面→體”自底而上的實(shí)體建模技術(shù)路線為指導(dǎo)，以某金礦為應(yīng)用實(shí)例，試圖利用ANSYS 軟件對(duì)復(fù)雜工程地質(zhì)體建立相應(yīng)數(shù)值模型，再通過數(shù)據(jù)轉(zhuǎn)化實(shí)現(xiàn) FLAC3D有限元模型的自動(dòng)生成.針對(duì)建模過程中所遇到的諸如數(shù)據(jù)獲取和模型簡化等問題提出了具體的解決方法，以供其它類似礦山建立應(yīng)力場計(jì)算模型時(shí)參考.

1 實(shí)體建模理論與方法

目前，三維地質(zhì)建模使用的方法基本是基于空間曲面插值擬合方法來實(shí)現(xiàn)，主要有趨勢面分析法、克里格法、形函數(shù)法和人工神經(jīng)網(wǎng)絡(luò)模擬法等[10].本模型主要在形函數(shù)和趨勢面分析法的理論基礎(chǔ)上，借助ANSYS自底向上的建模技術(shù)，建立三維地質(zhì)體模型.ANSYS各圖元之間的關(guān)系如下：

(低→高)

1.1 由關(guān)鍵點(diǎn)生成樣條曲線

以研究目的為前提，對(duì)模型進(jìn)行適當(dāng)?shù)暮喕?，結(jié)合具體的資料，確定模型數(shù)據(jù)采樣點(diǎn)(關(guān)鍵點(diǎn))的數(shù)量及三維坐標(biāo).在ANSYS實(shí)體建模中，所有線均用非均勻有理B樣條表示[11]，利用公式(1)對(duì)數(shù)據(jù)采樣點(diǎn)進(jìn)行插值計(jì)算，形成由關(guān)鍵點(diǎn)生成的樣條曲線.

(1)

Pi為數(shù)據(jù)采樣點(diǎn)；Ni,k(u)是k次規(guī)范B樣條基函數(shù)；wi為權(quán)因子(權(quán)重系數(shù))；若基函數(shù)的節(jié)點(diǎn)是非均勻的，則稱為非均勻有理B樣條.

1.2 由樣條曲線構(gòu)造樣條曲面

對(duì)于復(fù)雜金屬礦山地表模型的建立，主要運(yùn)用“蒙皮法”.蒙皮法傳統(tǒng)上稱為放樣，其任務(wù)就是擬合一張曲面(即“皮”)，使其通過一組有序的截面線.本文引用文獻(xiàn)[12-13]對(duì)樣條曲線進(jìn)行蒙皮操作，借助公式(2)生成一張雙三次B樣條曲面.具有(m+1)(Nmax+2)個(gè)控制點(diǎn)的雙三次B樣曲面的方程為：

(2)

式(2)中：di,j(i=0,1,…,Nmax+1)是曲面的控制點(diǎn)；Ni,3(u)是定義在同一節(jié)點(diǎn)矢量上的B樣條函數(shù)，Nj,3(v)是定義在另一個(gè)參數(shù)方向的節(jié)點(diǎn)矢量V上的B樣條函數(shù).

針對(duì)礦巖地質(zhì)界面擬合，主要采用趨勢面分析方法.本模型采用白潤才[14]提出的趨勢面分析法的改進(jìn)方法，由關(guān)鍵點(diǎn)插值生成曲面.如果平面上分布的某參數(shù)用Z表示，那么趨勢面分析可簡單看成尋求一個(gè)方程，按二次多項(xiàng)式擬合空間曲面為：

Z=f(X,Y)

(3)

P(X,Y)=C00+C10X+C20Y+C30X2+

C40XY+C50Y2

(4)

?Q/?Crs(r,s=0,1,2)

(5)

可得到6個(gè)未知數(shù)的6個(gè)方程組，求得式(5)的系數(shù)Crs.

式(3)～(5)中，(Xi,Yi)為巖體所在范圍某平面的數(shù)據(jù)點(diǎn)坐標(biāo)；Zi為點(diǎn)(Xi,Yi)上的估測點(diǎn)；Q為所有數(shù)據(jù)點(diǎn)的擬合值P(Xi,Yi)和估測值Zi的誤差的平方和；n為估點(diǎn)數(shù)；w(di)=w[(Xi-a)2+(Yi-b)2]為權(quán)值，是距離的函數(shù).

1.3 三維實(shí)體的生成

由已生成的關(guān)鍵點(diǎn)建立樣條曲線，再根據(jù)樣條曲線構(gòu)造相應(yīng)的曲面，最后由一系列連續(xù)封閉的曲面定義體.強(qiáng)大的布爾運(yùn)算工具可以實(shí)現(xiàn)實(shí)體之間加、減、分類、搭接和分割等復(fù)雜運(yùn)算，大大提高了建立復(fù)雜地質(zhì)體三維模型的效率.對(duì)于實(shí)體模型的網(wǎng)格劃分，ANSYS提供了功能強(qiáng)大的控制工具，還可對(duì)實(shí)體模型圖直接劃分網(wǎng)格[15].

2 采用ANSYS實(shí)體建模的步驟

對(duì)于整個(gè)有限元模型的建立，其主要的建模步驟如下：

a.建模方案的確定.在開始進(jìn)入ANSYS之前，首先收集資料并采集多元數(shù)據(jù)，選擇合適的單元類型，考慮適當(dāng)?shù)木W(wǎng)格密度，并確定單元的實(shí)常數(shù)、材料屬性和單元坐標(biāo)系.

b.模型簡化.在建模前，首先熟悉工程地質(zhì)特征和明確自己的建模對(duì)象，通過經(jīng)驗(yàn)定性分析影響礦體開挖的主要因素和次要因素，確定是否具備簡化條件，即所謂的模型簡化.

c.幾何模型的建立.進(jìn)入前處理(PREP 7)模塊，利用幾何數(shù)據(jù)，遵循從“點(diǎn)→線→面→體”自外而內(nèi)、自底而上的建模技術(shù)，生成符合幾何形狀的獨(dú)立實(shí)體模型，并借助布爾操作將各個(gè)獨(dú)立的實(shí)體模型適當(dāng)?shù)倪B接在一起.

d.網(wǎng)格剖分.考慮研究對(duì)象的影響程度，依次設(shè)置不同的網(wǎng)格劃分精度，對(duì)實(shí)體模型劃分網(wǎng)格，生成節(jié)點(diǎn)和單元.

e.接口與轉(zhuǎn)換.運(yùn)用ANSYS to FLAC程序,將d中的節(jié)點(diǎn)和單元導(dǎo)入到FLAC有限差分軟件，進(jìn)行計(jì)算.

3 建模過程中存在的難點(diǎn)及解決方法

ANSYS的建模是一個(gè)復(fù)雜的過程，計(jì)算模型的建立是應(yīng)力計(jì)算中最為困難的一步，模型最終將直接影響到計(jì)算分析結(jié)果的真實(shí)性和可靠性[16].

3．1 數(shù)據(jù)獲取

建模的關(guān)鍵取決于數(shù)據(jù)采樣點(diǎn)坐標(biāo)的準(zhǔn)確性，所以首先要解決的問題是如何獲取精確數(shù)據(jù).由于資料的詳盡程度不同，圖形上可能有部分?jǐn)?shù)據(jù)缺失，此時(shí)可依據(jù)周圍地形地貌的走勢，借助“線性內(nèi)插法”對(duì)數(shù)據(jù)采樣點(diǎn)的高程值進(jìn)行估算.在結(jié)合水平平面圖以及剖面圖讀取數(shù)據(jù)點(diǎn)坐標(biāo)時(shí)，要注意角度轉(zhuǎn)化問題.在地表建模過程中，參照現(xiàn)場實(shí)際情況，對(duì)描述露天采場特征的數(shù)據(jù)要進(jìn)行適當(dāng)?shù)男拚?，力求與工程實(shí)際開挖情況相符.

3．2 模型簡化

復(fù)雜礦山的三維地質(zhì)體建模最大的難點(diǎn)就在模型簡化.目前對(duì)斷層的三維表示與模擬研究還處于研究和探索階段，要模擬形態(tài)各異的金礦，也存在技術(shù)上的困難.鑒于模型主要用于模擬礦體開挖，從系統(tǒng)工程的角度出發(fā)，堅(jiān)持“重點(diǎn)研究局部，局部服從整體”的原則，可對(duì)礦體模型和斷裂帶模型進(jìn)行適當(dāng)?shù)暮喕?

a.三維實(shí)體的建模技術(shù).該模型主要用于應(yīng)力場計(jì)算，借助基于面的構(gòu)模方法的思路，將三維實(shí)體假想成不規(guī)則封閉曲面構(gòu)成的六面體模型.雖然礦體的形狀、產(chǎn)狀、品味和大小各異，在建模過程中可其簡化為六面體模型.根據(jù)礦體大致呈脈狀以及礦體品位可圈定擬建礦體模型的外輪廓線，構(gòu)造封閉曲面，從而建立三維實(shí)體模型.值得注意的是：第一，將圈定范圍內(nèi)不同品位的礦體與存在的夾層，均視為品味相同的礦體，盡量使礦體形狀由不規(guī)則變?yōu)橐?guī)則.第二，為了便于建模和對(duì)模型進(jìn)行修改，一般選取實(shí)體輪廓線的特征點(diǎn)為數(shù)據(jù)采樣點(diǎn)(關(guān)鍵點(diǎn))，宜盡量保持采樣點(diǎn)數(shù)量的統(tǒng)一、采樣點(diǎn)編號(hào)順序的統(tǒng)一.

b.三維曲面的建模技術(shù).斷層一般由斷層面來表示其產(chǎn)狀、斷層帶的寬度與延伸尺寸由控制點(diǎn)處得幾何特征來定義.由地質(zhì)構(gòu)造判斷斷裂破碎帶的穩(wěn)定性，確定其對(duì)數(shù)值計(jì)算的影響程度.如果該斷裂帶較為穩(wěn)定，屬于控礦斷裂帶且礦體伴隨斷裂帶產(chǎn)生，可假定該斷裂帶靠近礦體的一側(cè)與礦體大致重合，避免夾層對(duì)網(wǎng)格剖分產(chǎn)生不利影響.考慮到斷裂帶的控制長度、寬度變化以及影響程度，取適當(dāng)?shù)臄?shù)值作為斷裂帶的平均寬度.以假定重合的曲面為基準(zhǔn)面，利用ANSYS的VDRAG命令將基準(zhǔn)面向外平移，平移的距離即為平均寬度，最終生成斷裂帶模型.

除此之外，ANSYS軟件的建模方法較為靈活，可利用菜單創(chuàng)建有限元模型或利用命令流的方式建模.在進(jìn)行布爾運(yùn)算之前，要判斷模型是否滿足運(yùn)算的要求，特別是在邊界處.而且盡量減少布爾運(yùn)算的次數(shù)和充分利用模型的對(duì)稱性，使模型的幾何形式盡量簡單.

4 應(yīng)用實(shí)例

某金礦產(chǎn)于伊犁-中天山板塊北部的古生代中天山北緣活動(dòng)大陸邊緣的吐拉素?cái)嘞菖璧刂衃17],區(qū)域地層及構(gòu)造線方向均呈北西西向展布，是古生代的低硫型淺成低溫?zé)嵋航鸬V.該金礦中心點(diǎn)地理坐標(biāo)為：東經(jīng)81°36′30″，北緯44°13′45″，礦區(qū)面積0.348 km2.如圖1所示，北露天采礦場為山坡地形，整體呈一個(gè)長軸近于南北的橢圓形凹陷露天采坑.礦體在地表呈-向西突出的帶狀分布，呈脈狀產(chǎn)出，傾向東-北東，傾角60°～80°.在露天開采范圍內(nèi)，與采礦工程相關(guān)的地質(zhì)構(gòu)造有2組，即F2斷裂破碎帶和不整合接觸面.在礦區(qū)的北部F2產(chǎn)狀傾向?yàn)?00°，向南部傾向變?yōu)?0°，傾角為55°～85°，局部傾角為90°，大體與礦體延伸方向一致.結(jié)合上述工程地質(zhì)特性，筆者主要以地表建模為例進(jìn)行詳細(xì)介紹.

圖1 北露天采礦場構(gòu)造地質(zhì)圖

4.1 地表建模

為了建模方便，地表模型采用相對(duì)坐標(biāo)系XOY，取某金礦地形地質(zhì)及井上下對(duì)照?qǐng)D(1∶2000)，將圖上的大地坐標(biāo)(4899000,27548200)定為XOY平面的起點(diǎn)，起點(diǎn)坐標(biāo)為(200.00,0.00)，坐標(biāo)精度保留兩位小數(shù)，以米為單位，地表模型范圍為1 200 m×1 600 m.地表建模的思路：在地形圖上構(gòu)造一系列大小相等、規(guī)則的M×N柵格(GRID)，然后根據(jù)已知等高線和控制點(diǎn)的高程值，通過線性內(nèi)插的方法求得每個(gè)柵格點(diǎn)的高程值，從而得到各柵格點(diǎn)的三維坐標(biāo)，按照點(diǎn)成線、線成面的步驟，建立地表光滑曲面.

其中，幾何模型的建模過程如下：

a.在地形圖上建立20.00 m×20.00 m的規(guī)則柵格，將每個(gè)柵格四周的拐角點(diǎn)作為數(shù)據(jù)采集點(diǎn)，本圖采集點(diǎn)共4 941個(gè).

b.將數(shù)據(jù)采集點(diǎn)相鄰的等高線值作為已知值，通過數(shù)據(jù)采集點(diǎn)做兩條等高線的公垂線，測量公垂線以及采樣點(diǎn)與其中一條等高線的最短垂直距離，利用線性內(nèi)插的方法，估算各采樣點(diǎn)的高程值.最后結(jié)合所在地形圖的XOY平面，確定各點(diǎn)的X坐標(biāo)和Y坐標(biāo)，并按照ANSYS軟件指定格式編輯命令流.

c.依據(jù)文中提到的非均勻有理B樣條插值方法，由關(guān)鍵點(diǎn)生成樣條曲線，如圖2所示；借助“蒙皮法”將樣條曲線生成蒙皮似的光滑曲面，如圖3所示；由光滑曲面向下拖拉1 200 m生成體，如圖4所示；在800 m高程處建立工作面，利用布爾運(yùn)算切割各拖拉體后，再將其合并，如圖5所示.

d.由于北露天采礦場的地形較為復(fù)雜，建模方法是依據(jù)北露天采礦場各勘探線剖面圖，讀取各臺(tái)階拐點(diǎn)處的X坐標(biāo)和Z坐標(biāo)，再將該點(diǎn)投影至該金礦地形地質(zhì)及井上下對(duì)照?qǐng)D(1∶2000)中各勘探線上的相應(yīng)位置，讀取該點(diǎn)Y坐標(biāo).同樣利用多次線性內(nèi)插法，由各臺(tái)階拐點(diǎn)的坐標(biāo)估算北露天采礦場區(qū)域內(nèi)各刪格四周的拐角點(diǎn)三維坐標(biāo)，做法同b所示.

圖2 點(diǎn)成線

圖3 線成面

圖4 面成體

圖5 修改后的體

e.在建模過程中考慮不整合界線，并對(duì)兩組不同的巖性分別進(jìn)行物理力學(xué)性質(zhì)參數(shù)賦值.不整合界限模型主要的建模思路同地表模型一致，即以地層不整合界線的特征點(diǎn)作為關(guān)鍵點(diǎn)，按照點(diǎn)成線和線成面的原理生成邊界不整合界面，最后利用VSBA命令，對(duì)地層不整合界面與整個(gè)地表模型布爾運(yùn)算，從而形成含邊界不整合界面的地表模型，如圖6所示.

圖6 簡化后的邊界不整合界面模型

4.2 計(jì)算模型的生成

根據(jù)已有的地質(zhì)勘探、礦體資源儲(chǔ)量、采礦工藝和流程等資料，依照相同的原理，在地表模型建成的基礎(chǔ)上，分別建立斷裂構(gòu)造、礦體、塌落體、上覆巖梁(隔離層)、礦房和礦柱等三維地質(zhì)體模型.通過定義實(shí)體模型的各種屬性以及借助ANSYS軟件強(qiáng)大的網(wǎng)格劃分功能，將該實(shí)體模型轉(zhuǎn)化為有限元計(jì)算模型.運(yùn)用ANSYS to FLAC3D程序,將有限元計(jì)算模型中的節(jié)點(diǎn)和單元導(dǎo)入到FLAC3D有限差分軟件，進(jìn)行礦體開挖的數(shù)值模擬計(jì)算，如圖7所示.

圖7 有限元計(jì)算模型

5 結(jié) 語

以上研究的基于ANSYS的金屬礦山三維地質(zhì)建模方法，不僅克服了建立復(fù)雜計(jì)算模型的困難，而且操作簡單、易行，所要的基本數(shù)據(jù)僅為幾何圖形數(shù)據(jù)點(diǎn)，大大減少了建模所需的時(shí)間、精力，具體如下:

a.采用“點(diǎn)→線→面→體”自外而內(nèi)、自底而上的建模技術(shù)，依次建立地表、斷裂構(gòu)造、礦體、塌落體、上覆巖梁(隔離層)、礦房和礦柱等三維地質(zhì)體模型，再利用布爾運(yùn)算組裝各個(gè)獨(dú)立的實(shí)體模型.不僅保證各地質(zhì)體模型之間的獨(dú)立性和有限元模型的整體性，而且便于模型的修改和更新，提高了建模效率.

b.ANSYS軟件的前處理功能適用于建立龐大復(fù)雜的模型,而FLAC3D有限差分軟件計(jì)算模擬速度較快,充分結(jié)合兩種軟件的優(yōu)點(diǎn)，解決了復(fù)雜工程地質(zhì)體模型的數(shù)值計(jì)算難題.該方法已運(yùn)用于某金礦的實(shí)際生產(chǎn)中，實(shí)踐證明其對(duì)礦山企業(yè)的安全生產(chǎn)及采礦方法的優(yōu)化具有指導(dǎo)意義.

c.遵循“重點(diǎn)研究局部，局部服從整體”的簡化原則，以建模的目的為出發(fā)點(diǎn)，以實(shí)際的工程地質(zhì)情況為前提，對(duì)模型進(jìn)行簡化.但是在礦山三維地質(zhì)建模的過程中，至今還沒有形成一套完善的簡化規(guī)律，該問題仍待進(jìn)一步探討和研究.

參考文獻(xiàn)：

[1] Simon W Houlding. 3D Geosciences modeling：computer technique for geological characterization[M]．Hong Kong: Springer.1993.

[2] 侯恩科，吳立新．面向地質(zhì)建模的三維體元拓?fù)鋽?shù)據(jù)模型研究[J]．武漢大學(xué)學(xué)報(bào)：信息科學(xué)版，2002,10(5):467-472.

[3] 王占剛，曹代勇．基于改進(jìn)三棱柱模型的復(fù)雜地質(zhì)體3D建模方法[J]．中國煤田地質(zhì)，2004,16(1):3-6.

[4] 齊安文，吳立新．基于類三棱柱的三維地質(zhì)模擬與拓?fù)溲芯縖J]．礦山測量，2003(3):65-67.

[5] 侯遵譯，許振邦．建立在三維空間上的地學(xué)計(jì)算機(jī)圖形圖像系統(tǒng)[J]．中國數(shù)學(xué)地質(zhì)進(jìn)展，1995，25(7):168-179.

[6] 張菊明．空間地質(zhì)模型的設(shè)計(jì)和顯示[J]．中國數(shù)學(xué)地質(zhì)進(jìn)展，1995,25(4):189-195.

[7] 柴賀軍，黃地龍，黃潤秋，等．巖體結(jié)構(gòu)三維可視化及其工程應(yīng)用研究[J]．巖土工程學(xué)報(bào)，2001,23(2):217-220.

[8] 毛善軍，許友志，張海榮，等．空間地質(zhì)模型及其可視化系統(tǒng)[J]．中國數(shù)學(xué)地質(zhì)進(jìn)展，2003，25(7):168-179.

[9] 盧宏建，羅濤，甘德清．基于ANSYS的礦山三維動(dòng)態(tài)建模系統(tǒng)研究與開發(fā)[J]．金屬礦山，2009(1):106-107.

[10] 潘結(jié)南，孟召平，甘莉．礦山三維地質(zhì)建模與可視化研究[J]．煤田地質(zhì)與勘探，2005,33(1):16-18.

[11] 潘玉松．淺談使用ANSYS時(shí)的建模方法[J]．成都電子機(jī)械高等?？茖W(xué)校學(xué)報(bào)，2000(2):28-30.

[12] Woodward C D． Skinning techniques for interactive B-spline surface interpolation[J].Computer-Aided Design,1988,20(8).

[13] 張秀萍．基于點(diǎn)云數(shù)據(jù)復(fù)雜曲面的三維重構(gòu)[D]．新疆大學(xué)，2004.

[14] 白潤才，趙樹賢，劉淑讓．關(guān)于趨勢面分析法在礦巖界面擬合應(yīng)用中的改進(jìn)[J]．阜新礦業(yè)學(xué)院學(xué)報(bào)，1995,14(2):101-103.

[15] 廖秋林，曾錢幫．基于ANSYS平臺(tái)復(fù)雜地質(zhì)體FLAC3D模型的自動(dòng)生成[J]．巖土力學(xué)與工程學(xué)報(bào)，2005,24(6):1010-1013.

[16] 李元松，李新平，姜天華．板桁結(jié)合梁斜拉橋3D-ANSYS建模技術(shù)研究[J]．武漢工程大學(xué)學(xué)報(bào)，2008,30(4):65-68.

[17] 翟偉，孫曉明，蘇麗薇，等．新疆阿希金礦：古生代的低硫型淺成低溫?zé)嵋航鸬V床[J]．地學(xué)前緣(中國地質(zhì)大學(xué))，2010,17(2):266-285.