A Fractal Model for Predicting the Relative Permeability of Rough-Walled Fractures

The relative permeability and saturation relationships through fractures are fundamental for modeling multiphase ﬂow in underground geological fractured formations. In contrast to the traditional straight capillary model from porous media, the realistic ﬂow paths in rough-walled fractures are tortuous. In this study, a fractal relationship between relative permeability and saturation of rough-walled fractures is proposed associated with the fractal characteristics of tortuous parallel capillary plates, which can be generalized to several existing models. Based on the consideration that the aperture distribution of rough-walled fracture can be represented by Gaussian and lognormal distributions, aperture-based expressions between relative permeability and saturation are explicitly derived. The developed relationships are validated by the experimental observations on Gaussian distributed fractures and numerical results on lognormal distributed fractures, respectively.


Introduction
Multiphase flow through fractured media is an important process for a number of practical applications, such as seepage control for rock slope [1,2], petroleum reservoirs [3], nuclear waste disposal [4], and geologic storage of carbon dioxide [5]. In order to understand the influence of multiphase flow behavior in fractured media on implementation with these applications, numerical simulations have been applied in fractured media using a continuum method [6,7] or a discrete fracture network model [8,9]. In contrast to the macroscale continuum model, the flow behavior in the fractured rock is mainly governed by the presence of fracture, which can highly improve the permeability of rock masses [10][11][12]. For the fracture-dominated flow, the relative permeability of a single fracture is one of the key parameters to describe the multiphase flow behavior that controls the accuracy of modeling results during numerical simulations.
Some well-known models such as the Burdine model [31], Corey model [32], and Brooks and Corey model [33] induced from porous media were also borrowed in fitting the relative permeability properties of fractures, according to the conceptualization of rough-walled fractures as two-dimensional porous media. In these conceptual models from porous media, the hypothesis that the flow paths are straight is in contradiction with the actual tortuous flow paths. e tortuosity factor is generally quantified by a power function of function with exponent parameter m and m � 2 is widely used in Brooks and Corey model [33] while the fitted value m of 1.0 is obtained during Liu et al.'s curve fitting [34]. e physical meaning of the tortuosity factor is still not clear for rough-walled fractures.
When the flow paths are controlled by capillarity during quasistatic displacement, the relative permeability of fracture is substantial influenced by aperture distribution. Some studies [23,24,27,28] have been made to understand the relationships between relative permeability and aperture distribution for rough-walled fractures including numerical approaches for multiphase flow in single fractures, the sensitivity of spatial distribution between apertures on relative permeability predictions, and the determination of empirical parameters related to the residual saturation. However, explicit mathematical relationships between relative permeability and fracture-aperture distribution are still lacking.
A systematic study of relative permeability models with respect to flow tortuosity and its dependence on fractureaperture distribution is a challenging task. erefore, the primary objectives of this study are (1) to propose a generalized model for relative permeability in rough-walled fractures with consideration of the tortuosity fractal dimension, (2) to develop the relationships between relative permeability and aperture distribution of fractures, and (3) to assess the validity of the proposed model through experimental data and numerical results from the literature.

The Existing Relative Permeability Models
e relative permeability (k r ) versus saturation (S w ) relationship through rough-walled fractures has been extensively investigated due to their importance in engineering applications. On the observation of two-phase flow between water and kerosene through artificial parallel-plate fractures, Romm [13] presented the k r -S w relationship firstly as a simple linear function: is linear relationship is also applied to the oil phase, so that the sum of relative permeabilities for water and oil phases was equal to one known as X-curve. It is indicated that each phase kept flowing in its self-governed capillary paths without impediment from the other phase.
According to experimental results from the distilled water-nitrogen gas flow through a transparent replica of a natural rock fracture [14], air-water flow in horizontal artificial fractures [15,18], kerosene-water flow in a parallel glass plate fracture [16], nitrogen-water flow in smooth-and rough-walled glass fractures [21], and gas-water flow in an induced fracture from a cylindrical basalt core [20], strong phase interference at the intermediate saturation was found and the relative permeability did not linearly depend on the saturation. us, several simple models were established or directly borrowed from porous media to describe the nonlinear behavior of relative permeability.
Based on the porous media approach, the rock fractures are conceptualized as two-dimensional connected porous media. Burdine [31] proposed one of the first theoretical models for k r -S w relationship as k r � S e 0 dS we /P 2 where S we and S wr are the effective saturation and residual saturation of the wetting fluid, respectively, and P c is the capillary pressure. Nevertheless, all the capillary tubes in the Burdine model are assumed to be straight along the hydraulic gradient direction. In fact, the flow paths for porous media are quite tortuous. To take into consideration the tortuosity of the capillary tube, an additional tortuosity factor obtained from measured relative permeability data could be approximated as a power function of effective saturation. en, a new version of the modified Burdine model is given as (4) e formulation of the Burdine model includes three parameters S w , S wr , and P c . e relative permeability function is not explicitly expressed within the Burdine model. Once P c as a function of saturation is known, the k r -S w relationship can be derived from equation (4).
Based on the laboratory experiments for a number of soil samples, Corey [32] and Brooks and Corey [33] proposed an empirical relationship between water saturation and capillary pressure in a general form as where P b is the air entry pressure and λ is the pore size distribution index. In Corey [32] relation, λ is equal to two. Inserting equation (5) into equation (4) yields When λ � 2, equation (6) is equivalent to the Corey model: Brooks and Corey [33] have discussed that media with a wide range of pore size distribution should have small values of λ. However, media with a uniform pore size could have λ 2 Advances in Civil Engineering values close to infinity such as the fractured media; equation (6) can be simplified as After the observation on several flow structures (bubbles, unstable bubbles, film flow, etc.) for water-air flows in an artificial fracture [15], Fourar and Lenormand [29] developed a viscous coupling model accounting for the coupling between the two fluids flowing simultaneously, whose interface is assumed to be a plane. Based on the pipe-flow model, the wetting phase is in contact with the capillary wall and the nonwetting phase flows in between. e k r -S w function of the water phase is expressed as As indicated by Huo and Benson [22] who presented a detailed comparison between the experimental data and existing models of the k r -S w relationship, the X-curve could only describe the data from Romm [13]. e Corey model and Brooks and Corey model presented similar results and gave good representations for Chen and Horne's [21] test on randomly rough-walled fractures. e viscous coupling model provided good agreement with the data from Fourar et al. [15] experiment in horizontal artificial fractures and Chen and Horne's [21] on smooth-and homogeneously rough-walled fractures. Some other experimental data [14,20] have a considerable deviation from any of the existing models.
While substantial progress has been made to investigate the relative permeability in rock fractures, there is no general model to predict the k r -S w relationship of a rock fracture. We will develop a reasonable model in this paper and also demonstrate the relationship between the relative permeability and aperture distribution in rough fractures.

Fractal Model of Relative Permeability in Rough-Walled Fractures
A schematic of parallel capillary plates within a single roughwalled fracture is shown in Figure 1. For simplicity, the water-air two-phase flow through rough-walled fracture is focused on in this study. However, the mechanism and results can be extended to other multiphase flow systems. e flow rate q in a parallel capillary plate within a single fracture is given by cubic law: where ΔP is the pressure drop across the parallel capillary plate of length L t (along the flow direction) and aperture b, μ is the viscosity of water, and Δw is the width of the capillary plate.
Based on the Young-Laplace equation, the parallel capillary plate with aperture b is related to the magnitude of capillary pressure P c as where T s is surface tension and α is contact angle. For a given capillary pressure P c , a parallel capillary plate with an aperture equal to or less than (2T s cos α/|P c |) is generally filled by water and the rest by air. Integrating equation (10) from minimum aperture b min to the aperture b, the total flow rate Q of the entire fracture is given by where dN is the number of capillary plates corresponding to aperture b with capillary pressure P c . When the aperture distribution of fracture is defined as a probability density function f (b), inserting equation (10) into equation (12) leads to where V is the total void volume within a single fracture and the residual saturation S wr corresponds to the minimum aperture b min . As water flow is restricted to fracture aperture below the critical aperture b min , the water phase cannot flow in a continuous way and the volume fraction of the immobile water is defined as the residual saturation S wr . In the traditional approaches for developing a closedform relative permeability model in porous media, the circle capillary plates are both assumed to be straight. In reality, the flow paths are generally tortuous based on the experimental observations. e essential assumption of the traditional approaches is that the tortuosity factor is quantified by a power function and the exponent value cannot be determined rigorously in theory, which highly depends on experience from the experimental measurements.
As a result of the geometrical complication of the aperture distribution, the length of the disorder and irregular parallel capillary plate is similar to the measures of coastline, rough surface, mountains, lakes, and islands, which follow the fractal geometry rather than the Euclidean geometry [35]. According to the self-similarity of a fractal capillary plate, the length of a fractal capillary plate can be expressed as where ε is the length scale of measurement. For the flow through heterogeneous porous media, the tortuous length of the capillary tube is generally calculated by a fractal scaling relationship L t � ε 1− D T L D T 0 [36][37][38][39]. Similarly, we argue that the aperture b of parallel capillaries is analogous to the length scale ε. Hence, the length of a capillary plate with aperture b can be rewritten as

Advances in Civil Engineering
where D T is the fractal dimension of the tortuous parallel capillary plate, L 0 is the straight distance between inflow and outflow boundaries. For the fractal dimension of tortuous parallel capillary D T � 1, L t represents a straight flow path the same as the Burdine model. Inserting equations (3) and (16) into (13) yields For a given water-air flow system through a fixed fracture, ΔP, V, and L 0 are assumed to be constants. By the concept of relative permeability proposed, the relative permeability k r is given as where Q sa is the total flow rate under the saturated condition for the water phase.
Applying the Young-Laplace equation and substituting for b in equation (18) gives When the tortuosity fractal dimension approaches one with respect to a straight parallel capillary plate, equation (19) can be reduced to the Burdine model equation (2).
Combining equations (5) and (19) gives For D T /λ � 0 when λ ⟶ ∞ in which the unsaturated flow behavior in the rough-walled fracture is approximated to that in porous media with uniform pore size, equation (20) is reduced to the X-curve equation (1). For D T /λ � 1, equation (20) is reduced to the Brooks and Corey model equation (8). For D T /λ � 1.5, equation (20) is reduced to the Corey model equation (7). us, with one more parameter

The Fractal Relative Permeability for Aperture-Based Fractures
Based on the Young-Laplace equation (11), equation (19) can also be presented in the form of aperture b: For the tortuosity fractal dimension D T � 1, the numerator and denominator of equation (20) [30]. Once the probability density function f (b) of aperture distribution is known, the relative permeability can directly be determined from equation (21).
Measurements of aperture distribution from natural and artificial fractures have been performed through numerous measuring techniques including void casting [40], 3D laser/ stereotopometric scanning [41][42][43], and computed tomography [22,44,45]. e measured aperture distributions are supported to be either Gaussian or lognormal. For the sake of completeness, Gaussian and lognormal distribution are both employed to evaluate the relationship between the relative permeability and aperture distribution.

Normal Distribution-Based Fractures.
e Gaussian aperture distribution can be given by the following expression: where u and δ are the mean value and standard deviation of the aperture, respectively. Before presenting the relative permeability relationship with respect to the aperture distribution, it is useful to define a derived integration function: By defining t � (|b − u|/ � 2 √ δ) > 0, the following equation can be obtained: According to Newton's generalized binomial theorem [46], equation (24) can be rewritten as where c refers to Gamma function:

Advances in Civil Engineering 5
A combination of equations (21) and (22) yields and, using equation (26), the relative permeability expression is obtained as For D T is assigned to be zero in mathematics, the water saturation can be calculated as Relative permeability and saturation equations (28) and (29) based on Gaussian aperture distribution are functions of the parameters D T and b min . Figure 3 illustrates the k r -S w relationships for Gaussian aperturebased fractures with different tortuosity fractal dimensions and minimum apertures. For a Gaussian aperturebased fracture with low D T close to 2 when b min is kept at zero, the variation of the relative permeability k r versus S w is more pronounced and the relative permeability drops off rather rapidly with increasing D T . is is because, for a higher D T , a continuous flow path turns to be more tortuous and its component flow rate gets smaller according to the cubic law equation (10).
While the tortuosity fractal dimension D T remains unchanged, with the growth of b min corresponding to the larger residual saturation, the slope of the k r -S w curves is markedly elevated and the relative permeability decreases considerably. is may be a result of the fact that the number of capillary flow paths is fewer when much more void space is invalid to generate immobile capillary.

Evaluation of the New Relationship Based on Lognormal
Distribution. e lognormal aperture distribution can be given by the following expression: Note that u and δ are the mean value and standard deviation of the logarithm of aperture, respectively.
Similarly, the relative permeability based on the lognormal aperture distribution can be expressed as , the relative permeability and saturation can be derived as Relative permeability and saturation based on lognormal aperture distribution equations (32) and (33) are also a function of the parameters D T and b min . A relatively large range of the k r -S w relationships for lognormal aperturebased fractures with different tortuosity fractal dimensions and minimum apertures is depicted in Figure 4, indicating considerable sensitivity to the tortuosity fractal dimension and minimum aperture. Similar to the relative permeability curves shown in Figure 3, the relative permeability decreases with the increment of D T for a specified b min , while it reduces considerably with the increment of b min for a specified D T .

Experimental Validation against Laboratory
Observations. In order to demonstrate the validity of the proposed k r -S W relationships in the form of equations (28) and (29), the comprehensive experimental measurements for the nitrogen-water flow of rough-walled fractures reported by Chen and Horne [21] are applied to be compared with our predictions. e k r -S W data of two analog fractures (homogeneously and randomly rough-walled fractures) with Gaussian aperture distributions are both presented. e statistical parameters for homogeneously and randomly rough-walled fractures are u � 0.145 mm and δ � 0.03 mm, and u � 0.24 mm and δ � 0.05 mm, respectively. Figure 5 shows good agreement between the proposed model and the data from Chen and Horne [21]. Based on the least square method, the fitted parameters are D T � 1.75 and b min � 0.136 mm for homogeneously rough-walled fracture and D T � 1.9 and b min � 0.24 mm for randomly rough-walled fracture. e values for the correlation coefficient (R 2 ) of curve fitting are above 0.95. e predictions by other models including X-model, Brooks and Corey model, and Corey model are also plotted in Figure 4. e Brooks and Corey model gives better results than the Xmodel and Corey model. However, the residual saturation predicted using the Brooks and Corey model is specified as zero and diverges from the actual residual saturation (S wr � 0.25 and 0.39 for homogeneously and randomly rough-walled fractures, resp.) based on the experimental measurements.
Significantly, it is realized that the tortuosity fractal dimension D T of the homogeneously rough-walled fracture is less than that of randomly rough-walled fracture. is can be understood from the characteristics of the Gaussian aperture distribution, in which the aperture variation of the former does not have correlated spatial correlation whereas that of the latter does. erefore, the flow paths formed in the homogeneously rough-walled fracture are less tortuous corresponding to lower D T .

Comparison with Numerical
Solutions. At the present time, there are a few systematic studies on experimental data based on lognormal distributed fracture with which theoretical prediction for fracture relative permeability can be compared. However, Pruess and Tsang [24] have developed a numerical model based on percolation theory to numerically study the two-phase flow properties in fractures with lognormal aperture distributions and the relative permeability data for two typical fractures with the same lognormal aperture distribution (u � 0.0818 mm and δ � 0.0043 mm), and different anisotropy of spatial correlation are both calculated by simulating each phase flows separately in their filled void space. Based on a general simulator "MULKOM," the fracture aperture is discretized with 20 × 20 grid blocks and a pressure difference is assigned to the inflow and outflow boundaries while the residual boundaries are impermeable. Figure 6 shows matches between the simulated relative permeability-saturation relations and curves predicted using the proposed model, X model, Corey model, and Brooks and Corey model. For a given saturation, the X model generally overestimates the relative permeability and the Corey model generally underestimates the relative permeability, being in disagreement with the numerical results. e proposed model and Brooks and Corey model are overall closer to the measurements than the X model and Corey model predictions.

Advances in Civil Engineering
As indicated in Figure 6, some discrepancies between the analytical predictions and numerical results are apparent; for example, relative permeabilities are somewhat too low for large saturations, while, at low saturations, relative permeabilities are somewhat too high. Agreement is rather good at intermediate saturations. Note that the resultant D T values are identical with one, indicating that the flow paths in these two different spatial correlated fractures are close to straight without the influence of anisotropy of aperture distribution. is may be due to the insufficient statistical quantity of aperture elements as mentioned by Pruess and Tsang [24]; each fracture realization includes only 400 aperture elements and lacks the chance to generate tortuous flow paths in large probability.

Conclusions
e relative permeability in rough-walled fractures is a fundamental parameter in the multiphase flow through fractured media. A fractal model for estimating the relationship between relative permeability and saturation of rough-walled fractures is derived analytically based on the cubic law and the tortuous capillary model. e proposed model has considered the fractal characteristics of tortuous parallel capillary plates and can generalize the primary existing relationship in the literature. Based on the aperture-distributed dependence of relative permeability properties, the relative permeability expressions in the form of aperture with respect to Gaussian and lognormal distribution are both developed explicitly. e predictions of relative permeabilities from the proposed model are shown to be consistent with experimental observations for Gaussian distributed fractures and numerical data for lognormal distributed fracture. e fractal relative permeability model for rough-walled fracture is the terms of the tortuosity fractal dimension, determination of the tortuosity fractal dimension associated with aperture distribution is an interesting and challenging topic, and this work is in processing.

Data Availability
e data used in the present study may be available upon request from the authors.

Conflicts of Interest
e authors declare that they have no conflicts of interest regarding the publication of this paper.