馬長年,徐國元,,江文武,劉曉明
(1. 中南大學 資源與安全工程學院,長沙 430083;2. 華南理工大學 土木與交通學院,廣州 510641)
伴隨著人類對自然資源、能源需求的不斷增加,各類巖土工程、采礦工程規模越來越大,巖土工程結構和開挖過程越來越復雜,對所涉及結構和過程進行更為貼近實際的三維力學分析、優化結構體力學結構和開挖過程、保證巖土工程的穩定性和安全性有著重要意義。FLAC3D是美國Itasca Consulting Group Inc.開發的三維有限差分力學計算程序,是一種基于拉格朗日算法的數值計算方法。FLAC3D提供了7種材料本構模型,適用于模擬材料的大變形和扭曲轉動,已被廣泛應用于模擬計算巖土類工程地質材料的力學行為[1-5]。
FLAC3D提供了一組命令集和內嵌式解釋語言Fish,用于編寫模擬各類連續體力學過程程序,但對于復雜形態工程幾何體、多屬性、幾十萬至幾百萬單元、大量結構體和復雜開挖過程,由用戶手工編寫程序代碼對上述過程進行三維力學分析存在諸多困難,主要表現在:(1)模型的建立只能用數據文件來實現,不能直觀顯示,建模過程可視化程度低;(2)對于復雜的大型工程地質體建模,需要控制各個邊界節點的坐標數據,容易出錯,檢查也不方便[3-5];(3)對于大量單元、結構體對其特定屬性進行檢索、查詢、分類和賦值等操作困難;(4)手工編制上述大規模FLAC3D程序效率低、耗時多、出錯概率大、查錯糾錯困難,常導致程序運行不收斂、結果失真等異常。
本文所述 FLAC3D力學仿真代碼生成系統是對解決此類問題進行的研究與嘗試。
任何復雜開挖過程都是由一些具體的、特定的、可劃分成模塊或階段的過程構成的,通過每一個模塊或階段,開挖過程模擬可實現整個開挖過程的力學仿真。以下按4種不同方法對復雜開挖過程進行剖分。
(1)按照坐標(X、Y、Z及其組合)、平面(或曲面)相對位置(上下/前后/左右)、邊界(內/外)對模型進行剖分;
(2)按一定的開挖順序(串行或并行)對一組或多組單元進行剖分;
(3)按不同的結構功能可對開挖過程進行多種方式的剖分,如隧洞的開挖過程可按結構功能剖分為:上拱、側墻、底拱、核部、錨桿和襯砌等結構;
(4)按不同的巖性、開挖狀態(開挖/未開挖)、充填狀態、礦體/巖石、原巖/破碎帶等單元屬性進行剖分。
(1)按礦體/巖石屬性剖分模型
[4-6]的方法,針對某采用下向進路膠結充填采礦法礦山,本文利用Surpac軟件為該礦礦體和地表地形按照模型整體要求的范圍分別建立3 D和DTM幾何模型[6-7],所建模型如圖1所示。

圖1 礦山三維幾何模型Fig.1 3D model of mine
在此基礎上建立模型的Surpac塊體模型,擴展其屬性得到結構如表1所示的包含單元質心坐標數據的塊體單元數據庫。

表1 塊體單元數據庫基本結構Table1 The basic construction of zone table
可通過礦體塊體單元約束條件向塊體單元數據庫巖性屬性賦值,礦體內賦1,礦體外賦0。礦體塊體單元約束中的塊體尺寸必須與FLAC3D模型的單元尺寸保持一致。由Surpac塊體模型質心坐標數據可生成六面體單元的FLAC3D模型[4]。
(2)按開挖順序剖分模型
參照礦山實際生產情況,按礦體高程以 20 m劃分分段;在分段內沿礦體走向按整數勘探線每100 m劃分為一個盤區;進路規格與分段高保持一致,進路布置有沿走向和垂直走向2種方式。
(3)按結構體剖分模型
進路中布設的鋼筋(吊筋和水平筋)可由FLAC3D所提供的錨索結構體模擬,系統默認所有回采進路都布設鋼筋,所有礦體內單元和部分圍巖中單元將與一條或多條錨索結構體關聯。這樣模型單元被劃分為布設錨索和不布設錨索2類。
(1)對模型單元直接操作的意義
大型 FLAC3D力學模型通常是由幾萬至幾百萬甚至更多的單元組成的,由于地質工程體的原始應力狀態、材料屬性、模型幾何形態、開挖體幾何形態、開挖順序和步驟、工程結構體等的多種多樣,每一個單元都有可能具有自身獨特的屬性。對任一單元進行直接操作,主要涉及對特定單元進行定位、捕捉和重組的問題,如此才可能對各類復雜開挖過程進行有效的力學模擬與仿真。
(2)FLAC3D操作模型單元的方法[2]
FLAC3D對模型單元定位的方法歸納起來有 2種:RANGE命令和單元操作函數。
RANGE命令有一組子關鍵字,可按環、柱、組、方向、平面相對位置、模型類型、坐標、單元號等多種方式確定單元位置或單元范圍。其中按坐標和單元號操作單元是最為基礎的。
單元操作函數有:

(3)基于單元數據庫的模型單元操作方法
將模型所有單元的質心坐標、巖性、開挖標志、單元上施加的結構體等屬性(還可以增加其他屬性)通過Surpac形成塊體單元數據庫,其單元數據表基本結構如表1所示。
通過數據庫查詢語言可對單元一個或多個屬性進行組合檢索,將單元質心坐標傳遞給FLAC3D的單元捕捉函數,如gp_near(x,y,z),z_near(x,y,z)等,再利用函數z_id(p_z)得到單元在模型中的單元號。這樣可利用單元號對需要相同操作的單元通過重組進行成批操作和處理。單元的重組和批處理操作可利用FLAC3D中數組、函數和循環語句完成。
3.2.1 Fish語言和FLAC3D命令特點[2]
(1)Fish語言是內嵌于FLAC3D內部的語言,它可以讓服務用戶自己定義變量和函數;
(2)Fish共有4種數據類型:整型、浮點、字符串和指針;
(3)Fish函數由一對關鍵字Define與End定義;
(4)在Fish語言中所有的變量和函數都是全局的(與Basic相似),函數沒有參數;
(5)一條完整的Fish語句占用一行,無法續行。如果公式太長,需要用臨時變量將其分開;
(6)Fish包含條件分支、條件選擇、循環語句等程序控制語句;
(7)Fish預定義了一些函數和變量,使用此類變量和相關函數,可獲得所有節點、單元或全局空間坐標點在內存中的地址,并對其進行操作;
(8)函數可重定義,重定義后,原先的函數將不能再被調用,但原函數所形成的內存空間仍存在,其中的內存變量不變,數組不可重定義。
以上這些Fish語言和FLAC3D命令特點和功能是開發代碼生成系統所必須了解和注意的。
3.2.2 字典變量與FLAC3D關鍵字映射
為實現由 FLAC3D代碼生成系統自動生成FLAC3D程序,首先要建立應用程序中變量與FLAC3D語句關鍵字的映射關系,即翻譯字典。應用程序由 VC++開發,在此將專門用于處理這一映射關系的字符串變量稱為字典變量,數據類型為CString,字典變量映射格式如表2所示。這樣在應用程序中建立起字典變量與所有FLAC3D關鍵字之間一一對應的關系,這是代碼生成系統自動生成FLAC3D應用程序代碼的基本前提。

表2 字典變量與FLAC3D關鍵字映射表Table2 The mapping table between dictionary variables and FLAC3Dkeywords
3.2.3 傳遞程序變量值到FLAC3D變量、數組和函數
下面以布設2條錨索的程序代碼為例來說明如何將程序變量值傳遞到FLAC3D變量、數組和函數,并形成2條錨索布設的FLAC3D程序代碼。
VC++源程序代碼片段




至此,我們建立了代碼生成系統與FLAC3D關鍵字、變量、數組、函數和控制語句之間的對應關系,這是應用程序生成FLAC3D程序代碼的基礎。
(1)仿真代碼生成系統的輸入和輸出
仿真代碼生成系統的輸入是由 Surpac塊體模型導出的包含單元質心數據、單元屬性及單元約束狀態的Access數據庫。仿真代碼生成系統的輸出是按用戶的要求,由系統生成的FLAC3D程序代碼。
(2)仿真代碼生成系統的工作原理
仿真代碼生成系統程序通過 Ado訪問 Surpac塊體模型導出的塊體單元數據庫,利用數據庫所提供的函數可以方便地對其中的數據進行多種類型的檢索和操作。同時程序提供一組交互式表單與用戶進行交互,由用戶對一類特定巖體回采開挖過程所需步驟和參數進行選擇和確認。系統接收這些步驟和參數后生成*.txt格式的FLAC3D仿真程序代碼,同時修改塊體單元數據庫中相應單元屬性值。
(1)系統結構
①單元數據由數據庫結構組織;②應用程序由單文檔應用程序框架來實現;③屬性表單向導結構:應用程序將采用向導方式按生產過程劃分的步驟,一步一步引導用戶完成滿足要求的FLAC3D代碼。每一步采用屬性表單的方式,為用戶提供多個選項,以確定 FLAC3D力學仿真程序的各關鍵屬性和參數。
(3)功能設計
本文所述 FLAC3D仿真代碼生成系統由 VC++開發,系統具體功能是為模擬某采用下向充填法礦山的實際生產過程力學仿真而設計的,主要功能如圖2所示。

圖2 仿真代碼生成系統功能設計Fig.2 The function design of system of generating code
各頁面主要完成以下功能:
①主頁面為一單文檔應用程序,主要功能包含在主菜單“仿真代碼生成”菜單中,如圖 2(a)所示。打開“仿真代碼生成菜單”,程序將引導用戶逐步完成單元數據庫調用、分段水平劃分、盤區劃分及進路回采的 FLAC3D力學仿真代碼,并以 FLAC3D可識別的*.txt文件格式輸出到用戶指定位置。
②點擊“礦體分段劃分”按鈕,程序將要求用戶調入礦體單元數據庫。通過訪問塊體單元數據庫,程序可按分段高度從數據庫中提取礦山開采分段的高程數據,本步驟是在垂直方向對礦體進行劃分,用戶可選擇將開采分段,界面功能設計如圖2(b),該步驟同時向程序變量提供單元的Z坐標。
③本步驟是對上一步所選擇的開采分段進行盤區劃分,通常只對礦體部分進行劃分。盤區寬度默認為100 m,盤區編號從西到東。盤區是按生產實際情況劃分的,盤區劃分有利于實現不同開采順序的力學模擬。界面功能設計如圖2(c),該步驟向程序變量提供盤區內單元X、Y坐標。
④由用戶確定盤區內進路布置的方向,選擇待回采進路,確定充填進路內是否布設鋼筋等參數(在本實例中鋼筋的力學參數以現場實際使用值為默認值)。本步驟主要用來實現生產實際中上下分段(分層)進路交錯布置的仿真,在一個盤區內可模擬隔一采一、隔二采一等回采方式。界面功能設計如圖2(d),該步驟向程序變量傳遞待采進路所有單元質心坐標及鋼筋(錨索)的起點、終點單元的質心坐標,同時將單元數據庫中相應單元的開挖標志和充填標志置為1。
⑤用戶可指定仿真代碼的存儲位置,每個盤區內的進路可按一定回采順序生成一個*.txt文件。完成后程序將進入下一個選單,點擊回采新盤區可對新的盤區進行回采,點擊回采新分段可對新的分段進行回采,直到按用戶的方案回采完畢,點擊完成按鈕退出代碼生成器,其界面功能設計如圖 2(e)、(f)。
(1)基本類
單文檔 MFC應用程序都具有的 4個類:CCodeGenApp、CCodeGenDoc、CCodeGenView、CMainFrame。
(2)表單類
該類為程序功能的實現提供與用戶交互的界面(見表 3),用戶通過這些界面的選項將自己的意圖一步步傳遞給程序,程序從用戶獲得所需的參數,最終由程序自動地形成FLAC3D力學仿真代碼。
(3)功能類
CDBP類:該類是用來訪問數據庫和進行數據處理的類,CDBP通過接收 CProPageLevel、CProPageSection、CProPageMining、CProPageReturn和CProSheet等類傳遞來的消息和參數,對Surpac塊體單元數據庫進行訪問和檢索。

表3 類的種類與功能Table3 The type and function of class
應用本系統成功地生成了下向進路式充填采礦法礦山生產過程的FLAC3D力學仿真代碼,可實現垂直方向、水平方向的開采順序,不同進路布置方式,充填體內鋼筋布置等主要生產過程的力學仿真。程序中所采用的有關礦山原巖應力場、礦巖力學參數、開采工藝[8-11]基本與礦山實際保持一致。模型共有139626個單元,錨索結構體1256個,是1個較為復雜的三維力學模型,所生成程序代碼調入FLAC3D運行正常。礦體上部5個分段回采充填后,礦體附近第1主應力分布和鋼筋(錨索)受力狀態仿真結果見圖3,圖3(b)中錨索淺色表示受壓,深色表示受拉,粗細表示受力大小。

圖3 部分仿真結果圖Fig.3 The part of simulating results
(1)FLAC3D代碼生成系統從分析 FLAC3D命令集和Fish語言特點入手,著重解決了模型單元質心定位、單元捕捉和重組,應用程序變量到FLAC3D變量、函數、語句、數據結構之間映射等問題,為編制FLAC3D程序提供了一種新的方法,能夠大幅提高編制復雜開挖過程FLAC3D程序代碼的效率,降低代碼編制出錯率。
(2)GCS系統可對大量模型單元進行直接和有效的操作,具有良好的適應性和可擴展性。
(3)在建立塊體單元數據庫過程中利用了礦用地質軟件Surpac的許多功能,架起了地質三維幾何造型軟件與三維力學計算軟件之間的橋梁。
參考文獻
[1]謝和平,陳忠輝. 巖石力學[M]. 北京: 科學出版社,2004.
[2]Itasca Consulting Group. Fast Lagrangian analysis of continua in three dimensions[M]. America: Itasca Consulting Group.,1997.
[3]廖秋林,曾錢幫. 基于 ANSYS平臺復雜地質體FLAC3D模型的自動生成[J]. 巖石力學與工程學報,2005,24(6): 1010-1013.LIAO Qiu-lin,ZENG Qian-bang. Automatic model generation of complex geologic body with FLAC3Dbased on ANSYS platform[J]. Chinese Journal of Rock Mechanics and Engineering,2005,24(6): 1010-1013.
[4]羅周全,吳亞斌,劉曉明,等. 基于Surpac的復雜地質體 FLAC3D模型生成技術[J]. 巖土力學,2008,29(5):1334-1338.LUO Zhou-quan,WU Ya-bin,LIU Xiao-ming,et al.FLAC3Dmodeling for complex geologic body based on Surpac[J]. Rock and Soil Mechanics,2008,29(5): 1334-1338.
[5]何忠明,彭振斌,曹平,等. 雙層空區開挖頂板穩定性的 FLAC3D數值分析[J]. 中南大學學報(自然科學版),2009,40(4): 1066-1071.HE Zhong-ming,PENG Zhen-bin,CAO Ping,et a1.Numerical analysis for roof stability of double gob area after excavation by FLAC3D[J]. Journal of Central South University (Science and Technology),2009,40(4): 1066-1071.
[6]Surpac Software International國際軟件公司.SURPAC Vision軟件用戶使用手冊(第四版)[M]. [s. l.]: SURPAC Software International國際軟件公司,2000.
[7]劉曉明,羅周全,楊 彪,等. 復雜礦區三維地質可視化及數值模型構建[J]. 巖土力學,2010,31(12): 4006-4011.LIU Xiao-ming,LUO Zhou-quan,YANG Biao,et al.Numerical modeling and geological body visualization for complex mine [J]. Rock and Soil Mechanics,2010,31(12): 4006-4011.
[8]路德維格,宋恕夏,田永妥.中國-瑞典關于金川二礦區巖石力學研究報告[R]. 北京: 北京科技大學,1998.
[9]陳俊彥,N G W,庫克,等. 金川鎳礦二礦區采礦方法的巖石力學研究報告[R]. 長沙: 中南礦冶學院,1986.
[10]廖椿庭,施兆賢. 金川礦區原巖應力實測及在礦山設計中的應用[J]. 巖石力學與工程學報,1983,2(1): 103-112.LIAO Chun-ting,SHI Zhao-xian. In-situ stress measurements and their application to engineering design in the Jinchuan mine[J]. Chinese Journal of Rock Mechanics and Engineering,1983,2(1): 103-112.
[11]吳滿路,馬宇,瘳椿庭,等. 金川二礦深部1000 m中段地應力測量及應力狀態研究[J]. 巖石力學與工程學報,2008,27(增刊2): 3785-3790.WU Man-lu,MA Yu,LIAO Chun-ting,et al. Study of recent state of stress in depth 1000 m of Jinchuan mine[J].Chinese Journal of Rock Mechanics and Engineering,2008,27(Supp.2): 3785-3790.