**1. Introduction**

Darcy's law [1] establishes that the water flow in porous media is proportional to the hydraulic gradient; the proportionality coefficient is denoted hydraulic conductivity (K). The law, discovered in the context of water flow in saturated soils, has since been generalized to flow in unsaturated soil [2]. In saturated soils, conductivity is independent of water pressure, whereas in unsaturated soils it is a highly nonlinear function of pressure (ψ), or volumetric water content (θ) [3–8].

The saturated hydraulic conductivity, denoted KS, is at most a function of spatial coordinates. In unsaturated soils the hydraulic conductivity is a function of water pressure K(ψ) as well as spatial coordinates. In such cases the soil-water retention curve θ(ψ) is needed to relate volumetric water content to soil-water pressure. The two curves θ(ψ) and K(ψ) are known as the soil hydrodynamic characteristics and are important to the study of mass and transfers such as infiltration, drainage and evaporation, and groundwater recharge [2,9–12].

The aim of the present work is to establish relationships between the soil-water retention curve and hydraulic conductivity curve, using concepts of probability theory and fractal geometry in order to reduce the amount of unknown functions in the unsaturated soil zone. Soil here is considered as a set of Lebesgue measure different than zero. This does not consider sets which porosity is a unity, such as Menger sponge, which may not represent natural soil.

#### **2. Materials and Methods**

#### *2.1. A General Model of Hydraulic Conductivity*

A conceptual model for the hydraulic conductivity based on Poiseuille law of water flow in capillary tubes has been proposed at the literature [13–18]. The model has the general form:

$$\mathbf{K} = \mathbf{f} \mathbf{C}\_{\mathbf{f}} \int\_{\Omega} (\mathbf{R} / \mathbf{T})^2 \mathbf{d}w \tag{1}$$

where Ω represents the water flow area divided by total area of the exposed face [13]; f = ρ<sup>w</sup> g/η is the fluidity, ρ<sup>w</sup> is water density, η is the dynamical viscosity coefficient, g is the gravitational acceleration; R is the pore radius; the non-dimensional coefficient (Cf) take in account the irregular shape of the pore perimeter, for a circular pore Cf = 1/8; if R is taken as the hydraulic radius then Cf is called Koseny coefficient [19] with Cf = 1/2 for a circle; T is the tortuosity factor defined as T = dzf/dz ≥ 1, where z is the rectilinear path of water particles following macroscopic direction of the movement and zf is the actual path of water particles [20]; ω is the water flow effective area relative to total area of soil or partial effective areal porosity.

The effective area definition is established by Fuentes et al. [14–18] from the probabilistic idea of Childs and Collis–George [21] and fractal geometry concepts. After a perpendicular cut to macroscopic trajectory of water we obtain two faces, which are located at z and z + dz positions (Figure 1); the radii of pores of the z-face are denoted by r and those of the z + dz-face are denoted by ρ. A water particle in a pore of z-face can continue its trajectory by the same pore or by other pore of equal or different radius. The introduction of the joint probability of the two faces at intermediate point z + 1/2 dz allows the modeling of these possibilities.

**Figure 1.** A cross section in the perpendicular macroscopic direction of the water flow. In the z-face the different pore radius are represented by r and in the z + dz-face the different pore radius are represented by ρ, where dz is the pore size order.

We consider a completely saturated soil. Childs and Collis–George [21] assume that the pore size distribution f(r) is the same at both faces and that dθ(r) = f(r)dr is the water content in the pore interval containing r, (r − 1/2 dr, r + 1/2 dr), and dθ(ρ) = f(ρ)dρ is the water content in the pore interval containing ρ, (ρ − 1/2 dρ, ρ + 1/2 dρ), at the other face. The completely random joint probability of all the pores represented by these intervals is equal to the product of these probabilities, provided that dz is of the order of a pore size and less than a characteristic particle size. The product of dθ(r) and dθ(ρ) represents the flow effective area dω(r,ρ) = dθ(r)dθ(ρ) = f(r)dr f(ρ)dρ, which integration over whole pore domain gives the total effective flow area μ = φφ = φ2, where μ represents also a total effective areal porosity and φ the total volumetric porosity. In a parallel capillary system, a water

particle moves always in the same capillary and the effective area is equal to volumetric porosity dω(r) = dθ(r) = f(r)dr, which integration gives μ = φ [22].

In probabilistic terms, the Purcell model represents a complete correlation between the two faces, whereas the Childs and Collis–George model represents a complete decorrelation, these models represent the possible extreme behaviors. Mualem and Dagan [13] show that the Purcell model can be formally deduced from the Childs and Collis–George model if the flow effective area is defined as dω(r, ρ) = f(r)dr δ(ρ − r) dρ, where δ indicates the Dirac delta.

An intermediate approach between the Purcell and Childs and Collis–George models is proposed by Millington and Quirk [23]. These authors suppose that, if the areal porosity at each face is φs, then the effective areal porosity at intermediate point is μ = φsφ<sup>s</sup> where s = 1/2 represent the Purcell model and <sup>s</sup> <sup>=</sup> <sup>1</sup> the Childs and Collis–George model. Further, since <sup>φ</sup>2s <sup>≤</sup> <sup>φ</sup>, Millington and Quirk proposed to add the larger solid area (1 − φ) <sup>s</sup> to <sup>φ</sup>2s, in order to obtain the total area, that is (1 <sup>−</sup> <sup>φ</sup>) <sup>s</sup> + φ2s = 1. In most studies about soil structure, it is studied as a fractal object [16,17,24–28]. In this context, we'll offer an interpretation of the relationship (1 − φ) <sup>s</sup> + φ2s = 1.

In Euclidean geometry we have V <sup>∝</sup> L3 and A <sup>∝</sup> L2, where V, A, and L are the volume, area and length in an object, respectively; for example in a sphere V = 4/3πr3, A = 4πr2, L = r; since <sup>L</sup> <sup>∝</sup> <sup>V</sup>1/<sup>3</sup> we have A <sup>∝</sup> <sup>V</sup>2/3. Using the Mandelbrot [29] area–volume relationship in fractal geometry we have <sup>A</sup> <sup>∝</sup> <sup>V</sup>D/E, where D is the particle surface fractal dimension or the particle-pore interface fractal dimension and E = 3 is the Euclidean dimension of the space where the object is embedded.

In a cross section, the total area is the sum of the solid cross-sectional area and the pore cross-sectional area. If μ<sup>s</sup> = 1 − μ is the solid cross-sectional area relative to the total area (call it the areal solidity), and φ<sup>s</sup> = 1 − φ is the solid volume relative to soil volume (call it volumetric solidity), then μ<sup>s</sup> = φ<sup>s</sup> s, or 1 − μ = (1 − φ) s, where s = D/E is the fractal dimension relative to Euclidean dimension. Following the probabilistic idea, the areal and volumetric porosities relationship resulting is μ = φsφ<sup>s</sup> = φ2s; s = 1/2 corresponds to Purcell model and s = 1 to Childs and Collis–George model. Notice that the solid area takes the exponent s rather than 1-s; this is because μ<sup>s</sup> = φ<sup>s</sup> <sup>s</sup> is Hausdorff measure of the solid phase. The exponent 1-s is important in the calculus of the parallel body volume of a fractal [30].

The r-parallel body of a set F is defined by Pr(F) <sup>=</sup> <sup>x</sup> <sup>∈</sup><sup>E</sup> : <sup>x</sup> <sup>−</sup> <sup>y</sup> <sup>≤</sup> r, y <sup>∈</sup> <sup>F</sup> where <sup>E</sup> is the E-cartesian product of real numbers and represents soil (both pores and solids), F is the solids and <sup>E</sup> <sup>−</sup> <sup>F</sup> the pores. Clearly, F is contained in Pr(F). The volume of a parallel body is obtained as the product of the cover set's volume, cr<sup>E</sup> where c is a form coefficient (c = 1 if all covers are parallelepipeds), and the number of covers; therefore: volE(Pr) <sup>=</sup> NrcrE. Considering that Nr <sup>∝</sup> <sup>r</sup>−<sup>D</sup> when <sup>r</sup> <sup>→</sup> 0, we obtain volE(Pr) <sup>∝</sup> <sup>r</sup>E−D. It must be noted that the body parallel volume of the solids is not the same as the porous volume; it would, however, be the same if porosity tends to unity and if F were dense in E.

Since <sup>μ</sup><sup>s</sup> <sup>+</sup> <sup>μ</sup> <sup>=</sup> 1, the relationship between the relative fractal dimension and the total volumetric porosity is defined implicitly by the equation:

$$(1 - \phi)^{s} + \phi^{2s} = 1\tag{2}$$

where s = D/E. It can be shown that μ ≤φ, s→1/2 when φ→0, and s→1 when φ→1, in other words 1/2 < s < 1 when 0 < φ < 1.

From mass additive property we deduce the relationship ρ<sup>t</sup> = ρ<sup>s</sup> φ<sup>s</sup> + ρ<sup>v</sup> φ, where ρ<sup>t</sup> is the total density of soil, ρ<sup>s</sup> the density of solids, and ρ<sup>v</sup> the density of pores; if this last one is considered null, we get the classic φ = 1 − ρt/ρ<sup>s</sup> formula for estimating porosity, where ρ<sup>t</sup> becomes the total density of dry soil. The comparison with Equation (2) allows us to deduce that ρs/ρ<sup>t</sup> = φs−<sup>1</sup> <sup>s</sup> and ρv/ρ<sup>t</sup> = φ2s−1. The first could also be written as ρs/ρ<sup>t</sup> = φs−<sup>1</sup> <sup>s</sup> /φs, which is the quotient of Hausdorff measure and Lebesgue measure of solids. The second one can be interpreted as the quotient of the pores Hausdorff measure and the pores parallel body volume if we write ρv/ρ<sup>t</sup> = φs/φ1−s. It can be shown that <sup>φ</sup><sup>s</sup> <sup>≤</sup> <sup>φ</sup>1−<sup>s</sup> <sup>s</sup> / and <sup>φ</sup> <sup>≤</sup> <sup>φ</sup>1−s.

Clearly the relation s(φ) defined by Equation (2) cannot by applied, for instance, to Menger sponge, where <sup>φ</sup>= 1 and s <sup>=</sup> log20/3log3 - 0.9089; in other words, this relationship cannot be applied to any abstract sets where Lebesgue measure is zero (φ<sup>s</sup> = 0) and s < 1. Because relation s(φ) established a one-to-one relationship between porosity and the fractal dimension of solid pore interphase, it may not work for certain types of soils. However, their implications in the modeling of hydraulic conductivity are being investigated.

For the modeling of the hydraulic conductivity of unsaturated soils we accept the classic hypothesis that the water is contained in saturated pores with radius r, where 0 < r < R, for a given water content θ(R). Consequently, the effective area of flux, or partial area porosity, is the generalization of μ = φ<sup>s</sup> φs:

$$\mathbf{d}\,\omega(\mathbf{r},\boldsymbol{\uprho}) = \mathbf{d}\theta^{\mathbf{s}}(\mathbf{r})\mathbf{d}\theta^{\mathbf{s}}(\boldsymbol{\uprho})\tag{3}$$

As concerns the tortuosity factor, it was demonstrated in Fuentes et al. [17] that R is the pore radius measured perpendicularly to actual trajectory of the water particles (zf) and Rs is its projection in the macroscopic direction (z), consequently Rs/R <sup>=</sup> dz/dzf <sup>=</sup> 1/T. According to the probabilistic idea we have Rs <sup>∝</sup> Rs <sup>R</sup><sup>s</sup> <sup>=</sup> <sup>R</sup>2s, which we can write as an equality Rs/Rso <sup>=</sup> (R/Ro) 2s with Rso <sup>=</sup> Ro/To, where Ro is a reference radius and To the associated tortuosity factor. The pore radius-tortuosity factor relationship resulting is:

$$\text{T(R)} = \text{T}\_{\text{o}} \text{(R}\_{\text{o}}/\text{R)}^{\text{\textdegree}} \tag{4}$$

where 0 <sup>&</sup>lt; <sup>δ</sup> = 2s <sup>−</sup> <sup>1</sup> <sup>&</sup>lt; 1.

The general hydraulic conductivity model results of the introduction of the Equations (2) and (4) in the Equation (1):

$$\mathbf{K} = \mathbf{f} \frac{\mathbf{C}\_{\mathbf{f}}}{\mathbf{T}\_{\mathbf{o}}^{2} \mathbf{R}\_{\mathbf{o}}^{2(2s-1)}} \int\_{\Omega} \mathbf{R}^{4s}(\mathbf{r}, \boldsymbol{\rho}) \mathbf{d}\theta^{s}(\mathbf{r}) \mathbf{d}\theta^{s}(\boldsymbol{\rho}) \tag{5}$$

The saturated hydraulic conductivity (Ks) is obtained from Equation (5), replace Ω by the total pore domain ΩT.

#### *2.2. Classical Models of Hydraulic Conductivity*

Classic models reported in literature may be deduced from the proposed general model, provided that certain hypotheses are introduced. From Equation (3) we can deduce <sup>d</sup>ω(r, <sup>ρ</sup>) = <sup>s</sup>2[θ(r)]s−<sup>1</sup> [θ(ρ)]s−<sup>1</sup> <sup>d</sup>θ(r)dθ(ρ), with 0 <sup>&</sup>lt; <sup>r</sup> <sup>&</sup>lt; <sup>R</sup> and 0 <sup>&</sup>lt; <sup>ρ</sup> <sup>&</sup>lt; R. Assuming that the multiplicative function of the θ(r) and θ(ρ) differentials can be replaced by a medium value, which clearly depends on a superior limit, we have:

$$\mathbf{d}d\boldsymbol{\omega}(\mathbf{r},\boldsymbol{\rho};\mathbf{R}) = \left[\boldsymbol{\theta}(\mathbf{R})\right]^{2s-2} \mathbf{d}\boldsymbol{\theta}(\mathbf{r})\mathbf{d}\boldsymbol{\theta}(\boldsymbol{\rho})\tag{6}$$

where the s<sup>2</sup> term has been eliminated to satisfy:

$$\int\_{0}^{\mathbb{R}} \int\_{0}^{\mathbb{R}} \text{d}\omega(\mathbf{r}, \mathfrak{p}; \mathbb{R}) = \omega = \mathfrak{e}^{2s} \tag{7}$$

Likewise, tortuosity will depend only on the major radius, that of θ(R). We know that for a small R θ(R) = φ(R/Ro) <sup>λ</sup>, where Ro is the radius of a reference pore and λ > 0 is an index of pores (Brooks and Corey, 1964). From Equation (4) we obtain that:

$$\mathrm{T}(\mathbb{R}) = \mathrm{T}\_{\mathrm{o}} \Big[ \frac{\Phi}{\theta(\mathbb{R})} \Big]^{\mathrm{y}} \text{ with } \mathrm{y} = \frac{\delta}{\lambda} \tag{8}$$

In several articles [24,25] has been suggested that partial volumetric porosity is proportional to the volume of the parallel body, which means that <sup>θ</sup>(R)∝RE-D, where the fractal dimension D is estimated from soil-water retention curve, a result that's valid when porosity tends to unity. If that's the case, we get:

$$
\lambda = \mathcal{E} - \mathcal{D} \tag{9}
$$

Considering Equations (6) and (8), Equation (1) becomes:

$$\mathbf{K} = \mathbf{f} \frac{\mathbf{C}\_{\mathbf{f}}}{\mathbf{T}\_{\mathbf{o}}^{2}} \boldsymbol{\Phi}^{2s-2} \begin{bmatrix} \boldsymbol{\Phi} \\ \boldsymbol{\Phi} \end{bmatrix}^{\mathbb{P}} \int \left[ \mathbf{R}(\mathbf{r}, \boldsymbol{\rho}) \right]^{2} \mathbf{d}\boldsymbol{\theta}(\mathbf{r}) \mathbf{d}\boldsymbol{\theta}(\boldsymbol{\rho}) \text{ with } \mathbf{p} = (2s - 2) + (2\chi) \tag{10}$$

In the power p, the first addend represents the global effects of the correlation among pores, whereas the second represents the global effects caused by the tortuosity of flow trajectories.

For each relation between the effective radius R(r,ρ) and the radii r and ρ we can deduce a special model of the hydraulic conductivity. We will obtain four models corresponding to four relations R(r,ρ) utilized in models of the hydraulic conductivity reported at the literature: (i) 'small pore' model R(r,ρ) <sup>=</sup> min(r,ρ) used by Childs and Collis–George [21]; (ii) 'geometric pore' model <sup>R</sup>(r, <sup>ρ</sup>) = <sup>√</sup>r<sup>ρ</sup> used by Mualem [31]; (iii) 'neutral pore' model R(r,ρ) = either r or ρ corresponding to Burdine (1953) model; and iv) 'large pore' model R(r,ρ) = max(r,ρ) used by Fuentes [14]. Following the general indications of Brutsaert [32] for the integration of the Equation (10), the resulting special models are:

Small pore model: R(r, <sup>ρ</sup>) <sup>=</sup> min(r, <sup>ρ</sup>).

$$\mathbf{K}(\boldsymbol{\theta}) = \mathbf{f} \frac{2\mathbf{C}\_{\mathbf{f}}}{\mathbf{T}\_{\mathbf{o}}^{2}} \boldsymbol{\Phi}^{2s-2} \left[\frac{\boldsymbol{\Phi}}{\boldsymbol{\Phi}}\right]^{\mathbf{P}} \int\_{0}^{\boldsymbol{\Phi}} [\boldsymbol{\Phi} - \boldsymbol{\Phi}] \mathbf{r}^{2} \mathbf{d}\boldsymbol{\Phi} \tag{11}$$

This is the model of Childs and Collis–George [21], with a correction factor [θ/φ] p. Geometric pore model: R(r, <sup>ρ</sup>) = <sup>√</sup>rρ.

$$\mathbb{K}(\boldsymbol{\theta}) = \mathbf{f} \frac{\mathbf{C}\_{\text{f}}}{\mathbf{T}\_{\text{o}}^{2}} \boldsymbol{\Phi}^{2s-2} \left[\frac{\boldsymbol{\Theta}}{\boldsymbol{\Phi}}\right]^{\mathbb{P}} \left[\int\_{0}^{\boldsymbol{\Theta}} \mathbf{r} \,\mathrm{d}\,\boldsymbol{\theta}\right]^{2} \tag{12}$$

This model presents the structure of Mualem's model [31] with p = 1/2. Neutral pore model: R(r, <sup>ρ</sup>) <sup>=</sup> either r or <sup>ρ</sup>.

$$\mathbf{K}(\boldsymbol{\theta}) = \mathbf{f} \frac{\mathbf{C}\_{\mathbf{f}}}{\mathbf{T}\_{\mathbf{o}}^{2}} \boldsymbol{\Phi}^{2\mathbf{s}-1} \left[\frac{\boldsymbol{\Theta}}{\boldsymbol{\Phi}}\right]^{\mathbb{P}+1} \int\_{0}^{\boldsymbol{\Theta}} \mathbf{r}^{2} \mathbf{d} \,\mathbf{s} \tag{13}$$

This model presents the structure of Burdine's model [33] with p = 1. Large pore model: R(r, <sup>ρ</sup>) <sup>=</sup> max(r, <sup>ρ</sup>).

$$\mathbf{K}(\boldsymbol{\theta}) = \mathbf{f} \frac{2\mathbf{C}\_{\mathbf{f}}}{\mathbf{T}\_{\mathbf{o}}^{2}} \boldsymbol{\Phi}^{2s-2} \left[\frac{\boldsymbol{\Theta}}{\boldsymbol{\Phi}}\right] \int\_{0}^{\boldsymbol{\Theta}} \mathbf{r}^{2} \boldsymbol{\Phi} \,\mathrm{d}\boldsymbol{\Phi} \tag{14}$$

Childs and Collis–George [21] suppose that the resistance to the flow is determined by the small pore. Whereas Fuentes [14] proposes that, to deduce the opposite behavior, the conductance is determined by the pore of greater size [17]. Mualem [31] gives a little more weight to the large pore by proposing the geometric mean. The Burdine model is deduced when the pores sizes have the same weight.

The classic models are reported for relative hydraulic conductivity [Kr(θ)] and they are in function of the retention curve ψ(θ). These are obtained from Equations (11)–(14) with the rule Kr(θ) = K(θ)/K(θs), and the introduction of the Young–Laplace–Jurin law for the capillary rise phenomena:

$$
\psi = -\frac{\ell\_\mathrm{L}^2}{\mathrm{R}} \cos(\alpha\_\mathrm{c}) \tag{15}
$$

where the scale or capillary number (L) is defined by <sup>L</sup> <sup>=</sup> 2σ/ρwg [34], σ is the interfacial tension, <sup>L</sup> - 0.386 cm at 20 ◦C; α<sup>c</sup> is the contact angle formed between the air–water interface and the solid particles, assumed generally constant and equal to zero.

Classic models also consider the residual water content θr, defined by Brooks and Corey [35] as K(θr) = 0. This can be incorporated in the precedent models replacing θ by the effective water content θef = θ − θr, and φ by the effective volumetric porosity φef = φ − θr. The exponent s must be calculated by replacing φ with effective porosity (φef) in Equation (3); θ<sup>r</sup> is added to solid particles, i.e., φsef = 1 − φef = φ<sup>s</sup> + θr. In the classical conductivity models, the porosity φ is replaced by the volumetric water content to natural saturation θs, when entrapped air is considered.

The power p that appears in models (11)–(14) has been considered as an empiric parameter. This power (p <sup>=</sup> p1 <sup>+</sup> p2) is the result of the effects of correlation between pores, [θ(R)] p1 with p1 <sup>=</sup> 2s <sup>−</sup> 2, and tortuosity T2(R) <sup>=</sup> T2 <sup>o</sup>[φ/θ(R)] p2 with p2 <sup>=</sup> <sup>2</sup>(2s <sup>−</sup> <sup>1</sup>)/λ. To know the order of magnitude of the power p2 we assume that the particle surface fractal dimension is roughly equal to the fractal dimension estimated from the soil-water retention curve, i.e., from Equation (9) λ - 3(1 − s), hence p2 - 2(2s − 1)/3(1 − s). In consequence, the value of p may be estimated from porosity, some of which are shown in Table 1.

**Table 1.** Predicted values of the exponent p of the classical hydraulic conductivity models, which results from the effect of the pore correlation (p1) and the tortuosity factor (p2), for some values of the total volumetric porosity.


The approximate value of p ≈ 1/2 was obtained by Mualem [31] from the calibration of Equation (12) over the experimental data of 45 soils reported in different works, with total volumetric porosity in the range of 0.4 < φ < 0.7. According to Table 1, these soils may be represented by a soil with an average porosity of roughly 0.6.

#### **3. Results and Discussion**

#### *3.1. Some New Models of Hydraulic Conductivity*

We may obtain new hydraulic conductivity models from the general model established by Equation (5), without the introductions of those hypotheses established in Equations (6) and (8) to deduce the classic conductivity models, which may be restricted.

Each of the R(r,ρ) relationships mentioned above corresponds to a specific model of hydraulic conductivity. And again, following Brutsaert [32] for the integration of the Equation (5), the resulting special models are:

Small pore model: R(r, <sup>ρ</sup>) = min(r, <sup>ρ</sup>)

$$\mathbf{K}(\mathbf{R}) = \mathbf{f} \frac{2\mathbf{C}\_{\mathbf{f}}}{\mathbf{T}\_{\mathbf{o}}^{2}\mathbf{R}\_{\mathbf{o}}^{2(2s-1)}} \int\_{0}^{\mathbf{R}} [\boldsymbol{\Theta}^{s}(\mathbf{R}) - \boldsymbol{\Theta}^{s}(\mathbf{r})] \cdot \mathbf{r}^{4s} d\boldsymbol{\Theta}^{s}(\mathbf{r}) \tag{16}$$

Geometric pore model: R(r, <sup>ρ</sup>) <sup>=</sup> <sup>√</sup>rρ.

$$\mathbf{K}(\mathbf{R}) = \mathbf{f} \frac{\mathbf{C}\_{\mathbf{f}}}{\mathbf{T}\_{\mathbf{o}}^{2} \mathbf{R}\_{\mathbf{o}}^{2(2s-1)}} \left[ \int\_{0}^{\mathbf{R}} \mathbf{r}^{2s} \mathbf{d}\theta^{s}(\mathbf{r}) \right]^{2} \tag{17}$$

Neutral pore model: R(r, <sup>ρ</sup>) <sup>=</sup> r or R(r, <sup>ρ</sup>) <sup>=</sup> <sup>ρ</sup>.

$$\mathbf{K}(\mathbf{R}) = \mathbf{f} \frac{\mathbf{C}\_{\mathbf{f}}}{\mathbf{T}\_{\mathbf{o}}^{2} \mathbf{R}\_{\mathbf{o}}^{2(2s-1)}} \boldsymbol{\Theta}^{s}(\mathbf{R}) \left[ \int\_{0}^{\mathbb{R}} \mathbf{r}^{4s} \mathbf{d} \boldsymbol{\Theta}^{s}(\mathbf{r}) \right] \tag{18}$$

Large pore model: R(r, <sup>ρ</sup>) <sup>=</sup> max(r, <sup>ρ</sup>).

$$\mathbf{K}(\mathbf{R}) = \mathbf{f} \frac{\mathbf{2C\_f}}{\mathbf{T}\_o^2 \mathbf{R}\_o^{2(2s-1)}} \int\_0^\mathbf{R} \mathbf{r}^{4s} \theta^s(\mathbf{r}) \mathbf{d}\theta^s(\mathbf{r}) \tag{19}$$

Note that the Equations (17) and (18) can be generalized assuming <sup>R</sup>(r, <sup>ρ</sup>) <sup>=</sup> <sup>r</sup>αρ1−α, where 0 ≤ α ≤ 1; the Equation (17) follow with α = 1/2 and the Equation (18) with α = 0 or α = 1.

To obtain specific functions from new special models is necessary to provide the function θ(R). This can be obtained from the soil-water retention curve θ(ψ), relating the soil-water content to soil-water pressure, and the Laplace law, defined by the Equation (15).

From Equations (16)–(19) we can obtain the corresponding models to calculate the relative hydraulic conductivity from the retention curve:

$$\frac{\mathcal{K}(\Theta)}{\mathcal{K}\_{\mathsf{s}}} = \int\_{0}^{\Theta} (\Theta^{s} - \mathfrak{P}^{s}) \frac{\mathfrak{P}^{s-1}}{\left| \psi(\mathfrak{P}) \right|^{4s}} d\mathfrak{P} \Big/ \int\_{0}^{1} (1 - \mathfrak{P}^{s}) \frac{\mathfrak{P}^{s-1}}{\left| \psi(\mathfrak{P}) \right|^{4s}} d\mathfrak{P} \tag{20}$$

$$\frac{\mathcal{K}(\Theta)}{\mathcal{K}\_{\mathfrak{s}}} = \left[ \int\_0^{\Theta} \frac{\mathcal{S}^{s-1}}{\left| \psi(\mathfrak{s}) \right|^{2s}} d\mathfrak{s} \Big| \int\_0^1 \frac{\mathcal{S}^{s-1}}{\left| \psi(\mathfrak{s}) \right|^{2s}} d\mathfrak{s} \right]^2 \tag{21}$$

$$\frac{\mathcal{K}(\Theta)}{\mathcal{K}\_{8}} = \Theta^{8} \left[ \int\_{0}^{\Theta} \frac{\mathfrak{P}^{s-1}}{\left| \psi(\mathfrak{d}) \right|^{4s}} \mathrm{d}\mathfrak{d} \bigg/ \int\_{0}^{1} \frac{\mathfrak{P}^{s-1}}{\left| \psi(\mathfrak{d}) \right|^{4s}} \mathrm{d}\mathfrak{d} \right] \tag{22}$$

$$\frac{\mathcal{K}(\Theta)}{\mathcal{K}\_{\mathsf{K}}} = \int\_{0}^{\Theta} \frac{\mathfrak{s}^{2s-1}}{\left| \psi(\mathfrak{s}) \right|^{4s}} d\mathfrak{s} \Big| \int\_{0}^{1} \frac{\mathfrak{s}^{2s-1}}{\left| \psi(\mathfrak{s}) \right|^{4s}} d\mathfrak{s} \tag{23}$$

where <sup>Θ</sup> <sup>=</sup> (<sup>θ</sup> <sup>−</sup> <sup>θ</sup>r)/(θ<sup>s</sup> <sup>−</sup> <sup>θ</sup>r) is an effective degree of saturation.

#### *3.2. Applications*

#### 3.2.1. Brooks and Corey Equation

The equation proposed by Brooks and Corey [35] to represent the soil-water retention curve is:

$$
\Theta = (\psi\_{\rm cr}/\psi)^{\lambda} \tag{24}
$$

if ψ < ψcr, and Θ = 1 if ψcr ≤ ψ; where ψcr is a critical pressure and λ > 0 is an index of the pore distribution. The assumption λ = E − D [24,25] is not used here.

The introduction of the Equation (24) in the Equations (20)–(23) gives the same expression for the relative hydraulic conductivity:

$$\mathbb{K}(\Theta)/\mathbb{K}\_{\sf s} = \Theta^{2\sf s(2/\lambda + 1)}\tag{25}$$

The saturated hydraulic conductivity is given by:

$$\mathbf{K}\_{\mathfrak{s}} = \mathbf{f} \mathbf{C}\_{\mathfrak{f}} (\boldsymbol{\theta}\_{\mathfrak{s}} - \boldsymbol{\theta}\_{\mathfrak{r}})^{2s} (\mathbf{R}\_{\mathfrak{o}} / \mathbf{T}\_{\mathfrak{o}})^{2} \boldsymbol{\Lambda} \tag{26}$$

where we have defined Ro <sup>=</sup> <sup>λ</sup><sup>L</sup> λL/ ψcr . The factor Λ is different for each model: small pore <sup>Λ</sup><sup>s</sup> <sup>=</sup> 1/[2(2/<sup>λ</sup> <sup>+</sup> 1/2)(2/<sup>λ</sup> <sup>+</sup> <sup>1</sup>)];geometric poreΛ<sup>g</sup> <sup>=</sup> 1/(2/<sup>λ</sup> <sup>+</sup> <sup>1</sup>) 2 ; neutral poreΛ<sup>N</sup> <sup>=</sup> 1/[2(2/<sup>λ</sup> <sup>+</sup> 1/2)]; and large pore <sup>Λ</sup><sup>L</sup> <sup>=</sup> 1/(2/<sup>λ</sup> <sup>+</sup> <sup>1</sup>). The following inequalities are satisfied <sup>Λ</sup><sup>s</sup> <sup>&</sup>lt; <sup>Λ</sup><sup>g</sup> <sup>&</sup>lt; <sup>Λ</sup><sup>N</sup> <sup>&</sup>lt; <sup>Λ</sup>L, the equalities are given at extremes (<sup>λ</sup> <sup>→</sup> 0, <sup>∞</sup>). When <sup>λ</sup> <sup>→</sup> <sup>0</sup> we have: <sup>Λ</sup><sup>s</sup> <sup>=</sup> <sup>λ</sup>2/8 <sup>&</sup>lt; <sup>Λ</sup><sup>g</sup> <sup>=</sup> <sup>λ</sup>2/4 <sup>&</sup>lt; <sup>Λ</sup><sup>N</sup> <sup>=</sup> <sup>λ</sup>/4 <sup>&</sup>lt; <sup>Λ</sup><sup>L</sup> <sup>=</sup> <sup>λ</sup>/2.

The corresponding saturated hydraulic conductivity value satisfies the inequalities: Kss < Ksg < KsN < KsL.

#### 3.2.2. Generalized Power Function

One of the larger groups of models used to represent the soil–water retention curve is the following power function [36]:

$$
\psi = \psi\_{\rm d} \Theta^{-1/\lambda} \big( 1 - \Theta^{1/\rm m} \big)^{1/n} \tag{27}
$$

where ψ<sup>d</sup> is a pressure scale; m > 0, n > 0 and λ > 0 are three form parameters.

In Equation (27), we can note when Θ → 0 we obtain the Brooks and Corey [35] equation, and when λ = mn we obtain the van Genuchten equation [37]:

$$\Theta(\psi) = \left[1 + \left(\psi/\psi\_{\rm d}\right)^{n}\right]^{-m} \tag{28}$$

Introducing the Equation (27) in the Equations (20)–(23) we obtain, respectively:

$$\frac{\mathbf{K}(\Theta)}{\mathbf{K}\_{\mathbf{s}}} = \frac{\Theta^{\mathbf{s}} \mathbf{B}\_{\mathbf{l}} \big( \Theta^{1/\mathbf{m}}; 4\mathbf{s}\mathbf{m}/\lambda + \mathbf{s}\mathbf{m}, \ 1 - 4\mathbf{s}/\mathbf{n} \big) - \mathbf{B}\_{\mathbf{l}} \big( \Theta^{1/\mathbf{m}}; 4\mathbf{s}\mathbf{m}/\lambda + 2\mathbf{s}\mathbf{m}, \ 1 - 4\mathbf{s}/\mathbf{n} \big)}{\mathbf{B}(4\mathbf{s}\mathbf{m}/\lambda + \mathbf{s}\mathbf{m}, \ 1 - 4\mathbf{s}/\mathbf{n}) - \mathbf{B}(4\mathbf{s}\mathbf{m}/\lambda + 2\mathbf{s}\mathbf{m}, \ 1 - 4\mathbf{s}/\mathbf{n})} \tag{29}$$

$$\mathbf{K}(\Theta)/\mathbf{K}\_{\mathbf{s}} = \left[ \beta\_{\mathrm{l}} \big( \Theta^{1/\mathrm{m}}; 2\mathrm{s}\mathbf{m}/\lambda + \mathrm{s}\mathbf{m}, \ 1 - 2\mathrm{s}/\mathbf{n} \big) \right]^{2} \tag{30}$$

$$\mathcal{K}(\Theta)/\mathcal{K}\_{\sf s} = \Theta^{\sf s} \beta\_{\sf l} \Big( \Theta^{1/\sf m}; 4\sf{sm}/\lambda + \sf{sm}, \, 1 - 4\sf{s}/\bf n\right) \tag{31}$$

$$\mathbf{K}(\Theta)/\mathbf{K}\_{\mathbf{s}} = \beta\_{\mathrm{I}} \Big( \Theta^{1/\mathrm{m}}; 4\mathrm{s}\mathbf{m}/\lambda + 2\mathrm{s}\mathbf{m}, \,\mathbf{1} - 4\mathrm{s}/\mathbf{n} \Big) \tag{32}$$

where <sup>β</sup>I(x; p, q) <sup>=</sup> BI(x; p, q)/B(p, q), BI(x; p, q) is the incomplete beta function of variable x and parameters p > 0 and q > 0 and B(p,q) = BI(1:p,q) is the complete beta function.

We can obtain closed-form equations accepting the van Genuchten [37] idea consisting in to assign integral values to parameter p of the beta function and specially p = 1; this conduces to impose relationships between the form parameters of the soil–water retention curve. This idea is only applicable to models defined by Equations (30)–(32) because the unicity of these relationships. From Equation (29) we can obtain only closed-form equations of the first or second integral of the numerator, the result is an incomplete closed-form formula (semi closed-form) of the conductivity.

Small pore model:

$$\frac{\mathbf{K}(\Theta)}{\mathbf{K}\_{\mathbf{s}}} = \frac{\Theta^{\mathbf{s}} \mathbf{B}\_{\mathbf{l}} \Big(\Theta^{1/\mathbf{m}}; 1, 1 - 4\mathbf{s}/\mathbf{n}\Big) - \mathbf{B}\_{\mathbf{l}} \Big(\Theta^{1/\mathbf{m}}; 1 + \text{sm, } 1 - 4\mathbf{s}/\mathbf{n}\Big)}{\mathbf{B}(1, 1 - 4\mathbf{s}/\mathbf{n}) - \mathbf{B}(1 + \text{sm, } 1 - 4\mathbf{s}/\mathbf{n})}, \\ \boldsymbol{\lambda} = \frac{4\mathbf{s}\mathbf{m}}{1 - \text{sm}} \tag{33}$$

$$\frac{\mathbf{K}(\Theta)}{\mathbf{K}\_{\mathbf{s}}} = \frac{\Theta^{\mathbf{s}} \mathbf{B}\_{\mathbf{l}}^{\mathbf{s}} \mathbf{B}\_{\mathbf{l}}^{\mathbf{l}} \Theta^{\mathbf{l}/\mathbf{m}} ; 1 - \text{sm, } 1 - 4\mathbf{s}/\mathbf{n} \text{)} - \mathbf{B}\_{\mathbf{l}} \big( \Theta^{\mathbf{l}/\mathbf{m}} ; 1, \, 1 - 4\mathbf{s}/\mathbf{n} \big)}{\mathbf{B} (1 - \text{sm, } 1 - 4\mathbf{s}/\mathbf{n}) - \mathbf{B} (1, \, 1 - 4\mathbf{s}/\mathbf{n})}, \, \lambda = \frac{4\mathbf{s} \text{m}}{1 - 2\mathbf{s} \text{m}} \tag{34}$$

with n > 4s and BI <sup>Θ</sup>1/m; 1, 1 <sup>−</sup> 4s/n = (<sup>1</sup> <sup>−</sup> 4s/n) −1 1 − <sup>1</sup> <sup>−</sup> <sup>Θ</sup>1/m 1−4s/n Geometric pore model:

$$\mathbf{K}(\Theta)/\mathbf{K}\_{\mathbf{s}} = \left[ 1 - \left( 1 - \Theta^{1/\mathbf{m}} \right)^{1 - 2\mathbf{s}/\mathbf{n}} \right]^2, \lambda = \frac{2\mathbf{s}\mathbf{m}}{1 - \mathbf{s}\mathbf{m}} \tag{35}$$

.

with n > 2s.

Neutral pore model:

$$\mathcal{K}(\Theta)/\mathcal{K}\_{\mathbb{S}} = \Theta^{\mathbb{S}} \Big[ 1 - \left( 1 - \Theta^{1/\text{m}} \right)^{1 - 4\kappa/\text{n}} \Big], \ \lambda = \frac{4\text{sm}}{1 - \text{sm}} \tag{36}$$

with n > 4s.

Large pore model:

$$\mathbf{K}(\Theta)/\mathbf{K}\_{\mathbf{s}} = 1 - \left(1 - \Theta^{1/\mathbf{m}}\right)^{1 - 4\kappa/\mathbf{n}}, \ \lambda = \frac{4\mathbf{s}\mathbf{m}}{1 - 2\mathbf{s}\mathbf{m}}\tag{37}$$

with n > 4s.

We can note that soil–water retention curves induced by Equations (33) and (34) are equals to those induces by Equations (36) and (37), respectively.

The use of models (33)–(37) reduces the form parameters number of the soil–water retention curve defined by Equation (27): the three independent parameters {m, n, λ} are reduced to two parameters {m, n}.

The form parameters can even be reduced to one. If we assume λ = mn in the Equation (27) we obtain van Genuchten equation, Equation (28), which makes the function Θ(ψ) explicit where the form parameters {m, n} are still independent. If we accept the relationships between λ and m used to obtain Equations (33)–(37), the Equation (27) will have only one form parameter (m). The corresponding models of the conductivity associated to van Genuchten equation with a form parameter are the following:

Small pore model:

$$\frac{\mathcal{K}(\Theta)}{\mathcal{K}\_{\mathsf{s}}} = \frac{\Theta^{\mathsf{s}} \mathsf{B} \Big( \Theta^{1/\mathsf{m}}; 1, \text{ sm} \Big) - \mathsf{B} \Big( \Theta^{1/\mathsf{m}}; 1 + \text{sm}, \text{ sm} \Big)}{\mathsf{B} (1, \text{ sm}) - \mathsf{B} (1 + \text{sm}, \text{ sm})}, \\ 0 < \text{sm} = 1 - 4\mathsf{s}/\mathsf{n} < 1 \tag{38}$$

$$\frac{\mathbf{K}(\Theta)}{\mathbf{K}\_{\mathbf{s}}} = \frac{\Theta^{\mathbf{s}} \mathbf{B}\_{\mathbf{l}} \Big( \Theta^{1/\mathbf{m}}; 1-\text{sm, sm} \Big) - \mathbf{B}\_{\mathbf{l}} \Big( \Theta^{1/\mathbf{m}}; 1,\text{sm} \Big)}{\mathbf{B}(1-\text{sm, sm}) - \mathbf{B}(1,\text{sm})}, \ 0 < 2\text{sm} = 1 - 4\text{s}/\text{m} < 1\tag{39}$$

Geometric pore model:

$$\mathrm{K}(\Theta)/\mathrm{K}\_{\mathrm{s}} = \left[1 - \left(1 - \Theta^{1/\mathrm{m}}\right)^{\mathrm{sm}}\right]^2, \ 0 < \mathrm{sm} = 1 - 2\mathrm{s}/\mathrm{n} < 1\tag{40}$$

Neutral pore model:

$$\mathcal{K}(\Theta)/\mathcal{K}\_{\sf s} = \Theta^{\sf s} \Big[ 1 - \left( 1 - \Theta^{1/\sf m} \right)^{\rm sm} \Big], \ 0 < \text{sm} = 1 - 4\sf s /\bf n < 1 \tag{41}$$

Large pore model:

$$\mathcal{K}(\Theta)/\mathcal{K}\_{\mathbb{S}} = 1 - \left(1 - \Theta^{1/\text{m}}\right)^{2\text{sm}}, \ 0 < 2\text{sm} = 1 - 4\text{s/n} < 1\tag{42}$$

A first evaluation of the predictive capacity of the relative hydraulic conductivity models defined by Equations (38)–(42), on fifty soils of the GRIZZLY database reported by Haverkamp et al. [38], was presented by Fuentes et al. [18] with acceptable results.

The classical version of the models defined by Equations (21) and (22) corresponds to the Mualem [31] and Burdine [33] models that are defined in Equations (12) and (13). The use of these models has been proved by different authors [37]. We will show the capacity of prediction of the model that corresponds to the hypothesis of the large pore defined in Equation (23) and in particular the close-form equation of the relative hydraulic conductivity defined in Equation (42).

We use three of the five soils that van Genuchten [37] analyzes, holding the values θ<sup>r</sup> and θ<sup>s</sup> that the author use and the parameter s is estimated from θs. In Table 2 we present some properties of the three soils. The parameters ψ<sup>d</sup> and m obtained by least squares method corresponding to three different models are presented in Table 3.


**Table 2.** Some physical properties of the three analyzed soils.

**Table 3.** Parameters values of the soil water retention curve of the three analyzed soils corresponding to three different models, Equations (40)–(42).


Figures 2–4 present the fitted soil water retention curves, and the predicted relative hydraulic conductivity curves by the geometric pore, neutral pore and large pore models for the three studied soils. We can observe that the predictions are good enough in these three soils. In addition, in Figure 5 we illustrate the prediction capability of the small pore model using the relationships between m and n provides by the neutral pore and large pore models.

The analysis performed on the classical models and the comparison between experimental and predicted relative hydraulic conductivity with the four special new models allows us to show that the different models may be used to estimate relative hydraulic conductivity.

**Figure 2.** Observed and calculated curves of the soil hydrodynamic characteristics of hygiene sandstone, Equations (40)–(42).

**Figure 3.** Observed and calculated curves of the soil hydrodynamic characteristics of Touchet silt loam G.E.3, Equations (40)–(42).

**Figure 4.** Observed and calculated curves of the soil hydrodynamic characteristics of silt loam G.E.3, Equations (40)–(42).

**Figure 5.** Observed and calculated curves of the relative hydraulic conductivity with small pore model, Equations (38) and (39).

#### **4. Conclusions**

The proposed hydraulic conductivity model has been deduced using classic probabilistic theory and fractal geometry concepts in order to approach the flow effective area and the tortuosity of flow trajectories for each pore radius. The introduction of four definitions used at the literature of the pore effective radius that is involved in the general model has permitted to establish four new models to predict the relative hydraulic conductivity.

Some additional considerations related to the definitions of flow effective area and the tortuosity factor have allow us to deduce four classical models that are extensively used in different studies. In particular, we have given some interpretations of its empirical parameters in the fractal geometry context. The power function proposed by Brooks and Corey (1964) to represent the water retention curve leads to a power function to represent the relative hydraulic conductivity. The difference between the new models consists in the prediction of the relative hydraulic conductivity.

We have studied a general power function with three form parameters to represent the water retention curve reported by Braddock et al. [36], which presents as a special case the equation of Brooks and Corey [35] with a form parameter and the van Genuchten equation [37] with two form parameters. With this function we have gotten close-form equations of hydraulic conductivity with two and one form parameter. The general power function of the soil water retention curve and the close-form equations of hydraulic conductivity, obtained through special fractal models, may be used for the study of mass and energy transferences through soil, as infiltration, drainage, evaporation and ground water recharge.

**Author Contributions:** Formal analysis, C.F., C.C. and F.B.; investigation, C.F., C.C. and F.B.; software, C.F., C.C. and F.B.; C.F., C.C. and F.B.; writing—original draft, C.F. and F.B.; writing—review and editing, C.C. All authors have read and agreed to the published version of the manuscript.

**Funding:** This research received no external funding.

**Conflicts of Interest:** The authors declare no conflict of interest.
