土體網(wǎng)狀裂隙分形插值法生成技術(shù)

王媛 ，馮迪 ，陳尚星

(1.河海大學(xué) 巖土力學(xué)與堤壩工程教育部重點實驗室，江蘇 南京，210098；2.河海大學(xué) 隧道與城市軌道工程研究所，江蘇 南京，210098)

由于構(gòu)造運動、自然氣候和人類活動等作用，自然界的土壤表面可能會發(fā)育著大量的裂隙。裂隙的存在對土體的力學(xué)性質(zhì)[1?2]和滲透特性[3]產(chǎn)生了重要的影響，可能會導(dǎo)致土體的工程性質(zhì)變差，裂土地區(qū)的各類工程建設(shè)常因此而遭受破壞，造成巨大的經(jīng)濟損失。此外，裂隙土體中水分的流動和溶質(zhì)的運移在農(nóng)業(yè)工程[4]和環(huán)境工程[5]中也得到廣泛關(guān)注。土體裂隙網(wǎng)絡(luò)的研究[6?9]是開展裂隙土體工作的基礎(chǔ)，只有建立了與實際情況近似的網(wǎng)絡(luò)，其力學(xué)和滲流分析結(jié)果才是可靠的。許多學(xué)者借鑒巖體裂隙網(wǎng)絡(luò)的生成方法，將裂隙簡化為線段，通過蒙特卡洛方法來生成土體的裂隙網(wǎng)絡(luò)[10?11]。然而由此得到的裂隙網(wǎng)絡(luò)與實際網(wǎng)絡(luò)形態(tài)差別比較大，其主要原因是土體中裂隙形態(tài)與巖體中裂隙不同，這種裂隙形態(tài)的不同會引起其性質(zhì)的不同，例如連通率、滲透張量以及表征單元體積等。巖體中的裂隙多為 1組(或幾組)沿某優(yōu)勢方向的裂隙的隨機分布，而土體中的裂隙多為互相聯(lián)通的多邊形，如果用傳統(tǒng)的生成巖體裂隙的方法，無法模擬土體中裂隙的形態(tài)。本文作者基于分形插值理論，建立了土體裂隙網(wǎng)絡(luò)的生成方法，并編制了可視化程序，得到了與實際裂隙網(wǎng)絡(luò)形態(tài)較為接近的二維裂隙網(wǎng)絡(luò)，為研究裂隙土體的力學(xué)性質(zhì)和滲流特性建立了數(shù)值模型。

1 裂隙網(wǎng)絡(luò)幾何特征的調(diào)查與統(tǒng)計

1.1 土體裂隙網(wǎng)絡(luò)的現(xiàn)場調(diào)查

采用改進的測窗法對裂隙網(wǎng)絡(luò)進行測量與調(diào)查[12]。首先，利用數(shù)碼相機進行裂隙圖片的拍攝；然后對拍攝的裂隙圖片進行處理，得到二元化的裂隙灰度圖片，最后將處理過的裂隙圖片導(dǎo)入AutoCAD中，對裂隙網(wǎng)絡(luò)進行素描，得到裂隙空間分布圖像，并獲取裂隙網(wǎng)絡(luò)的信息。具體步驟如下：

(1) 選擇好場地以后，對拍攝區(qū)域進行選擇，要求拍攝區(qū)域平整不存在凸起和凹陷，同時拍攝區(qū)域內(nèi)的土體顏色要均一、少有雜色物質(zhì)。在拍攝之前要對拍攝區(qū)域進行清理，去除裂隙中的充填雜物，使裂隙充分暴露出來；去除土體表面上的其他雜物，剔除它們可能對影像造成的干擾。在拍攝區(qū)域邊緣放置標度尺，作為確定分析區(qū)域?qū)嶋H尺寸的依據(jù)和比例換算的標準。拍攝時鏡頭主光軸與拍攝面保持垂直，這樣能夠真實反映裂隙的形態(tài)，提高計算準確性。

(2) 將拍攝好的照片轉(zhuǎn)儲到計算機上，根據(jù)圖像中的標度尺計算確定圖像與實物的比例尺。剪掉圖像中的標度尺部分，挑選圖像質(zhì)量較好并且沒有干擾的區(qū)域作為計算區(qū)域，計算區(qū)域應(yīng)該盡量處于照片的中央位置。將圖像調(diào)整為灰度模式。對圖像進行增益除噪處理，使裂隙與其他圖像信息最大程度地分離。設(shè)置適當?shù)拈y值，將灰度圖轉(zhuǎn)換為只存在純黑和純白 2種顏色的二元位圖，即實際裂隙網(wǎng)絡(luò)圖(如圖1所示)，以此作為土體裂隙幾何學(xué)特征研究的素材。

(3) 將處理過的土體裂隙圖像導(dǎo)入到AutoCAD軟件中，對照裂隙圖像，用線段素描得到裂隙網(wǎng)絡(luò)幾何特征圖(如圖 2所示)。裂隙網(wǎng)絡(luò)圖是進行裂隙分組、統(tǒng)計、計算和裂隙網(wǎng)絡(luò)幾何特征分析的基礎(chǔ)。根據(jù)圖像與實物的比例尺，就可以獲取裂隙幾何參數(shù)的真實信息，本文中分析區(qū)域的實際尺寸(長×寬)為 275 mm×220 mm。

圖1 實際裂隙網(wǎng)絡(luò)圖Fig.1 Graph of real crack network

圖2 裂隙網(wǎng)絡(luò)幾何特征圖Fig.2 Geometric features of crack network

1.2 裂隙網(wǎng)絡(luò)的幾何特征統(tǒng)計

通過調(diào)查發(fā)現(xiàn)裂隙土的裂隙大多首尾相接，最終形成由四邊形、五邊形、六邊形及七邊形組成的網(wǎng)狀裂隙網(wǎng)絡(luò)。為了能夠在一定程度上反映裂隙的發(fā)育，需要對裂隙網(wǎng)絡(luò)的幾何特征進行概化。定義節(jié)點為裂隙與其他裂隙相交的點(如圖 2中圓圈所示)，兩節(jié)點間即為一條裂隙；拐點為裂隙延展方向發(fā)生明顯改變的點(如圖 2中三角形所示)。整個裂隙網(wǎng)絡(luò)的幾何結(jié)構(gòu)由節(jié)點、拐點和裂隙段構(gòu)成，可以通過以下參數(shù)來進行描述：

(1) 節(jié)點密度：單位面積的節(jié)點個數(shù)。

(2) 節(jié)點的度：節(jié)點與其他節(jié)點連接的個數(shù)，它反映了裂隙網(wǎng)絡(luò)節(jié)點間的連接程度。如圖 3所示，A點的度為1，B點的度為2，C點的度為3。

(3) 拐點密度：每條裂隙的平均拐點個數(shù)。

(4) 拐點變幅：每條裂隙拐點相對于節(jié)點連線的垂直距離(如圖 4所示)，它反映了每條裂隙的細部結(jié)構(gòu)形態(tài)。

圖3 節(jié)點的度Fig.3 Degree of node

圖4 拐點變幅Fig.4 Amplitude of inflection

(5) 裂隙寬度：是指裂隙張開的寬度，是描述裂隙結(jié)構(gòu)特征的一個重要指標。

裂隙網(wǎng)絡(luò)節(jié)點度的統(tǒng)計結(jié)果見表 1。可見：裂隙網(wǎng)絡(luò)共有節(jié)點41個，其中36個節(jié)點的度為3。裂隙網(wǎng)絡(luò)拐點統(tǒng)計結(jié)果見表 2。可見：裂隙網(wǎng)絡(luò)共有拐點72個，裂隙共有49條，拐點個數(shù)為2的裂隙共有25條，其次為拐點個數(shù)為0和1的裂隙，分別為11條和9條。拐點的最大變幅為17.24 mm，最小變幅為0.06 mm，平均變幅為3.64 mm。裂隙寬度最大為4.12 mm，最小為1.09 mm，平均為2.57 mm。

表1 裂隙網(wǎng)絡(luò)節(jié)點統(tǒng)計結(jié)果Table 1 Statistical results of crack network nodes

表2 裂隙網(wǎng)絡(luò)拐點統(tǒng)計結(jié)果Table 2 Statistical results of crack network inflections

2 分形插值法生成土體的裂隙網(wǎng)絡(luò)

裂隙土體中的網(wǎng)狀裂隙形態(tài)復(fù)雜，裂隙延展方向及隙寬變化較大，如用直線段模擬，就會與實際形態(tài)產(chǎn)生較大的差別。為了更接近實際裂隙網(wǎng)絡(luò)，需用另一種方法來進行模擬。易順民等[13?14]在研究膨脹土裂隙時，發(fā)現(xiàn)裂隙面的起伏形態(tài)也具有統(tǒng)計意義上的自相似性符合分形特征，依據(jù)這一結(jié)果，本文建立了一種利用分形插值來生成土體裂隙網(wǎng)絡(luò)的方法。

2.1 迭代函數(shù)系與分形插值

迭代函數(shù)系(IFS)是一種仿射映射變換方法，其實質(zhì)是通過仿射映射變換，將整體的形態(tài)變換到局部，應(yīng)用這一特性，人們模擬重構(gòu)出了大自然中許多具有自相似形態(tài)的圖案[15]。

并滿足如下條件：

則可以確定以下方程組：

對于每個仿射變換ωi有ai，ci，di，ei，fi5 個參數(shù)，如果要確定這5個參數(shù)，就必須增加1個方程或固定其中一個參數(shù)。在這里選取參數(shù)d代表變換中垂直方向變化的參數(shù)，稱其為垂直比例因子。令d為任意一個取定的實數(shù)，聯(lián)解方程可得：

2.2 土體裂隙網(wǎng)絡(luò)的生成方法

土體裂隙網(wǎng)絡(luò)的生成中采用的參數(shù)是節(jié)點、拐點數(shù)和分維數(shù)。我們通過統(tǒng)計實際土體單位面積的平均節(jié)點數(shù)來生成模擬區(qū)內(nèi)裂隙的節(jié)點數(shù)，然后根據(jù)節(jié)點之間位置關(guān)系生成拐點。裂隙土的裂隙符合分形特征，根據(jù)這一結(jié)論，裂隙土裂隙的形態(tài)應(yīng)該具有自相似特征，我們可以根據(jù)節(jié)點、拐點的位置特征，插值出節(jié)點、拐點間的裂隙形態(tài)特征，獲得的裂隙和實際裂隙具有相同的分形特征。

分形插值法生成土體裂隙網(wǎng)絡(luò)的過程如下：

(1) 生成域和分析域。考慮到邊界周圍裂隙的節(jié)點可能會在域外，實際裂隙網(wǎng)絡(luò)生成域應(yīng)比分析域大一些，然后用分析區(qū)域邊界去切割生成的裂隙網(wǎng)絡(luò)，就可以得到分析域的裂隙網(wǎng)絡(luò)圖。

(2) 節(jié)點的生成。根據(jù)節(jié)點密度，計算生成域內(nèi)節(jié)點數(shù)，按均勻分布隨機生成節(jié)點。

(3) 確定節(jié)點間的連接關(guān)系。節(jié)點間的連接需要遵從一定的關(guān)系，即裂隙之間的拓撲關(guān)系。通過對實際裂隙網(wǎng)絡(luò)觀察發(fā)現(xiàn)：節(jié)點總是與周圍距離相近的幾個節(jié)點相連。節(jié)點的度k表示節(jié)點與其他節(jié)點的連接數(shù)，因此本文假定每個節(jié)點與其距離最近的k個節(jié)點相連，由此可以確定整個區(qū)域中所有節(jié)點間的連接關(guān)系。

(4) 拐點生成。統(tǒng)計節(jié)點間的拐點數(shù)，以平均值作為模擬時節(jié)點間生成的拐點數(shù)。量測每條裂隙拐點的變幅，然后對整個區(qū)域的拐點變幅值進行平均，將其作為裂隙模擬時的變幅值。設(shè)兩節(jié)點的坐標分別為(x0,y0), (x1,y1)，則第i個拐點的坐標(xi,yi)為：

式中：m為拐點數(shù)；為拐點變幅平均值；α為該裂隙與x軸的夾角。

(5) 裂隙寬度。測量每條裂隙的寬度，然后對整個區(qū)域的裂隙寬度進行統(tǒng)計，根據(jù)下式確定第i條裂隙的寬度bi：

式中:bmax為最大隙寬；bmin為最小隙寬；random(bmax?bmin)表示在 0～(bmax?bmin)之間的均勻隨機數(shù)。

(6) 分形插值。根據(jù)節(jié)點坐標和節(jié)點間的拐點坐標，由式(5)獲得迭代函數(shù)的參數(shù)，通過下式進行迭代計算：

獲得插值點的坐標，經(jīng)過L次插值迭代后，插值點足夠多便連成線，表現(xiàn)出裂隙的形態(tài)變化和粗細變化。

2.3 裂隙網(wǎng)絡(luò)的生成結(jié)果與分析

2.3.1 參數(shù)的選取

裂隙網(wǎng)絡(luò)生成域為344 mm×275 mm(長×寬)，分析域為275 mm×220 mm(長×寬)，根據(jù)實際裂隙網(wǎng)絡(luò)幾何參數(shù)的統(tǒng)計結(jié)果可知：相同節(jié)點密度下，生成域內(nèi)節(jié)點數(shù)n為51，節(jié)點的度k平均值為3，節(jié)點間的拐點數(shù)m均值為2，拐點變幅h均值為3.64 mm，最大隙寬為4.12 mm，最小隙寬為1.09 mm，插值迭代次數(shù)L為6。改變節(jié)點坐標隨機數(shù)產(chǎn)生種子可以得到不同的裂隙網(wǎng)絡(luò)樣本，如圖5所示。

圖5 不同的裂隙網(wǎng)絡(luò)樣本Fig.5 Different samples of crack network

2.3.2 節(jié)點數(shù)n的影響

生成域內(nèi)節(jié)點數(shù)n對分形插值網(wǎng)絡(luò)結(jié)果的影響如圖6所示。由分維數(shù)計算方法[14]可以得到裂隙網(wǎng)絡(luò)的維數(shù)，從圖6可以看出：隨著n的增大，裂隙網(wǎng)絡(luò)的裂隙密度和分維數(shù)也在不斷的增大。

圖6 參數(shù)n對分形插值網(wǎng)絡(luò)結(jié)果的影響Fig.6 Effect of n on crack network generated by fractal interpolation

2.3.3 垂直比例因子d的影響

圖7 參數(shù)d對分形插值網(wǎng)絡(luò)結(jié)果的影響Fig.7 Effect of d on crack network generated by fractal interpolation

垂直比例因子d對分形插值網(wǎng)絡(luò)結(jié)果的影響如圖7所示，從圖7可以看出：d雖然是個自由參數(shù)，但其對分形插值的結(jié)果影響很大。當d=0時，分形插值退化為線性插值，此時得到的裂隙網(wǎng)絡(luò)不能反映土體裂隙的細部形態(tài)；隨著d的不斷增大，土體裂隙的細部起伏形態(tài)表現(xiàn)的越來越明顯，裂隙網(wǎng)絡(luò)的分維數(shù)逐漸增大；但是，當d=0.6時，裂隙已呈現(xiàn)出不光滑形態(tài)。因此，在裂隙網(wǎng)絡(luò)生成過程中可以通過改變垂直比例因子d來控制土體裂隙網(wǎng)絡(luò)的精細程度。

2.3.4 關(guān)于裂隙網(wǎng)絡(luò)生成結(jié)果的討論

在裂隙網(wǎng)絡(luò)生成程序中嵌入分維數(shù)計算程序，當?shù)玫降姆中尾逯稻W(wǎng)絡(luò),分維數(shù)等于實際裂隙網(wǎng)絡(luò)圖像的分維數(shù)時，插值結(jié)束，此時得到的分形插值網(wǎng)絡(luò)如圖8所示。

圖8 分形插值網(wǎng)絡(luò)(d=0.281，分維數(shù)1.466)Fig.8 Crack network generated by fractal interpolation

目前，土體裂隙網(wǎng)絡(luò)的生成主要是通過蒙特卡洛方法來生成[9]，按照裂隙網(wǎng)絡(luò)幾何參數(shù)統(tǒng)計結(jié)果，得到裂隙網(wǎng)絡(luò)如圖9所示。該方法主要存在2個缺點：一是不能反映土體裂隙的細部起伏形態(tài)；二是與土體裂隙多為互相聯(lián)通多邊形的實際形態(tài)不符。將圖8和圖9與土體實際的裂隙網(wǎng)絡(luò)(圖1)對比可以發(fā)現(xiàn)：分形插值法得到的裂隙網(wǎng)絡(luò)顯然更符合實際形態(tài)。

圖9 蒙特卡洛法生成的裂隙網(wǎng)絡(luò)Fig.9 Crack network generated by Monte-Carlo method

3 土體裂隙網(wǎng)絡(luò)的表征單元體積

由于土體中裂隙的存在，使得原來的連續(xù)介質(zhì)假定變得不合理，需要用非連續(xù)介質(zhì)的方法進行研究。在現(xiàn)階段的研究中較為成熟并行之有效的方法是等效連續(xù)介質(zhì)模型，采用等效連續(xù)介質(zhì)模型可以將由于裂隙引起的力學(xué)或滲流的集中效應(yīng)等效平均到整個土體中，將裂隙土體模擬為連續(xù)的土體，然后應(yīng)用經(jīng)典的連續(xù)介質(zhì)理論進行研究分析。在應(yīng)用等效連續(xù)介質(zhì)進行分析時存在有效性問題，即研究對象必須存在表征單元體積(REV)[18]。表征單元體積(REV)是指當試件的體積大于某一特定值V0時，由試驗獲得的參數(shù)不在隨試件的體積變化而變化，該體積就稱為表征單元體積(如圖10所示)。為此，本文引入表征單元體積的概念來對生成的裂隙網(wǎng)絡(luò)質(zhì)量做進一步地評價和討論。

圖10 表征單元體積示意圖Fig.10 Schematic diagram of REV

在多孔介質(zhì)中，一般根據(jù)空隙率來定義表征單元體積[19]。對于裂隙土體，其不連續(xù)性主要是由于裂隙引起的，求表征單元體積的過程就是將裂隙影響效應(yīng)均一化的過程，因此可以用裂隙率來定義表征單元體積。當裂隙率不隨深度變化而變化時，體裂隙率就等于面裂隙率，面裂隙率定義為：

式中：nc為面裂隙率；Ac為裂隙面積；A為裂隙土體面積。

用一系列不同尺寸的窗口(如圖11所示)去切割裂隙網(wǎng)絡(luò)，計算出每個尺寸窗口對應(yīng)的裂隙率，繪制裂隙率與計算窗口面積之間的關(guān)系曲線圖。當裂隙率變化連續(xù)小于誤差范圍且無明顯變化趨勢時，即認為裂隙率穩(wěn)定，確定不出現(xiàn)大的波動的那個最小計算窗口面積，即為表征單元體積所對應(yīng)的面積。

選取裂隙網(wǎng)絡(luò)中間的220 mm×220 mm(長×寬)區(qū)域作為表征單元體積的計算區(qū)域，用步長為 13.75 mm×13.75 mm(長×寬)的16個窗口去切割裂隙網(wǎng)絡(luò)，分形插值網(wǎng)絡(luò)與實際裂隙網(wǎng)絡(luò)的裂隙率隨計算窗口的變化如圖12所示。

圖11 一系列不同尺寸的窗口Fig.11 A series of windows with different scales

圖12 裂隙率隨計算窗口變化圖Fig.12 Variation of crack ratio with different computing windows

從圖12可以看出：實際裂隙網(wǎng)絡(luò)的表征單元面積大約為第12個計算窗口，即165 mm×165 mm(長×寬)；模擬得到的分形插值裂隙網(wǎng)絡(luò)的表征單元體積大約為第10個計算窗口，即137.5 mm×137.5 mm(長×寬)；從第12個計算窗口開始，兩者的裂隙率趨于相同。結(jié)果表明：生成網(wǎng)絡(luò)與實際網(wǎng)絡(luò)的表征單元體積相差不大，兩者的整體分布特征較為接近。當然，兩者之間還是存在一定的誤差，說明用分形插值法所得到的裂隙網(wǎng)絡(luò)還需要進一步的優(yōu)化，使其更加符合實際裂隙發(fā)育的情況。用拐點數(shù)、拐點分布特征及分維值控制生成裂隙還存在一些缺點，以后若從裂隙形態(tài)的頻譜特征去統(tǒng)計和模擬裂隙會改善單條裂隙的相似性，同時建立裂隙間的拓撲關(guān)系函數(shù)，會使整個裂隙網(wǎng)絡(luò)更趨于合理。

4 結(jié)論

(1) 采用改進的測窗法對裂隙網(wǎng)絡(luò)進行調(diào)查較為簡便，統(tǒng)計結(jié)果表明土體裂隙網(wǎng)絡(luò)中節(jié)點的度以3居多，裂隙間的拐點個數(shù)一般為2，其次為0和1。

(2) 分形插值法生成的裂隙網(wǎng)絡(luò)能很好地反映出土體裂隙的細部變化特征，基本符合其變化形態(tài)。裂隙網(wǎng)絡(luò)的裂隙密度和分維數(shù)隨著節(jié)點密度的增大而增大；垂直比例因子d對結(jié)果的影響很大，隨著d的不斷增大，土體裂隙的細部起伏形態(tài)表現(xiàn)得越來越明顯，且裂隙網(wǎng)絡(luò)分維數(shù)逐漸增大。

(3) 分形插值生成的裂隙網(wǎng)絡(luò)與土體實際裂隙網(wǎng)絡(luò)的表征單元體積和穩(wěn)定裂隙率相差不大，表明兩者的整體分布特征相似。

(4) 基于分形插值理論的土體裂隙網(wǎng)絡(luò)生成技術(shù)切實可行，生成的二維裂隙網(wǎng)絡(luò)與土體實際裂隙網(wǎng)絡(luò)形態(tài)較為接近。不可忽略裂隙網(wǎng)絡(luò)的隨機性和裂隙深度對于土體性質(zhì)的影響，如何全面考慮參數(shù)的隨機性并建立土體三維裂隙網(wǎng)絡(luò)模型需要進一步的研究。

[1] 胡卸文, 李群豐, 趙澤三, 等.裂隙性粘土的力學(xué)特性[J].巖土工程學(xué)報, 1994, 16(4): 81?88.HU Xiewen, LI Qunfeng, ZHAO Zesan, et al.Mechanical properties of fissured clay[J].Chinese Journal of Geotechnical Engineering, 1994, 16(4): 81?88.

[2] 包承鋼.非飽和土的性狀及膨脹土邊坡穩(wěn)定問題[J].巖土工程學(xué)報, 2004, 26(1): 1?15.BAO Chenggang.Behavior of unsaturated soil and stability of expansive soil slope[J].Chinese Journal of Geotechnical Engineering, 2004, 26(1): 1?15.

[3] 馬佳.裂土優(yōu)勢流與邊坡穩(wěn)定分析方法[D].武漢: 中國科學(xué)院武漢巖土力學(xué)研究所, 2007: 6?9.MA Jia.Preferential flow and stability analysis method for fissure clay slopes[D].Wuhan: Institute of Rock and Soil Mechanics.The Chinese Academy of Sciences, 2007: 6?9.

[4] Tuong T P, Cabangon R J, Wopereis M C S.Quantifying flow processes during land soaking of cracked rice soils[J].Soil Science Society of America Journal, 1996, 60(3): 872?879.

[5] Alonso E E, Springman S M, Ng C W W.Monitoring large-scale tests for nuclear waste disposal[J].Geotechnical and Geological Engineering, 2008, 26(6): 817?826.

[6] Velde B.Structure of surface cracks in soil and muds[J].Geoderma, 1999, 93: 101?124.

[7] 盧再華, 陳正漢, 蒲毅彬.膨脹土干濕循環(huán)脹縮裂隙演化的CT試驗研究[J].巖土力學(xué), 2002, 23(4): 417?422.LU Zaihua, CHEN Zhenghan PU Yibin.A CT study on the crack evolution of expansive soil during drying and wetting cycles[J].Rock and Soil Mechanics, 2002, 23(4): 417?422.

[8] Bear J U, Kent T F, Anderson S H.Image analysis and fractal geometry to characterize soil desiccation cracks[J].Geoderma,2009, 154: 153?163.

[9] Li J H, Zhang L M.Geometric parameters and REV of a crack network in soil[J].Computers and Geotechnics, 2010, 37:466?475.

[10] 袁俊平.非飽和膨脹土的裂隙概化模型與邊坡穩(wěn)定研究[D].南京: 河海大學(xué)巖土工程研究所, 2003: 22?25.YUAN Junping.Generalized Model of Fissures Distribution &Slope Stability Analysis for Unsaturated Expansive Soils[D].Nanjing: Hohai University.Institute of Geotechnical, 2003:22?25.

[11] Li J H, Zhang L M, Wang Y.Permeability tensor and representative elementary volume of saturated cracked soil[J].Canadian Geotechnical Journal, 2009, 46(8): 928?942.

[12] 李亞萍.甘肅北山花崗巖裂隙幾何特征研究[D].北京: 中國地震局地質(zhì)研究所, 2005: 6?9.LI Yaping.Geometrical Characteristics of Fractures in Granite in the Beishan Area, Gansu Province[D].Beijing: China Earthquake Administration.Institute of Geology, 2005: 6?9.

[13] 易順民, 黎志恒, 張延中.膨脹土裂隙結(jié)構(gòu)的分形特征及其意義[J].巖土工程學(xué)報, 1999, 21(3): 294?298.YI Shunmin, LI Zhiheng, ZHANG Yanzhong.The fractal characteristics of fractures in expansion soil and its significance[J].Chinese Journal of Geotechnical Engineering,1999, 21(3): 294?298.

[14] 陳尚星.基于分形理論的土體裂隙網(wǎng)絡(luò)研究[D].南京: 河海大學(xué)巖土工程研究所, 2006: 44?52.CHEN Shangxing.Reseach on crack network of soil by fractal dimension theory[D].Nanjing: Hohai University.Institute of Geotechnical, 2006: 44?52.

[15] 楊杰.分形插值及其分形維數(shù)研究[J].武漢工業(yè)學(xué)院學(xué)報,2006, 25(1): 9?11.YANG Jie.Study on fractal interpolated and its fractal dimension[J].Journal of Wuhan polytechnic University, 2006,25(1): 9?11.

[16] 王興元, 梁慶永, 沈立新.真彩色IFS吸引子的計算機構(gòu)造[J].大連理工大學(xué)學(xué)報, 2005, 45(1): 142?146.WANG Xingyuan, LIANG Qingyong, SHEN Lixin.True color IFS attractors constructed on computer[J].Journal of Dalian University of Technology, 2005, 45(1): 142?146.

[17] 陳慧琴.分形插值及分形維數(shù)的圖解法[J].江西科學(xué), 2010,28(2): 167?169.CHEN Huiqin.Study on fractal interpolation and its fractal dimension by graphic method solving[J].Jiangxi Science, 2010,28(2): 167?169.

[18] 向文飛, 周創(chuàng)兵.裂隙巖體表征單元體研究進展[J].巖石力學(xué)與工程學(xué)報, 2005, 24(增2): 5686?5692.XIANG Wenfei, ZHOU Chuangbing.The advances in investigation of representative elementary volume for fractures rock mass[J].Chinese Journal of Rock Mechanics and Engineering, 2005, 24(supp2): 5686?5692.

[19] Bear J.Dynamics of fluids in porous media[M].New York:Elsevier, 2007: 15?27.