Pore Space Connectivity and the Transport Properties of Rocks
Connectivité de l’espace poreux et propriétés de transport des roches
^{1}
Earth, Atmospheric and Planetary Sciences Department, Massachusetts Institute of Technology, Cambridge, Massachusetts – USA
^{2}
State Key Laboratory of Oil and Gas Reservoir Geology and Exploitation, Southwest Petroleum University, Chengdu – China
email: yvb@mit.edu
^{*} Corresponding author
Pore connectivity is likely one of the most important factors affecting the permeability of reservoir rocks. Furthermore, connectivity effects are not restricted to materials approaching a percolation transition but can continuously and gradually occur in rocks undergoing geological processes such as mechanical and chemical diagenesis. In this study, we compiled sets of published measurements of porosity, permeability and formation factor, performed in samples of unconsolidated granular aggregates, in which connectivity does not change, and in two other materials, sintered glass beads and Fontainebleau sandstone, in which connectivity does change. We compared these data to the predictions of a KozenyCarman model of permeability, which does not account for variations in connectivity, and to those of Bernabé et al. (2010, 2011) model, which does [Bernabé Y., Li M., Maineult A. (2010) Permeability and pore connectivity: a new model based on network simulations, J. Geophys. Res. 115, B10203; Bernabé Y., Zamora M., Li M., Maineult A., Tang Y.B. (2011) Pore connectivity, permeability and electrical formation factor: a new model and comparison to experimental data, J. Geophys. Res. 116, B11204]. Both models agreed equally well with experimental data obtained in unconsolidated granular media. But, in the other materials, especially in the low porosity samples that had undergone the greatest amount of sintering or diagenesis, only Bernabé et al. model matched the experimental data satisfactorily. In comparison, predictions of the KozenyCarman model differed by orders of magnitude. The advantage of the Bernabé et al. model was its ability to account for a continuous, gradual reduction in pore connectivity during sintering or diagenesis. Although we can only speculate at this juncture about the mechanisms responsible for the connectivity reduction, we propose two possible mechanisms, likely to be active at different stages of sintering and diagenesis, and thus allowing the gradual evolution observed experimentally.
Résumé
La connectivité porale est probablement l’un des plus importants facteurs affectant la perméabilité des roches réservoir. En outre, les effets de connectivité ne sont pas limités aux seuls matériaux approchant un seuil de percolation mais peuvent avoir lieu progressivement dans des roches soumises à des processus tels que la diagénèse mécanique ou chimique. Dans cette étude, nous avons compilé plusieurs jeux, préalablement publiés, de mesures de porosité, perméabilité et facteur de formation, obtenues sur des agrégats granulaires nonconsolidés, pour lesquels la connectivité ne change pas, ainsi que dans deux autres matériaux, billes de verre frittées et grès de Fontainebleau, pour lesquels elle varie. Nous avons comparé ces données aux prédictions d’un modèle de perméabilité de type KozenyCarman, qui ne tient pas compte des variations de connectivité, et à celles du modèle de Bernabé et al. (2010, 2011), qui les inclut [Bernabé Y., Li M., Maineult A. (2010) Permeability and pore connectivity: a new model based on network simulations, J. Geophys. Res.115, B10203; Bernabé Y., Zamora M., Li M., Maineult A., Tang Y.B. (2011) Pore connectivity, permeability and electrical formation factor: a new model and comparison to experimental data, J. Geophys. Res.116, B11204]. Les deux modèles sont apparus également en bon accord avec les données expérimentales pour les agrégats granulaires nonconsolidés. En revanche, pour les deux autres matériaux, en particulier dans les échantillons les moins poreux, qui ont été le plus longtemps soumis au frittage ou à la diagénèse, les données expérimentales n’ont pu être expliquées que par le seul modèle de Bernabé et al. En comparaison, les prédictions du modèle de type KozenyCarman étaient erronées de plusieurs ordres de grandeur. La valeur du modèle de Bernabé et al. réside dans sa capacité à reproduire la réduction continue et progressive de la connectivité ayant lieu lors du frittage ou de la diagénèse. Bien que nous ne puissions, à ce stade, que spéculer sur les mécanismes responsables de la réduction de connectivité, nous proposons deux mécanismes possibles, susceptibles d’être actifs à des étapes différentes du frittage et de la diagénèse, et donc de permettre l’évolution continue observée expérimentalement.
© Y. Bernabé et al., published by IFP Energies nouvelles, 2016
This is an Open Access article distributed under the terms of the Creative Commons Attribution License (http://creativecommons.org/licenses/by/4.0), which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.
Introduction
Permeability is the most important rock property controlling oil and gas production rates. There is great interest in the oil and gas industry for devising reliable and accurate models of permeability that could be applied to core, welllog and geophysical data. Dimensional analysis shows that permeability is the product of a length scale squared by a dimensionless, scaleinvariant quantity (Berryman, 1992ab, 1993; Bernabé et al., 2010), where “scaleinvariant” refers to invariance under application to the porous medium of a transform such that the distance between any pair of points in the medium is scaled by a constant factor. The choice of a lengthscale parameter is arbitrary (Bernabé et al., 2010). The only qualities required of an effective lengthscale parameter are that it should be unambiguously defined, and that it should be measurable in rock samples, preferably through a variety of independent techniques. One suitable lengthscale is the hydraulic radius r_{H}, here defined as twice the volume to surface ratio of the pore space (r_{H} = 2V_{p}/A_{p}, where V_{p} and A_{p} are the pore volume and internal wetted surface area, respectively). The hydraulic radius can be measured by several methods, including picnometry, the BET gasadsorption technique, or microstructural image analysis. The form of the scaleinvariant term is more difficult to determine. Fluid flow through porous media is affected to various degrees by a vast number of scaleinvariant factors (e.g., normalized moments of the pore size and porelength distributions, pore aspect ratio, relative porewall roughness, pore connectivity, and so forth). Therefore, identifying the major factors and their functional form in a permeability model is a challenging task.
Our main goals here are to show that pore connectivity is one of the most important scaleinvariant factors and, furthermore, that connectivity effects are not limited to media going through a percolation transition but can gradually occur in rocks undergoing common geological processes such as mechanical and chemical diagenesis, microfracturing during uplift, and so forth. These continuous changes in connectivity can alter permeability in ways that are not predicted using widely held models such as the equivalent channel model of Paterson (1983) and Walsh and Brace (1984) (PWB), which do not explicitly include pore connectivity. Our approach will be to analyze experimental datasets in porous materials subject to processes expected to affect differently the length scale r_{H} and the pore connectivity. As a result, these materials display very different interrelationships among porosity (ϕ), permeability (k) and formation factor (F), which can then be used to identify and quantify the effect of the connectivity variations. In particular, we will show that the experimental data can be better explained using the connectivity models such as Bernabé et al. (2010, 2011).
The paper is organized as follows. In Section 1, we briefly present the permeability and formation factor models used to interpret the experimental data. The materials selected for this study and the associated datasets are described in Section 2. The data are analyzed in Section 3 using the models of Section 1. The results are discussed in Section 4 and our conclusions presented in at the end.
1 Permeability and formation factor models
The PWB model belongs to the KozenyCarman class of models. According to PWB, the pore space can be represented by a single characteristic channel having the appropriate dimensions, shape and orientation. When simple channel shapes such as circular cylinders or flat slits are considered, it is easy to derive the following equations:(1) (2)and(3)where D is a geometric form factor equal to 1/8 for a circular cylinder and 1/12 for a flat slit, and τ denotes the tortuosity, i.e., the ratio of the channel length to the nominal length of the corresponding porous medium. Inspection of Equations (13) shows that the PWB model includes three scaleinvariant parameters, D, ϕ and τ, corresponding to the shape, volume fraction and orientation of the pores. If a strict geometrical definition is used for D and τ, they should both be expected to display a very limited range of variation (i.e., between 1/8 to 1/12 for D and around square root of 2 for τ). These ranges have often been found inadequate, and tortuosity, in particular, is frequently assumed to account implicitly for a host of effects not directly included in PWB (Clennell, 1997, for a thorough review of tortuosity). However, this practice transforms tortuosity into an illdefined parameter, which cannot be measured independently. Equation (2) can be replaced or supplemented by Archie’s law F = ϕ^{−m} (Revil and Cathles, 1999; Glover et al., 2006), in which the cementation exponent m is interpreted as an indicator of the pore space “connectedness” (Glover, 2009).
Bernabé et al. (2010, 2011) recently proposed an alternative model, in which connectivity is explicitly included. The model assumes that porescale heterogeneity and connectivity are the major scaleinvariant factors controlling k and F. For the sake of simplicity, porescale heterogeneity is restricted to crosssectional pore size variability and is measured by the standard deviation (σ) of the pore radius distribution normalized to the mean radius. Additional heterogeneity parameters may be needed for rocks with very complex pore shapes (e.g., carbonates) and/or a bimodal pore size distribution (e.g., presence of a microporous phase). Bernabé et al. (2010, 2011) notion of connectivity is similar to Glover’s “connectedness”, i.e., “general availability of pathways for transport” (Glover, 2009). According to this definition, a natural measure of pore connectivity is the mean coordination number (z), i.e., the average number of channels connected to a nodal pore (Bernabé and Maineult, 2015, for a discussion of z as pore connectivity measure). Bernabé et al. model is summarized in Equations (46):(4) (5)and(6)where the exponents (β, γ, α = β/γ) are functions of the porescale heterogeneity measure σ while the prefactors (C_{k}, C_{F}, C = C_{k}C_{F}^{−α}) depend on σ and ε, where ε refers to the mean aspect ratio of the approximately elliptical pore crosssections. Approximate functional forms for the prefactors and exponents were determined from the results of network simulations (Bernabé et al., 2010, 2011). They are repeated here for convenience in Table 1. The four scaleinvariant parameters used in the model, i.e., σ, z, ε and the normalized pore length l/r_{H}, can, in principle, all be measured in two and threedimensional images of the pore space of rock samples.
Approximate expressions of the exponents and prefactors of Equations (46) as functions of σ and ε (from Bernabé et al., 2011)
Before moving on to the rest of the study, we note that the classic definition of the formation factor, i.e., the ratio of the sample resistivity to that of the saturating brine is formally incorrect because an electrical double layer usually forms at the mineral/water interface, giving rise to anomalous surface conduction that may significantly affect the rock conductivity (Revil and Glover, 1997, 1998). Neglecting surface conduction leads to an underestimation of the formation factor that is approximately:(7)
where F* is the apparent formation factor, Σ_{S} the specific surface conductivity and σ_{br} the brine conductivity. In clayfree sandstones saturated with a highsalinity aqueous solution the ratio Σ_{S}/r_{H}σ_{br} is much smaller than 1 and surface conduction can be neglected. This condition is assumed to hold in the rest of the paper.
2 Materials: transport properties and hydraulic radius
2.1 Unconsolidated Glass Beads and Sand
If properly prepared, random packs of monodisperse spherical grains should have statistically equivalent structures with identical scaleinvariant parameters, independent of the radius of the spheres. The PWB and Bernabé et al. models both predict that, in such sphere packs, permeability scales as r_{H}^{2} while the formation factor remains constant. The same should be approximately true of unconsolidated packs of wellsorted, slightly nonspherical particles. To test this we considered the experimental data collected on unconsolidated packs of monodisperse spherical glass beads (Glover and Walker, 2009; Glover and Déry, 2010) and wellsorted, rounded sands (Biella and Tabacco, 1981; Biella et al., 1983). These datasets include values of porosity, grain radius, permeability and formation factor (note that F was measured in the high salinity limit, with negligible surface conduction). As can be seen in Figure 1, the porosity ϕ and formation factor F of unconsolidated, wellsorted granular media are nearly constant, equal to 0.39 and 4.3, respectively, while the permeability k varies over 6 orders of magnitude from about 0.1 to 10 000.10^{−12} m^{2} (or Darcy). These values of ϕ and F correspond to a cementation exponent of 1.5 in good accordance with theoretical predictions (Glover et al., 2000). A minor trend is visible in the data with high permeabilities corresponding to low porosities and high formation factors, and vice versa (Fig. 1a, c). This effect is more pronounced in sands than in glass beads and may reflect the existence of a weak correlation of grain size to grain asphericity (Bernabé et al., 2011) and/or packing structure. In the papers cited above, the experimental errors were on the order of few percent, i.e., substantially lower than the sampletosample fluctuations occurring during the preparation of the bead and sand packs. Sampletosample variability was even more prevalent in sintered glass beads and Fontainebleau sandstone as shown in Figure 1.
Figure 1. Compilation of published experimental data, a) k versus ϕ, b) F versus ϕ and c) k versus F, in unconsolidated monodisperse glass beads (Glover and Walker, 2009; Glover and Déry, 2010; black dots), wellsorted, rounded sands (Biella and Tabacco, 1981; Biella et al., 1983; blue dots), sintered glass beads (Wong et al., 1984; Guyon et al., 1987; Li et al., 1995; Blair et al., 1996; purple, red and orange squares corresponding to different grain sizes as indicated in the inset) and Fontainebleau sandstone (Doyen, 1988; Fredrich et al., 1993; Zamora, unpublished data; black diamonds; Revil et al., 2014; grey diamonds). The colored lines are inserted to highlight powerlaw trends in the various materials (note that the three families of sintered glass obey fairly similar power laws). The values of the exponents are given when necessary. 
The hydraulic radius of a pack of perfectly monodisperse spheres is given by:(8)
where R_{G} denotes the sphere radius (see the derivation in Appendix A). Note that Equation (8) considers a sphere pack of infinite extent and does not include the end effect associated with the walls of the container that must always be used in practice. Obviously, this end effect is negligible when the grain size is significantly lower than the container size, a condition that appears to be met in the studies mentioned above, except possibly for the millimetersize glass beads of Glover and Walker (2009), who used relatively small containers (2.5 cm diameter). Based on the small experimental errors reported for ϕ and R_{G}, the relative uncertainty δr_{H}/r_{H} was estimated to be lower than ±5% (Fig. 2).
Figure 2. Influence of the hydraulic radius on a) permeability and b) formation factor. The diagrams use the same symbol and line conventions as in Figure 1. Three horizontal error bars (dark purple) help visualize the relative uncertainty δr_{H}/r_{H} determined for the various materials (i.e., ±10, ±15 and ±20% for the unconsolidated sands, sintered glass beads and Fontainebleau sandstone, respectively). 
In the case of sand packs, the sand grains are not spherical and Equation (8) must be modified as:(9)
where the grain radius R_{G} is now defined as the radius of the sphere with an equal volume V_{G} = 4πR_{G}^{3}/3 and the asphericity factor f_{G} > 1 is equal to the ratio of the grain surface area A_{G} to that of the equivalent sphere, f_{G} = A_{G}/4πR_{G}^{2}. In order to estimate f_{G}, we applied two previously proposed formulas linking f_{G} and ϕ for packs of ellipsoidal grains (Coelho et al., 1997) and found f_{G} values increasing from 1 for the coarsest sands to about 1.3 for the finest ones (this increase of ellipticity with decreasing grain size was also noted in Bernabé et al., 2011). However, in addition to being slightly elongated or flattened (roughly corresponding to prolate or oblate ellipsoids), sand grains have irregular, uneven surfaces that further increase A_{G} and therefore f_{G}. Since the Biella and Tabacco (1981) and Biella et al. (1983) sands were described as rounded, we only increased f_{G} by a limited amount. Our final f_{G} estimates ranged from about 1.5 for the coarsest sands to 1.9 for the finest ones. Owing to the relative imprecision of f_{G}, we evaluated δr_{H}/r_{H} around ±10% (Fig. 2).
In loglog plots of permeability and formation factor versus the estimated hydraulic radius (Fig. 2), we observe power law dependences of k on r_{H} with exponents equal to 1.91 and 1.96 for glass beads and sand, respectively (notice that the exponent increases to 1.98 when Glover and Walker (2009) three coarsest glass bead packs are omitted to avoid the end effect mentioned in the previous paragraph). Conversely, F displays very weak variations associated with a possible minor correlation of grain size with grain shape and/or packing structure. In summary, wellsorted, unconsolidated granular media represent an excellent approximation of materials possessing a variable length scale r_{H} and constant scaleinvariant factors (with the possible exception of some minor factors that are not investigated in the present study).
2.2 Sintered Glass Beads
Sintering of granular materials is a densification process driven by the excess free energy associated with the surface of the grains. Sintering is activated by raising temperature to a sufficiently high value. Sintered glass beads have often been used as rock proxies to study the effect of porosity on the physical properties of rocks. In particular, Wong et al. (1984), Guyon et al. (1987), Li et al. (1995) and Blair et al. (1996) reported measurements of k, F and ϕ in sintered monodisperse glass beads. They used three main types of glass beads with radii of around 25, 50 and 100 μm (materials with grain radii as low as 90 μm and as high as 160 μm are included in the last group). Despite large fluctuations, the three categories display very similar powerlaw trends of permeability versus ϕ and F (k α ϕ^{4.3} and k α F^{−2.1}; Fig. 1a, c) and of the formation factor with respect to porosity (F α ϕ^{−2}, i.e., the classic Archie’s relation; Fig. 1b).
During sintering solid matter is transferred over very short distances and the total volume of solid does not change, even though the density of the aggregate increases. Densification of the granular assembly is caused entirely by redistribution of the solid phase and decrease of the pore volume. This process produces an augmentation of the number of grain contacts per grain very similar to the one experimentally observed in plastically compressed bronze powders (Fischmeister et al., 1978). Using these properties of sintering, we developed a method to estimate the hydraulic radius of monodisperse sintered glass beads from porosity measurements and the initial radius of the spheres (Appendix A for more details). Here, we estimated the relative uncertainty δr_{H}/r_{H} to be ±15% (Fig. 2).
We observe that permeability, formation factor and porosity follow welldefined power laws of the hydraulic radius, k α r_{H}^{3.7} (Fig. 2a), F α r_{H}^{−1.7} (Fig. 2b) and ϕ α r_{H}^{0.9} (not shown). The observed exponents are substantially different from those predicted by the PWB and Bernabé et al. models (Eq. 16) when the scaleinvariant factors are assumed constant. These differences imply that, during sintering, variations in hydraulic radius were coupled with changes in the major scaleinvariant factors, in particular pore connectivity.
2.3 Fontainebleau Sandstone
Fontainebleau sandstone is a 99% pure quartz arenite from the Paris basin, France. It is relatively well sorted with a mean grain radius of about 120 μm. One characteristic that makes Fontainebleau sandstone an excellent material for testing rock physics models, is that the formation includes a very wide range of porosities (from as low as 0.03 to about 0.25). In a seminal study, Bourbié and Zinszner (1985) measured permeability, porosity and the acoustic properties of a large set of Fontainebleau sandstone specimens. Moreover, they performed a detailed microstructure investigation using a variety of techniques, including the preparation of pore casts. Similarly, Doyen (1988) and Fredrich et al. (1993) combined quantitative microstructure observations and laboratory measurements of k, F and ϕ. They used twodimensional images of the pore structure to estimate the hydraulic radius (Fredrich et al., 1993), the mean pore coordination number z and the mean pore radius <r> (Doyen, 1988). The statistical distributions of porosity, pore coordination number z, throat radius r and pore length l were measured in threedimensional images of several Fontainebleau sandstone samples by Lindquist et al. (2000). Furthermore, we complemented the datasets cited above with the laboratory measurements of ϕ, k and F (including surface conduction correction) of Revil et al. (2014) and Zamora (unpublished data).
As reported by Bourbié and Zinszner (1985), Fontainebleau sandstone displays a very well defined relationship between porosity and permeability (Fig. 1a). It is then reasonable to assume that other transport properties and microstructural quantities may also be accurately expressed as functions of porosity, allowing the individually incomplete datasets of Doyen (1988), Fredrich et al. (1993), Revil et al. (2014) and Zamora (unpublished data) to be blended together into a more complete one (for example, see the approach used by Bernabé et al., 2010). However, visual inspection shows that the formation factor values reported in Revil et al. (2014) are significantly lower and more scattered than the other data sets (as is particularly evident for the samples with ϕ < 0.12; Fig. 1b, c). This discrepancy is probably not caused by surface conduction since the presumably uncorrected formation factors of Doyen (1988) and Fredrich et al. (1993) are actually higher than the corrected ones, while the reverse would be expected if surface conduction played a major role. Moreover, note that Doyen (1988) and Fredrich et al. (1993) used highly saline NaCl solutions, which cause negligible surface conduction in clayfree rocks such as Fontainebleau sandstone.
In order to acknowledge the differences in F, we divided the data into two groups, Fontainebleau 1 containing Doyen (1988), Fredrich et al. (1993) and Zamora (unpublished data) data, and Fontainebleau 2 including only the Revil et al. (2014) data. The two groups of (k, ϕ) data points appear well superposed on a loglog scatterplot (Fig. 1a) and display the classic downwardcurved trend originally reported by Bourbié and Zinszner (1985). We also note that the Fontainebleau sandstone data approximately coincide with those of the coarsest sintered glass beads, as should be expected since they nearly have the same grain size. The formation factor of the Fontainebleau 1 group is consistent with the classic Archie’s relation, F α ϕ^{−2}, while a smaller cementation exponent of about 1.8 is observed for the Fontainebleau 2 group (Fig. 1b). The difference between the Fontainebleau 1 and 2 groups is also visible in Figure 1c, where permeability is represented as a function of formation factor.
Bernabé et al. (2010) constructed a functional transform to calculate r_{H} as a function of ϕ in Fontainebleau sandstone based on Fredrich et al. (1993) results. This transform was established using only 4 data points and, therefore, may have limited applicability. Here, we tried to construct a more robust r_{H}(ϕ) transform (Tab. 2) by assuming that the average pore and/or throat radius values reported by Doyen (1988) and Lindquist et al. (2000) could be used as hydraulic radius proxies. Notice that the relative difference of the two transforms is modest, between 2 and 8%, and that the new one has a slightly greater range of variation, corresponding to a relative uncertainty δr_{H}/r_{H} of about ±20%.
The transforms used to complete the Fontainebleau sandstone datasets of Doyen (1988), Fredrich et al. (1993), Revil et al. (2014) and Zamora (unpublished data)
We observe that permeability declines much more rapidly with decreasing r_{H} in Fontainebleau sandstone than in the sintered glass beads (Fig. 2a) whereas the formation factor increases more sharply (Fig. 2b, in which the difference between Fontainebleau 1 and Fontainebleau 2 is visible). The implication is that the evolution of the pore structure during diagenesis of Fontainebleau sandstone (mostly quartz cementation) differs from that during sintering.
3 Data analysis and modeling
The procedure for testing the PWB model is as follows. According to a strict geometrical interpretation of PWB, we assign fixed constant values to the input parameters, D = 0.1 and τ^{2} = 2. We then calculate the PWB predictions of k and F using Equations (13) and compare the results to the experimental data. The model is deemed acceptable if the calculated k and F values fall within some intervals ±δk and ±δF, respectively, of the corresponding experimental values. Since D and τ^{2} are not allowed to vary, the model tolerances are given by δk/k ≈ 2 δr_{H}/r_{H} + δϕ/ϕ and δF/F ≈ δϕ/ϕ, where the δϕ/ϕ term is small compared to δr_{H}/r_{H}.
We specifically chose this procedure to test whether or not a model exclusively sensitive to variations in the characteristic length scale r_{H} could explain the experimental data. If a different goal were pursued, flexibility could be increased by separately considering the hydraulic and electrical tortuosities, τ_{h}^{2} = Dϕr_{H}^{2}/k and τ_{e}^{2} = ϕF. In this second approach, the effects of scaleinvariant factors such as pore connectivity, are implicitly incorporated into τ_{h}^{2} and τ_{e}^{2}. Note that, in order to be consistent with Equations (1) and (2), Equation (3) requires formal equality of τ_{h}^{2} and τ_{e}^{2} but it does not exclude large variations of tortuosity as a function of ϕ.
Similarly, in order to test the Bernabé et al. (2011) model, we need to assign input values to the mean pore coordination number, z, the normalized pore length, l/r_{H}, the heterogeneity measure, σ, and the crosssectional aspect ratio, ε. For this purpose, we prefer to use actual measurements if they exist (a rare case). When direct experimental data are not available, the best option is to infer the values of the desired parameters from measurements of other related quantities. The inference scheme may be based on a theoretical model (e.g., the estimation of r_{H} for unconsolidated glass beads using Eq. 8) or on a purely empirical relationship (e.g., the r_{H}(ϕ) transform established for Fontainebleau sandstone). When even indirect information is unavailable, the last option is to identify plausible, reasonably narrow ranges of variation of the input parameters and test the models using values falling within these ranges.
Examination of the functions C_{k}(σ, ε), C_{F}(σ, ε), β(σ) and γ(σ) in Equations (4) and (5) shows that the Bernabé et al. (2011) model is more sensitive to σ and z than l/r_{H} and ε. Furthermore, the sensitivity of the model strongly depends on the magnitude of σ and z. For example, imposing a change δz = ±0.1 (or, equivalently, δσ = ±0.05) results in a relative variation of calculated k as low as ±1% for σ = 0.05 and z = 10 and as high as ±60% for σ = 0.95 and z = 2 (moderately smaller relative variations are obtained for F in the same conditions).
3.1 Unconsolidated Glass Beads and Sands
It is rather obvious that unconsolidated monodisperse bead packs and wellsorted sands have a wellconnected pore space, for which the appropriate value of z is on the order of 6, and a relatively homogeneous pore structure with σ between 0.3 and 0.5 (Bernabé et al., 2011, and references therein). It is also clear that the crosssections of the pores, although not circular, are not elongated in a particular direction; thus, ε should be taken to be unity (Bernabé et al., 2011). Finally, the average pore length, l, must be related to the grain radius. Bernabé et al. (2011) estimated l ≈ 1.2 R_{G} for a Face Centered Cubic (FCC) array of identical spheres, but random sphere packs are less dense than FCC and a slightly greater value may therefore be more realistic (here, we used l = 1.3 R_{G}). In sand packs, the situation is different. We must account for a relatively broad grain size distribution and recognize that, in the average, the pore length should be more strongly controlled by the fine particles than be the coarse ones (see the analysis of binary sand mixtures in Bernabé et al., 2011). Wellsorted sands have grain radius distributions with a standard deviation lower than about 1/3 of the mean grain radius and hence contain a significant fraction of grains with a radius as low as half of the mean grain radius. Thus, values of l as low as 0.65 R_{G} can reasonably be assumed (here, we used l = 0.85 R_{G}). For convenience, the assigned values of the input parameters are summarized in Table 3.
Summary of the values assigned to the input parameters of the Bernabé et al.’s models
Owing to moderately low δr_{H}/r_{H}, high z and low σ, the δk/k and δF/F tolerances for both models were found to be relatively small, i.e., δk/k < ±30% and δF/F < ±10%. Despite these relatively tight tolerances, both the PWB and Bernabé et al. models fitted the experimental k and F data quite well for both glass bead packs and sands (Fig. 3). We note, however, a slight tendency of the models (especially, when expressed in Eq. 3 and 6) to overestimate the permeability of the glassbead packs (Fig. 3a, c). It is also worth emphasizing that a strictly geometrical definition of tortuosity (i.e., τ^{2} = 2) is sufficient to ensure a good fit of the PWB model.
Figure 3. Comparison of calculated a) k from Equations (1, 4), b) F from Equations (2, 5) and c) k from Equations (3, 6) to the observed values of permeability and formation factor in the unconsolidated glass beads and sands. The meaning of the different symbols is explained in the inset, where k(z) refers to Bernabé et al. Equation (4), F(z) to Equation (5) and k(F) to Equation (6). The vertical error bars (dark purple) represent the tolerances estimated for the Bernabé et al. model, δk/k ≈ ±30% and δF/F ≈ ±10%. 
3.2 Sintered Glass Beads
Clearly, sintering is bound to produce significant structural changes in the glass bead packs. We expect pore connectivity to decrease and porescale heterogeneity to increase during sintering. Also the reduction of the intergranular distances should lead to a simultaneous decrease of the pore length. One of our main assumptions in this paper is that structural changes accompanying processes such as sintering or diagenesis, should be gradual and that the model parameters should be expressible as simple, continuous functions of porosity. Here, for the sake of simplicity, we limited our analysis of z, σ and l to linear functions of ϕ. At very low levels of sintering (i.e., high porosity, ϕ = 0.4), z, σ and l should evidently take the values mentioned in Section 3.1 for unconsolidated glass bead packs (i.e., z = 6, σ = 0.45 and l = 1.3 R_{G}). To determine the linear functions z(ϕ),σ (ϕ) and l(ϕ) we only need to find the adequate values at the low porosity endmember. Having no microstructure information on the sintered glass beads at low porosities, we proceeded by trial and error and found that z = 4, σ = 0.7 and l = 0.4 R_{G} at ϕ = 0.02, yielded satisfactory results. The assigned values of the input parameters are summarized in Table 3.
The main result is that, in sintered glass beads, the Bernabé et al. Equations (4, 5) worked significantly better than the PWB Equations (1, 2) (Fig. 4a, b). As expected from the results described in Section 3.1, the two models gave almost identical, wellmatching results at porosities near 0.4 but PWB produced increasingly strong overestimation of k and underestimation of F with decreasing porosity. The δk/k and δF/F tolerances calculated for the PWB model were clearly too small to explain the observed discrepancies of up to one order of magnitude at very low ϕ. However, the tolerances for both models were nearly sufficient to account for the experimental sampletosample variability. Surprisingly, we found identically good fits with both the PWB Equation (3) and Bernabé et al. Equation (6) (Fig. 4c). The success of Equation (3) is remarkable since this equation was derived from Equations (1, 2), which both appear to be invalid in sintered glass beads. This outcome gives support to the classic idea of using electrical conductivity measurements to estimate the “connectedness” of the pore space and incorporate it into a permeability model (Revil and Cathles, 1999; Glover et al., 2006; Glover, 2009).
Figure 4. Comparison of calculated a) k from Equations (1, 4), b) F from Equations (2, 5) and c) k from Equations (3, 6) to the observed values of permeability and formation factor in the sintered glass beads. The meaning of the different symbols is explained in the inset, where k(z) refers to Bernabé et al. Equation (4), F(z) to Equation (5) and k(F) to Equation (6). The vertical error bars (dark purple) represent minimal estimates of the model tolerances (i.e., δk/k ≈ ±30% and δF/F ≈ ±10% for the PWB model as well as the Bernabé et al. model applied to high porosity samples; δk/k ≈ ±50% and δF/F ≈ ±20% for the Bernabé et al. model in low porosity samples). 
3.3 Fontainebleau Sandstone
Based on the reported microstructural data (Doyen, 1988; Fredrich et al., 1993; Lindquist et al., 2000), Bernabé et al. (2010, 2011) obtained approximate expressions of the pore length l and pore coordination number z as functions of porosity (Tab. 2). These transforms describe a large reduction of z from about 6 to 1.6 (very near the percolation threshold) and a substantial increase of l/R_{G} from 1.2 to 2.2, over a wide range of porosities (i.e., from 0.2 to 0.04). Bernabé et al. (2010, 2011) argued that pore space images (Bourbié and Zinszner, 1985; Doyen, 1988) suggest that the pores in Fontainebleau sandstone tend to have elongated crosssections and that the crosssectional aspect ratio ε decreases with decreasing porosity. The fact that the highest porosities reported in Doyen (1988), Fredrich et al. (1993), Revil et al. (2014) and Zamora (unpublished data) were much lower than 0.4, the porosity of unconsolidated sands, suggests that the corresponding rock samples had already experienced a significant amount of diagenesis and were therefore likely to have pore aspect ratios lower than unity. For the sake of simplicity, we assumed that ε was a linear function of ϕ, decreasing from 0.31 to 0.03 in the porosity range mentioned above. The pore crosssection area and pore length data reported by Lindquist et al. (2000) approximately satisfied highly skewed exponential distributions, which have equal standard deviations and means, suggesting that the heterogeneity measure σ should also be relatively high. By trial and error, we found that a constant value σ = 0.75 gave satisfactory results. The assigned values of the input parameters are summarized in Table 3.
The results obtained for the Fontainebleau sandstone datasets are quite similar to those for the sintered glass beads. The PWB Equations (13) yielded adequate results only for the samples with highest porosities but overestimated k by up to three orders of magnitude and underestimated F by one and half at the lowest porosities (Fig. 5). The observed discrepancies were much greater than the δk/k and δF/F tolerances calculated for the PWB model. In contrast, Bernabé et al. Equation (4) produced a very good fit with data for both Fontainebleau 1 and Fontainebleau 2 (Fig. 5a). Equation (5) was also satisfactory, although the estimated δF/F tolerance could not offset the large overestimation of F observed for the low porosity samples of the Fontainebleau 2 datasets (Fig. 5b). Equation (6), while still relatively acceptable, displayed a visibly poorer fit quality than Equation (4), particularly for the Fontainebleau 2 dataset (Fig. 5c). The data scatter associated with both the PWB Equation (3) and Bernabé et al. Equation (6) was much larger in Fontainebleau sandstone (even if we only consider the Fontainebleau 1 dataset) than in the sintered glass beads (compare Fig. 4c and 5c), suggesting that additional, relatively important scaleinvariant factors were probably omitted in the models.
Figure 5. Comparison of calculated a) k from Equations (1–4), b) F from Equations (2, 5) and c) k from Equations (3, 6) to the observed values of permeability and formation factor in the Fontainebleau sandstone. The meaning of the different symbols is explained in the inset, where k(z) refers to Bernabé et al. Equation (4), F(z) to Equation (5) and k(F) to Equation (6) and where the digits 1 and 2 indicate the Fontainebleau 1 and 2 datasets. The difference between the two datasets is particularly visible in the diagrams b) and c). The vertical error bars (dark purple) represent minimal estimates of the δk/k and δF/F model tolerances (i.e., δk/k ≈ ±40% and δF/F ≈ ±10% for the PWB model; δk/k ≈ ±60% and δF/F ≈ ±20% for the Bernabé et al. model in high porosity samples and δk/k ≈ ±90% and δF/F ≈ ±30% in low porosity samples). 
4 Discussion
The main observation reported here is that the PWB model correctly predicts the transport properties of unconsolidated, wellsorted granular media (Fig. 3) but becomes gradually more inexact when pore connectivity is altered by processes such as sintering or cementationcontrolled diagenesis (Fig. 4, 5). Our task in this discussion is to identify the origin of the PWB breakdown.
The granular materials investigated by Biella and Tabacco (1981), Biella et al. (1983), Glover and Walker (2009) and Glover and Déry (2010) cover nearly three orders of magnitude in grain size and nevertheless display very narrow ranges of porosity (0.37 to 0.41, except a couple of outliers as high as 0.45) and formation factor (4 to 5). This insensitivity to grain size demonstrates that these materials have essentially identical pore structures, characterized by nearly constant scaleinvariant factors. Their permeability, therefore, must almost exclusively depend on r_{H}, the length scale considered here. As discussed in Section 2.1, Equation (8) allows an accurate estimation of r_{H} for glass bead packs. For the sands, Equation (9) is satisfactory, although the asphericity factor f_{G} is somewhat uncertain. Moreover, good fits were generally obtained with both the PWB and Bernabé et al. models using realistic values of the input parameters. Accordingly, we can conclude that these data supports both the PWB and Bernabé et al. models with constant scaleinvariant parameters, summarized in the following power laws, k α r_{H}^{2} and F α r_{H}^{0}.
In sintered glass beads and Fontainebleau sandstone, on the other hand, the predictions of Equations (1, 2) (denoted k_{PWB} and F_{PWB}, respectively) become strongly erroneous (Fig. 4, 5). As mentioned earlier, another way to look at Equations (1, 2) is to calculate the hydraulic and electrical tortuosities τ_{h}^{2} and τ_{e}^{2}. Both parameters are found to increase with decreasing porosity (Fig. 6). One interesting observation is that τ_{h}^{2} and τ_{e}^{2} have comparable values in sintered glass beads (consistent with the validity of Eq. 3 demonstrated by Fig. 4c), whereas τ_{h}^{2} increases much faster than τ_{e}^{2} in Fontainebleau sandstone (Fig. 5c). Furthermore, the rate of increase of τ_{e}^{2} is about the same for both materials (Fig. 6). As a quantitative summary, we note that, at porosities of about 0.02 to 0.04, the ratio k_{PWB}/k reached values as low as 1/10 and 1/1 000 in sintered glass beads and Fontainebleau sandstone, respectively, while F_{PWB}/F approached 15 in both materials. What causes these discrepancies? PWB has three easily identified unrealistic assumptions, (a) the pores have idealized shapes, namely, straight cylinders or flat slits, (b) the pore space is fully connected and pore connectivity cannot change, (c) fluctuations in pore size, length and other geometric characteristics are negligible.
Figure 6. The hydraulic and electrical tortuosities τ_{h}^{2} and τ_{e}^{2}versus porosity for the sintered glass beads and Fontainebleau sandstone. The meaning of the different symbols is explained in the inset, where H and E stand for hydraulic and electrical. For the sake of visibility, the Fontainebleau 2 dataset is not included in this diagram but the associated trends are comparable to the Fontainebleau 1 trends represented here. 
We will first discuss hypothesis (a). The very small range of variation of the form factor D mentioned in Section 1 suggests that the crosssectional shape of the pores does not have a substantial effect on fluid flow (indeed, square and equilateral triangular conduits have form factors of 1/9 and 1/9.6, respectively). Although complex crosssectional shapes such as stars with narrow spokes may produce a greater effect (Yale, 1984), unrealistic crosssectional pore shapes can be ruled out as the main cause of the breakdown of PWB in sintered glass beads and Fontainebleau sandstone.
Alternatively, we note that channels in porous rocks have constricted shapes with smaller crosssections (usually known as throats) in their middle parts than at the ends. A gradual narrowing of the constrictions during sintering and diagenesis will increasingly impede fluid flow and may explain the breakdown of PWB. In order to estimate this effect, we calculated the hydraulic and electrical conductances of a frustumshaped pipe (i.e., truncated cone) and compared them to those of an equivalent straight cylindrical pipe (Appendix B).
In order to achieve identical values of porosity and hydraulic radius, the frustum and cylindrical pipes must have identical volumes and wetted surface areas, yielding the following equations:(10)
where r_{1}, r_{2} > r_{1} and l_{c} denote the endface radii and length of the frustum pipe, and, r_{H} and l refer to the estimated pore radius and pore length discussed in Sections 2 and 3. Although intuitively appealing, assuming l_{c} = l is not practical since the system of Equations (10) does not necessarily admit positive, real solutions for r_{1} and r_{2} besides the obvious but irrelevant r_{1} = r_{2} = r_{H}. Alternatively, we can complement the system with either Equations (B5) or (B7) of Appendix B, in which g_{H} is assumed equal to the ratio k/k_{PWB} and g_{E} to F_{PWB}/F. We then solve both expanded systems for r_{1}, r_{2} and l_{c}, and check whether or not the two sets of solutions are consistent with each other.
For the sintered glass beads, the lowest porosity (ϕ ≈ 0.02) occurred in the coarsest group (R_{G} = 100 μm), corresponding to r_{H} ≈ 1.7 μm, l ≈ 39 μm, g_{H} ≈ 1/10 and g_{E} ≈ 1/15. The system with Equation (B5) and that with (B7) produced significantly different values for the constriction radius, i.e., r_{1} ≈ 0.6 μm and 0.1 μm, respectively, although the single values r_{2} ≈ 2.5 μm and l_{c} ≈ 45 μm were obtained in both cases. The failure to find consistent solutions for r_{1} is due to the fact that g_{H} and g_{E} have comparable magnitudes although a given constriction has a greater effect on fluid flow than electrical conduction. Thus, tightening of pore constrictions is not a plausible cause of the breakdown of PWB in sintered glass beads.
For the Fontainebleau sandstone, at low porosities (ϕ ≈ 0.04) we have r_{H} ≈ 14 μm, l ≈ 250 μm, g_{H} ≈ 1/1 000 and g_{E} ≈ 1/15, which yield consistent solutions, i.e., r_{1} ≈ 0.9 μm, r_{2} ≈ 21 μm and l_{c} ≈ 320 μm, for the two systems. Nevertheless, pore constrictions can be ruled out as major PWB breakdown factors in Fontainebleau sandstone also because r_{1} ≈ 0.9 μm is much smaller than Lindquist et al. (2000) reported mean throat radii (from 23 μm at a porosity of 0.22 to 18 μm at 0.075). Notice that, given the very strong throat size heterogeneity observed in Fontainebleau sandstone, we cannot conclude that very narrow constrictions do not exist in this rock and that throat tightening does not happen at all. Our conclusion is merely that throat tightening cannot alone explain the breakdown of the PWB model.
We are left with hypotheses (b) and (c), which are addressed in Bernabé et al. (2010, 2011). We saw in Sections 3.2 and 3.3 that sets of reasonable values of z, l/r_{H}, σ and ε producing good fits with the experimental data, can be found for both the sintered glass beads and Fontainebleau sandstone. These results suggest that sintering and cementationcontrolled diagenesis of unconsolidated granular media lead to reduced pore connectivity and increased pore scale heterogeneity. In Fontainebleau sandstone, diagenesis is also associated with flattening of the pores.
One important question is: how does reduction of connectivity occur? As proposed by Zhu et al. (1995, 1999), one possible mechanism is pinching off of conduits, which transforms connecting conduits into pairs of deadends and strings of isolated pores. The specific pinching off mechanism invoked by Zhu et al. (1999) is tube ovulation, a surface tensioncontrolled process that occurs increasingly rapidly with decreasing tube radius. Pinching off should, therefore, be primarily active during the late stages of sintering and diagenesis when appropriately narrow throats are present.
We can also expect rearrangement of grains to take place during the early stages of sintering and diagenesis, owing to the sufficiently large pore space available for grain motion and to the presence of intergranular mechanical forces (gravity and surface tension in the case of sintering and overburden pressure for diagenesis). Such a grain rearrangement event is schematically illustrated in Figure 7a. The large nodal pore in the upper cartoon is split into two interconnected smaller cavities, each one acting as a new nodal pore and connected to a lower number of conduits than the original one. Hence, splitting of nodal pores is associated with a decrease of the pore coordination number and, therefore, of connectivity (Glover and Walker, 2009, also attributed connectivity loss to grain rearrangements). Notice that, in the schematic example of Figure 7a, the coordination number reduction is not associated with a change in topology. In graph theory, the appropriate topological invariant is the first Betti number (or genus), ζ = N_{b} – N_{n} + N_{c}, where N_{b} denotes the number of branches, N_{n} the number of nodes and N_{c} the number of separate connected clusters (in Fig. 7a, N_{c} = 1). Indeed, ζ is left unchanged by the grain rearrangement event represented in Figure 7a whereas the mean coordination number z decreases (see also the example of Fig. 7b, where branchnode pairs are added at diverse locations without changing ζ).
Figure 7. a) Schematic representation of the node splitting mechanism. The cartoon is twodimensional but it is meant to illustrate a threedimensional arrangement of sintered grains and pores between them. The large blue dots indicate the center of the nodes and the blue arrows symbolize the (out of plane) channels emanating from them. The large red arrow indicates the evolution associated with sintering. b) Example of graph transformations that leave the graph topological invariant (genus) unchanged while reducing the mean coordination number. In this cartoon, extra node/branch combinations (represented in red) are inserted at different places in the graph without altering the genus. By definition, nodes must always be present at the extremities of a branch (the dotted segments are intended to symbolize the unspecified continuation of the graph). 
As argued by Bernabé and Maineult (2015), pore connectivity is not a pure topological property. In network simulations, the dissipation of energy associated with fluid flow occurs in the branches while the nodes are just locations were mass conservation is enforced. If the resistance of a particular branch is negligibly small, this branch can be removed and the two nodes connected to it merged together without any change in the numerical results (Bernabé and Maineult, 2015). Numerical merging of nodes was indeed found to be necessary in threedimensional microstructure studies using skeletonization algorithms, which tend to define nodes only at triple junctions (z = 3) and blind ends (z = 1). These algorithms thus identify many hydraulically irrelevant branches and require implementation of corrective steps to eliminate them. For example, Lindquist et al. (2000) devised a rule to merge nodes connected by branches shorter than some threshold value. Similarly, Petford et al. (2001) measured z by considering only the branches possessing a significant constriction (in other words, by counting the throats) and by merging the redundant nodes. The Bernabé et al. (2010, 2011) model, of course, assumes hydraulically representative values of z (for example, they found that the z(ϕ) transform based on Doyen’s, 1988, data gave better results that the one from Lindquist et al., 2000).
In summary, the processes referred to above as pinching off and node splitting could explain the regular reduction in pore connectivity inferred from the sintered glass beads and Fontainebleau sandstone experimental data. Although they both reduce the mean pore coordination number, these two processes affect the pore space topology differently. Pinching off of a conduit adds two nodes and a branch and hence decreases ζ by one. In the example of Figure 7a, node splitting leaves ζ constant, but, if a node splits into more than two connected nodes and if the new nodes form a loop, the genus ζ will increase. Therefore, determining the genus of the skeletonized pore space of rock samples may provide useful insight into the rock genetic processes, although we do not expect ζ (unlike the coordination number z) to be formally related to the rock transport properties.
Conclusions
Our main conclusions are:

one of the most important scaleinvariant factors controlling permeability is pore connectivity;

a recent permeability model that explicitly includes the effect of connectivity and takes connectivity variations into account, performs much better than an older, widely used one that does not;.

the application of the Bernabé et al. model indicates that pore connectivity was continuously and gradually reduced during sintering and cementationcontrolled diagenesis;

examination of the sintered glass beads and Fontainebleau sandstone experimental data also suggests that poreplugging processes (e.g., pinching off) should be active only during the late stages of sintering and diagenesis, and should be complemented by other processes during the early stages. For example, grain rearrangements occurring in highly porous samples may have lead to the division of large nodal pores into smaller, less connected ones, thus causing the pore coordination number to decrease.
Acknowledgments
This work was partially funded by the Open Fund PLN1501 of the State Key Laboratory of Oil and Gas Reservoir Geology and Exploitation (Southwest Petroleum University, Chengdu, China) and by the U.S. Department of Energy under grant DEFG0997ER14760.
References
 Bernabé Y., Li M., Maineult A. (2010) Permeability and pore connectivity: a new model based on network simulations, J. Geophys. Res. 115, B10203. [CrossRef] [Google Scholar]
 Bernabé Y., Zamora M., Li M., Maineult A., Tang Y.B. (2011) Pore connectivity, permeability and electrical formation factor: a new model and comparison to experimental data, J. Geophys. Res. 116, B11204. [Google Scholar]
 Bernabé Y., Maineult A.(2015) Physics of porous media: fluid flow through porous media, Treatise on Geophysics, 2nd ed., Vol. 11, Resources in the NearSurface Earth, Schubert G. (ed.) , pp. 19–41. [Google Scholar]
 Berryman J.G. (1992a) Effective stress for transport properties of inhomogeneous porous rock, J. Geophys. Res.  Sol. Earth 97, 17409–17424. [Google Scholar]
 Berryman J.G. (1992b) Exact effectivestress rules in rock mechanics, Phys. Rev. A 46, 6, 3307–3311. [CrossRef] [Google Scholar]
 Berryman J.G. (1993) Effectivestress rules for porefluid transport in rocks containing two minerals, Int. J. Rock Mech. Min. Sci. & Geomech. Abstr. 30, 1165–1168. [Google Scholar]
 Biella G., Lozej A., Tabacco I. (1983) Experimental study of some hydrogeophysical properties of unconsolidated porous media, Ground Water 21, 741–751. [CrossRef] [Google Scholar]
 Biella G., Tabacco I. (1981) The influence of grain size on the relations between resistivity, porosity and permeability in unconsolidated formations, Boll. Geof. Teor. Appl. 23, 43–58. [Google Scholar]
 Blair S.C., Berge P.A., Berryman J.G. (1996) Using twopoint correlation functions to characterize microgeometry and estimate permeabilities of sandstones and porous glass, J. Geophys. Res.  Sol. Earth 101, 20359–20375. [CrossRef] [Google Scholar]
 Bourbié T., Zinszner B. (1985) Hydraulic and acoustic properties as a function of porosity in Fontainebleau sandstone, J. Geophys. Res.  Sol. Earth 90, 11524–11532. [Google Scholar]
 Clennell M.B. (1997) Tortuosity: a guide through the maze, Geol. Soc. London, Spec. Pub. 122, 299–344. [CrossRef] [Google Scholar]
 Coelho D., Thovert J.F., Adler P.M. (1997) Geometrical and transport properties of random packings of spheres and aspherical particles, Phys. Rev. E 55, 1959–1978. [CrossRef] [Google Scholar]
 Doyen P.M. (1988) Permeability, conductivity, and pore geometry of sandstone, J. Geophys. Res.  Sol. Earth 93, 7729–7740. [CrossRef] [Google Scholar]
 Fischmeister H.F., Arzt E., Olsson L.R. (1978) Particle deformation and sliding during compaction of spherical powders: a study by quantitative metallography, Powder Metallurgy 4, 179–187. [CrossRef] [Google Scholar]
 Fredrich J.T., Greaves K.H., Martin J.W. (1993) Pore geometry and transport properties of Fontainebleau sandstone, Int. J. Rock Mech. Min. Sci. & Geomech. Abstr. 30, 691–697. [CrossRef] [Google Scholar]
 Glover P.W.J., Hole M.J., Pous J. (2000) A modified Archie’s law for two conducting phases, Earth Planet Sci. Lett. 180, 369–383. [CrossRef] [Google Scholar]
 Glover P.W.J., Zadjali I.I., Frew K.A. (2006) Permeability prediction from MICP and NMR data using an electrokinetic approach, Geophys. 71, F49–F60. [CrossRef] [Google Scholar]
 Glover P.W.J. (2009) What is the cementation exponent? A new interpretation, The Leading Edge 28, 82–85. [CrossRef] [Google Scholar]
 Glover P.W.J., Walker E. (2009) Grainsize to effective poresize transformation derived from electrokinetic theory, Geophys. 74, E17–E29. [CrossRef] [Google Scholar]
 Glover P.W.J., Déry N. (2010) Streaming potential coupling coefficient of quartz glass bead packs: dependence on grain diameter, pore size, and pore throat radius, Geophys. 75, F225–F241. [CrossRef] [Google Scholar]
 Guyon E., Oger L., Plona T.J. (1987) Transport properties in sintered porous media composed of two particle sizes, J. Phys. D: Appl. Phys. 20, 1637–1644. [CrossRef] [Google Scholar]
 Li S.X., Pengra D.B., Wong P.z. (1995) Onsager’s reciprocal relation and the hydraulic permeability of porous media, Phys. Rev. E 51, 5748–5751. [CrossRef] [Google Scholar]
 Lindquist W.B., Venkatarangan A., Dunsmuir J., Wong T.f. (2000) Pore and throat size distributions measured from synchrotron Xray tomographic images of Fontainebleau sandstones, J. Geophys. Res.  Sol. Earth 105, 21509–21527. [Google Scholar]
 Paterson M.S. (1983) The equivalent channel model for permeability and resistivity in fluidsaturated rocks – a reappraisal, Mech. Mater. 2, 345–352. [CrossRef] [Google Scholar]
 Petford N., Davidson G., Miller J.A. (2001) Investigation of the petrophysical properties of a porous sandstone sample using confocal scanning laser microscopy, Pet. Geosci. 7, 99–105. [CrossRef] [Google Scholar]
 Revil A., Cathles L.M. (1999) Permeability of shaly sands, Water Resour. Res. 35, 651–662. [Google Scholar]
 Revil A., Glover P.W.J. (1997) Theory of ionic surface electrical conduction in porous media, Phys. Rev. B 55, 1757–1773. [CrossRef] [Google Scholar]
 Revil A., Glover P.W.J. (1998) Nature of surface electrical conductivity in natural sands, sandstones and clays, Geophys. Res. Lett. 25, 691–694. [CrossRef] [Google Scholar]
 Revil A., Kessouri P., TorresVerdin C. (2014) Electrical conductivity, induced polarization, and permeability of the Fontainebleau sandstone, Geophys. 79, D301–D318. [CrossRef] [Google Scholar]
 Roberts J.N., Schwartz L.M. (1985) Grain consolidation and electrical conductivity in porous media, Phys. Rev. B 31, 5990–5997. [CrossRef] [Google Scholar]
 Schwartz L.M., Kimminau S. (1987) Analysis of electrical conduction in the grain consolidation model, Geophys. 52, 1402–1411. [CrossRef] [Google Scholar]
 Walsh J.B., Brace W.F. (1984) The effect of pressure on porosity and the transport properties of rock, J. Geophys. Res.  Sol. Earth 89, 9425–9431. [Google Scholar]
 Wong P.z., Koplik J., Tomanic J. (1984) Conductivity and permeability of rocks, Phys. Rev. B 30, 6606–6614. [CrossRef] [Google Scholar]
 Yale D.P. (1984) Network modeling of flow, storage and deformation in porous rocks, Ph.D. Thesis, Stanford University. [Google Scholar]
 Zhu W., David C., Wong T.f. (1995) Network modeling of permeability evolution during cementation and hot isostatic pressing, J. Geophys. Res. 100, 15451–15464. [CrossRef] [Google Scholar]
 Zhu W., Evans B., Bernabé Y. (1999) Densification and permeability reduction in hotpressed calcite: a kinetic model, J. Geophys. Res. 104, 25501–25511. [CrossRef] [Google Scholar]
Appendix
Appendix A: Hydraulic Radius
Unconsolidated packing of monodisperse spheres: We consider a volume V_{0} containing n spheres. The porosity and hydraulic radius are given by:(A1)and(A2)
We can calculate V_{0} from Equation (A1) and plug it into (A2). Since V_{0} is proportional to the number of spheres, n cancels out and Equation (7) is obtained. The same derivation applies to sand grains, except that the surface of the grains must be written n f_{G} 4π R_{G}^{2}, where f_{G} is the asphericity factor.
Sintered glass beads: After sintering the glass beads volume V_{0} is reduced to V_{snt}. Since the volume of solid does not change during sintering, the porosity of the sintered glass beads can be written as:(A3)
Similarly, the hydraulic radius is given by:(A4)
where A_{G} is the average surface area of a single sintered sphere. A sintered sphere can approximately be represented as a truncated sphere with z spherical caps removed (Roberts and Schwartz, 1985; Schwartz and Kimminau, 1987), where z is the average number of grain contacts per sphere, which can be estimated as z = 6 ϕ^{−1/5} based on the experimental data of Fischmeister et al. (1978) on plastically compressed packings of monodisperse bronze spheres. Because the grain contacts are not all formed at the same time during sintering, we must evaluate the number dz of new grain contacts forming when the porosity decreases from ϕ* to ϕ* − dϕ* (where * denotes a running porosity variable, between the initial and final values, ϕ_{0} and ϕ < ϕ_{0}, respectively) and the depth h* reached by the corresponding spherical caps at the end of the sintering process, i.e., when ϕ* approaches ϕ. By differentiation we obtain:(A5)
The linear dimension of the volume of the sintered spheres must be proportional to the distance between contacting spheres, implying:(A6)
where only the new grain contacts formed at ϕ* are considered. Remembering that the surface of a spherical cap of depth h* is 2πR_{G}h*, we obtain A_{G} equal to:(A7)
Appendix B: Hydraulic and Electrical Conductance of Constricted Pipes
We consider a frustumshaped pipe of length l_{c}. We call r_{1} and r_{2} > r_{1} the radii of the narrow and wide bases, respectively, and x the coordinate along the pipe axis. We assume that fluid flow through this pipe is laminar and obeys Reynold’s equations, i.e., the locally expressed Poiseuille law is valid for any value of x, a condition most likely met in long pipes. The flow rate q is constant along the pipe and given by:(B1)
where p(x) is the fluid pressure at x, r(x) = r_{1} + x(r_{2}−r_{1})/l_{c} the crosssectional radius and η the fluid viscosity. We can thus write:(B2)
After double integration of Equation (B2), we find:(B3)
where Δp is the absolute value of the apparent pressure difference between the pipe faces. A similar expression can be written for a straight, cylindrical pipe of radius r = r_{H} and length l:(B4)
Thus, the ratio g_{H} of the frustum pipe hydraulic conductance to that of a cylindrical pipe is:(B5)
The electrical resistance R_{c} of a frustum pipe saturated with a solution of resistivity ρ_{s} is:(B6)
Likewise, the resistance of a cylindrical pipe is equal to lρ_{s}/(πr_{H}^{2}) and the ratio g_{E} of the electrical conductances is:(B7)
Cite this article as: Y. Bernabé, M. Li, Y.B. Tang and B. Evans (2016). Pore Space Connectivity and the Transport Properties of Rocks, Oil Gas Sci. Technol 71, 50.
All Tables
Approximate expressions of the exponents and prefactors of Equations (46) as functions of σ and ε (from Bernabé et al., 2011)
The transforms used to complete the Fontainebleau sandstone datasets of Doyen (1988), Fredrich et al. (1993), Revil et al. (2014) and Zamora (unpublished data)
Summary of the values assigned to the input parameters of the Bernabé et al.’s models
All Figures
Figure 1. Compilation of published experimental data, a) k versus ϕ, b) F versus ϕ and c) k versus F, in unconsolidated monodisperse glass beads (Glover and Walker, 2009; Glover and Déry, 2010; black dots), wellsorted, rounded sands (Biella and Tabacco, 1981; Biella et al., 1983; blue dots), sintered glass beads (Wong et al., 1984; Guyon et al., 1987; Li et al., 1995; Blair et al., 1996; purple, red and orange squares corresponding to different grain sizes as indicated in the inset) and Fontainebleau sandstone (Doyen, 1988; Fredrich et al., 1993; Zamora, unpublished data; black diamonds; Revil et al., 2014; grey diamonds). The colored lines are inserted to highlight powerlaw trends in the various materials (note that the three families of sintered glass obey fairly similar power laws). The values of the exponents are given when necessary. 

In the text 
Figure 2. Influence of the hydraulic radius on a) permeability and b) formation factor. The diagrams use the same symbol and line conventions as in Figure 1. Three horizontal error bars (dark purple) help visualize the relative uncertainty δr_{H}/r_{H} determined for the various materials (i.e., ±10, ±15 and ±20% for the unconsolidated sands, sintered glass beads and Fontainebleau sandstone, respectively). 

In the text 
Figure 3. Comparison of calculated a) k from Equations (1, 4), b) F from Equations (2, 5) and c) k from Equations (3, 6) to the observed values of permeability and formation factor in the unconsolidated glass beads and sands. The meaning of the different symbols is explained in the inset, where k(z) refers to Bernabé et al. Equation (4), F(z) to Equation (5) and k(F) to Equation (6). The vertical error bars (dark purple) represent the tolerances estimated for the Bernabé et al. model, δk/k ≈ ±30% and δF/F ≈ ±10%. 

In the text 
Figure 4. Comparison of calculated a) k from Equations (1, 4), b) F from Equations (2, 5) and c) k from Equations (3, 6) to the observed values of permeability and formation factor in the sintered glass beads. The meaning of the different symbols is explained in the inset, where k(z) refers to Bernabé et al. Equation (4), F(z) to Equation (5) and k(F) to Equation (6). The vertical error bars (dark purple) represent minimal estimates of the model tolerances (i.e., δk/k ≈ ±30% and δF/F ≈ ±10% for the PWB model as well as the Bernabé et al. model applied to high porosity samples; δk/k ≈ ±50% and δF/F ≈ ±20% for the Bernabé et al. model in low porosity samples). 

In the text 
Figure 5. Comparison of calculated a) k from Equations (1–4), b) F from Equations (2, 5) and c) k from Equations (3, 6) to the observed values of permeability and formation factor in the Fontainebleau sandstone. The meaning of the different symbols is explained in the inset, where k(z) refers to Bernabé et al. Equation (4), F(z) to Equation (5) and k(F) to Equation (6) and where the digits 1 and 2 indicate the Fontainebleau 1 and 2 datasets. The difference between the two datasets is particularly visible in the diagrams b) and c). The vertical error bars (dark purple) represent minimal estimates of the δk/k and δF/F model tolerances (i.e., δk/k ≈ ±40% and δF/F ≈ ±10% for the PWB model; δk/k ≈ ±60% and δF/F ≈ ±20% for the Bernabé et al. model in high porosity samples and δk/k ≈ ±90% and δF/F ≈ ±30% in low porosity samples). 

In the text 
Figure 6. The hydraulic and electrical tortuosities τ_{h}^{2} and τ_{e}^{2}versus porosity for the sintered glass beads and Fontainebleau sandstone. The meaning of the different symbols is explained in the inset, where H and E stand for hydraulic and electrical. For the sake of visibility, the Fontainebleau 2 dataset is not included in this diagram but the associated trends are comparable to the Fontainebleau 1 trends represented here. 

In the text 
Figure 7. a) Schematic representation of the node splitting mechanism. The cartoon is twodimensional but it is meant to illustrate a threedimensional arrangement of sintered grains and pores between them. The large blue dots indicate the center of the nodes and the blue arrows symbolize the (out of plane) channels emanating from them. The large red arrow indicates the evolution associated with sintering. b) Example of graph transformations that leave the graph topological invariant (genus) unchanged while reducing the mean coordination number. In this cartoon, extra node/branch combinations (represented in red) are inserted at different places in the graph without altering the genus. By definition, nodes must always be present at the extremities of a branch (the dotted segments are intended to symbolize the unspecified continuation of the graph). 

In the text 