# A Complex Network Model for Analysis of Fractured Rock Permeability

Authors: Xinmin Zhu; Guannan Liu; Feng Gao; Dayu Ye; Jixi Luo
Journal: Advances in Civil Engineering
 (2020)

---

## Abstract

Fractured rocks exist widely in nature. The fracture network is an effective storage space and main seepage channel of low-permeability oil and gas reservoirs, which controls the seepage system of low-permeability oil and gas reservoirs. The connection characteristics of fracture networks are complex and evolve dynamically with time. The rise of complex network research can provide reliable analysis for the relationship between network structures and network behaviors. In this work, the fracture network is considered as a hierarchical network with self-similarity, and complex network theory is applied to analyze the permeability of fractured rocks. According to the power-law relationship of degree distribution of network nodes, the number of nodes is corresponding to the number of network edges and a new power-law distribution relationship of edges with degree of nodes is proposed. Eventually, the permeability model of fractured rocks is derived and it is found that permeability of fractured rocks is a function of degree of maximum nodekmax, self-similarity index γ, power index dk, and other structural parameters. Compared with the existing numerical simulations, the validity of the model is verified. By calculating the influence of model parameters on the permeability, the following results are obtained: (1) fracture porosity is directly proportional to permeability; (2) fracture surface density is linearly increasing with permeability; (3) power index is inversely proportional to permeability; and (4) permeability is exponentially increasing with the maximum degree of a node.


---

## Body

## 1. Introduction

The networks of rock fracture are formed by structural deformation and physical diagenesis [1]. On the inner surface of rock stratum, the scale of naturally formed fracture network is expanding and redistributing randomly with different fracture development degrees, which is always difficult to identify. However, researchers usually use the dip angle and azimuth to determine the spatial orientation of fractures. The structural characteristics are consistent with the two directional attributes of structural geology: tendency and trend. In geological objects, the complex trace analysis is used to calculate the obliquity estimation of the three-dimensional data body, i.e., steering cube [2], and then gets the dip angle and azimuth information of each data point. The permeability of fractured reservoirs is very low, and the fracture network controls the fluid flow [3]. Hence, it has an important influence on oil or gas exploitation [4] and geothermal energy extraction [5].In recent years, researchers around the world have studied permeability characteristics of fracture networks and put forward corresponding models [6–8]. Snow [9] in his study established the parallel plate model and obtained a tensor analytical formula of permeability of fracture network. Koudina et al. [10] studied the permeability of fracture network in three-dimensional space by means of numerical simulation. The fracture network was composed of polygonal shapes and the flow of fluid in each fracture satisfied Darcy’s law, while comparing it with Snow’s model. Xia [11] established the dynamic model of permeability and opening of fracture network under different confining pressures. Van Stappen et al. [12] also connected the seepage model with fracture opening by determining the relationship between fracture permeability and confining pressure. Li et al. [13, 14] broke away from the traditional practice of thinking fractured reservoirs as dual media and established a percolation model with equivalent continuum suitable for low-permeability fractured shale reservoirs. De Dreuzy et al. [15] studied the permeability of randomly generated two-dimensional fracture network by numerical and theoretical methods and compared it with natural fractures to verify the accuracy of the model. Klimczak et al. [16] used the parallel plate model to obtain the permeability formula of a single crack under conditions that the fracture length and opening satisfied the power-law relationship and verified the accuracy of the model through numerical simulations. Wei et al. [17] derived a forecasting model of permeability using the electrokinetic relationship between fluid flow and current in microfractures and analyzed the influence of connectivity between fractures on permeability. Li [18] proposed a new model considering fracture connectivity according to the hydraulic fracture morphology of raw coal, the “matchstick” seepage model and the cubic law. However, the above models do not quantitatively relate the permeability of fracture network with porosity, surface density, and microstructural parameters of fracture network, such as fracture connectivity, openness, the dip angle, and azimuth.The randomly distributed fracture networks in rocks have been shown to have statistical self-similarity, which is a basic feature of fractal. Interested readers may consult [19–26] for details. Watanabe and Takahashi [5] used the fractal theory to study the permeability of fracture networks and the extraction of heat in dry hot rocks, but they did not put forward a permeability expression with micro parameters. Yu et al. [27] based on the study of seepage characteristics of porous media in fracture networks by using fractal methods put forward an explicit expression with micro parameters, such as the structure of fracture network and porosity, and then gave the scaling relationship between permeability and the structure of fracture network. Li et al. [28] established the mathematical model of equivalent permeability tensors in fractured reservoirs, based on fracture statistics, the simulation technique of fracture network, and equivalent flow assumption, and then obtained the equivalent permeability tensor of fractured media by using boundary element method. Jafari and Babadagli [29] obtained fractal permeability expressions of random fractures by using multiple regression analysis based on logging data but their empirical relationship contained many empirical constants. Recently, Miao et al. [6] obtained the analytical expression of fracture network permeability according to the basic fractal theory. This model quantitatively connected the fracture length, aperture, the fracture dip angle, and fracture azimuth with permeability of fractured rocks, which did not include any empirical constant.Most of the above models initiate from the statistical parameters of isolated fractures and macroscopic homogenization. The connectivity of fracture networks is not considered, particularly the influence of the connectivity of a small number of local fractures (maximum degree) on the overall permeability. Starting from the topological structure of fracture network and based on the complex network theory, this paper establishes the network permeability model of fractured rocks and probes into the internal mechanism of the influence of structure parameters of fracture network on permeability, including fracture porosity∅M, fracture density D, power index dk, and the maximum node degree kmax.


## 2. Self-Similarity of Complex Networks

### 2.1. Degree Distribution of Hierarchical Networks

In order to illustrate the modularity, local clustering, and scale-free topological characteristics of many complex network systems, it is necessary to assume that the modules generate a hierarchical network in some iterative way [30]. Recent studies show that [31] some of the topology modules are well organized hierarchically in the network. The hierarchical network seems to have a very conspicuous feature; that is, the local is similar to the whole in a sense, i.e., the self-similarity. Hierarchical network integrates scale-free topology with internal module structure. Song et al. [32] further reveal that self-similarity and degree distribution of scale-free hold true at all coarse-grained stages of the network by adopting renormalization procedure and the degree distribution Pk of the renormalized network is invariant under renormalization. The power-law relationship can be expressed as follows [32]:(1)Pk∝k−γ,kmin≤k≤kmax,where P represents the total number of nodes in the network with degree k, k represents the number of other nodes connected to a node, and γ represents the self-similarity index with the range of 1–3, which is transformed by the exponential formula [32].


### 2.2. Basic Features of Fractals

Most of the trace length of the fracture satisfies the power-law (scale-free) distribution [33, 34]. The fractal power-law distribution refers to the fact that the fracture length in nature is random and disordered, showing the characteristics of similarity and fractal. The power-law expression is [35](2)N≥l∝l−Df,where Df is the fractal dimension of fracture length, l is the track length of fracture, and N is the total number of fractures. This is the basic expression of fractal scaling law and basis of box counting method.


### 2.3. Power-Law Expression of Fracture Complex Networks

Covariant analogy is also known as mathematical similarity analogy. Power-law relations (1) and (2) have obviously similar functional relations and equation (1) multiplied by k can be analogous to equation (2). The number of edges of complex network is associated with the number of edges of fracture network. The following expression is obtained:(3)l−Df=k1−γ.The relationship betweenγ in power-law distribution formula (1) and Df in power-law expression (2) is as follows [32]:(4)γ=1+Dfdk,where the power exponent dk is 1.5 times of Df [32]. Substitute equation (4) into equation (3) to get a proportional relationship:(5)l=k1/dk.The parallel plate model is usually used to represent the effective aperture of the fracture and the relationship between crack length and effective aperture has also been studied by a large number of researchers [36, 37]. This relationship is given by(6)a=βln,where β is the proportionality coefficient, which is related to the mechanical properties of the medium around the fracture in the range of 10−3∼10−1 [16]. a is the effective aperture of the fracture and n is the power exponent.When the power exponentn = 1, the fracture network has the characteristics of self-similarity and fractal [37]. So, for the fracture network with self-similarity [16], equation (6) can be rewritten as(7)a=βl.Equation (1) can be rewritten as(8)Mk=kPk=αk1−γ,where M represents the number of network edges, and α is the proportionality coefficient.Differentiating equation (8), we can get the number of edges whose node degrees are in the range k to k+dk:(9)−dMk=−kdPk=αγk−γdk,wherein the negative sign indicates that the number of edges of a complex network decreases with the increase of node degree, which is in line with the actual situation and −dMk>0.The probability density of an edge with node degreesk is expressed as(10)−dMkMt=αMtγk−γdk,where Mt represents the total number of edges the network has, and fk=α/Mtγk−γ is the probability density function of the edge with node degrees k, which satisfies the normalization principle:(11)∫kminkmaxfkdk=∫kminkmaxαMtγk−γdk=1.Thus, it can be obtained that(12)αMtγγ−11kminγ−11−kminkmaxγ−1=1.Evidently, whenkmin≪kmax, equation (12) can be expressed as(13)α=Mtkminγ−11−1γ,generally, kmin/kmax≤10−2 can be taken and complex networks in nature usually meet this requirement.Yu [38, 39] studied the power-law relation of fractal distribution of pores in porous media. Likewise, Majumdar and Bhushan gave the cumulative size distribution of islands on the Earth’s surface [40]:(14a)N>s=smaxsD/2,where N is the total number of islands with the largest area smax greater than s, and D is the fractal dimension for the size distribution of islands. Equation (14a) indicates that there is the largest island on the Earth’s surface; in addition, Majumdar and Bhushan [40] used this power-law formula to describe the contact points on engineering surfaces, where the maximum point area smax=gλmax2, a point area s=gλ2 with λ being the diameter of a point and g being a geometry coefficient.Since self-similarity is one of the basic characteristics of fractal, the self-similarity of porous media with fractures needs to satisfy a certain power-law relationship [41]. Hence, equation (14a) is used to describe islands on the Earth’s surface and points on the engineering surface can be extended to describe the size distribution of nodes on the surface of a fracture network. In the complex network theory, the characteristic size of a single node includes out-degree and in-degree [30]:(14b)N>s=komaxkimaxkokiγ/2,where komaxkimax represents the maximum node size with komax and kimax, respectively, being the maximum out-degree and maximum in-degree, koki is a node size with the out-degree and in-degree being ko and ki, respectively. When the direction of the degree is ignored, equation (14b) can be simplified as(14c)N>s=kmaxkmaxkkγ/2.From equation (14c), the cumulative number of nodes whose degrees are greater than k can be expressed as(14d)N>k=kmaxkγ,where N is the cumulative number of nodes in a fracture network.From equation (14d), the total number of nodes in a complex fracture network is obtained:(15)Nt=kmaxkminγ.Because the contribution of one edge to degree is 2, the average degree of complex networks is(16)kav=2MtNt.Inserting equation (15) into equation (16),(17)Mt=12kmaxkminγkav.Inserting equation (17) into equation (13) to get the proportionality coefficient,(18)α=12kmaxγkmin1−1γkav.Then, we insert equation (18) into equation (9) which gives(19)−dMk=−kdPk=12γ−1kmaxγkminkavk−γdk.Equation (19) is an important power-law distribution relation of edges with certain node degree in complex networks. Furthermore, by the same logic, the average degree of complex networks can be obtained as(20)kav=γ2γ−1kmin1−kminkmax1−γ.


### 2.4. Surface Porosity of Fracture Networks

Self-similarity is closely related to the fractal. Yu and Li [42] deduced the relationship between porosity and fractal dimension in porous media based on the fractal theory:(21)∅m=λminλmaxdE−Dp,where λmin and λmax are, respectively, the minimum pore diameter and the maximum pore diameter. Dp is the fractal dimension of pores. dE is the Euclidean dimension: in two dimensions, dE = 2; in three dimensions, dE = 3.Equation (21) is appropriate not only for precise fractal geometry but also for statistical fractal geometry. As long as the pores of porous media fall within the self-similar range of λmin∼λmax, forming a fractal set, equation (21) holds accurate regardless of the shape of the pores. Therefore, in a hierarchical complex network, it is embedded into the matrix as a fracture, forming a network model with fracture properties. The edges of the complex network, that is, the fractures with a node, satisfy the above equation (21) within the self-similar range of kmin∼kmax and are independent of the shape of the node. It can be rewritten as(22)γ=dE−ln∅Mlnkmin/kmax,where ∅M is the effective porosity of fractures in the rock and kmin and kmax are the minimum and maximum of the nodes, respectively.On the cross section of the representative elementary volume, the surface porosity of the fracture network is defined as [6](23)∅M=APMAM,where AM represents the cross-sectional area of the representative elementary volume in which the fracture network is located, and APM represents the total area of fracture pores on this area.According to equations (5), (7), and (19), we can get the total cross section area of the fracture [6]:(24)APM=−∫kminkmaxa⋅l⋅dMk=12γ−1βkavkminkmax2/dk+12/dk−γ+11−kminkmax2/dk−γ+1.Inserting equation (22) into equation (24),(25)APM=12γ−1βkavkminkmax2/dk+12/dk−γ+11−∅M1/dk,where the porosity ∅M is used in the two-dimensional space of equation (22), i.e., dE = 2.


## 2.1. Degree Distribution of Hierarchical Networks

In order to illustrate the modularity, local clustering, and scale-free topological characteristics of many complex network systems, it is necessary to assume that the modules generate a hierarchical network in some iterative way [30]. Recent studies show that [31] some of the topology modules are well organized hierarchically in the network. The hierarchical network seems to have a very conspicuous feature; that is, the local is similar to the whole in a sense, i.e., the self-similarity. Hierarchical network integrates scale-free topology with internal module structure. Song et al. [32] further reveal that self-similarity and degree distribution of scale-free hold true at all coarse-grained stages of the network by adopting renormalization procedure and the degree distribution Pk of the renormalized network is invariant under renormalization. The power-law relationship can be expressed as follows [32]:(1)Pk∝k−γ,kmin≤k≤kmax,where P represents the total number of nodes in the network with degree k, k represents the number of other nodes connected to a node, and γ represents the self-similarity index with the range of 1–3, which is transformed by the exponential formula [32].


## 2.2. Basic Features of Fractals

Most of the trace length of the fracture satisfies the power-law (scale-free) distribution [33, 34]. The fractal power-law distribution refers to the fact that the fracture length in nature is random and disordered, showing the characteristics of similarity and fractal. The power-law expression is [35](2)N≥l∝l−Df,where Df is the fractal dimension of fracture length, l is the track length of fracture, and N is the total number of fractures. This is the basic expression of fractal scaling law and basis of box counting method.


## 2.3. Power-Law Expression of Fracture Complex Networks

Covariant analogy is also known as mathematical similarity analogy. Power-law relations (1) and (2) have obviously similar functional relations and equation (1) multiplied by k can be analogous to equation (2). The number of edges of complex network is associated with the number of edges of fracture network. The following expression is obtained:(3)l−Df=k1−γ.The relationship betweenγ in power-law distribution formula (1) and Df in power-law expression (2) is as follows [32]:(4)γ=1+Dfdk,where the power exponent dk is 1.5 times of Df [32]. Substitute equation (4) into equation (3) to get a proportional relationship:(5)l=k1/dk.The parallel plate model is usually used to represent the effective aperture of the fracture and the relationship between crack length and effective aperture has also been studied by a large number of researchers [36, 37]. This relationship is given by(6)a=βln,where β is the proportionality coefficient, which is related to the mechanical properties of the medium around the fracture in the range of 10−3∼10−1 [16]. a is the effective aperture of the fracture and n is the power exponent.When the power exponentn = 1, the fracture network has the characteristics of self-similarity and fractal [37]. So, for the fracture network with self-similarity [16], equation (6) can be rewritten as(7)a=βl.Equation (1) can be rewritten as(8)Mk=kPk=αk1−γ,where M represents the number of network edges, and α is the proportionality coefficient.Differentiating equation (8), we can get the number of edges whose node degrees are in the range k to k+dk:(9)−dMk=−kdPk=αγk−γdk,wherein the negative sign indicates that the number of edges of a complex network decreases with the increase of node degree, which is in line with the actual situation and −dMk>0.The probability density of an edge with node degreesk is expressed as(10)−dMkMt=αMtγk−γdk,where Mt represents the total number of edges the network has, and fk=α/Mtγk−γ is the probability density function of the edge with node degrees k, which satisfies the normalization principle:(11)∫kminkmaxfkdk=∫kminkmaxαMtγk−γdk=1.Thus, it can be obtained that(12)αMtγγ−11kminγ−11−kminkmaxγ−1=1.Evidently, whenkmin≪kmax, equation (12) can be expressed as(13)α=Mtkminγ−11−1γ,generally, kmin/kmax≤10−2 can be taken and complex networks in nature usually meet this requirement.Yu [38, 39] studied the power-law relation of fractal distribution of pores in porous media. Likewise, Majumdar and Bhushan gave the cumulative size distribution of islands on the Earth’s surface [40]:(14a)N>s=smaxsD/2,where N is the total number of islands with the largest area smax greater than s, and D is the fractal dimension for the size distribution of islands. Equation (14a) indicates that there is the largest island on the Earth’s surface; in addition, Majumdar and Bhushan [40] used this power-law formula to describe the contact points on engineering surfaces, where the maximum point area smax=gλmax2, a point area s=gλ2 with λ being the diameter of a point and g being a geometry coefficient.Since self-similarity is one of the basic characteristics of fractal, the self-similarity of porous media with fractures needs to satisfy a certain power-law relationship [41]. Hence, equation (14a) is used to describe islands on the Earth’s surface and points on the engineering surface can be extended to describe the size distribution of nodes on the surface of a fracture network. In the complex network theory, the characteristic size of a single node includes out-degree and in-degree [30]:(14b)N>s=komaxkimaxkokiγ/2,where komaxkimax represents the maximum node size with komax and kimax, respectively, being the maximum out-degree and maximum in-degree, koki is a node size with the out-degree and in-degree being ko and ki, respectively. When the direction of the degree is ignored, equation (14b) can be simplified as(14c)N>s=kmaxkmaxkkγ/2.From equation (14c), the cumulative number of nodes whose degrees are greater than k can be expressed as(14d)N>k=kmaxkγ,where N is the cumulative number of nodes in a fracture network.From equation (14d), the total number of nodes in a complex fracture network is obtained:(15)Nt=kmaxkminγ.Because the contribution of one edge to degree is 2, the average degree of complex networks is(16)kav=2MtNt.Inserting equation (15) into equation (16),(17)Mt=12kmaxkminγkav.Inserting equation (17) into equation (13) to get the proportionality coefficient,(18)α=12kmaxγkmin1−1γkav.Then, we insert equation (18) into equation (9) which gives(19)−dMk=−kdPk=12γ−1kmaxγkminkavk−γdk.Equation (19) is an important power-law distribution relation of edges with certain node degree in complex networks. Furthermore, by the same logic, the average degree of complex networks can be obtained as(20)kav=γ2γ−1kmin1−kminkmax1−γ.


## 2.4. Surface Porosity of Fracture Networks

Self-similarity is closely related to the fractal. Yu and Li [42] deduced the relationship between porosity and fractal dimension in porous media based on the fractal theory:(21)∅m=λminλmaxdE−Dp,where λmin and λmax are, respectively, the minimum pore diameter and the maximum pore diameter. Dp is the fractal dimension of pores. dE is the Euclidean dimension: in two dimensions, dE = 2; in three dimensions, dE = 3.Equation (21) is appropriate not only for precise fractal geometry but also for statistical fractal geometry. As long as the pores of porous media fall within the self-similar range of λmin∼λmax, forming a fractal set, equation (21) holds accurate regardless of the shape of the pores. Therefore, in a hierarchical complex network, it is embedded into the matrix as a fracture, forming a network model with fracture properties. The edges of the complex network, that is, the fractures with a node, satisfy the above equation (21) within the self-similar range of kmin∼kmax and are independent of the shape of the node. It can be rewritten as(22)γ=dE−ln∅Mlnkmin/kmax,where ∅M is the effective porosity of fractures in the rock and kmin and kmax are the minimum and maximum of the nodes, respectively.On the cross section of the representative elementary volume, the surface porosity of the fracture network is defined as [6](23)∅M=APMAM,where AM represents the cross-sectional area of the representative elementary volume in which the fracture network is located, and APM represents the total area of fracture pores on this area.According to equations (5), (7), and (19), we can get the total cross section area of the fracture [6]:(24)APM=−∫kminkmaxa⋅l⋅dMk=12γ−1βkavkminkmax2/dk+12/dk−γ+11−kminkmax2/dk−γ+1.Inserting equation (22) into equation (24),(25)APM=12γ−1βkavkminkmax2/dk+12/dk−γ+11−∅M1/dk,where the porosity ∅M is used in the two-dimensional space of equation (22), i.e., dE = 2.


## 3. Relationship between Surface Density and the Self-Similarity Index

According to equation (5), the total length of fractures on the cross section of the representative elementary volume is as follows:(26)Ltotal=−∫kminkmaxl⋅dMk=12γ−1kavkminkmax1/dk+11/dk−γ+11−kminkmax1/dk−γ+1.Inserting equation (22) into equation (26),(27)Ltotal=12γ−1kavkminkmax1/dk+11/dk−γ+11−∅M1/dk−γ+1/2/dk−γ+1.In the two-dimensional fracture network, the surface density refers to the density of the cross section in a unit cell of fractures (not a single fracture), which is defined by [43](28)D=LtotalAM,where D is the surface density of fractures, and Ltotal is the total length of all fissures on the cross section of the representative elementary volume body, which is related to the complex network model.Equations (23), (25), and (27) are inserted into equation (28) to get the surface density of the fracture, i.e.,(29)D=2/dk−γ+11−∅M1/dk−γ+1/2/dk−γ+1∅M1/dk−γ+11−∅M1/dkβkmax1/dk.Equation (29) shows that the surface density in the two-dimensional complex fracture network is a function of the fracture porosity ∅M, self-similarity index γ, power index dk related to fractional dimension, proportion coefficient β, and the degree kmax of the largest node.In order to study the relationship between the surface density of these four fracture networks and the self-similarity index, the prediction results of the surface density for the complex fracture network are compared with the four random fracture networks generated by Zhang and Sanderson [43] through the numerical method of self-avoiding walking. In their simulation, the critical fractal dimensions lie in a narrow range from 1.22 to 1.38 (average 1.30) for those critical clusters with variations in the lower limit of length from 0.005 to 1.5 m, in the dispersion angle of fracture direction from 0 to 50°, and in exponents from 1.2 to 1.8. Therefore, during calculation, the minimum degree of the node is taken to be 1 and the maximum degree is taken to be 300 according to the same ratio coefficient. Meanwhile, the average power index is 2 and the average self-similarity index is 1.65 through equation (4), and the average porosity ∅M is calculated through equation (22). From Figure 1, it can be observed that the predicted results are in good agreement with numerical simulations. Meanwhile, Figure 1 shows that the surface density of fracture network increases with the increase of fracture self-similarity index.Figure 1
A comparison between theoretical model prediction and numerical simulation of surface density of complex fracture network [43] β=0.011. Groups A–D represent the four-component fracture network generated in the numerical simulation.Figure2 shows the relationship between fracture surface density and porosity when the maximum degree of node kmax is 300 and β is 0.006. It can be observed from Figure 2 that the surface density of fracture network increases with the increase of fracture porosity. This is because the larger is the porosity, the greater will be the pore area of the fracture network. Under certain conditions of β, the longer is the total length of the fracture, the stronger the connectivity will be as mentioned. This result is consistent with the simulations by Miao et al. [6]. It can be explained that the change of numerical values will not affect the general trend between them.Figure 2
The density of fracture surface varies with the porosity of fracture network.


## 4. Complex Network Model of Permeability of Fractured Rocks

Generally, the production of low-permeability reservoirs often depends on the seepage system of fracture network. When there are differences of temperature and pressure in the system, there will be fluid flow or heat transfer between the fracture networks. In these processes, however, the laws of mass, momentum, and energy transfer among fluids are very complicated. Moreover, the geometric aspects of fractures cannot be determined relatively, including density and surface roughness. For fractured reservoirs, we can proceed from macroheterogeneity, because the degree distribution, length, aperture, and orientation of the fracture are often random and disordered. Complex network can provide an effective method for representing irregular objects.The topological model of complex networks normally considers the positional relationship between nodes but not their shape and size. Therefore, network space determined by the azimuth and dip angle has an important influence on the seepage characteristics of fracture network. Nevertheless, the spatial orientation of fractures is usually random and the number of fractures in space is so large that it is almost impossible to express the orientation of each fracture precisely [44]. Generally, the statistical method in the field of engineering has been adopted to show the location of fracture network, which is to take the average value of fracture dip angle and fracture azimuth [45], and this is often used in petroleum engineering, shale gas exploitation, and geothermal energy extraction. Therefore, in this paper, we assume that the average dip angle of the complex fracture network is θ and the average azimuth of the fracture is α, as shown in Figure 3.Figure 3
The average orientation of fractures in the three-dimensional space, the plane of the coordinate axis is the horizontal plane and the direction of water flow is along thex-axis. The included angle between the fracture direction and y-axis is α, that is, the azimuth of the fracture. The θ angle between the fracture plane and the horizontal plane is the dip angle of the fracture.The cubic law of single fracture is based on the model of parallel plate, which becomes the basic theory of network seepage of fractured rocks and this is usually considered simple and effective. The flow rate along the flow direction through a fracture can be described by the famous cubic law [46, 47]:(30)ql=a3l12μΔpL0,where L0 denotes the length of the representative elementary volume, a denotes the fracture aperture, l denotes fracture trace length, Δp denotes the pressure drop across a fracture along flow direction, and μ denotes dynamic viscosity coefficient of the fluid.If the spatial orientation of fracture is considered, the flow rate of a single fracture can be expressed as [6, 48](31)ql=a3l1−cos2αsin2θ12μΔpL0.The total flow rate of fluid through a set of complex fracture networks can be obtained by integrating equation (31) from minimum degree to maximum degree in a unit cross section; i.e.,(32)Q=−∫kminkmaxqldMk=β312μγ−121−cos2αsin2θ4/dk−γ+1ΔpL0kavkminkmax4/dk+11−kminkmax4/dk−γ+1.In general,kmin≪kmax. According to equation (4) and [35], since 1<γ<2.3 in the two-dimensional plane, and kmin/kmax4/dk−γ+1≪1, consequently, equation (32) can be simplified as(33)Q=β312μγ−121−cos2αsin2θ4/dk−γ+1ΔpL0kavkminkmax4/dk+1.It can be seen from equation (33) that the total flow rate of fluid in the complex fracture network is related to the index of self-similarity γ, fractional-dimension related power index dk, fracture azimuth α, and fracture dip angle θ and the flow rate is very sensitive to the maximum degree kmax of nodes.Darcy’s law for Newtonian fluid flow in porous media is given by [6](34)Q=KAMμΔpL0.The permeability of complex fracture network can be obtained by inserting equation (33) into equation (34):(35)K=β324AMγ−11−cos2αsin2θkav4/dk−γ+1kminkmax4/dk+1.By inserting equation (28) and equation (29) into equation (35), the permeability of fracture network can be expressed by the surface density of fracture:(36)K=β3⋅D121/dk−γ+14/dk−γ+1kmax3/dk1−cos2αsin2θ1−∅M1/dk−γ+1/2/dk−γ+1.Equation (36) suggests that permeability is a function of self-similarity index γ, power index dk, structural parameters (maximum degree of a node kmax, fracture surface density D, fracture azimuth α, and fracture dip angle θ), and fracture porosity ∅M in a medium formed by a complex fracture network. Equation (36) further reveals that permeability is strongly dependent on the maximum degree kmax of the node. The higher is the node degree, the stronger is the connectivity of fracture network. The fluid capacity increases with increase of the flow path, leading to higher permeability. Therefore, this model has more advantages than the traditional model and can better explain the influence of node failure on fluid flow in the fracture network.


## 5. Results and Discussion

Jafari and Babadagli [49] analyzed 22 different fracture networks in nature. The digitized fracture patterns were exported to commercial fracture modeling software (FRACA) to calculate their equivalent fracture network permeability. A 3D model with a grid block size of 100 m × 100 m × 10 m was constructed. Each digitized 2D fracture pattern (i.e., the digitized mapped fracture traces from outcrops) was imported into the 3D model in such a way that all fractures were considered to be vertically touching the top and the bottom of the layer, wherein the maximum fracture length is 2 m and the dip angle of the fracture is 0°. Therefore, in the calculation, the minimum degree of the node is 1, and the maximum degree is 6. Furthermore, since the model of parallel plate mainly depends on the effective aperture of a single fracture, the actual tortuosity of the fracture is not considered by using this simplified model. Via equation (4) and equation (20), the average power index and the average degree of the node are calculated. All the structural parameters used in theoretical calculations are listed in Table 1. Figure 4 shows that the predicted values of our model are in good agreement with the results of numerical simulations.Table 1
Structural parameters of fractured rock media.
Structural parametersValueDescription of parameterskmin1Minimum degreekmax6Maximum degreedk2.6Average power indexkav4Average degree of fracture networkθ (°)0Average dip angle of fracture networkα (°)—Azimuth average value of fracture networkβ0.001Ratio of fracture aperture to lengthAM (m2)10000Cross-sectional area of the 3D modelFigure 4
Comparison of permeability between model prediction and numerical simulation results of complex fracture network [49] (β = 0.001).We discuss the influence of model parameters on permeability. From equation (36), it is observed that the parameters that play a decisive role mainly include fracture porosity ∅M, fracture dip angle θ, fracture surface density D, power index dk, and maximum node degree kmax. Figure 5 shows the relationship between permeability and fracture porosity of the complex network model at different dip angles. In the calculation, the maximum degree kmax = 398 of the fracture node (at β = 0.006) is taken. It can be seen from Figure 5 that the permeability of fracture network increases with the increase of fracture porosity. In addition, with the same porosity, the larger is the fracture dip angle, the smaller is the permeability of fracture network. This is because the flow resistance of fluid increases with the increase of the fracture dip angle.Figure 5
The relationship of fracture permeability with porosity at different dip angles.Figure6 shows the relationship between permeability and fracture surface density in the complex network model. In this calculation, the maximum degree kmax = 398 of fracture node, fracture dip angle θ = 45°, fracture azimuth α = 0°, and β = 0.006 are taken. It can be seen from Figure 6 that the permeability of fracture network increases with the increase of fracture surface density. This is because as the density of fracture surface increases, the porosity of fracture also increases. Hence, the permeability of fracture network increases.Figure 6
The relationship between fracture permeability and fracture surface density.Figure7 shows the relationship between permeability of complex network model and the power index. In the calculation, the minimum and maximum degrees of fracture network nodes as 1 and 6 are taken, respectively. Equations (4) and (20) are used to calculate the average self-similarity index γav = 1.67 and the average degree of node kav = 4 and we take the dip angle of fracture θ = 0° and β = 0.001. It can be seen from Figure 7 that the permeability of fracture network decreases slowly with the increase of power index. Miao et al. [6] have verified that the permeability of the fractal fracture network model increases slowly with the increase of fractal dimension. By considering equation (4), that is, the internal correlation between scale-free property of complex networks and fractal scaling law, it can be concluded that there will be a competitive relationship between the inhibition of seepage flow by power index and the promotion of seepage flow by fractal dimension. Henceforth, it leads to the discontinuous phase that is not always occurring.Figure 7
The relationship between fracture permeability and power index.Figure8 shows the relationship between the permeability of complex network model and the maximum node degree. When a node in a network has multiple edges connected to it, the number of edges is the degree of the node, regardless of its direction. In the calculation, we take the dip angle (from 0° to 180°) of fracture θ = 45°, azimuth of the fracture α = 0°, β = 0.006 with the range of 0.001∼0.1 [16], and average surface density D = 10 (m/m2). It can be seen from Figure 8 that permeability of fracture network increases sharply with the increase of the maximum degree of nodes. Since connectivity of the entire fracture network is strongly dependent on maximum degree of a node, it is equivalent to the connection hub of entire complex network. When a small number of edges are removed from the network, the overall connectivity of the network will not be greatly affected. Thus, the complex network has a high robustness to the node destruction. At the same time, if a node with the maximum degree is deliberately attacked, the entire network will paralyze quickly and the fluid can only flow through a few paths. This is also the vulnerability of complex network to deliberate attacks on nodes.Figure 8
The relationship between fracture permeability and maximum node degree.


## 6. Conclusion

This paper applies the complex network theory and topological model to fractured rocks, while describing the fracture network as a hierarchical network with self-similarity. Meanwhile, the fracture network model of surface density is obtained based on the power-law distribution relation of network edges. Then, the permeability model of fractured rocks is deduced in accordance with the famous cubic law, Darcy’s law, and complex network theory.Compared with the existing numerical simulations, the predicted results show that the above models are accurate. Besides, the effect of structural parameters on the permeability of fractured media is also discussed. The permeability of fracture networks increases with the increases of porosity and surface density. The permeability of fracture networks increases exponentially with the increase of the maximum node degree and its power exponent is3/dk.



---

*Source: 8824082-2020-10-30.xml*