基于損傷理論的巖石細(xì)觀尺度微裂隙量化分析

基于損傷理論的巖石細(xì)觀尺度微裂隙量化分析

1巖石細(xì)觀結(jié)構(gòu)巖石是由各種礦物結(jié)晶、空隙和結(jié)構(gòu)物組成的綜合體。經(jīng)過億萬年的地質(zhì)演變和多次復(fù)雜的構(gòu)造運(yùn)動(dòng)，使巖石含有不同階次的隨機(jī)分布的微觀孔隙和裂紋。巖石受載后的宏觀斷裂失穩(wěn)和破壞與其變形時(shí)內(nèi)部微裂紋的分布以及微裂紋的產(chǎn)生、擴(kuò)展和貫通密切相關(guān)。巖石是一種多尺度的材料，M.E.Kassner等認(rèn)為只有建立多尺度的力學(xué)模型才能真正了解微裂隙起裂的機(jī)制，并應(yīng)在細(xì)觀尺度上對(duì)脆延性材料作進(jìn)一步研究。因此，本文從巖石的掃描電鏡(SEM)圖像研究入手，獲得了巖石微裂隙的細(xì)觀結(jié)構(gòu)信息，并從細(xì)觀角度出發(fā)，認(rèn)為巖石在細(xì)觀上同樣為一結(jié)構(gòu)體；同時(shí)，即為巖石細(xì)觀結(jié)構(gòu)體引入巖石細(xì)觀結(jié)構(gòu)體量化信息運(yùn)用于巖體損傷理論中，基于G.Swoboda損傷理論得到巖石細(xì)觀結(jié)構(gòu)體的損傷本構(gòu)模型，從而將量化巖石細(xì)觀結(jié)構(gòu)特征信息運(yùn)用于損傷理論之中。2巖石開裂程度表征巖石細(xì)觀結(jié)構(gòu)(rockmesostructure)量化試驗(yàn)研究就是研究巖石體內(nèi)的微裂隙和微孔洞的三維形狀參數(shù)?；谖⒘严兜拈L度、方位角、間距和面積等形狀參數(shù)，得到巖石的初始細(xì)觀結(jié)構(gòu)損傷張量，從而將細(xì)觀結(jié)構(gòu)損傷張量引入到損傷理論得到巖石的本構(gòu)模型。K.I.Kanatani于1984年通過對(duì)微結(jié)構(gòu)體的參數(shù)——微裂隙密度來表征巖石的開裂程度。D.F.Howarth和J.C.Rowlands于1987年通過對(duì)巖石切片圖像分析，建立了量化巖石微細(xì)顆粒參數(shù)(形狀、方位角和相對(duì)比例等)的方法，并與巖樣的力學(xué)強(qiáng)度進(jìn)行了比較，得到兩者的相互關(guān)系。D.H.Campbell和J.S.Galehouse于1991年討論了利用光學(xué)顯微鏡進(jìn)行巖石、磚和水泥等材料微觀量化的可行性，并且與X射線方法進(jìn)行了比較和研究。R.H.CWong等于1996年研究了初始微裂隙密度和巖石微顆粒尺寸對(duì)于YuenLong大理石巖樣單軸壓縮強(qiáng)度的影響。Y.H.Hatzor等于1997年研究了白云石的細(xì)觀結(jié)構(gòu)與微裂隙起裂的初始應(yīng)力和試樣最終強(qiáng)度之間的關(guān)系，發(fā)現(xiàn)巖石微結(jié)構(gòu)對(duì)其強(qiáng)度極限的影響非常大。X.Y.Wu等于2000年利用光學(xué)顯微鏡和SEM研究了DarleyDale砂巖在壓縮破壞中的各向異性損傷的微觀力學(xué)演化過程，得到了微裂隙密度與應(yīng)變之間的關(guān)系。B.Menendez等于2001年探討了運(yùn)用激光共焦點(diǎn)掃描顯微鏡來研究巖石中微裂隙和微孔洞網(wǎng)絡(luò)的方法。2.1試樣制備和檢測巖石試樣采自四川錦屏II級(jí)水電站工程現(xiàn)場的大理巖，白色、塊狀、粒狀變晶結(jié)構(gòu)，均勻性好。為了保證試樣性質(zhì)的均一性，嚴(yán)格密集采樣，巖樣嚴(yán)格按照國際巖石力學(xué)學(xué)會(huì)(ISRM)試驗(yàn)規(guī)程來加工，加工成3組尺寸為φ50mm×100mm的圓柱體標(biāo)準(zhǔn)試件。其中2組進(jìn)行單軸壓縮試驗(yàn)和常規(guī)三軸壓縮試驗(yàn)，另外1組按照以下介紹的方法進(jìn)一步加工成3個(gè)切片進(jìn)行巖石細(xì)觀結(jié)構(gòu)量化試驗(yàn)。試驗(yàn)切片的制作流程見圖1，具體試樣加工步驟如下：(1)首先從加工好的大理巖試樣中選取1塊，另外2塊將用來做單軸壓縮和常規(guī)三軸壓縮試驗(yàn)。(2)從試樣頂部切割一高大約為10mm左右的圓柱，見圖1步驟aㄢ(3)在得到小部分圓柱上，切出一個(gè)大致尺寸為30mm×30mm×10mm的切片，其為XY切片(垂直于加載方向)，見圖1步驟b,d,eㄢ(4)在剩余的大部分圓柱上，豎向切出一個(gè)足夠大小的部分，切成大致尺寸為30mm×30mm×10mm，為YZ切片，見圖1步驟c,f,h,iㄢ(5)從剩余的部分切出一個(gè)垂直于YZ切片的XZ切片，見圖1步驟g,jㄢ最后對(duì)切片進(jìn)行相應(yīng)標(biāo)號(hào)，并在切片的左上角標(biāo)記，以確保三者之間的空間關(guān)系，即觀測到恰好是試樣內(nèi)的一個(gè)結(jié)構(gòu)立方體的3個(gè)面。在進(jìn)行掃描電鏡圖像采集前，需要對(duì)切片進(jìn)行清洗和干燥以及鍍金，采用的步驟為：(1)在乙醇中超聲波處理5min以清除表面雜質(zhì)。(2)在100?C下干燥12hㄢ(3)進(jìn)行鍍金。加工好的大理石切片見圖2ㄢ2.2圖像放大個(gè)數(shù)采用HitachiS–3000h型掃描電鏡對(duì)相應(yīng)的切片進(jìn)行觀測，采集到的圖像的分辨率為640像素×480像素。圖3為一張放大300倍的大理巖細(xì)觀結(jié)構(gòu)SEM圖像，圖中一個(gè)像素的距離代表0.6667μm。圖像的放大倍數(shù)，決定了一張圖像上能夠顯示的微裂隙數(shù)量。如果放大倍數(shù)太大，則可能出現(xiàn)一張圖像上只能顯示微裂隙的一部分這種情況；相反，如果放大倍數(shù)太小，則會(huì)由于微裂隙在圖像上顯示的太小，在圖像處理過程中，丟失圖像信息，造成誤差過大，所以在參考F.E.Mohamed研究成果的基礎(chǔ)上，作者進(jìn)行相應(yīng)的嘗試，最后決定采用的放大倍數(shù)為300倍。圖像的采集數(shù)量決定了是否能獲取到足夠多的微裂隙信息，進(jìn)而形成微裂隙的信息統(tǒng)計(jì)分布。在試驗(yàn)中，每個(gè)面采集最少200張以上圖片。在切片的表面上獲取圖片，從上而下，避免圖像的重復(fù)獲取。3巖石細(xì)觀結(jié)構(gòu)的量化分析近年來，隨著計(jì)算機(jī)技術(shù)和數(shù)字圖像處理的發(fā)展使得通過分析巖石切片細(xì)觀圖像得到巖石細(xì)觀結(jié)構(gòu)參數(shù)量化信息成為可能。同時(shí)數(shù)字圖像處理技術(shù)的發(fā)展使得從切片微觀圖像上獲取數(shù)據(jù)收集更方便、更準(zhǔn)確、統(tǒng)計(jì)上更有代表性。3.1區(qū)域生長——局部分割過程圖像分割是將圖像細(xì)分為構(gòu)成它的子區(qū)域或?qū)ο?。分割的程度取決于要解決的問題以及如何來識(shí)別SEM圖像中的微裂隙部分。圖像分割算法的一般是基于亮度值的兩個(gè)基本特性之一：不連續(xù)性和相似性。第1個(gè)性質(zhì)的應(yīng)用途徑是基于亮度的不連續(xù)變化分割圖像，比如圖像的邊緣；第2個(gè)的主要應(yīng)用途徑是依據(jù)事先制訂的準(zhǔn)則將圖像分割為相似的區(qū)域，如門限處理、區(qū)域生長、區(qū)域分離和聚合都是這類方法。由于巖石細(xì)觀SEM圖像的特征，決定了采用諸如邊緣檢測算子如Robison,Roberts,Sobel,Prewitt,Krisch和Log等，無法很好地達(dá)到從圖像中分割出微裂隙的目的，所以采用后者即采用相似性的原則運(yùn)用區(qū)域生長來實(shí)現(xiàn)圖像分割的目的。區(qū)域生長是一種事先根據(jù)定義的準(zhǔn)則將像素或子區(qū)域聚合成更大區(qū)域的過程?；痉椒ㄊ且砸唤M“種子”點(diǎn)開始將與種子性質(zhì)相近(諸如灰度級(jí)或顏色的特定范圍)的相鄰像素附加到生長區(qū)域的每個(gè)種子上。通常根據(jù)所解決的問題的性質(zhì)選擇一個(gè)或多個(gè)起點(diǎn)。當(dāng)一個(gè)先驗(yàn)信息無效時(shí)，這一過程將對(duì)每個(gè)像素計(jì)算相同的特性集，最終這個(gè)特性集在生長過程中將像素歸屬某個(gè)區(qū)域。如果這些計(jì)算結(jié)果呈現(xiàn)不同簇的值，則那些由于自身性質(zhì)而處于這些簇中心附近的像素可以作為種子。區(qū)域生長的另一個(gè)問題是用公式描繪一個(gè)終止規(guī)則。基本上，在沒有像素滿足加入某個(gè)區(qū)域條件時(shí)，區(qū)域生長就會(huì)停止。因此，區(qū)域生長的2個(gè)基本原則為：(1)任何像素和種子之間的灰度級(jí)絕對(duì)差必須小于某一個(gè)值β，這個(gè)數(shù)值是根據(jù)某一算法得到的，足以區(qū)別圖上2個(gè)不同的區(qū)域，如果β超過了其合理大小的值，則分割出來的微裂隙將偏大；反之，則分割出來的微裂隙將偏小。(2)要添加入某一個(gè)區(qū)域的像素必須與此區(qū)域至少一個(gè)像素是八連通的。如果某個(gè)像素被發(fā)現(xiàn)于多于一個(gè)區(qū)域相聯(lián)系，就將這些區(qū)域合并起來。形象地說，生長的過程中必須有“根”。而β值可以通過多種算法獲取，如雙峰法、迭代方法(最優(yōu)方法)、大津法(OTSU法)、灰度拉伸法(改進(jìn)的OTSU法)以及Krish算子法，具體參見劉爽的研究成果。3.2人工智能以及模糊聚類獲取種子點(diǎn)，算法充選取種子點(diǎn)是通過在微裂隙區(qū)域人工點(diǎn)擊一定種子點(diǎn)，當(dāng)然在今后可以通過其他一些算法結(jié)合人工智能以及模糊聚類也可以實(shí)現(xiàn)智能獲取種子點(diǎn)，本文不作考慮；圖4為設(shè)置的區(qū)域生長算法的種子點(diǎn)。(2)區(qū)域生長運(yùn)算符使用MATLAB進(jìn)行編程運(yùn)算(見圖5)，得到微裂隙二值化后的圖像(見圖5(a))，5(b)為區(qū)域生長算法識(shí)別結(jié)果。(3)邊緣有小球質(zhì)性識(shí)別出來的微裂隙圖像的圖像中微裂隙(圖6中白色部分)的邊緣存在不光滑以及有一些毛刺，在微裂隙中可能有一些小的孔洞，所以有必要采用形態(tài)學(xué)的閉操作以及孔洞填充進(jìn)行處理，處理后的結(jié)果見圖7ㄢ3.3值化圖像的識(shí)別在識(shí)別出來的二值化微裂隙圖像(見圖7)上得到微裂隙細(xì)觀結(jié)構(gòu)參數(shù)，MATLAB中有Bwlabel和Regionprops兩個(gè)函數(shù)可以對(duì)二值化圖像后顆粒或孔隙進(jìn)行標(biāo)識(shí)和統(tǒng)計(jì)各參數(shù)值。可以獲取的參數(shù)值為顆?；蚩紫兜目偯娣e、孔隙率、總周長、面積、周長、形狀系數(shù)、粒徑或孔徑、扁圓度和方位角。圖中一個(gè)像素代表的長度為0.6667μm，從計(jì)算機(jī)識(shí)別出來的各量像素單位再轉(zhuǎn)化為微米單位。從圖7識(shí)別出來的微裂隙細(xì)觀結(jié)構(gòu)參數(shù)如表1所示。微裂隙的平面間距的表達(dá)式為式中：L為切片上沿某一方向剖線的長度，NLC為長為L的剖線通過的微裂隙數(shù)量。4巖石薄視結(jié)構(gòu)信息分析對(duì)所有的3個(gè)切片上細(xì)觀結(jié)構(gòu)圖片進(jìn)行了處理，得到了大量的微裂隙的細(xì)觀結(jié)構(gòu)參數(shù)數(shù)據(jù)，對(duì)這些數(shù)據(jù)進(jìn)行相應(yīng)的統(tǒng)計(jì)分析。4.1密度函數(shù)的估計(jì)圖8為XY面微裂隙方位角直方圖和理論概率分布函數(shù)，從圖8中可以看出，微裂隙的方位角服從正態(tài)分布，其正態(tài)分布的密度函數(shù)為采用點(diǎn)估計(jì)的方法算出μ和σ的無偏估計(jì)量的估計(jì)值為得到XY,YZ,XZ面的微裂隙方位角的統(tǒng)計(jì)參數(shù)的無偏估計(jì)量見表2。將表2中相應(yīng)數(shù)值代入式(2)得到理論概率分布函數(shù)曲線與實(shí)際數(shù)據(jù)直方圖的比較(見圖8～10)，其分別為XY,YZ,XZ面的微裂隙方位角數(shù)據(jù)的直方圖和代入統(tǒng)計(jì)參數(shù)的無偏估計(jì)量的理論概率分布函數(shù)曲線。4.2密度函數(shù)的統(tǒng)計(jì)描述圖11為XY面微裂隙長度直方圖和理論分布密度函數(shù)曲線。由圖11可看出，微裂隙長度分布基本是服從對(duì)數(shù)正態(tài)分布，其對(duì)數(shù)正態(tài)分布的密度函數(shù)為同樣得到XY,YZ,XZ面的微裂隙長度的統(tǒng)計(jì)參數(shù)的無偏估計(jì)量見表2。同樣將表2中相應(yīng)數(shù)值代入式(4)得到理論概率分布函數(shù)曲線與實(shí)際數(shù)據(jù)直方圖的比較(見圖11～13)，其分別為XY,YZ,XZ面的微裂隙長度數(shù)據(jù)的直方圖和理論概率分布函數(shù)曲線。4.3微裂隙間距的統(tǒng)計(jì)參數(shù)由圖14可以看出，XY面的微裂隙間距分布為指數(shù)分布，其密度函數(shù)為同樣得到XY,YZ,XZ面的微裂隙間距的統(tǒng)計(jì)參數(shù)的無偏估計(jì)量見表2。將表2中相應(yīng)數(shù)值代入式(5)得到理論概率分布函數(shù)曲線與實(shí)際數(shù)據(jù)直方圖的比較(見圖14～16)，其分別為XY,YZ,XZ面的微裂隙間距數(shù)據(jù)的直方圖和理論概率分布函數(shù)曲線。5巖石細(xì)觀結(jié)構(gòu)體的模擬微裂隙的方位角服從正態(tài)分布，長度服從對(duì)數(shù)正態(tài)分布，間距服從指數(shù)分布。模擬巖石細(xì)觀結(jié)構(gòu)體是通過模擬符合以上3種統(tǒng)計(jì)分布的微裂隙來實(shí)現(xiàn)的，程序采用MATLAB編制，MATLAB中的隨機(jī)變量模擬是采用逆變換方法的得到，即所有的分布都是從均勻分布中得到。5.1xy面的微裂隙解析微裂隙的中心坐標(biāo)是通過模擬微裂隙的間距得到的，將表2中微裂隙間距分布參數(shù)輸入，編程模擬得到一定量隨機(jī)微裂隙間距，其為一矩陣S(a,b)。采用式(6)計(jì)算得到各模擬微裂隙的中心坐標(biāo)，即式中：Cx(a,b)，Cy(a,b)為(a,b)處的微裂隙中心。模擬得到XY面的微裂隙中心見圖17ㄢ在得到微裂隙中心以后，同樣，按照微裂隙方位角的分布參數(shù)和長度的分布參數(shù)模擬得到相應(yīng)的微裂隙方位角和長度，微裂隙兩端點(diǎn)的坐標(biāo)分別為式中：為第i條微裂隙的起點(diǎn)坐標(biāo)；為第i條微裂隙的終點(diǎn)坐標(biāo)；為第i條微裂隙的中點(diǎn)坐標(biāo)；iL為第i條微裂隙的長度；θi為第i條微裂隙的方位角。采用式(6)分別計(jì)算巖石細(xì)觀結(jié)構(gòu)的3個(gè)表面上微裂隙的起點(diǎn)和終點(diǎn)坐標(biāo)，運(yùn)用程序繪制XY,YZ和XZ三個(gè)面上的微裂隙，分別見圖18～20。最后，模擬細(xì)觀結(jié)構(gòu)體的三維顯示見圖21ㄢ5.2微裂隙的尺寸大小由于微裂隙間距的分布情況直接決定了微裂隙的密集程度，也即反映了微裂隙單位面積的數(shù)量，所以認(rèn)為巖石細(xì)觀REV的大小由微裂隙的間距決定，模擬得到在同一間距分布而不同數(shù)量的微裂隙具有的尺寸大小，得到面積與單位面積微裂隙數(shù)量之間的關(guān)系。圖22為以XY面微裂隙間距分布參數(shù)得到的以上兩者關(guān)系圖，從圖22中可看出，當(dāng)面積為40×106μm2時(shí)，單位面積的微裂隙數(shù)趨于穩(wěn)定，因此得到細(xì)觀結(jié)構(gòu)體即為REV，即可得到REV各個(gè)面上的微裂隙數(shù)量、REV的體積以及統(tǒng)計(jì)分析得到微裂隙方位角、長度和間距各分布函數(shù)的數(shù)學(xué)期望，見表3ㄢ6巖石細(xì)觀損傷的本構(gòu)模型6.1微裂隙損傷張量的計(jì)算T.Kyoya等基于Murakami-Ohno蠕變損傷理論給出了節(jié)理巖體的損傷基本理論，仿照巖體損傷張量計(jì)算方法，計(jì)算得到巖石細(xì)觀結(jié)構(gòu)初始損傷張量。設(shè)微裂隙結(jié)構(gòu)體在其體積V內(nèi)的裂隙平均面積為a，平均微裂隙數(shù)為N，則損傷張量可以表示為式中：n為微裂隙的法向矢量，n={n1,n2,n3}T，且有同時(shí)，參數(shù)λ滿足：另外，不同面上的方位角之間滿足：式中：Ni(i=1,2)為第in面上平均微裂隙數(shù)。對(duì)于存在多組的微裂隙，以下為計(jì)算存在多組微裂隙的損傷張量的方法，首先依據(jù)以上理論求出每組的損傷張量，然后依據(jù)能量等效原理，即由于全局損傷張量引起的能量消散等于各個(gè)損傷張量引起的能量消散的總合。計(jì)算全局損傷張量的表達(dá)式為式中：?g為全局損傷張量，?n為第n組微裂隙的損傷張量，N為微裂隙組別的數(shù)量，I為二階單位張量。將表3中細(xì)觀結(jié)構(gòu)量化參數(shù)數(shù)據(jù)代入式(8a)和(13)得到初始損傷張量矩陣為6.2g.swoboda損傷理論6.2.1連續(xù)介質(zhì)力學(xué)把細(xì)觀損傷信息引入剛度矩陣是利用連續(xù)介質(zhì)力學(xué)方法求解的關(guān)鍵。G.Swoboda和Q.Yang引入了一種計(jì)算損傷彈性系數(shù)張量的方法。這種方法所得到的損傷彈性系數(shù)張量符合對(duì)稱、正定等連續(xù)介質(zhì)力學(xué)條件。其計(jì)算公式為式中：λ和μ為Lamé常數(shù)；C1,C2和3C均為?ij的不變量的函數(shù)，且1C,C2和C3滿足：由于1C,C2和3C為常數(shù)，而且Φij為非負(fù)性，且有6.2.2主應(yīng)變損傷張量的確定在節(jié)理巖體中，一般認(rèn)為損傷與應(yīng)變方向相關(guān)聯(lián)。為了考慮裂縫的開裂與閉合以及由于裂縫中的夾雜材料或表面粗糙不平而引起的不同的裂縫表面接觸情況對(duì)應(yīng)力傳遞的影響，基于主應(yīng)變的有效損傷張量可定義為式中：pv和εv為第v個(gè)的主應(yīng)力的主矢量和主值；H(·)是修正了的階躍函數(shù)，其表達(dá)式為6.2.3．非負(fù)材料一階段的損傷張量的面中,且對(duì)于等效狀態(tài)，基于共軛力的損傷演化法則對(duì)于任何應(yīng)力狀態(tài)都是適合的，且有其中，的推導(dǎo)在有關(guān)研究中有詳細(xì)介紹，這里不作介紹，可表示為其中，式中：和ωv分別為損傷張量?的主方向和主值，且關(guān)聯(lián)的損傷面為式中：β為材料常數(shù)，g為非負(fù)材料常數(shù)。6.2.4本結(jié)構(gòu)方程本構(gòu)方程可以表示為6.3woda損傷理論由于在G.Swoboda中有一些參數(shù)需要事先確定，而其中采用三軸試驗(yàn)可以確定一部分參數(shù)，其他的參數(shù)可