石 露,李小春,白 冰
(中國科學院武漢巖土力學研究所 巖土力學與工程國家重點實驗室,武漢 430071 )
自然界的巖體經(jīng)歷了多次反復的地質作用,經(jīng)受過變形破壞,是由一定巖性的巖塊和結構面組成的地質體[1]。結構面的存在是巖體的力學性質不同于完整巖石材料的根本原因。在組成巖體的兩大基本單元中,巖塊和結構面所處的地位不同,從強度上來比較,結構面強度遠小于巖塊的強度,為巖體中的薄弱部位,常常是控制巖體強度和變形的關鍵因素。從結構上來講,巖塊的形態(tài)、規(guī)模等特征由結構面所控制,完全取決于結構面的組合、交切與圍限的情況。因此,在巖體的結構中,通常結構面占據(jù)著更為重要的位置[2]。
關鍵塊體理論是由 Goodman和石根華[3]于1982年創(chuàng)立,具有完善的理論體系,包括有限性定理和可移動性定理,目前在分析地下硐室、邊坡和基礎穩(wěn)定性中得到廣泛的應用。塊體理論創(chuàng)立之初即尋找工程在自重條件下的不穩(wěn)定的關鍵塊體,并假定結構面無限大。鑒于工程中,除層面或斷層在研究區(qū)域內可視為無限平面,多數(shù)遍布節(jié)理均是有限長的,由此后續(xù)發(fā)展了基于三維隨機網(wǎng)絡的巖體關鍵塊體搜索[4-6],可對塊體的形態(tài)、大小及空間分布特征等進行統(tǒng)計分析。目前,塊體理論已發(fā)展成全空間內的塊體識別[7-10],不僅可用于關鍵塊體分析,也可用于三維巖體不連續(xù)變形分析。
雖然關鍵塊體理論在巖質邊坡、硐室的穩(wěn)定性評價中起到重要的作用,但目前尚無通用軟件、程序實現(xiàn)困難以及節(jié)理眾多造成的程序運行緩慢,嚴重阻礙該方法的進一步普及。
基于區(qū)域三維隨機節(jié)理網(wǎng)絡的關鍵塊體搜索方法主要有3種:
(1)基于組成自由面上閉合環(huán)路節(jié)理的圍交判斷的關鍵塊體搜索[2];
(2)基于有向性和封閉性原理,根據(jù)節(jié)理上的閉合回路進行全空間內的塊體拓撲搜索[8-9];
(3)基于節(jié)理區(qū)域無限切割以及子塊體拼裝的全空間塊體搜索[10]。
上述方法中,方法(1)完全基于關鍵塊體識別出發(fā),認為關鍵塊體的圍閉節(jié)理必然在自由面上形成閉合環(huán)路,而后根據(jù)組成閉合環(huán)路節(jié)理的交切判斷是否形成塊體,但塊體表面較多時,節(jié)理的空間交切非常復雜,且不能適應復雜的研究區(qū)域;方法(2)和(3)均能夠勝任全空間的塊體識別,且算法本身具有非常好的健壯性,能夠研究復雜的區(qū)域以及實現(xiàn)凹面體的識別,而方法(3)完全由節(jié)理切割凸形子塊體,非常易于編程實現(xiàn)。但全空間上的塊體識別無疑需要龐大的計算量。此外,全空間搜索算法也不易根據(jù)關鍵塊體搜索特點、工程特點對搜索過程進行適當簡化,避免極少數(shù)偶然出現(xiàn)的塊體所帶來巨大的計算量。
本文首先采用方法(1)的閉合環(huán)路搜索,然后根據(jù)方法(3)利用環(huán)路組成節(jié)理及相關后緣節(jié)理進行空間無限切割再拼裝來識別該環(huán)路是否對應關鍵塊體,方法針對性強,能夠很好地適應人為規(guī)定的判別條件,如只搜索楔形體以及有無后緣切割等,且能夠實現(xiàn)凹面體關鍵塊體的搜索,此外也繼承了方法(3)編程實現(xiàn)簡單的特點。
首先對本文的關鍵塊體搜索方法的基本步驟介紹如下:
(1)將研究區(qū)域劃分成凸形區(qū)域組合,如圖 1所示的多平臺邊坡和弧頂型硐室,這樣節(jié)理的切割完全轉化成平面切割凸形面和體的計算幾何問題,避免凹形區(qū)域切割帶來的麻煩。
(2)根據(jù)現(xiàn)場實測統(tǒng)計獲得的結構面產(chǎn)狀、長度與間距的統(tǒng)計特征,采用計算機偽隨機生成三維結構面網(wǎng)絡,并假定結構面為圓盤狀;確定的結構面如 III級以上的結構面可根據(jù)其實際空間形態(tài)添置進去。此外,三維結構面中心在區(qū)域中的生成是采用文獻[2]提出的歸一化方法,可以避免常規(guī)均勻分布帶來的節(jié)理間距不呈指數(shù)分布的情況。
(3)將一個環(huán)路的組成節(jié)理作為無限平面切割研究區(qū)域,判斷該環(huán)路對應的空間多面體的每個面是否真正存在,如不存在再進行相關的后緣節(jié)理切割。如果是非凸的環(huán)路,或者研究區(qū)域非凸還需要進行塊體的拼裝,如圖2所示的一個多平臺邊坡,其中一塊體跨平臺,切割后需要進行拼裝得到最終塊體。
(4)最后,對于第(3)步找到的塊體,計算出體積與滑動方向。對于單面滑動和雙面滑動計算其安全系數(shù)[11],如果安全系數(shù)小于1則此塊體即為關鍵塊體,而對于墜落式塊體則直接判定為關鍵塊體。
具體的流程如圖3所示。接下來將對關鍵塊體實現(xiàn)過程中的閉合環(huán)路搜索和節(jié)理無限切割這兩個關鍵環(huán)節(jié)進行詳細闡述。
圖1 研究區(qū)域的凸形分解Fig.1 Convex decomposition of study region
圖2 非凸塊體的拼裝Fig.2 Assembling of non-convex blocks
圖3 關鍵塊體識別流程Fig.3 Flowchart of key block identification
正如前文所述,任一個可滑動的塊體,肯定對應著自由面上跡線圍成的一個閉合環(huán)路,且閉合環(huán)路的形狀也基本上決定了關鍵塊體的形狀。另外,根據(jù)實際工程情況,如若只考慮重力作用,可根據(jù)環(huán)路最低點以及與該點共棱的塊體內部點的高低位置關系,首先排除一些不能滑出的塊體,也可以非常方便地根據(jù)工程經(jīng)驗對搜索過程作必要的簡化,如工程上經(jīng)常搜索的楔形體對應的環(huán)路為三角形。因此,搜索塊體之前則需要找到自由面上的閉合環(huán)路。一般情況下,對于面積較小的環(huán)路可不做工程考慮。下面則重點介紹閉合環(huán)路的篩選過程。
對于在自由面上可能構成閉合環(huán)路的跡線,必須與兩個或兩個以上的跡線包括自由邊界(如邊坡平臺的眉線)相交,這是一個必要條件,否則該跡線肯定不能構成閉合回路,顯然應被剔除在該閉合環(huán)路的第1類節(jié)理集合中,該過程即為文獻[6]所述的“砍樹”(tree cut)。上述過程是個循環(huán)過程,如圖4所示。圖4(a)為自由面的初始跡線,可以看出與跡線1和2相交的跡線數(shù)目分別為0和1,不具備構成閉合環(huán)路的必要條件,因此,應當被剔除,而當跡線2被剔除后,與跡線3相交的跡線數(shù)目則降為了1條,同樣也需要剔除,依次循環(huán)下去,直至所有的跡線都有兩條或兩條以上相交的跡線而終止,最終如圖圖4(b)所示。以上過程,前提是必須已知各跡線間的相交關系。
圖4 無效跡線剔除示意圖Fig.4 Schematic diagram of rejecting invalid traces
由于節(jié)理數(shù)目眾多,不可能無限制地生成復雜的閉合環(huán)路,在搜索過程中構成閉合環(huán)路的跡線數(shù)目需要根據(jù)工程經(jīng)驗以及實際情況設一個上限,以降低問題的復雜性。事實上,根據(jù)大量的工程實踐,邊坡的關鍵塊體多以四面體的楔形體為主,即三角形的閉合環(huán)路,而洞室圍巖的可動塊體多為四面體、五面體或六面體,有時也會有七面體,但多于七面的可能性較小[2]。
對于邊數(shù)較少的閉合環(huán)路搜索比較容易實現(xiàn),例如與一跡線相交的兩條跡線如果也相交的話,則可形成封閉的三角形環(huán)路,而如果此兩條跡線有一個共同的相交跡線,則可形成一個封閉的四邊形回路。然而,隨著環(huán)路邊數(shù)的增加,以上這種尋找閉合環(huán)路的編程思路則較難實現(xiàn)。為此,需要一個通用的環(huán)路搜索方法。采用邏輯矩陣(元素只能存儲0和1,分別代表不相交與相交)表示跡線相交的信息,是一種圖的存儲結構,如圖5所示的跡線相交圖,其相交信息A見式(1)。
圖5 跡線示意圖Fig.5 Schematic diagram of traces
為了能夠順利進行閉合環(huán)路搜索,還需設計一個堆棧管理數(shù)據(jù)進出,其基本的規(guī)則就是越晚得到的數(shù)據(jù)越早輸出[12],堆棧的目的則是存儲邏輯矩陣的元素的位置信息,該堆棧的兩個基本函數(shù)為壓入(PUSH)和彈出(POP)函數(shù),另外還有重復檢測函數(shù)以防止搜索過程中原路返回的現(xiàn)象發(fā)生。環(huán)路搜索的程序流程見圖6,圖中R表示堆棧中存儲行信息的向量,N為堆棧中數(shù)據(jù)的組數(shù),從廣義上說,此屬于圖遍歷方法中的深度優(yōu)先搜索方法。假定搜索閉合環(huán)路邊數(shù)上限為m,依次從一對角元素(i,i)開始,將該對角元素位置壓入堆棧中形成初始堆棧,而后進行條件循環(huán)。在每次循環(huán)中彈出一個位置數(shù)據(jù) (j0,k0),如果 j0=i且 A(j0,k0)=1那么堆棧中行位置元素加上i即形成一個環(huán)路。接下來,找到該位置同列中的下一個非0元素位置 (j1,k0),將這一數(shù)據(jù)壓入到堆棧中,判斷 j1和堆棧中數(shù)據(jù)的重復性,如果不重復且堆棧個數(shù)據(jù)未數(shù)達到m-1,壓入(i,j1),如果達到m-1,則不壓入(i,j1),但如果A(i,j1)= 1,同樣也形成一個閉合環(huán)路。之后,繼續(xù)進行重復彈棧操作。當棧中元素為空時,壓入對角元素(i+1,i+1),如此直至對角元素循環(huán)結束。
圖5所示的跡線相交圖,經(jīng)過上述的算法可以得到的閉合環(huán)路有[2,4,1]、[2,4,5,3,1]、[2,5,3,1]、[2,5,4,1]、[3,5,2,4,1]、[3,5,4,1]和[4,5,2]這 7 個結果,這其中[2,5,4,1]這種交叉型的互連不能夠成閉合環(huán)路需要對其排除。一般情況下,研究塊體的閉合環(huán)路為凸性。否則造成其安全系數(shù)難以求解,即使如圖2(d)所示的非凸塊體,其自由面上的閉合環(huán)路仍為凸的。因此,根據(jù)環(huán)路的對應交點是否均在每個環(huán)路節(jié)理的同側可排除掉交叉型和非凸性的回路,即排除了[2,5,4,1]和[3,5,2,4,1],從而剩下 5個閉合回路。
圖6 跡線閉合回路搜索流程圖Fig.6 Flowchart of searching closed loop
已知了自由面上的閉合環(huán)路,則以閉合環(huán)路為單位進行關鍵塊體搜索。
如圖2所示,對于一個閉合回路,確定其可能形成的最大塊體,在該過程中必須有規(guī)則的數(shù)據(jù)鏈表和拓撲意義上的數(shù)據(jù)結構。一個完整的多面體結構為一樹狀的結構(如圖7所示),其層次為體-面-環(huán)-棱-點,除點含有坐標表示幾何信息外,其余均為拓撲信息。
圖7 多面體的邏輯表示Fig.7 Logic representation of polyhedron
多面體的存儲采用類半邊數(shù)據(jù)結構[13],但不使用“半邊”層,而是在環(huán)到邊的拓撲信息中采用一個邏輯數(shù)記錄該邊棱的方向(如圖8所示),共有5個雙向的鏈表來依次表示體、面、環(huán)、邊以及點這5層幾何或拓撲信息,對于任一層上的某節(jié)點,均能通過指針訪問其上下層拓撲關系,也可以訪問該層內其前一個與后一個指針,如采用FORTRAN95定義環(huán)(LOOP)結構體如下
TYPE :: loop_type
integer :: id
type(includeE_type),pointer :: includeE
type(face_type),pointer :: belongF_1 !面指針
type(face_type),pointer :: belongF_2 !面指針
type(loop_type),pointer :: prev
type(loop_type),pointer :: next
END TYPE loop_type
這其中有
TYPE :: includeE_type
type(edge_type),pointer :: edge !邊指針
logical :: dir !邊方向
type(includeE_type),pointer :: prev
type(includeE_type),pointer :: next
END TYPE includeE_type
可以看出,一個環(huán)節(jié)點即含有了環(huán)到邊的拓撲信息,亦含有環(huán)屬于面的拓撲信息,其他各層節(jié)點亦是如此。另外,includeF、includeL和includeE分別表示體、面和環(huán)的組成,亦是一個雙向鏈表,includeF節(jié)點包含著指向面的指針和面的外法向信息;includeL節(jié)點則是指向環(huán)的指針,且該鏈表的頭節(jié)點定為“外部”邊界,其余節(jié)點均為面上的“孔”;includeE的節(jié)點如上FORTRAN代碼定義所示,包括指向邊的指針和判別邊方向的邏輯變量。
圖8 多面體的層次結構與數(shù)據(jù)結構Fig.8 Hierarchy and data structure of polyhedron
有了上述完整地表征多面體幾何和拓撲信息的數(shù)據(jù)結構,則很方便進行下面的切割與拼裝的計算幾何操作。
規(guī)定空間內平面的正法向為:當平面法向與 z軸夾角小于5°時,該法向與 z軸正方向夾角為銳角;否則,當其與x軸夾角小于45°,那么該法向與x軸正方向夾角為銳角,如仍不成立,則以與y軸正方向夾角為銳角的方向為正。已知空間平面P方程為
式中:(a,b,c)為該平面的單位正法向向量;d為系數(shù)。對于任一點 v(X,Y,Z),則點與平面的位置關系L為
式中:ε為個小正數(shù)。L=0、1和-1分別表示在面上、面的正側和負側。
若邊、平面以及多面體的所有頂點和平面的位置關系中既存在1亦存在-1,則很容易判斷其與該平面相交。圖9為多面體與節(jié)理的相交示意圖,由式(3)可以判斷出節(jié)理B與C所在的平面肯定與多面體不相交。此外,對于節(jié)理D與多面體是沒有任何關系的,則應避免節(jié)理切割該塊體,但當塊體數(shù)量比較龐大時,如果嚴格地判斷塊體和節(jié)理是否交切的自身計算量亦會非常大,因此,需要較為經(jīng)濟的方法判斷兩者是否呈切割關系。文獻[10]通過節(jié)理與多面體的范圍箱是否有重疊來決定兩者是否進行相交計算。范圍箱指的是一個長、寬、高分別平行于空間坐標軸的長方體,使節(jié)理或多面體恰好落在該長方體內。對于多面體其范圍箱可根據(jù)所有頂點的坐標最大、最小值求出;而對于節(jié)理則需根據(jù)其在各面上的投影求得。有了范圍箱后,只需進行大小判斷即可斷定范圍箱是否有重疊,從而決定是否進行平面與多面體的切割計算。有了該步判斷,圖9的節(jié)理D則不會進行切割計算。
圖9 多面體的與節(jié)理交切關系示意圖Fig.9 Schematic diagram of cutting relationship between polyhedron and joints
下面分步驟介紹平面切割塊體:
(1)平面切割邊。該子程序將返回一個頂點指針。當邊與平面相交時,在頂點鏈表的尾部增加一個節(jié)點,并更新邊到環(huán)以及點到邊的拓撲信息,將返回指針指向該新頂點;若平面恰與邊交與邊的終端點,無需更新拓撲信息,直接返回終端點指針;若以上都不滿足則返回空指針。
(2)平面切割環(huán)。該子程序返回一個邊的指針,當環(huán)的某一邊在平面上時,直接返回該邊的指針。圖10為一平面切割環(huán)E12→E23→E34→E41,對環(huán)所有的邊進行第(1)步操作后,環(huán)變成E15→E52→E23→E36→E64→E41,并根據(jù)第(1)步返回的兩個頂點指針在邊鏈表尾部生成一個新邊 E56(或E65)。接下來,配置兩個環(huán)指針,對環(huán)所有的邊進行循環(huán),并以一邏輯變量為指示,若為真將該邊原封不動地放入環(huán)一,反之則放入環(huán)二,若當前邊的終點為5或6的一點時,將對應E56或E65放入當前環(huán)的邊鏈表中,并對邏輯變量取反,根據(jù)該循環(huán)可得到兩個環(huán) E15→E56→E64→E41和 E52→E23→E36→E65,用兩個新環(huán)中一個覆蓋原環(huán),另一個添加到環(huán)鏈的尾部,最后再將新邊E56作為指針返回。
圖10 平面切割環(huán)Fig.10 Loops cutting by plane
(3)平面切割體。該程序返回一個面的指針。圖 11是一體被平面切割示例,首先是體的每個面(環(huán))進行第(2)步操作后,根據(jù)返回的所有邊指針,生成一個新的環(huán)和面;配置兩個體指針,體一和體二分面位于平面的正負方向,將所有的面按照與平面的位置關系分別放置至兩個體中,并保持外法向不變,再將新面添加到兩個體中,其外法向分別對應著該平面的負法向和正法向;最后將平面正方向的體替換原體,而將負方向上的體放置在體鏈表的尾部。
圖11 平面切割體Fig.11 Solids cutting by plane
有了切割后的子塊體,可以將節(jié)理恢復成圓盤狀,從而判斷子塊體是否有效,如該子塊體含有研究區(qū)域邊界面,則同樣也是無效的。
首先,由形成自由環(huán)路的節(jié)理切割研究區(qū)域,含有閉合環(huán)路的塊體則稱為主塊體。若主塊體是有效的,則無需進行后緣切割;否則,則采用與閉合環(huán)路相交的第2類節(jié)理(與第1類節(jié)理相交的節(jié)理)對主塊體進行切割,通常為了減少復雜程度,不考慮第3類節(jié)理(與第2類節(jié)理相交的節(jié)理)。有了一系列將主塊體切割后的子塊體,則可進行子塊體的拼裝。
文獻[10]給出了通用的合并方法,拼裝的基礎是共面刪除邊、共邊刪除點的拓撲信息改變過程,在此不做過多贅述。主塊體被切割后的子塊體進行拼裝后,若不含有研究區(qū)域邊界,則此塊體則為該閉合環(huán)路對應的塊體。
本節(jié)以某掛幫礦開采的礦柱為例搜索一頂柱的關鍵塊體。
該區(qū)域的節(jié)理層面信息如表1所示,此外區(qū)域的層面產(chǎn)狀為345°∠73°,層面的間距為1.21 m。頂柱的尺寸為15 m×40 m×12 m,頂柱的幾何模型及三維節(jié)理網(wǎng)絡如圖12所示。由于頂柱的塊體基本為墜落式,因此,不需要根據(jù)摩擦角和凝聚力計算其安全系數(shù)。設H為塊體的頂點距自由面的最大距離,根據(jù)上述關鍵塊體搜索過程獲得H>0.2 m關鍵塊體如圖13所示。
該頂柱的關鍵塊體以四面體(106個)為主,五面體(31個)次之,亦有少量的六面體(20個),無七面體,主要是節(jié)理1、2、4和層面相交生成四面體與五面體,節(jié)理3近水平從而對塊體起到后緣切割作用。頂柱的關鍵塊體最大體積為2.739 m3,開采過程可表現(xiàn)出小規(guī)模的塌方。
表1 節(jié)理產(chǎn)狀統(tǒng)計Table1 Statistical parameters of joint occurrence
圖12 頂柱的三維節(jié)理網(wǎng)絡Fig.12 3D joint network in the region of top pillar
圖13 頂柱關鍵塊體的分布圖Fig.13 Distribution of key blocks in top pillar
本文針對關鍵塊體搜索的特點,將研究區(qū)域分解成凸子區(qū)域,搜索自由面上的閉合環(huán)路,然后利用環(huán)路組成節(jié)理及后緣相交節(jié)理進行無限切割和拼裝來識別該環(huán)路是否對應關鍵塊體,得到如下結論:
(1)將自由面節(jié)理跡線相交信息用圖表示,利用堆棧進行圖的遍歷操作,可以適應任意邊數(shù)目的閉合環(huán)路搜索。
(2)采用類半邊數(shù)據(jù)結構存儲空間多面體,除點鏈含有幾何信息外,其余均為拓撲信息,拓撲結構明確,在節(jié)理平面切割邊、環(huán)和面時能夠自由的更新切割后的拓撲關系。
(3)結合自由面上的閉合環(huán)路搜索與節(jié)理平面無限切割的關鍵塊體搜索方法,避免了采用環(huán)路節(jié)理圍閉在節(jié)理較多時復雜的判斷以及全空間塊體搜索塊體龐大的計算量,方法針對性強,較好地適應人為規(guī)定的判別條件,可使關鍵塊體搜索達到效率與效果的平衡。
[1]蔡美峰,何滿潮,劉東燕. 巖石力學與工程[M]. 北京:科學出版社,2002.
[2]賈洪彪,唐輝明,劉佑榮,等. 巖體結構面三維網(wǎng)絡模擬與工程應用[M]. 北京: 科學出版社,2008.
[3]GOODMAN R E,SHI G H. Geology and rock slope stability—Application of the key block concept for rock slopes[C]//Proceedings of 3rdInternational Conference on Stability in Surface Mining. New York: AIME/SME,1982.
[4]SHI G H,GOODMAN R E. The key blocks of unrolled joint traces in developed maps of tunnel walls[J].International Journal for Numerical Analytical Methods in Geomechanics,1989,13(2): 131-158.
[5]鄔愛清,任放,郭玉. 節(jié)理巖體開挖面上塊體隨機分布及錨固方式研究[J]. 長江科學院院報,1991,8(4): 27-34.WU Ai-qing,REN Fang,GUO Yu. A study of block random distribution and anchor design for joint rock excavation[J]. Journal of Yangtze River Scientific Research Institute,1991,8(4): 27-34.
[6]張奇華. 塊體理論的應用基礎研究與軟件開發(fā)[D]. 武漢: 武漢大學,2004.
[7]LU Jun. Systematic identification of polyhedral rock blocks with arbitrary joints and faults[J]. Computers and Geotechnics,2002,29(1): 49-72.
[8]SHI Gen-hua. Producing joint polygons,cutting joint blocks and finding key blocks from general free surfaces[J]. Chinese Journal of Rock Mechanics and Engineering,2006,25(11): 2161-2170.石根華. 一般自由面上的多面節(jié)理生成、節(jié)理塊切割與關鍵塊搜尋方法[J]. 巖石力學與工程學報,2006,25(11): 2161-2170.
[9]張奇華,鄔愛清. 邊坡及洞室?guī)r體的全空間塊體拓撲搜索研究[J]. 巖石力學與工程學報,2008,27(10): 2072-2078.ZHANG Qi-hua,WU Ai-qing. Study of spatial block topological identification of slope and cavern rock mass[J]. Chinese Journal of Rock Mechanics and Engineering,2008,27(10): 2072-2078.
[10]YU Q,OHNISHI Y,XUE G,et al. A generalized procedure to identify three-dimensional rock blocks around complex excavations[J]. International Journal for Numerical and Analytical Methods in Geomechanics,2009,33(3): 355-375.
[11]于青春,薛果夫,陳德基. 裂隙巖體一般塊體理論[M].北京: 中國水利水電出版社,2007.
[12]劉大有,唐海鷹,孫舒楊,等. 數(shù)據(jù)結構[M]. 北京: 高等教育出版社,2001.
[13]伏玉琛,周洞汝. 計算機圖形學——原理、方法與應用[M]. 武漢: 華中科技大學出版社,2003.