關(guān)鍵詞：離散裂隙網(wǎng)絡(luò)模型;等效連續(xù)介質(zhì)模型;哈爾克山;離散裂隙-等效連續(xù)介質(zhì)耦合模擬;水文地質(zhì)doi：10.13278/j.cnki.jjuese.20240011 中圖分類號(hào)：TU45；P641.135 文獻(xiàn)標(biāo)志碼：A

Abstract：Halk Mountains are situated in the northern Tianshan range. Their complex geological structure and intense topographic dissection have resulted in well-developed fractures and wel-developed fissures within river valleys. This area significantly impacts the overall water resources and ecological environment of the northern Tianshan piedmont. To evaluate groundwater infiltration capacity in this region and explore the permeability characteristics of bedrock fractures in Halk Mountains area，this study proposes a coupled simulation method of discrete fracture-equivalent continuous medium based on the theoretical foundations of discrete fracture network （DFN） model and equivalent continuous medium （ECM） model. Firstly it analyzes the spatial characteristics of regional fissures by integrating the hydrogeological background and natural fissure network characteristics，then performs discrete fissure network simulation according to the natural fissure network characteristics to establish a discrete fissure network model consistent with the regional actual geological conditions，and finally achieves fissuregeological unit coupling based on the equivalent medium model theory to calculate the equivalent permeability of the fissure network.The results indicate that：Fractures in the study area are influenced by Halk Mountain fault structure， with an overall strike predominantly NWW and NEE； The fissure widths follow a log-normal distribution，and the spatial distribution of fissures adheres to a Fisher distribution;The equivalent radius of the fissures ranges from 312.9-780.6m ， the fissure width spans 1.43-25.19mm ， the permeability varies from 2.3-66.7mD ，and the fissure density ranges from 0.10- （204號(hào) 5.25m^-1 ; The equivalent permeability values are 0.050-1.877mD in the xx direction， 0.050-1.947mD （204號(hào) in the yy direction，and 0.100-2.825mD in the zz direction，with the permeability of the entire region being less than 3mD. . The bedrock fissure network in Halk Mountains area exhibits relatively low permeability， with fault zones and river valley fissure development zones serving as the primary pathways for groundwater flow in the mountainous region. Comparing the simulation results with previously conducted hydrogeological test results shows good agreement，validating the reliability of the coupled discrete fracture-equivalent continum simulation method.

Key words： discrete fracture network model; equivalent continuous medium model; Halk Mountain；discrete fracture-equivalent continuous medium coupled simulation； hydrogeology

0 引言

隨著社會(huì)的發(fā)展和進(jìn)步，諸如石油和天然氣開采、二氧化碳和核廢料的地質(zhì)封存、地下開挖等地質(zhì)工程建設(shè)越來(lái)越多，對(duì)地下水環(huán)境的保護(hù)越來(lái)越重要。裂隙巖體是地下水滲流遷移的主要通道，因此對(duì)于裂隙巖體滲透性的準(zhǔn)確描述和定量估計(jì)逐漸成為研究熱點(diǎn)[1-3]。宋曉晨等[4]認(rèn)為在地質(zhì)演變的漫長(zhǎng)過(guò)程中，巖體受到構(gòu)造運(yùn)動(dòng)、卸荷作用、風(fēng)化作用等影響，巖石內(nèi)部和表面形成了各種類型的不連續(xù)結(jié)構(gòu)面，形成裂隙巖體。對(duì)于裂隙巖體而言，裂隙的過(guò)流能力會(huì)比多孔基質(zhì)的過(guò)流能力高出數(shù)個(gè)數(shù)量級(jí)，這意味著裂隙及其網(wǎng)絡(luò)可以強(qiáng)烈影響甚至支配巖體的滲透能力和流動(dòng)路徑[5-7]。Zhao等[8]在研究中發(fā)現(xiàn)裂隙對(duì)流體流動(dòng)的影響比孔隙更顯著，裂隙網(wǎng)絡(luò)的幾何特征和連通性對(duì)巖體的滲透特性具有重要影響。Zhang等[9]發(fā)現(xiàn)由于裂隙在巖體中發(fā)育極不規(guī)則，裂隙網(wǎng)絡(luò)的幾何構(gòu)成非常復(fù)雜，由此帶來(lái)的非均質(zhì)性異常明顯，使得地下流體的滲流過(guò)程具有很強(qiáng)的不連續(xù)性和復(fù)雜性。因此，對(duì)于裂隙介質(zhì)中水流的運(yùn)動(dòng)路徑和滲透性的刻畫一直都是地下工程及地下環(huán)境保護(hù)的重難點(diǎn)，裂隙巖體滲透參數(shù)的確定對(duì)工程實(shí)踐及開展地下環(huán)境保護(hù)都有著巨大的參考價(jià)值，同時(shí)裂隙巖體流場(chǎng)的準(zhǔn)確預(yù)測(cè)在確保工程實(shí)踐的成功方面起著重要作用[10]

目前對(duì)于研究裂隙介質(zhì)中流體流動(dòng)和溶質(zhì)運(yùn)移的模擬方法主要分為兩大類：等效連續(xù)介質(zhì)模型（ECM）和離散裂隙網(wǎng)絡(luò)模型（DFN）。

等效連續(xù)介質(zhì)模型最早由Barenblatt等[11]提出，之后Snow[12]、Long 等[13]、Oda[14]、杜延齡等[15]、周志芳等[16]、張有天[17]和田開銘等[18]國(guó)內(nèi)外學(xué)者都對(duì)此模型進(jìn)行了理論研究。Nian等[19]基于雙連續(xù)體模型分析了裂隙巖質(zhì)邊坡降雨入滲邊界條件與穩(wěn)定性。Frehner 等[20]、Ritz 等[21]、Sladek等[22]、Starnoni等[23]則分別采用有限元、有限差分、有限體積、邊界元與無(wú)網(wǎng)格法等連續(xù)介質(zhì)模型，對(duì)巖體的動(dòng)力學(xué)、滲流行為與穩(wěn)定性等進(jìn)行了研究，并在工程實(shí)踐中取得了良好的應(yīng)用效果。李博南等[24]提出了一種基于頻變AVO理論與Chapman等效介質(zhì)模型的同步反演方法，用于反演裂縫密度和流體時(shí)間尺度因子。該方法以等效連續(xù)介質(zhì)模型作為理論基礎(chǔ)，該模型將裂隙介質(zhì)近似為多孔介質(zhì)，其核心在于能夠充分利用連續(xù)介質(zhì)領(lǐng)域中成熟的理論體系展開分析研究。其不僅擁有堅(jiān)實(shí)的理論支撐，還積累了豐富的案例，在解決問(wèn)題方面展現(xiàn)出可靠性和有效性[25-27]。然而，等效連續(xù)介質(zhì)模型也存在局限性，其僅適用于裂隙分布廣泛且密集的裂隙巖體，對(duì)于孤立裂隙較多或裂隙連通性差的巖體則不適用。此外，該模型在應(yīng)用中存在兩個(gè)難點(diǎn)，一是模型的有效性難以判斷，二是巖體裂隙等效滲透張量的確定較為復(fù)雜。

離散裂隙網(wǎng)絡(luò)模型最早由Louis提出，之后Koyama 等[28]、Min等[29]、Nordqvist等[30]基于該模型，針對(duì)二維和三維滲流問(wèn)題，對(duì)裂隙巖體的滲透性能展開了深入研究。張奇華等[31]則基于單純形積分理論，利用全空間塊體搜索技術(shù)構(gòu)建了離散裂隙網(wǎng)絡(luò)模型，實(shí)現(xiàn)了三維任意結(jié)構(gòu)面的網(wǎng)絡(luò)滲流計(jì)算。張雨晴等[32]以渤?；◢弾r地區(qū)為例，建立了離散裂隙網(wǎng)絡(luò)模型，對(duì)其儲(chǔ)存裂縫進(jìn)行成因分析，構(gòu)建了強(qiáng)風(fēng)化帶裂縫分布趨勢(shì)模型。Li等[33基于離散裂隙網(wǎng)絡(luò)模型，提出了一種研究裂隙巖體中地下水流動(dòng)與耦合行為的方法，發(fā)現(xiàn)裂隙中的水分對(duì)流會(huì)顯著影響溶質(zhì)濃度的分布，表明離散裂隙分布在加速物質(zhì)運(yùn)移過(guò)程中起到了重要作用。Karimzade等[34]利用代碼構(gòu)建了離散裂隙網(wǎng)絡(luò)模型，對(duì)地下工程實(shí)際涌水量進(jìn)行了模擬預(yù)測(cè)，預(yù)測(cè)結(jié)果的準(zhǔn)確性在伊朗某地下工程中得到了驗(yàn)證。自前，離散裂隙網(wǎng)絡(luò)模型發(fā)展迅速，但這種模型因忽略巖石基質(zhì)的透水性，僅適用于裂隙網(wǎng)絡(luò)貫通良好的場(chǎng)景。同時(shí)，其無(wú)法有效反映孤立裂隙對(duì)整體滲透性的影響，有一定局限性。

哈爾克山地處天山北部，受復(fù)雜地質(zhì)結(jié)構(gòu)和強(qiáng)烈地形切割作用影響，區(qū)域內(nèi)斷裂發(fā)育密集，河谷裂隙廣泛分布。該區(qū)域?qū)μ焐奖甭吹乃Y源及生態(tài)環(huán)境具有重要影響，因此開展哈爾克山地區(qū)地下水滲透能力評(píng)價(jià)研究意義重大。深入探究該區(qū)裂隙巖體的滲透能力，對(duì)有效保護(hù)和合理利用自然資源、促進(jìn)天山北麓生態(tài)保護(hù)具有重要作用。本文以哈爾克山為研究對(duì)象，通過(guò)定量分析研究區(qū)裂隙空間展布特征與發(fā)育規(guī)律，采用離散裂隙-等效介質(zhì)混合模型方法，構(gòu)建離散裂隙網(wǎng)絡(luò)模型，并在此基礎(chǔ)上，計(jì)算區(qū)域內(nèi)裂隙等效滲透張量。研究成果旨在為評(píng)價(jià)該區(qū)地下水滲流能力提供依據(jù)，并為未來(lái)哈爾克山區(qū)域裂隙水研究提供新的技術(shù)方法和理論支持。

1 地質(zhì)概況

研究區(qū)位于新疆哈爾克山區(qū)域內(nèi)。哈爾克山呈近東西向條帶狀展布，地貌上為侵蝕中一高山區(qū)。其南北兩側(cè)邊界為NW向和NEE向的逆沖斷裂，北界為昭蘇盆地南緣斷裂，南界為哈爾克山斷裂。那拉提斷裂以北為那拉提隆起。該隆起帶呈NEE向，遭受強(qiáng)烈擠壓，主要由中深變質(zhì)巖系和花崗巖組成。其高程 1 700～4 000m ，相對(duì)高差大于 1500m ，地形切割劇烈，溝谷發(fā)育。研究區(qū)邊界北部以昭蘇盆地南緣斷裂帶基巖為界，西、南邊界為阿克牙孜河分水嶺，東部邊界為庫(kù)克蘇河。研究區(qū)地層巖性主要為古元古界木扎爾特巖群斜長(zhǎng)片麻巖，石炭系大哈拉軍山組中段安山玢巖、英安斑巖夾流紋斑巖，石炭系阿克沙組灰?guī)r等。侵入巖為華力西中期侵入巖，巖性主要為花崗巖（圖1），巖石堅(jiān)硬，礦物成分主要為石英、斜長(zhǎng)石、鉀長(zhǎng)石、黑云母等。

研究區(qū)地表水由南向北徑流，物理風(fēng)化作用強(qiáng)烈導(dǎo)致地表巖體破碎，巖屑堆積物厚度巨大。在地表較大坡降條件下，即使是穩(wěn)定融化的雪水也能迅速下滲進(jìn)入堆積層和基巖裂隙，形成基巖裂隙水，僅少量形成地表徑流。該地段潛水位很高，含水帶厚度不均勻且起伏較大；同時(shí)受冰雪季節(jié)性消融及降水影響，地表水水量變化及地下水水位變幅較大。

.晚更新世洪積層； Qh^al .全新世沖積層； HtM .古元古界木扎爾特巖群； C₁d² .大哈拉軍山組中段； C_1-2a .阿克沙組； N_1-2k .康村組； C₂δ .晚石炭世閃長(zhǎng)巖； C₁γβ 早石炭世黑云母花崗巖； C₁ηγ .早石炭世二長(zhǎng)花崗巖； C₁δo .早石炭世石英閃長(zhǎng)巖； C₁γδσ .早石炭世英云花崗 巖； Pξγ ，二疊紀(jì)正長(zhǎng)花崗巖； Pβμ .二疊紀(jì)輝綠巖； 志留紀(jì)英云閃長(zhǎng)巖； S₃γ ，晚志留世破碎花崗巖； sinγ .志留紀(jì)二長(zhǎng)花崗巖； ChA_.³ .阿克 蘇巖群白云母石英片巖。

圖1 研究區(qū)區(qū)域地質(zhì)圖

2 模型理論

2.1 離散裂隙網(wǎng)絡(luò)模型

離散裂隙網(wǎng)絡(luò)模型將每個(gè)裂隙抽象為特定幾何形態(tài)的薄盤，通常假定為圓形或正多邊形。每條裂隙由其中心點(diǎn)位置、直徑、產(chǎn)狀和開度定義。模型由下述性質(zhì)定義：圓盤中心點(diǎn)構(gòu)成一個(gè)三維泊松點(diǎn)過(guò)程；圓盤直徑相互獨(dú)立，具有相同的分布；圓盤產(chǎn)狀相互獨(dú)立，具有相同的分布；直徑和產(chǎn)狀相互獨(dú)立；裂隙開度相互獨(dú)立，具有相同的分布[4]。離散裂隙的生成模型包括分形模型（fractalmodel）、最近距離模型（nearestneighbormodel）與Baecher模型，經(jīng)過(guò)前人[35]多種假設(shè)和研究后認(rèn)為，Baecher模型將裂隙假定為圓盤最為可靠。離散裂隙網(wǎng)絡(luò)模型是基于已有的裂隙進(jìn)行概率統(tǒng)計(jì)后，根據(jù)得到的概率分布參數(shù)進(jìn)行確定的隨機(jī)模擬，從而生成離散裂隙網(wǎng)絡(luò)。

2.1.1 裂隙密度

裂隙密度是衡量區(qū)域裂隙發(fā)育程度的參數(shù)，其不僅反映裂隙發(fā)育的密集程度，而且控制著巖體的完整性。一般情況下，裂隙面發(fā)育越密集，巖體的強(qiáng)度越低，完整性越差，巖塊越小，巖體力學(xué)性質(zhì)越差，滲透性越強(qiáng)。根據(jù)尺度的不同，裂縫密度又分為線密度 （P₁₀ ）、面密度 （P₂₀） ）、體密度 （P₃₀） 、累積密度等表示方式[36]。其中：線密度定義為在裂縫法線方向上單位長(zhǎng)度的裂縫數(shù)量；累積密度能夠直觀反映裂縫密度的垂向分布，便于進(jìn)行明確的垂向分區(qū)。本文采用線密度進(jìn)行密度約束，采用累積密度進(jìn)行井段分區(qū)。

2.1.2 裂隙直徑

裂隙大小是決定巖體滲流、力學(xué)等性質(zhì)的重大因素，目前現(xiàn)有的技術(shù)水平無(wú)法直接獲取裂隙大小。在Baecher模型中一般用裂隙直徑表征裂隙大小，其獲取手段是根據(jù)地表出露的跡長(zhǎng)結(jié)合區(qū)域?qū)嶋H大小、地質(zhì)背景進(jìn)行推測(cè)分析。前人[37計(jì)算裂隙直徑的基本研究思路為：首先假設(shè)裂隙圓盤的大小服從某種概率分布；再通過(guò)對(duì)裂隙跡長(zhǎng)的觀測(cè)值進(jìn)行擬合分析，進(jìn)而確定該分布特征的幾何參數(shù)。一般認(rèn)為，跡長(zhǎng)和裂隙直徑符合對(duì)數(shù)正態(tài)分布[4]

2.1.3 裂隙開度

裂隙開度一般指裂隙兩壁間的垂直距離，是地下流體流動(dòng)的主要通道，開度的大小一般決定了流體流動(dòng)速度的大小。裂隙開度是當(dāng)前研究的熱點(diǎn)課題之一。眾多學(xué)者研究[37]表明，裂隙開度的分布特征通常服從對(duì)數(shù)正態(tài)分布或冪律分布。但是在天然裂隙系統(tǒng)中，很難描述每個(gè)裂隙開度的分布。因此，通常用最大開度 （d_max ）或平均開度 （d_ave=πd_max/4） ）來(lái)代表單裂縫的開度。Olson[38]在研究中使用了一個(gè)非線性冪律關(guān)系：

式中： L 為裂隙長(zhǎng)度； α 為冪律指數(shù)； λ 為與圍巖力學(xué)參數(shù)相關(guān)的指數(shù)，表示為

式中： K_lc 為裂隙的韌性； u 為巖石的泊松比； E 為巖 石楊氏模量。Olson[3]進(jìn)一步指出，密集裂縫間的 力學(xué)作用會(huì)導(dǎo)致指數(shù) λ 隨裂隙間距而變化。

大量的統(tǒng)計(jì)測(cè)量表明， 0.22lt;λlt;0.69 ，均值為0.50。而裂隙的滲透率則通過(guò)立方定律[37計(jì)算：

式中， K_f 為裂縫的滲透率。

2.1.4 裂隙空間分布特征

在實(shí)際地質(zhì)情況中，裂隙巖體中流體的流動(dòng)方向一般由裂縫方位控制，因此在構(gòu)建地質(zhì)裂隙模型的時(shí)候，需要著重注意裂隙方位的刻畫，使其更好地符合實(shí)際的裂縫分布情況，從而使擬合效果更具有代表性和真實(shí)性。目前，F(xiàn)isher分布是最常見的描述裂縫方位分布的概率密度函數(shù)，單變量Fisher分布是球面上的正態(tài)（高斯）分布，當(dāng)裂縫方向受巖石性質(zhì)、應(yīng)力條件等多個(gè)不相關(guān)變量控制時(shí)，單變量Fisher分布在理論上是合理的[39]

Fisher分布由概率密度函數(shù)定義：

f（θ^′，φ^′）=

式中： κ 為裂隙分布參數(shù)，表示裂隙濃聚度； θ^′ 為裂隙方位角； φ^′ 為偏離平均裂縫軸的角度。

Fisher分布是單峰的，且圍繞 φ^′ 軸對(duì)稱。增加κ 會(huì)使分布更加集中在 φ^′ 軸周圍。

2.2 等效連續(xù)介質(zhì)模型

等效連續(xù)介質(zhì)模型將裂隙巖體概化成具有連續(xù)性的滲透等效介質(zhì)。該模型未考慮單個(gè)裂隙結(jié)構(gòu)面在流體滲流過(guò)程中的具體作用機(jī)制，而是將其影響統(tǒng)一整合到介質(zhì)的等效滲透特性之中。等效連續(xù)介質(zhì)模型一般將研究區(qū)域剖分成若干有限單元體。首先根據(jù)裂隙網(wǎng)絡(luò)構(gòu)建合適的REV（representativeelementaryvolume）單元體；然后采用有限差分法將單元體內(nèi)裂隙耦合進(jìn)單元體，計(jì)算單元體等效滲透張量，通過(guò)有限單元體構(gòu)建的連續(xù)體模型，從而刻畫區(qū)域的宏觀滲透性。目前常用的裂隙耦合方法是基于 Oda^[40] 提出來(lái)的方法。Oda裂隙張量是在滲流區(qū)域內(nèi)任取一體元，地下水僅沿裂隙流動(dòng)，忽略巖塊滲流，將裂隙概化為半徑為 r 的圓盤，裂隙開度為b ，裂隙面法向單位矢量為 ，某一裂隙可用 （n，2r ，b ）來(lái)表示，水流沿 （n，2r，b） 裂隙流動(dòng)可概化為理想化的厚度為 η（η=b ）的板內(nèi)層流問(wèn)題，裂隙的滲流遵循立方定律。含有一定數(shù)量裂隙單元體的滲透張量可用下式表示：

k_ij=_μ（A_kkδ_ij-A_ij）;

$E ＼left（ n ， 2 r ， b ＼right） ＼mathrm { d } ＼Omega ＼mathrm { d } （ 2 r ） ＼mathrm { d } b ＼ 。$

式中： k_ij 為裂隙滲透張量， i，j 為三維空間坐標(biāo)索引，取值1，2，3對(duì)應(yīng) x，y，z 方向； A_ij 為裂隙幾何張量； A_kk 為張量 A 的跡，三維空間中 A_kk=A₁₁+A₂₂+ A₃₃：δ_ij 為克羅內(nèi)克（Kronecker）符號(hào)； b_m 為開度最大值； 2r_m 為直徑最大值； 為空間角度域; 為單位向量 的分量； 為裂隙分布的概率密度函數(shù)； μ 為無(wú)量綱常數(shù)，如果測(cè)定裂隙連通性好，則 μ 取 1/12 。

3 模型構(gòu)建

3.1 區(qū)域裂隙網(wǎng)絡(luò)生成

本文基于研究區(qū)前期測(cè)井?dāng)?shù)據(jù)及14d的野外踏勘數(shù)據(jù)，對(duì)已有的裂隙數(shù)據(jù)進(jìn)行分析處理，構(gòu)建區(qū)域離散裂隙網(wǎng)絡(luò)。先基于六組測(cè)井裂隙數(shù)據(jù)，繪制了井中裂隙分布圖（圖2）、井中裂隙走向玫瑰圖（圖3）和井中裂隙等值線圖（圖4）。再根據(jù)井中裂隙分布特征得到了測(cè)井裂隙線密度分布（表1）及測(cè)井裂隙空間Fisher分布模型分析結(jié)果（表2）。根據(jù)前期地質(zhì)勘查及實(shí)際野外踏勘數(shù)據(jù)，假定研究區(qū)裂隙隙寬服從對(duì)數(shù)正態(tài)分布，則裂隙直徑由隙寬與直徑關(guān)系式推導(dǎo)計(jì)算。

圖2井中裂縫分布圖

圖3井中裂縫走向玫瑰圖 （30^°） （20

Fig.3 Rosechartoffracturestrikeofthewell （30^° 0

研究區(qū)域在地質(zhì)形成過(guò)程中受到塔里木和伊犁兩個(gè)古陸塊的對(duì)擠形成了一系列斷裂帶，其中主要為昭蘇盆地南緣斷裂帶及喀拉古溫泉斷裂帶。根據(jù)地質(zhì)測(cè)繪及遙感解譯對(duì)研究區(qū)域斷裂進(jìn)行分析刻畫，初步判定研究區(qū)域內(nèi)發(fā)育26條中小型斷裂，斷裂帶內(nèi)巖石破碎，充填大量的膠結(jié)物，孔隙發(fā)育，其滲透能力遠(yuǎn)遠(yuǎn)大于巖體裂隙；因此斷裂帶是研究區(qū)域不可忽視的滲流通道，在區(qū)域滲透率研究中必須考慮斷層的影響。同時(shí)研究區(qū)域內(nèi)還發(fā)育三條中大型河谷，河谷區(qū)域因卸荷作用形成大量裂隙，滲透率顯著高于正常巖體。這使得河流入滲補(bǔ)給成為區(qū)域地下水流系統(tǒng)的主要補(bǔ)給方式。因此，在數(shù)值模型構(gòu)建過(guò)程中，也必須考慮河谷區(qū)域這個(gè)關(guān)鍵因素。通過(guò)坐標(biāo)投影和層面概化的方式，將斷裂與河谷導(dǎo)入數(shù)值模型中（圖5）。依據(jù)區(qū)域裂隙特征分析獲得的裂隙參數(shù)，結(jié)合測(cè)井與踏勘數(shù)據(jù)分布位置及研究區(qū)域地質(zhì)構(gòu)造的發(fā)育情況，在數(shù)據(jù)分布較為集中的海拔 1370m 至地表的垂向范圍開展數(shù)值模擬。

圖4井中裂縫等值線圖

Fig.4 Fracture contour map of the well

表1測(cè)井裂隙線密度分布

Table 1 Linear density distribution of logging fractures

本文采用Baecher模型構(gòu)建區(qū)域裂隙網(wǎng)絡(luò)進(jìn)行數(shù)值模擬。建模過(guò)程對(duì)研究區(qū)域裂隙網(wǎng)絡(luò)作如下假定：1）所有裂隙均可以用圓盤近似表達(dá)，裂隙直徑服從對(duì)數(shù)正態(tài)分布，裂隙隙寬與裂隙長(zhǎng)度服從相同分布且具有數(shù)學(xué)關(guān)系；2）此次建模過(guò)程不考慮裂隙之間相互截?cái)噙^(guò)程，即裂隙間截?cái)喔怕试O(shè)置為0；3）研究區(qū)域具有相同的地層巖性，區(qū)域裂隙網(wǎng)絡(luò)具有一致性；4）單個(gè)裂隙內(nèi)為光滑平行面，認(rèn)為根據(jù)隙寬與長(zhǎng)度關(guān)系得到的隙寬是真實(shí)隙寬，滲透率計(jì)算過(guò)程中不考慮實(shí)際地質(zhì)的表面摩擦系數(shù)及隙寬誤差，即計(jì)算的絕對(duì)滲透率為實(shí)際滲透率。

3.2 區(qū)域等效介質(zhì)

本文裂隙模擬垂向范圍為海拔 1370m 至地表區(qū)域，考慮到區(qū)域地形起伏較大，結(jié)合裂隙直徑大小及計(jì)算機(jī)算力，在垂向上將整個(gè) 1 370m 至地表的空間均勻剖分成5層。在水平方向考慮研究區(qū)域面積及裂隙分布情況，將網(wǎng)格剖分成 100m×100m 的規(guī)則網(wǎng)格，即將研究區(qū)域在空間上剖分成529690個(gè)網(wǎng)格。首先根據(jù)Oda裂隙耦合方法，利用有限差分計(jì)算將網(wǎng)格體內(nèi)裂隙耦合嵌入網(wǎng)格體，同時(shí)考慮巖體基質(zhì)滲流作用，對(duì)整體的網(wǎng)格體賦予基質(zhì)滲透率；再根據(jù)研究區(qū)域地層巖性特點(diǎn)，結(jié)合前期的巖心流動(dòng)實(shí)驗(yàn)，設(shè)置模型基質(zhì)滲透率為 z 方向 0.1mD x，y 方向 0.05mD ;最后根據(jù)前期室內(nèi)實(shí)驗(yàn)結(jié)果可知，區(qū)域內(nèi)基巖孔隙度在 6% 以下，主要集中在1%～4% 之間，結(jié)合區(qū)域地質(zhì)背景，設(shè)置模型基質(zhì)孔隙度為 4% 。

表2測(cè)井Fisher分布模型分析結(jié)果

Table2Analysisresults of well logging Fisher distribution model

圖5裂隙模擬區(qū)地質(zhì)示意圖

Fig.5 Geological diagramoffracture simulation area

研究區(qū)域地質(zhì)條件復(fù)雜，發(fā)育昭蘇盆地南緣斷裂帶及喀拉溫古泉斷裂帶兩條大型斷裂帶，并由此衍生出26條中小型斷裂帶。根據(jù)前期開展的地質(zhì)調(diào)查和現(xiàn)場(chǎng)物探數(shù)據(jù)可知，研究區(qū)斷裂帶充填較好，破碎程度一般，計(jì)算認(rèn)為研究區(qū)斷裂帶區(qū)域的平均滲透率為 230.71mD ，區(qū)內(nèi)斷裂帶平均寬度為160m;因此在本次模擬中對(duì)于距離斷裂中心位置 80m 的網(wǎng)格滲透率設(shè)置為 230.71mD 。區(qū)域發(fā)育有三個(gè)較大的河谷裂隙區(qū)，根據(jù)前期水文地質(zhì)調(diào)查，設(shè)定河谷裂隙區(qū)的平均滲透率為 420.5mD ，影響范圍為300m ，即在數(shù)值模型中將距河谷層面 150m 網(wǎng)格的滲透率設(shè)置為 420.5mD 。結(jié)合裂隙網(wǎng)絡(luò)耦合的等效滲透率，疊加計(jì)算得到研究區(qū)域綜合滲透率。

4模型成果及檢驗(yàn)

4.1 模型成果

4.1.1 區(qū)域裂隙網(wǎng)絡(luò)結(jié)果

根據(jù)3.1小節(jié)對(duì)區(qū)域測(cè)井的裂隙數(shù)據(jù)進(jìn)行分析處理得到模型裂隙參數(shù)，以及不同走向及其對(duì)應(yīng)的線密度數(shù)據(jù)，共同構(gòu)建了10個(gè)離散裂隙集，與區(qū)域斷裂帶、河谷裂隙區(qū)共同組成了研究區(qū)域的裂隙網(wǎng)絡(luò)（圖6）。對(duì)整個(gè)區(qū)域裂隙的傾角和走向進(jìn)行統(tǒng)計(jì)，繪制裂隙傾角直方圖和裂隙走向直方圖（圖7）。從傾角直方圖（圖7a）可以看出，模型區(qū)裂隙主要為高角度裂隙，主要原因在于裂隙發(fā)育受到區(qū)域斷裂帶的控制。該研究區(qū)域內(nèi)的26條斷裂帶幾乎貫穿整個(gè)模擬區(qū)，且其傾角在 70^°～88^° 之間，屬于近垂直斷層的高角度斷層。此外，模型裂隙的主要走向?yàn)镹WW向和NEE向（圖7b），這是由研究區(qū)構(gòu)造時(shí)期的NW向和NEE向逆沖斷層影響所導(dǎo)致。模型模擬的裂隙網(wǎng)絡(luò)等效半徑為 312.9～780.6m ，裂隙隙寬為 1.43～25.19mm ，而裂隙網(wǎng)絡(luò)滲透率則在 2.3～ 66.7mD 之間。研究區(qū)域裂隙網(wǎng)絡(luò)隙寬大多分布在8～18mm 之間，滲透率集中在 6.3～25.0mD 之間。裂隙網(wǎng)絡(luò)整體滲透率較低，主要原因是研究區(qū)地質(zhì)構(gòu)造復(fù)雜，裂隙之間相互截?cái)嗖⒈舜擞绊?，致使裂隙網(wǎng)絡(luò)整體呈現(xiàn)低滲特征。結(jié)合3.2小節(jié)對(duì)斷裂帶及河谷裂隙區(qū)的分析可知，模型中斷裂滲透率為230.71mD ，影響范圍為 160m ；河谷裂隙區(qū)滲透率為 420.5mD ，影響范圍為 300m 。

圖6研究區(qū)裂隙網(wǎng)絡(luò)模型示意圖

4.1.2 區(qū)域等效介質(zhì)結(jié)果

法，計(jì)算得到研究區(qū)的單元體滲透率與孔隙度。滲透率耦合結(jié)果為 xx、yy、zz、xy、xz、yz 共6個(gè)方向的滲透率，其中： xy?xz?yz 方向的滲流能力極差，滲透率均小于 0.001mD ，可忽略不計(jì)； xx 方向滲透率為 0.050～1.877mD 圖 8a）sy 方向滲透率為 0.050～1.947mD （圖8b）； zz 方向滲透率為 （圖8c），等效介質(zhì)孔隙度為1.00%～13.62% （圖8d）。結(jié)果表明，模擬區(qū)域整體滲透性較低，孔隙發(fā)育程度較弱。同時(shí)，考慮區(qū)域斷裂帶的影響，對(duì)受斷裂影響網(wǎng)格進(jìn)行單獨(dú)賦值（圖9）。計(jì)算所得區(qū)域裂隙密度為 0.10～5.25m^-1 ，表明區(qū)域裂隙發(fā)育程度較低。區(qū)域綜合滲透率的結(jié)果表明，研究區(qū)域巖體滲流能力較弱，整體小于 3mD 屬于低滲區(qū)域。研究區(qū)斷裂破碎帶和河谷裂隙區(qū)的滲透能力顯著強(qiáng)于基巖裂隙區(qū)，是區(qū)域內(nèi)典型的高滲透區(qū)域，也是研究區(qū)地下水滲流的主要通道。高滲通道不僅是地下水滲流的主要通道，也是地下水環(huán)境保護(hù)中需要著重考慮的區(qū)域，在實(shí)際應(yīng)用中，需采取有效措施阻隔污染物通過(guò)高滲通道侵入地下水系統(tǒng)，避免引發(fā)區(qū)域性地下水污染。

4.2 模型檢驗(yàn)

前期水文地質(zhì)勘查工作在模擬區(qū)內(nèi)開展抽水試驗(yàn)12次、壓水試驗(yàn)301段。其中抽水試驗(yàn)選擇各單一巖性鉆孔段開展，壓水試驗(yàn)在裂隙發(fā)育及基巖地區(qū)均開展實(shí)驗(yàn)。通過(guò)抽水試驗(yàn)和壓水試驗(yàn)的數(shù)據(jù)，大致計(jì)算出各種巖性區(qū)域的滲透率平均值。根據(jù)地質(zhì)背景分析，該地區(qū)內(nèi)主要巖石類型包括玻屑晶屑凝灰?guī)r、黑云母二長(zhǎng)花崗巖、花崗巖、輝長(zhǎng)巖、閃長(zhǎng)巖和斜長(zhǎng)片麻巖，它們的滲透率分布在 0.72～3.5mD 之間。根據(jù)模擬區(qū)巖性分布情況，對(duì)研究區(qū)內(nèi)對(duì)應(yīng)巖性地層滲透率進(jìn)行初步劃分，結(jié)果顯示其主要集中在 0.8～3.2mD 之間（圖10）。為了進(jìn)一步了解裂隙特性，在LZK5、LZK7和LZK35等裂隙發(fā)育程度較高的地區(qū)開展了裂隙壓水試驗(yàn)。通過(guò)對(duì)滲透率數(shù)據(jù)和裂隙網(wǎng)絡(luò)模型進(jìn)行對(duì)比分析可知，鉆井裂隙壓水試驗(yàn)測(cè)得的滲透率范圍為 1.07～88.65mD ，而裂隙網(wǎng)絡(luò)模型給出的結(jié)果則為 2.35～66.98mD （圖11）。結(jié)果表明，數(shù)值模擬結(jié)果與水文地質(zhì)試驗(yàn)結(jié)果近乎相同，表明本文構(gòu)建的數(shù)值模型能較好地反映實(shí)際區(qū)域的水文地質(zhì)特征，模型模擬結(jié)果具有較高的準(zhǔn)確性和可靠性。由于數(shù)值模擬未能全面考慮到裂隙與巖石基質(zhì)間相互作用對(duì)滲流場(chǎng)產(chǎn)生的影響，僅考慮了裂隙流動(dòng)與基質(zhì)滲透簡(jiǎn)單疊加的情況，因此模擬結(jié)果的滲透率與實(shí)際現(xiàn)場(chǎng)的試驗(yàn)結(jié)果存在偏差。

由3.2小節(jié)區(qū)域地質(zhì)網(wǎng)格等效滲透系數(shù)的方

圖7裂隙網(wǎng)絡(luò)屬性統(tǒng)計(jì)直方圖

Fig.7 Statistical histogram of fracture network attributes

圖8等效模型模擬結(jié)果示意圖

圖9 斷裂河谷滲透率示意圖

Fig.9Permeabilitydiagramofthefractured valley

圖10抽水試驗(yàn)巖性滲透率與模擬結(jié)果對(duì)比

圖11壓水試驗(yàn)鉆孔裂隙段滲透率與模擬結(jié)果對(duì)比 Fig. 11 Comparison between permeability and simulation results of borehole fracture section in water pressure test

5 結(jié)論與建議

1）受復(fù)雜地質(zhì)結(jié)構(gòu)和強(qiáng)烈地形切割作用的影響，哈爾克山斷裂發(fā)育密集，河谷裂隙廣泛分布。區(qū)域內(nèi)裂隙網(wǎng)絡(luò)受到哈爾克山構(gòu)造成因的影響，整體裂隙走向主要為NWW向和NEE向。

2）哈爾克山地區(qū)地層巖性主要為花崗巖和變質(zhì)巖，裂隙發(fā)育程度較低，張開度較小，致使整體區(qū)域基巖滲透能力較低。因此，區(qū)域內(nèi)發(fā)育的斷裂帶及河谷裂隙成為區(qū)域地下水主要的補(bǔ)給通道與滲流通道，也是地下水遷移的主要路徑。

3）本文采用了離散裂隙-等效連續(xù)介質(zhì)耦合方法進(jìn)行模擬，該方法能夠較為準(zhǔn)確地模擬出裂隙網(wǎng)絡(luò)，并計(jì)算出區(qū)域的等效滲透系數(shù)。模擬結(jié)果與實(shí)際水文地質(zhì)試驗(yàn)參數(shù)相近，證明本文提出的離散裂隙-等效連續(xù)介質(zhì)耦合模擬方法具有可行性。本文方法不僅考慮了裂隙在滲流中的優(yōu)勢(shì)作用，還能簡(jiǎn)化后期滲流過(guò)程的計(jì)算量，從而更高效地分析計(jì)算大尺度區(qū)域的滲透率。

4）受軟件計(jì)算能力的限制，本文數(shù)值模型構(gòu)建的離散裂隙網(wǎng)絡(luò)主要選取對(duì)區(qū)域滲透率有顯著影響的大裂隙，未考慮區(qū)域內(nèi)的中小裂隙（直徑 lt;100 m）。裂隙模型精細(xì)度不足且未考慮基質(zhì)與裂隙滲流場(chǎng)的相互作用，僅對(duì)基質(zhì)與裂隙的滲流能力進(jìn)行了簡(jiǎn)化處理。在未來(lái)研究中應(yīng)進(jìn)一步完善不足，使模型更加符合實(shí)際情況。

致謝：北京遠(yuǎn)望景盛科技發(fā)展有限公司提供了Fracman軟件及技術(shù)支持，特此致謝。

參考文獻(xiàn)（References）：

[1]Huang N，Jiang Y J，Li B，et al. A Numerical Method for Simulating Fluid Flow Through 3 -D Fracture Networks[J]. Journal of Natural Gas Science and Engineering，2016，33：1271-1281.

[2] 劉宗賓，李超，路研，等.基于孔隙結(jié)構(gòu)表征的低滲透砂 巖流體賦存特征及滲透率評(píng)價(jià)[J].吉林大學(xué)學(xué)報(bào)（地 球科學(xué)版），2024，54（4）：1124－1136. Liu Zongbin，Li Chao，Lu Yan，et al. Fluid Occurrence State and Permeability Evaluation of Low-Permeability Sandstone Based on Pore Structure Characterization [J].Journal of Jilin University （Earth Science Edition），2024，54（4）：1124-1136.

[3]袁益龍，唐佳偉，李鵬，等.熱-力耦合作用下花崗巖粗 糙裂隙滲透率演化規(guī)律［J].吉林大學(xué)學(xué)報(bào)（地球科學(xué) 版），2025，55（1）：179-187. Yuan Yilong，Tang Jiawei，Li Peng，et al.Permeability Evolution of Granite with Rough Fracture Under Thermal-Mechanical Coupling[J]. Journal of Jilin University（Earth Science Edition），2025，55（1） ：179 - 187.

[4]宋曉晨，徐衛(wèi)亞.裂隙巖體滲流模擬的三維離散裂隙網(wǎng) 絡(luò)數(shù)值模型：I：裂隙網(wǎng)絡(luò)的隨機(jī)生成[J].巖石力學(xué)與工 程學(xué)報(bào)，2004，23（12）：2015-2020. Song Xiaochen，Xu Weiya. Numerical Model of ThreeDimensional Discrete Fracture Network for Seepage in Fractured Rocks：I ：Generation of Fracture Network LJ」：Uninese Journal oI Kock Iviecnanicsand Engineering，2004，23（12）：2015-2020.

[5]王庶丞，唐華風(fēng)，鄒明倬，等.渤海灣盆地萊州灣凹陷南 部斜坡帶中生界火山巖裂縫特征及分布規(guī)律[J].吉林 大學(xué)學(xué)報(bào)（地球科學(xué)版），2024，54（6）：1846－1859. Wang Shucheng，Tang Huafeng，Zou Mingzhuo，et al. CharacteristicsandDistributionofFracturesin Mesozoic Volcanic Rocks in Southern Slope of Laizhou Bay Depression，Bohai Bay Basin[J].Journal of Jilin University（Earth Science Edition），2024，54（6）： 1846-1859.

[6]宋清，孫盼科，相金元，等.鄂爾多斯盆地蘇里格東二區(qū) 走滑斷裂結(jié)構(gòu)特征及其對(duì)流體分布的控制作用[J].現(xiàn) 代地質(zhì)，2023，37（5）：1110－1122. Song Qing，Sun Panke，Xiang Jinyuan，et al. Structural Characteristics of Strike-Slip Fault and Its Control on Fluid Distribution in the East Block of the Sulige Gas Field，Ordos Basin[J]. Geoscience，2023，37（5） ：1110 - 1122.

[7]路芳，高明星，周書賢，等.阿爾金東段斷裂帶流域地貌 特征及其構(gòu)造活動(dòng)[J].現(xiàn)代地質(zhì)，2023，37（5）：1100 - 1109. Lu Fang，Gao Mingxing， Zhou Shuxian，et al. Fluvial Geomorphology of the Eastern Altyn Tagh and Its Tectonic Activity Implications[J]. Geoscience，2023，37 （5）：1100 -1109.

[8]Zhao Y S，F(xiàn)eng Z C，Lv Z X，et al. Percolation Laws of a Fractal Fracture-pore Double Medium[J]. FractalsComplex Geometry Patterns and Scaling in Nature and Society，2016，24（4）：1650053.

[9]Zhang X，Sanderson D J. Numerical Study of Critical BehaviourofDeformationandPermeabilityof Fractured Rock Masses[J]. Marine and Petroleum Geology，1998，15（6）：535- 548.

[10]Wang L C，Zheng L Z，Wang T J，et al.Development ofDepth-Averaged Nonlinear Flow Modelfor Fractures Considering Heterogeneity and Significant Inertial Effects [J]. Engineering Geology，2022，309： 106849.

[11]Barenblatt G I， Zheltov I P，Kochina I N. Basic Concepts in the Theory of Seepage of Homogeneous Liquids in Fissured Rocks[J]. Journal of Applied Mathematics and Mechanics，1960，24（5）：1286- 1303.

[12]Snow D T. Anisotropie Permeability of Fractured Media[J].Water Resources Research，1969，5（6）： 1273-1289.

[13] Long J C S，Remer J S， Wilson C R，et al. Porous Media Equivalents for Networks of Discontinuous Fractures[J]. Water Resources Research，1982，18 （3）：645-658.

[14]Oda M. An Equivalent Continuum Model for Coupled Stress and Fluid Flow Analysis in Jointed Rock Masses[J]. Water Resources Research，1986，22（13） ： 1845 -1856.

[15］杜延齡，許國(guó)安.復(fù)雜巖基三維滲流有限元分析研究 [J].水利學(xué)報(bào)，1991，35（7）：19－26，33. Du Yanling，Xu Guoan. Multipoint Quasi-Coincidence Method of Inverse Analysis for Arch Dams[J]. Journal of Hydraulic Engineering，1991，35 （7） ：19 26，33.

[16]周志芳，錢孝星，顧長(zhǎng)存.三峽工程壩基揚(yáng)壓力三維有 限元分析[J].河海大學(xué)學(xué)報(bào)，1993，21（5）：94－99. Zhou Zhifang，Qian Xiaoxing，Gu Changcun. 3- Dimensional Finite Element Analysis of Uplift in the Dam Foundation of Three Gorges Project[J]. Journal of Hohai Univesity，1993，21（5）：94-99.

[17]張有天.裂隙巖體滲流數(shù)學(xué)模型研究現(xiàn)況[J].人民長(zhǎng) 江，1991，22（3）：1-10. Zhang Youtian. Research Status of Mathematical Model of Seepage in Fractured Rock Mass [J]. Yangtze River，1991，22（3）：1-10.

[18]田開銘，萬(wàn)力.各向異性裂隙介質(zhì)滲透性的研究與評(píng) 價(jià)[M].北京：學(xué)苑出版社，1989. Tian Kaiming，Wan Li. Research and Evaluation of Permeability of Anisotropic Fractured Media[M]. Beijing： Academy Press，1989.

[19]Nian G Q，Chen Z G，Bao M，et al. Rainfall Infiltration Boundary Conditions and Stability of a FracturedRock Slope Based on a Dual-Continuum Model[J]. Hydrogeology Journal，2022，30（3） ：829 - 847.

[20]Frehner M，Schmalholz S M，Saenger E H，et al. Comparison of Finite Difference and Finite Element Methods for Simulating Two-Dimensional Scattering of Elastic Waves [J]. Physics of the Earth and Planetary Interiors，2008，171（1） ：112-121.

[21]Ritz E，Mutlu O，Pollard D D. Integrating Complementarityintothe2D Displacement Discontinuity Boundary Element Method to Model Faults andFractureswith Frictional Contact Properties[J]. Computers amp; Geosciences，2012，45： 304-312.

[22]Sladek J，Sladek V，Schanz M. A Meshless Method for Axisymmetric Problems in Continuously Computers and Geotechnics，2014，62：100 -109.

[23]Starnoni M，Pokrajac D，Neilson J E. Computation of Fluid Flow and Pore-Space Properties Estimation on Micro-CT Images of Rock Samples[J]. Computers amp; Geosciences，2017，106：118-129.

[24]李博南，劉財(cái)，郭智奇.基于等效介質(zhì)模型和頻變 AVO反演的裂縫儲(chǔ)層參數(shù)估算方法[J].吉林大學(xué)學(xué) 報(bào)（地球科學(xué)版），2017，47（1）：234－244. Li Bonan，Liu Cai，Guo Zhiqi. Estimation of Fractured Reservoir Parameters Based on Equivalent MediaModelandFrequncy-DependentAVO Inversion[J]. Journal of Jilin University （Earth Science Edition），2017，47 （1） ：234- 244.

[25] Zhou L，Su X P，Lu Y Y，et al. A New Three Dimensional NumericalModelBasedonthe Equivalent Continuum Method to Simulate Hydraulic Fracture Propagation in an Underground Coal Mine [J]. Rock Mechanics and Rock Engineering，2019，52 （8）：2871-2887.

[26] Zhang K，Zhang X M，Zhang L M，et al. Inversion of Fractures Based on Equivalent Continuous Medium Model of Fractured Reservoirs [J]. Journal of Petroleum Science and Engineering，2017，151： 496 - 506.

[27]王笑雨，張可霓，李毅.裂隙巖體滲流及污染物遷移模 型的數(shù)值模擬研究[J].北京師范大學(xué)學(xué)報(bào)（自然科學(xué) 版），2015，51（5）：527-532. Wang Xiaoyu， Zhang Keni， Li Yi. Modeling Fluid Flow and Tracer Transport in Fractured Rock by the Effective Continum Method[J]. Journal of Beijing Normal University （Natural Science），2O15，51（5）： 527 -532.

[28]Koyama T F，Jing L R. Effects of Model Scale and Particle Size on Micro-Mechanical Propertiesand Failure Processes of Rocks：A Particle Mechanics Approach[J].Engineering Analysis with Boundary Elements，2007，31（5） ：458- 472.

[29] Min K B，Rutqvist J，Tsang C F，et al. StressDependent Permeability of Fractured Rock Masses ： A Numerical Study[J]. International Journal of Rock Mechanics and Mining Sciences，2004，41（7） ： 1191 1210.

[30] Nordqvist A W，Tsang Y W，Tsang C F，et al.A Variable Aperture Fracture Network Model for Flow and Transport in Fractured Rocks [J]. Water Resources Research，1992，28（6） ：1703-1713. 法[J].巖石力學(xué)與工程學(xué)報(bào)，2010，29（4）：720－730. ZhangQihua， WuAiqing. Three-Dimensional Arbitrary Fracture Network Seepage Model and Its Solution[J]. Chinese Journal of Rock Mechanics and Engineering，2010，29（4）：720-730.

[32]張雨晴，王暉，范廷恩，等.花崗巖潛山儲(chǔ)層裂縫建模 表征方法：以渤?；◢弾r潛山 A油田為例.[J].吉林 大學(xué)學(xué)報(bào)（地球科學(xué)版），2016，46（5）：1312－1320. Zhang Yuqing，Wang Hui，F(xiàn)an Ting'en，et al.Granite Buried Hill Reservoir Characterization and Modeling ： Taking Offshore a Oilfield in Bohai Bay as an Example[J]. Journal of Jilin University （Earth Science Edition），2016，46（5）：1312-1320.

[33] Li X X，Li D Q，Xu Y，et al. A DFN Based 3D NumericalApproachforModelingCoupled Groundwater Flow and Solute Transport in Fractured Rock Mass[J]. International Journal of Heat and Mass Transfer，2020，149：119179.

[34] Karimzade E，Sharifzadeh M，Zarei H R，et al. Prediction ofWater Inflow into Underground Excavations in Fractured Rocks Using a 3D Discrete Fracture Network （DFN） Model[J]. Arabian Journal of Geosciences，2017，10（9） ：206.

[35]Baecher G B，Lanney N A. Trace Length Biases in Joint Surveys[C]//Proc of 19th U. S. Symposium on Rock Mechanics，AIME. New York：[s. n.]，1978： 56-65.

[36]Gl?ser D，F(xiàn)lemisch B，Class H，et al. Frackit：A Framework for Stochastic Fracture Network Generation and Analysis[J]. The Journal of Open Source Software，2020，5（56）： 2291.

[37]Ghanbarian B，Perfect E，Liu H H.A Geometrical Aperture -Width Relationship for Rock Fractures[J]. Fractals，2019，27（1）：1940002.

[38］Olson J E. Sublinear Scaling of Fracture Aperture Versus Length：An Exception or the Rule？ [J]. Journal of Geophysical Research： Solid Earth，2003， 108（B9） ：2413.

[39]Gottron D，Henk A. Upscaling of Fractured Rock Mass Properties： An Example Comparing Discrete Fracture Network （DFN） Modeling and Empirical RelationsBasedon EngineeringRockMass Classifications[J]. Engineering Geology，2021， 294： 106382.

[40]Oda M. Permeability Tensor for Discontinuous Rock Masses[J]. Geotechnique，1985，35（4） ：483- 495.