用于 FLAC3D力學(xué)分析的礦床塊段模型優(yōu)化建模技術(shù)研究

寇向宇,賈明濤,王李管,吳 霞

(1.中南大學(xué)資源與安全工程學(xué)院,湖南 長沙 410083)

(2.中南大學(xué)數(shù)字礦山研究中心,湖南 長沙 410083)

0 引言

FLAC3D是一種基于連續(xù)介質(zhì)快速拉格朗日差分法的三維空間分析軟件,在采礦與巖土工程中廣泛應(yīng)用。但 FLAC3D的前處理能力不強(qiáng),在建立計(jì)算模型時(shí)主要采用鍵入數(shù)據(jù)/命令行文件方式,雖然可以采用內(nèi)置 Fish語言編程實(shí)現(xiàn)復(fù)雜模型的建立,但其獨(dú)特的源代碼表達(dá)方式、加之實(shí)際采礦工程介質(zhì)和工程形體的復(fù)雜性,由 FLAC3D直接構(gòu)建復(fù)雜模型十分困難[1～2]。

三維礦床建模軟件可以利用地質(zhì)鉆孔、平剖面圖或其他測量數(shù)據(jù),生成復(fù)雜地質(zhì)和工程對象界(表)面的 T IN模型,并在此模型控制下,將其離散為可以逼近形體表面的規(guī)則六面體單元,從而用于對地質(zhì)體屬性(如品位)等的估值。

侯恩科等[3～4]介紹了將礦床建模軟件與FLAC3D軟件進(jìn)行耦合的應(yīng)用實(shí)例,但由于六面體單元在 2類軟件中的應(yīng)用目的不同,關(guān)于如何在塊段建模過程中對模型進(jìn)行優(yōu)化以使之滿足力學(xué)分析要求,在這些文獻(xiàn)中并沒有詳細(xì)討論。

本文從 FLAC3D網(wǎng)格劃分的原則出發(fā),根據(jù)礦床建模軟件 SURPAC構(gòu)建塊段模型的方法,對用于力學(xué)分析的塊段模型構(gòu)建的步驟、原則等進(jìn)行研究,使建立的模型得到優(yōu)化,使其能夠直接作為 FLAC3D分析的前處理模型進(jìn)行模擬計(jì)算。

1 FLAC3D網(wǎng)格劃分原則及 SURPAC塊段模型構(gòu)建方法

1.1 FLAC3D網(wǎng)格劃分原則

有限元網(wǎng)格劃分是進(jìn)行數(shù)值模擬分析至關(guān)重要的一步,它直接影響著后續(xù)數(shù)值計(jì)算分析結(jié)果的精確性。網(wǎng)格劃分涉及單元的形狀及其網(wǎng)格的密度、數(shù)量、質(zhì)量,單元階次,網(wǎng)格分界面和分界點(diǎn),位移協(xié)調(diào)性,網(wǎng)格布局,節(jié)點(diǎn)和單元編號等[5]。

網(wǎng)格數(shù)量的多少將影響計(jì)算結(jié)果的精度和計(jì)算規(guī)模的大小。網(wǎng)格疏密是指在結(jié)構(gòu)不同部位采用大小不同的網(wǎng)格,這是為了適應(yīng)計(jì)算數(shù)據(jù)的分布特點(diǎn)。在計(jì)算數(shù)據(jù)變化梯度較大的部位(如應(yīng)力集中處),為了較好地反映數(shù)據(jù)變化規(guī)律,需要采用比較密集的網(wǎng)格。而在計(jì)算數(shù)據(jù)變化梯度較小的部位,為減小模型規(guī)模,則應(yīng)劃分相對稀疏的網(wǎng)格。選用高階單元可提高計(jì)算精度,因?yàn)楦唠A單元的曲線或曲面邊界能夠更好地逼近結(jié)構(gòu)的曲線和曲面邊界,且高次插值函數(shù)可更高精度地逼近復(fù)雜場函數(shù),所以當(dāng)結(jié)構(gòu)形狀不規(guī)則、應(yīng)力分布或變形很復(fù)雜時(shí)可以選用高階單元。但高階單元的節(jié)點(diǎn)數(shù)較多,在網(wǎng)格數(shù)量相同的情況下由高階單元組成的模型規(guī)模要大得多,因此在使用時(shí)應(yīng)權(quán)衡考慮計(jì)算精度和時(shí)間。網(wǎng)格質(zhì)量好壞將影響計(jì)算精度。結(jié)構(gòu)中的一些特殊界面和特殊點(diǎn)應(yīng)分為網(wǎng)格邊界或節(jié)點(diǎn)以便定義材料特性、物理特性、載荷和位移約束條件。當(dāng)結(jié)構(gòu)形狀對稱時(shí),其網(wǎng)格也應(yīng)劃分對稱網(wǎng)格,以使模型表現(xiàn)出相應(yīng)的對稱特性(如集中質(zhì)矩陣對稱)。不對稱布局會引起一定誤差。節(jié)點(diǎn)和單元的編號影響結(jié)構(gòu)總剛矩陣的帶寬和波前數(shù),因而影響計(jì)算時(shí)間和存儲容量的大小,因此合理的編號有利于提高計(jì)算速度。

1.2 SURPAC塊段模型構(gòu)建方法

塊段模型是一種分割礦體、模擬礦體內(nèi)部特征的有效數(shù)據(jù)模型,就是對目標(biāo)地質(zhì)體進(jìn)行三維柵格化處理[6],把礦體劃分為一系列小長方體單元,其方法是,首先將研究的范圍形成最小包絡(luò)長方體,將其定義為原型；然后根據(jù)地質(zhì)勘探網(wǎng)度、采礦方法、地質(zhì)條件以及地質(zhì)統(tǒng)計(jì)學(xué)等方面對塊段要求來確定單元塊尺寸,以此對原型進(jìn)行三維柵格化；然后對地質(zhì)體進(jìn)行塊段劃分時(shí),應(yīng)選用最合適的、且最具代表性的礦體平剖面圖建立地質(zhì)體表面模型(包括礦體、地層和斷面等模型)對實(shí)體模型進(jìn)行邊界約束,在邊界處進(jìn)一步細(xì)分,以逼近地質(zhì)體的空間形態(tài)；每個(gè)小長方體單元的內(nèi)部屬性相同,且該屬性代表了礦體該位置的內(nèi)部特征,這樣所有小長方體單元的屬性變化規(guī)律就是整個(gè)礦體的特征變化規(guī)律。因此,在建立好塊段模型之后,就可以進(jìn)行相應(yīng)的統(tǒng)計(jì)分析工作了。

2 用于力學(xué)分析的 SURPAC塊段模型優(yōu)化建模技術(shù)

為使 S URPAC建立的塊段模型能夠直接用于FLAC3D力學(xué)分析,保證力學(xué)分析的效率和精度,根據(jù)前文分析,塊段模型建立時(shí)應(yīng)遵從以下步驟和方法。

2.1 模型分區(qū)

根據(jù)網(wǎng)格劃分原則,為了適應(yīng)計(jì)算數(shù)據(jù)的分布特點(diǎn),較好地反映數(shù)據(jù)變化規(guī)律,應(yīng)力可能集中的地方應(yīng)密,其他地方可稀疏一些。網(wǎng)格數(shù)量也要考慮精度和計(jì)算規(guī)模而確定。因此,考慮增加網(wǎng)格的經(jīng)濟(jì)性和分析數(shù)據(jù)的類型,按照不同的劃分方法和力學(xué)參數(shù)建立原始塊段模型,并根據(jù)三維模型進(jìn)行約束,最后,利用編制的接口程序,將塊段模型轉(zhuǎn)換為 FLAC3D分析前處理模型,然后對比不同網(wǎng)格劃分方法的模擬結(jié)果,以確定一個(gè)合適的值。這就是為滿足疏密、密度等的要求,在建模時(shí),可根據(jù)實(shí)際情況分成相關(guān)區(qū)段,外圍區(qū)段單元大、往里逐漸減小,只對最里面有形狀變化的區(qū)段細(xì)分。這樣既保證計(jì)算精度,又保證單元數(shù)量不會太多,且虛擬節(jié)點(diǎn)不會很多,找不全的問題,進(jìn)而提高了分析效率,實(shí)現(xiàn)了網(wǎng)格劃分及其優(yōu)化。

2.2 邊界單元細(xì)分

邊界單元細(xì)分一是為逼近幾何邊界、二是達(dá)到疏密網(wǎng)格的目的,從而反映應(yīng)力集中區(qū)的應(yīng)力變化。在通過 SURPAC軟件輸出網(wǎng)格時(shí),要進(jìn)行塊段模型的建立,在建立的過程中,塊段模型選取地質(zhì)體線框模型進(jìn)行約束,塊段模型的單元塊在邊界處將自動按照一定的準(zhǔn)則發(fā)生細(xì)分(變塊技術(shù))。同時(shí)在建模時(shí)細(xì)分次數(shù)也可超過 2次以上,導(dǎo)入 FLAC3D中使一個(gè)單元跟相鄰單元的邊長之比小于 1∶1或 1∶2,這樣對復(fù)雜形態(tài)邊界的逼近能力更高,效果更好。另外,在細(xì)分后導(dǎo)入 FLAC3D必須利用Attach FAce或 Gen merge命令找出并處理這些虛節(jié)點(diǎn)。因?yàn)楫?dāng)模型較大時(shí),如細(xì)分次數(shù)過多,會導(dǎo)致接觸面虛節(jié)點(diǎn)數(shù)量多,FLAC3D可能找不全這些點(diǎn),導(dǎo)致計(jì)算結(jié)果出錯(cuò)。

2.3 結(jié)果驗(yàn)證

為了驗(yàn)證這種網(wǎng)格劃分實(shí)現(xiàn)前處理模型優(yōu)化方法的正確性,根據(jù)二維彈性力學(xué)的原理[7],設(shè)計(jì)了一個(gè)邊長 48 m×48 m的板中有一個(gè)半徑為 3 m圓孔,在左右兩邊受均布拉力,大小為 1 MPa,然后根據(jù)對模型不同的細(xì)分方法得出疏密不同的單元體網(wǎng)格,進(jìn)行分析,再與理論值作比較,最終得出最優(yōu)的網(wǎng)格劃分方法,進(jìn)而也驗(yàn)證了這種優(yōu)化方法在力學(xué)分析中的重要性和正確性。所以以此方法做出理論應(yīng)力分布圖(11200Zones and 22402 Grid),遠(yuǎn)近疏密一致的應(yīng)力分布圖 (9988zones and 15448 grid points),疏密程度不一致的應(yīng)力分布圖 (9256zones and 18544 grid points),比較應(yīng)力的主要分布,顯示出不同網(wǎng)格劃分對數(shù)值模擬效果的影響,如圖 1～圖 3所示。圖中顯示,理論上的應(yīng)力應(yīng)該集中在中心的兩側(cè),而在不同位置的不同疏密程度模擬出來的效果不同,并不是網(wǎng)格越密效果越好,也不是網(wǎng)格越疏效果越好,只有在應(yīng)力相對應(yīng)該集中的地方密一些,其他地方疏一些,圖 2在兩側(cè)時(shí)應(yīng)力出現(xiàn)了發(fā)散,而圖 3模擬的效果最接近于理論值,效果最好,因此得出了網(wǎng)格劃分的最優(yōu)選擇,圖 4為網(wǎng)格劃分的顯示,進(jìn)而也得到了前處理模型的優(yōu)化。

圖 1 理論應(yīng)力分布圖

圖 2 疏密程度一致應(yīng)力分布圖

圖 3 疏密程度不一致應(yīng)力分布圖

圖 4 最優(yōu)網(wǎng)格劃分效果圖

3 工程應(yīng)用

3.1 工程概況

云南某銅礦經(jīng)過數(shù)十年的露天及地下大規(guī)模工業(yè)開采,井下采空區(qū)越來越多、越來越大,且許多早期、淺部空區(qū)已經(jīng)發(fā)生垮塌,多個(gè)中段上下相鄰或同一水平相鄰的空區(qū)發(fā)生貫通,甚至部分已經(jīng)貫通地表,對地表各種工業(yè)和社會活動,以及生態(tài)環(huán)境造成一定的威脅,存在較大的隱患。如果沒有有效的采空區(qū)穩(wěn)定性及地表移動規(guī)律的評價(jià)研究,則采空區(qū)垮落不僅加大礦石回采的貧化與損失,也必然對礦井及礦區(qū)相關(guān)工農(nóng)業(yè)的安全生產(chǎn)造成嚴(yán)重的威脅。為保證礦山有序的正常安全生產(chǎn)以及后續(xù)開采的安全,對該礦采空區(qū)穩(wěn)定性及由此引發(fā)的地表移動問題進(jìn)行研究具有非常重要的現(xiàn)實(shí)意義。因此,運(yùn)用FLAC3D數(shù)值方法模擬了影響采場穩(wěn)定性的空區(qū)圍巖應(yīng)力及地表位移變化規(guī)律,并在模型建立時(shí)以不同網(wǎng)格劃分方法進(jìn)行對比,實(shí)現(xiàn)前處理模型的優(yōu)化并進(jìn)行數(shù)值模擬,進(jìn)而有效地提高了采空區(qū)穩(wěn)定性分析和評價(jià)結(jié)果的準(zhǔn)確性和可靠性。

3.2 力學(xué)模擬及其分析

在 SURPAC中,在地表三維模型的約束下,建立原始塊段模型。其次,在各中段空區(qū)三維表面模型的約束下,建立相應(yīng)的空區(qū)塊段模型,然后,把所有塊段模型合并,形成一個(gè)統(tǒng)一的模型,該模型內(nèi)單元塊數(shù)目為 638762塊。最后,利用編制的接口程序,將 SURPAC塊段模型轉(zhuǎn)換為 FLAC3D分析前處理模型,該模型節(jié)點(diǎn)數(shù)為 810780個(gè)。

在前處理模型中,將各采空區(qū)所包含的單元設(shè)置為不同的組 (Group),賦力學(xué)參數(shù)、初始地應(yīng)力、邊界條件,就可以在分析過程中,通過指定某個(gè)組的移出,以達(dá)到模擬開挖,考察在采空區(qū)形成后,周邊圍巖應(yīng)力及位移分布情況的目的,最后得出結(jié)果。本文選取礦區(qū)某中段的空區(qū)作為實(shí)例進(jìn)行分析比較,按照相關(guān)力學(xué)參數(shù)進(jìn)行賦值,做出空區(qū)圍巖應(yīng)力分布圖和地表位移分布圖(圖 5、圖 6)。

圖 5 中段開采后最小主應(yīng)力空間分布圖

圖 6 中段開采后地表位移圖

圖 5為第一步開挖后模型內(nèi)的最小主應(yīng)力分布區(qū)及其大小。由圖 5可見,在采空區(qū)周邊圍巖中出現(xiàn)了拉應(yīng)力,多在0.3MPa以下,少數(shù)單元拉應(yīng)力較大,接近但小于 0.8 MPa。由圖 6可以看出,地表位移量整體上不是很大,在 1.5 cm以內(nèi)；最大位移發(fā)生在露天坑坑底中心位置,隨著離坑底中心漸遠(yuǎn),整體位移逐漸減??；然后按采場模擬開挖順序?qū)Ω髦卸芜M(jìn)行開挖,直至結(jié)束,并且分析各中段開挖后模型內(nèi)的最小主應(yīng)力分布區(qū)及其大小,以及采場中各應(yīng)力的變化情況以及開采后的地表位移變化,實(shí)現(xiàn)對該礦區(qū)采空區(qū)穩(wěn)定性的準(zhǔn)確模擬,以制定措施對井下空區(qū)進(jìn)行處理,排除礦區(qū)生產(chǎn)的安全隱患,保障礦區(qū)的正常生活和作業(yè)。

4 結(jié)論

為了實(shí)現(xiàn)數(shù)值模擬的效果更加準(zhǔn)確,研究了有限元網(wǎng)格劃分原則和 SURPAC構(gòu)建塊段模型的方法,然后建立了最優(yōu)的 FLAC3D前處理模型,并構(gòu)建了一個(gè)實(shí)例,依據(jù)彈性力學(xué)的原理,判定其有效性和準(zhǔn)確性,同時(shí)也驗(yàn)證此種前處理模型優(yōu)化方法的正確性,從而實(shí)現(xiàn)了 FLAC3D的前處理模型優(yōu)化。以此在云南某銅礦進(jìn)行設(shè)計(jì)實(shí)踐,通過按照礦區(qū)實(shí)際開采順序進(jìn)行模擬,做出采空區(qū)的應(yīng)力分布和地表的位移移動情況,判斷空區(qū)的穩(wěn)定性和地表生產(chǎn)生活的安全性,并制定一系列措施進(jìn)行治理,排除了礦區(qū)生產(chǎn)的安全隱患,保障了礦區(qū)的正常生活和安全作業(yè),在實(shí)際中驗(yàn)證了此模擬的有效性和精確性,同時(shí)也驗(yàn)證了該前處理模型優(yōu)化方法的可行性和必要性。

[1] 鄧紅衛(wèi),朱和玲,周科平,等.基于 FLAC3D數(shù)值模擬的前后處理優(yōu)化研究 [J].礦業(yè)研究與開發(fā),2008,28 (2)：60-62.

[2] 胡斌,張倬元,黃潤秋,等.FLAC3D前處理程序的開發(fā)及仿真效果檢驗(yàn) [J].巖石力學(xué)與工程學(xué)報(bào),2002, 21(9)：1387-1391.

[3] 劉科偉,李夕兵,宮鳳強(qiáng),等.基于 CALS及 Surpac-FLAC3D耦合技術(shù)的復(fù)雜空區(qū)穩(wěn)定性分析[J].巖石力學(xué)與工程學(xué)報(bào),2008,27(9)：1924-1931.

[4] 侯恩科,吳立新,李建民,等.三維地學(xué)模擬與數(shù)值模擬的耦合方法研究[J].煤炭學(xué)報(bào),2002,27(4)：388 -392.

[5] 杜平安.有限元網(wǎng)格劃分的基本原則[J].機(jī)械設(shè)計(jì)與制造,2000：34-36.

[6] 馮興隆,王李管,畢林,等.基于三維模擬技術(shù)的礦巖可崩性評價(jià)[J].煤炭學(xué)報(bào),2008,33(9)：971-976.

[7] 徐芝綸.彈性力學(xué)簡明教程 [M].人民教育出版社, 1983：200-203.