1 Introduction

A slightly compressible fluid injected into a saturated rock is partly transported through the connected pore space following the prevailing pore-pressure gradient and partly stored owing to its compressibility and that of the pore space. Thus, to describe the complete hydraulic behavior of a material, it is necessary to determine a transport property, as well as a storage property. Early experiments relied on the determination of permeability investigating steady-state flow at constant pressure gradients (e.g., Bernabé 1987; Darcy 1856). Since these methods are less suitable for low-permeability rocks, it has now been practice for many decades to determine hydraulic properties of rocks from pressure transients recorded during field and laboratory tests, such as pulse tests (Bauer et al. 1995; Brace et al. 1968; Bredehoeft and Papadopulos 1980; Hsieh et al. 1981; Neuzil 1982; Selvadurai and Carnaffan 1997; Zeynaly-Andabily and Rahman 1995), slug tests (e.g., Cooper et al. 1967; Butler 1997), or experiments with a constant flow rate (Lin 1977; Trimmer et al. 1980; Song et al. 2004). A basic disadvantage of these methods is their requirement for constant boundary conditions of the investigated rock volume during the experiment, regarding in particular temperature, posing a challenge for temperature control in the laboratory (e.g., Bernabé 1987). Initially developed in physical chemistry (Turner 1958, 1959) and then applied to rock samples (Stewart et al. 1961; Kranz et al. 1990; Fischer 1992), the axial oscillatory pore-pressure method is insensitive to changes in boundary conditions for the samples and permits determination of permeability k and specific storage capacity s in a single experiment when the period is appropriately chosen. As an extension, Song and Renner (2007) introduced the axial oscillatory pore-flow method. In the current laboratory experiments, we also employ a radial oscillatory pore-flow method analogous to periodic-pumping tests in wells (Cardiff et al. 2013; Rabinovich et al. 2015; Rasmussen et al. 2003; Renner and Messar 2006).

For the inversion of pressure and/or flow-rate transients to hydraulic parameters, all currently used laboratory methods rely on the assumption that the medium investigated is homogeneous and isotropic. Therefore, they actually provide equivalent or effective hydraulic properties, i.e., the properties of a hypothetical homogeneous and isotropic medium that would give the same response to the imposed hydraulic perturbation as the observed one. Frequency dependence of effective permeability, gained from applying the axial pore-pressure method at a range of frequencies, was previously interpreted to indicate sample heterogeneity (Bernabé et al. 2004; Fischer et al. 1992; Song and Renner 2006a). Frequency dependence is, however, also related to a transition from laminar to turbulent flow or to nonlinear pressure sensitivity of sample properties. Also, the role of anisotropy has not been analyzed in detail.

The relation between pore-space geometry and effective properties has frequently been studied with a range of methods, especially for Fontainebleau sandstone that exhibits a natural variability in porosity from about 3 to 30 % (e.g., Bourbie and Zinszner 1985; David and Darot 1989; Doyen 1988; Fredrich et al. 1993; Song and Renner 2006b, 2007). Recently, David et al. (2018b) performed a measurement exercise for Grimsel granodiorite with a rather low permeability of \(\approx\) 10\(^{-18}\) m\(^2\) involving 24 laboratories using three different methods to measure permeability (steady-state flow method, transient pulse method, and axial oscillatory pore-pressure method) complemented by in-depth microstructural analyses David et al. (2018a). This laboratory inter-comparison was, however, performed on different samples rather than in the form of a round robin; thus, the results may be affected by sample-to-sample variability, be it natural or related to preparation.

The objective of this study is to test to what extent information about the structure of a rock’s pore space can be obtained by combining the different oscillatory methods, and varying Terzaghi’s effective pressure and the oscillation period and thus the nominal penetration depth. We selected two rocks, Westerly granite (WG) and Wilkeson sandstone (WS), covering connected porosity from 1 to 10 % and representing pore spaces from micro-fracture dominated to characteristic for granular media. Three of the four investigated methods, the axial pore-pressure and flow method and the steady-state flow or Darcy method, were applied on one sample of the two rocks, and thus, the method comparison is not affected by sample-to-sample variability. Only the radial flow method required preparation of separate samples.

2 Methods

2.1 Material

Wilkeson sandstone (WS) forms part of the Pudged Group of the Mid Miocene (Gard 1968). The used block originates from an outcrop in Pierce County, Washington (USA). The light gray to brown, compact, poorly sorted sandstone has an average composition of 50 % quartz and 35 % feldspar, and approximately 10 % porosity (Ahrens et al. 2018; Arndt and Kent 1980; Duda 2011; Moen 1967). The pore space of WS consists of equant pores between grains but also inter- and intragranular microcracks (Ahrens et al. 2018). The grains exhibit an average size of about 0.2 mm as determined by the line-intercept method. Plagioclase grains are partly altered to sericite and saussurite. Some quartz grains contain fluid inclusions and show undulatory extinction as well as recrystallization features. Muscovite grains are partly kinked (Ahrens et al. 2018). Ultrasonic p-wave velocity measurements on the block, from which samples were cored, indicate an anisotropy of about 10 %.

Westerly Granite (WG) originates from the southeast of Rhode Island (USA). The granite is composed of 41 % plagioclase, 25.2 % quartz, 26.3 % potassium feldspar, 6.3 % biotite, and 1.2 % other constituents (Moore 1993). The samples used here have an average grain size of about 0.4 mm, determined with the line-intercept method. Our measurements of ultrasonic p-wave velocity indicate isotropy. The porosity of WG is dominated by microcracks, as typically acquired by a magmatite during cooling and decompression after solidification due to anisotropic elastic and thermal properties of the components.

2.2 Preparation and Physical Characterization

Cylindrical samples with a diameter of 30 mm were diamond-drilled from cubic blocks of the two rocks with an edge length of about 250 mm. Samples were sawn to a length of approximately 78 mm for axial and 45 to 48 mm for radial tests. The end faces were ground square within a tolerance of 0.03 mm to reach lengths of 75 mm and 40 to 45 mm for the axial and radial measurements, respectively. For the latter, a hole with a diameter of 3 mm was drilled into the center of the sample in axial direction. Drilling, sawing, and grinding were carried out with water as the cooling medium. All prepared samples were characterized by basic physical properties, such as bulk density, matrix density, total and connected porosity, and P- and S-wave velocity.

Bulk density \(\rho _{\mathrm{{geo}}}\) of specimens was determined from their masses after drying at 60 \(^\circ\)C and geometrical volume determination. Matrix density \(\rho _{\mathrm{{m}}}\) was gained from pycnometer measurements on powders produced by crushing and grinding pieces of the rocks. Total porosity is calculated from the densities according to \(\phi _{\mathrm{{tot}}} = 1 - (\rho _{\mathrm{{geo}}}/\rho _{\mathrm{{m}}})\). Connected porosity \(\phi _{\mathrm{{con}}}\) was determined from the difference in mass before and after evacuation and saturation with distilled water (for more details see Duda and Renner 2013). Calculated total porosity is consistently larger than connected porosity, but differences are barely significant considering the uncertainties of the two independent constraints on pore volume (see “Appendix A.1”).

2.3 Experimental Setup and Procedure

The experimental setup consists of a pore-pressure system, a confining-pressure system, and a computer-assisted control unit. Pore pressure \({p _{\mathrm {p}}}\) and confining pressure \({p _{\mathrm {c}}}\) are controlled by two identical intensifiers via servohydraulic valves. Two rubber tubes surround the cylindrical sample, saturated with distilled water, to separate the two fluid systems for pore pressure (distilled water) and confining pressure (oil, for details see Ahrens et al. 2018). In axial flow experiments (Fig. 1a), the sample connects two reservoirs in the pore-pressure system, downstream and upstream, where the pressure perturbation is exerted, both equipped with a pressure transducer. Porous metal disks (SIKAR S) with a thickness of 3 mm, a porosity of 31 %, a permeability larger than 10\(^{-12}\) m\(^2\), and pore diameters in the range of 7 to 14 µm cover both ends of the specimen to serve two purposes. The porous disks ensure that the pore fluid is uniformly distributed over the sample ends and that the confining pressure is transferred to the sample in axial direction. In radial flow experiments (Fig. 1b), the axial hole is only connected to an upstream reservoir. The two sample ends are equipped with a self-adhesive Teflon foil to prevent fluid penetration between the bottom plug and the sample, and the top spacer and the sample (see “Appendix A.2”).

The rubber tubes used for jacketing have to be stretched substantially to fit on the assembly and thus seat tightly on sample and end plugs, even when no confining pressure is applied, preventing any water loss during installation. Capillary forces tend to hinder fluid to drip out of the fluid-filled capillary tubes of the end-plugs in the vessel before they are connected to the external parts of the pore-pressure system. Furthermore, we minimized evaporation or other fluid losses by following a standardized routine of rapid assemblage. Once the assembly was placed in the vessel and all pore-pressure tubing was connected until an immediate pressure response to even small volumometer movements demonstrated the absence of trapped air in the tubing. The employed protocol minimizes, if not totally removes, trapped air in the system and since the solubility of air in water is substantial at the elevated pore pressures applied in the experiments, single-phase flow in the entire system, including the sample, is reached at run conditions even if minute amounts of air were initially trapped somewhere in the system at ambient pressure.

The assembly was first subjected to a confining pressure of 10 MPa. Subsequently, confining pressure and pore pressure were gradually increased to 20 MPa and 10 MPa, respectively. When the movement of the pore-pressure intensifier piston ceased after the pressurization indicating equilibration of pore pressure in the entire system, the bypass valve between the two reservoirs was closed for axial flow experiments. Finally, harmonic pressure changes in the upstream reservoir were applied with oscillation periods of 5, 10, 20, 50, 100, and 333 s. During a test sequence, confining pressure was increased to 40, 70, and 110 MPa, i.e., a sequence consisted of measurements at mean Terzaghi’s effective pressures, the difference between confining pressure or mean stress and pore-fluid pressure here denoted as \(\Delta p\), of 10, 30, 60, and 100 MPa. Permeability was also determined using the Darcy method (Darcy 1856) at confining pressures of 11, 31, 61, and 101 MPa. The upstream pressure was controlled to a constant value of 1 MPa, while the downstream reservoir remained at ambient pressure.

2.4 Determination of Hydraulic Parameters

The oscillatory methods rely on an evaluation of amplitude ratio \(\delta\) and phase shift \(\varphi\) between two harmonic signals for the determination of the effective hydraulic properties, permeability, and specific storage capacity. For the axial methods, the evaluation is based on an analytic solution of the one-dimensional diffusion equation for a homogeneous, isotropic, porous, and finite medium (e.g., Turner 1958; Stewart et al. 1961; Bennion and Goss 1977; Kranz et al. 1990; Bernabé et al. 2006; Song and Renner 2007). The radial method relies on an analytic solution of the diffusion equation for radial symmetry with a no-flow boundary at a finite distance \(r=R\) from the borehole with radius \(r_{\mathrm {i}}\) (Renner and Messar 2006) that gives phase shift and amplitude ratio between flow rate and pressure at the borehole wall as

$$\begin{aligned} {\varphi _{\mathrm{{rQp}}}}= & {} \arg \left( \hat{\eta } \frac{I_{1}(\hat{\eta }R)K_{1}(\hat{\eta }r_{\mathrm {i}})-I_{1}(\hat{\eta }r_{\mathrm {i}})K_{1}(\hat{\eta }R)}{K_{1}(\hat{\eta }R)I_{0}(\hat{\eta }r_{\mathrm {i}})+K_{0}(\hat{\eta }r_{\mathrm {i}})I_{1}(\hat{\eta }R)} \right) ~{\mathrm {and}} \end{aligned}$$
(1)
$$\begin{aligned} {\delta _{\mathrm{{rQp}}}}= & {} \frac{{2\pi r_{\mathrm {i}} k h}}{{\mu }} \left| \hat{\eta } \frac{I_{1}(\hat{\eta }R)K_{1}(\hat{\eta }r_{\mathrm {i}})-I_{1}(\hat{\eta }r_{\mathrm {i}})K_{1}(\hat{\eta }R)}{K_{1}(\hat{\eta }R)I_{0}(\hat{\eta }r_{\mathrm {i}})+K_{0}(\hat{\eta }r_{\mathrm {i}})I_{1}(\hat{\eta }R} \right| \end{aligned}$$
(2)

with \(\hat{\eta }=(i\omega /D)^{0.5}\), fluid viscosity \(\mu\), permeability k, sample height h, angular frequency \(\omega\), and \(I_{0}\), \(I_{1}\), \(K_{0}\) and \(K_{1}\) denoting the modified Bessel functions of the first and second kind of zero and first order, respectively.

Fourier analysis (FFT) was performed to determine amplitude ratio \(\delta\) and phase shift \(\varphi\) between recorded pressure and flow signals using a sliding window of three periods and a step width of one period over a data set of 10 periods recorded for each step in confining pressure. Fluid flows into or out of the sample when pressure in the upstream reservoir increases or decreases resulting in characteristic phase shift and amplitude ratio between upstream-pore pressure and involved fluid volume, evaluated in the case of pore-flow methods. The recorded fluid-volume signal, however, is not completely harmonic due to technical limitations. When the pressure-intensifier piston changes its direction of travel, an immediate change in the upstream pressure occurs without detectable movement of the piston, likely related to slight reordering of the sealing package and tilting of the piston. We address this effect as hysteresis and correct for it (for details of the procedure see “appendix A.3”) before the FFT is carried out.

It is not possible to analytically invert the solutions of the one-dimensional diffusion equation for hydraulic parameters. Instead, dimensionless permeability \(\eta\) and specific storage capacity \(\xi\) corresponding to the observed amplitude ratio and phase shift are determined using a grid search. Hydraulic parameters are calculated from the found parameters according to the relations \(\eta =AT/(\pi \mu h S_{\mathrm {d}})\) and \(\xi =Ahs/S_{\mathrm {d}}\) (see Bernabé et al. 2006), with the storage capacity of the downstream reservoir \(S_{\mathrm {d}}\) and the cross-sectional area of the sample A. It is also not possible to analytically invert the solutions of the radial diffusion equation for hydraulic parameters. One has to perform numerical calculations (details in Renner and Messar 2006) for two intermediate parameters, the argument \(\hat{\eta }\) and the modulus of the Bessel functions involved in (1, 2).

In addition to permeability k and specific storage capacity s, we report hydraulic diffusivity \(D=k/(s\mu )\) and apparent porosity \(\phi _{\mathrm{{ap}}}\), defined as the ratio between specific storage capacity and fluid compressibility c\(_{\mathrm {f}}\)

$$\begin{aligned} \phi _{\mathrm{{ap}}}=\frac{s}{{{c_{\mathrm{{f}}}}}} = {\phi _{\mathrm{{con}}}}\left( {1 + \frac{{{c_{\mathrm{{pp}}}}}}{{{c_{\mathrm{{f}}}}}}} \right) , \end{aligned}$$
(3)

where we follow the notation of Zimmerman et al. (1986) for the various compressibilities. This description reflects the two contributions to specific storage capacity s, the compressibility of the pore-fluid \({c _{\mathrm {f}}}\) and the pore compressibility in response to pore-pressure variations, \({c _{\mathrm {pp}}}\). For a homogeneous sample, the apparent porosity should exceed the connected porosity because the compressibility \(c_{\mathrm{{pp}}}\) in (3) is always positive.

Fig. 1
figure 1

Setup of the oscillatory pore pressure and pore-flow tests. For the axial tests, porous metal plates are installed at the sample’s end faces; in the radial tests, the downstream reservoir is separated from the system with an impermeable spacer and Teflon foil, the upstream is equipped with an impermeable spacer and Teflon foil both with a hole leading the fluid into the samples borehole

Full size image
Fig. 2
figure 2

Dependence of permeability on Terzaghi’s effective pressure for a Wilkeson sandstone and b Westerly granite for all tested periods. Results from previous studies of Duda (2011) (D11), Interlab-test (Int, see acknowledgement), David et al. (1994) (D94) and Brace (1965) (B65) using axial pore-pressure method (axpp), axial pore-pressure method (axpp), constant-flow rate tests, and hydraulic pulse tests, respectively, are shown for comparison. The legends in (a) and (b) apply to both plots. (axpp: axial pore-pressure method, axVp: axial pore-flow method, rQp: radial pore-flow method)

Full size image
Fig. 3
figure 3

Dependence of permeability and specific storage capacity on oscillation period for Wilkeson sandstone (a, b) and Westerly granite (c, d) at a Terzaghi’s effective pressure of 60 MPa determined by the four employed methods. The legend applies to all subplots. Error bars reflect the error-propagation analysis (see A.4). Where not visible, the error bars do not exceed symbol size. (axpp: axial pore-pressure method, axVp: axial pore-flow method, rQp: radial pore-flow method)

Full size image
Fig. 4
figure 4

Permeability moduli as a function of period applied in the oscillatory tests (axpp: axial pore-pressure method, axVp: axial pore-flow method, rQp: radial pore-flow method) for a Wilkeson sandstone and b Westerly granite in comparison with results of the Darcy tests (black lines). Uncertainties of the performed linear regression are represented by the vertical error bars (oscillatory methods) and the dashed horizontal lines (Darcy method). Where not visible, the error bars do not exceed symbol size

Full size image

The standard deviations of amplitude ratio and phase shift gained from the sliding-window analyses serve as uncertainty measures that we propagate to uncertainties of the dimensionless parameters. The uncertainty of all hydraulic parameters reflects a full error propagation of the relevant uncertainties (see “Appendix A.4” for details of the error analysis).

As previous studies showed, the resolution of the storage capacity for the axial pore-pressure method is limited for certain experimental conditions (Bernabé et al. 2006; Rutter and Faulkner 1996). The limits of the solution space are defined by the limits of the downstream reservoir size \({S_{\mathrm{{d}}}} = 0\) and \({S_{\mathrm{{d}}}} \rightarrow \infty\), to which limits of the dimensionless specific storage capacity of \(\xi \rightarrow \infty\) and \(\xi = 0\) correspond. In the range of these limits and for vanishing phase shift, \(\eta\) and \(\xi\)-isolines lie arbitrarily close to each other forming a “tail” of the solution space, where the resolution is low (Song and Renner 2007). For the axial pore-flow method, the solution space bends over for phase shifts smaller than \(-\)0.125 periods and the solutions are no longer unambiguous. Results of the axial pore-pressure method with a phase shift of more than \(-\)0.25 periods and of the axial pore-flow method with a phase shift smaller than \(-\)0.125 periods are thus not taken into account for further processing due to their uncertainty (see “Appendix” A.3). Typically, a specific oscillation period does not simultaneously lead to poor resolution of the two methods. The signal characteristics observed for the Westerly granite samples correspond to a small dimensionless storage capacity in both solution spaces, indicating that the pore space of Westerly granite is only 0.01 to 0.1 % of the size of the entire measuring system. For such an unfavorable ratio, a high uncertainty of hydraulic parameters results, in particular for specific storage capacity, but for periods that avoid the critical regions of the solution spaces (Song and Renner 2007).

Hydraulic properties gained from the axial methods are here presented as a function of Terzaghi’s effective pressure, i.e., \({\Delta } p\)=\({p _{\mathrm {c}}}\)-\({p _{\mathrm {p}}}\). For a given confining pressure, the stress states, to which the solid cylinders and the thick-walled hollow cylinders are subjected in axial testing and radial testing, respectively, differ by the stress concentration due to the central bore. At the borehole wall, an anisotropic stress state prevails with principal normal stresses \(\sigma _{\mathrm {ax}}=p_{\mathrm {c}}\) (axial), \({p _{\mathrm {p}}}\) (radial), and \(2p _{\mathrm {c}}\) (tangential), corresponding to a mean stress of \({\sigma _{\mathrm {m}}}\)=\({p _{\mathrm {c}}}\)+\({p _{\mathrm {p}}}/3\). We replace confining pressure in the calculation of Terzaghi’s effective pressure for solid samples by the mean stress for the calculation of the “pressure” difference for radial-flow symmetry, i.e., \({\Delta } p\)=\({p _{\mathrm {c}}}\)-\({2 p _{\mathrm {p}}/3}\).

2.5 Mechanical Properties Related to Hydraulic and Ultrasonic Measurements

For the full description of the volumetric mechanical response of an isotropic, homogeneous, linear poro-elastic medium, three mechanical parameters are needed out of the many commonly used ones, for example, the drained (also dry or skeleton) bulk modulus \(K_{\mathrm {d}}\), the undrained bulk modulus \(K_{\mathrm {ud}}\), the average mineral modulus \(K_{\mathrm {s}}\), and the fluid modulus \(K_{\mathrm {f}} (= c_{\mathrm {f}}^{-1})\). The ultrasonic measurements on dry and saturated samples nominally constrain \(K_{\mathrm {d}}\) and \(K_{\mathrm {ud}}\) at ambient conditions, respectively. The drained modulus at elevated pressure can be determined from specific storage capacity s obtained from the hydraulic testing according to

$$\begin{aligned} \frac{1}{{{K_{\mathrm {d}}}}} = s - {\phi _{{\mathrm {con}}}}{c_{\mathrm {f}}} + \left( {1 + {\phi _{{\mathrm {con}}}}} \right) \frac{1}{K_{\mathrm {s}}} \end{aligned}$$
(4)

(e.g., Zimmerman et al. 1986; Kümpel 1991). For the evaluation of equation (4), we estimated the average mineral modulus \(K_{\mathrm {s}}\) from the rocks’ composition and literature values for the constituting minerals (Gebrande 1982) to \({K_{\mathrm{{s,WG}}}} = 57.0\pm 6.7\;{\mathrm {GPa}}\) for WG and \({K_{\mathrm {s,WS}}} = 33.0\pm 4.3\;{\mathrm {GPa}}\) for WS. Water compressibility at the test conditions, i.e., room temperature and 10 MPa, was determined with the software FLUIDCAL (Wagner and Pruß 2002) to \({c_{\mathrm {f}}} = 4.474 \cdot {10^{ - 10}}\;{\mathrm {Pa}}^{-1}\).

We express the decrease in permeability with Terzaghi’s effective pressure by the permeability modulus \(K_{\mathrm {k}}\) defined as

$$\begin{aligned} K_{\mathrm {k}} =\frac{\mathrm {d}{\Delta } p}{\mathrm {d} \ln k} \end{aligned}$$
(5)

(e.g., Riepe et al. 1983; Debschutz et al. 1989; Rice 1992; David et al. 1994; Luo and Feng 2013). A small modulus corresponds to a strong pressure dependence and vice versa.

3 Results

3.1 Material Dependence and Effect of Terzaghi’s Effective Pressure

Independent of method, the results reflect that the investigated rocks represent what one would loosely address as a moderately permeable sandstone (WS) and a low-permeability granite (WG). The hydraulic properties of WS range from \(k_{\mathrm{{WS}}}={10^{ - 17}}\) to \({10^{ - 16}}\) \({\mathrm {m}}^2\) and \(s_{\mathrm {WS}}={10^{ - 11}}\) to \({10^{ - 9}}\) Pa\(^{-1}\) and those of WG from \(k_{\mathrm {WG}}={10^{ - 19}}\) to \({10^{ - 17}}\) \({\mathrm {m}}^2\) and \(s_{\mathrm {WG}}={10^{ - 12}}\) to \({10^{ - 10}}\) Pa\(^{-1}\) (Figs. 2 and 3). Corresponding hydraulic diffusivity values range from \(D_{\mathrm {WS}}={5 \cdot 10^{ - 5}}\) to \(D_{\mathrm {WS}}={9 \cdot 10^{ - 4}}\) \({\mathrm {m}}^2/{\mathrm {s}}\) for WS and \(D_{\mathrm {WG}}={1 \cdot 10^{ - 5}}\) to \(D_{\mathrm {WG}}={4 \cdot 10^{ - 3}}\) \({\mathrm {m}}^2/{\mathrm {s}}\) for WG. The quoted ranges encompass the variability with experimental parameters and the used method. All methods give a negative correlation between permeability and Terzaghi’s effective pressure (Fig. 2). Permeability moduli, determined after (5), of WS are larger than those of WG regardless of the used method (Fig. 4), i.e., the sandstone exhibits a lower sensitivity of permeability to Terzaghi’s effective pressure than the granite.

3.2 Effect of Method and Oscillation Period

Permeability and specific storage capacity show method-dependent differences for the two tested rock types to different degrees (Fig. 5). Of all determined hydraulic parameters, permeability values for Wilkeson sandstone exhibit the least variability with applied method. Permeability values gained by the axial pore-pressure method exhibit close agreement with the results of the Darcy tests for both rocks (Figs. 2, 3 and 5). The results of the two flow methods (axVp, rQp) for permeability and specific storage capacity typically exceed those from the axial pore-pressure method, with the radial pore-flow method yielding the largest values. The difference, however, remain rather small for WS, i.e., within an order of magnitude for permeability and up to one and a half orders of magnitude for storage capacity , while for WG \(k_{\mathrm {axpp}}\) and \(k_{\mathrm {rQp}}\) differ by up to two orders of magnitude, and \(s_{\mathrm {axpp}}\) and \(s_{\mathrm {rQp}}\) even by three orders of magnitude (Figs. 3 and 5). Hydraulic diffusivity values agree well for Westerly granite but less so for Wilkeson sandstone for the two axial methods; the radial method tends to give smaller values than the other two (Fig. 5).

Permeability and specific storage capacity of WS are almost insensitive to period for all methods, with the exception of \(s_{\mathrm {rQp}}\) (Fig. 3a, b). In contrast, the permeability of WG decreases significantly with increase in period. The specific storage capacity gained by the axial pore-pressure method increases with increase in period while for the flow methods it remains constant (Fig. 3c, d).

The central borehole of the radially tested samples is associated with a stress concentration and an anisotropic stress state at the borehole wall that—in conjunction with the pressure dependence of the hydraulic parameters—could potentially contribute to differences between the investigated methods. Yet, even for the granite samples, exhibiting a stronger pressure sensitivity than the sandstone samples (Figs. 2 and 4), the modest difference in mean stress at the “upstream” side between samples with a central bore and solid samples, tested at identical confining and fluid pressure, amounts to less than about 3 MPa, insufficient to explain the difference between permeability gained by the axial pore-flow and the radial pore-flow method when invoking the determined permeability moduli (Fig. 4). Moreover, the observed difference between results from axial and radial flow is opposite to the effect of pressure difference.

Fig. 5
figure 5

Permeability (a, b), specific storage capacity (c, d) and hydraulic diffusivity (e, f) of Wilkeson sandstone and Westerly granite determined by axial and radial pore-flow methods (axVp and rQp) in relation to axial pore-pressure method (axpp). The green lines give axial pore-pressure method (axpp) in relation to itself (solid) and the corresponding errors (dotted)

Full size image
Fig. 6
figure 6

Drained bulk modulus of a Wilkeson sandstone and b Westerly granite determined from ultrasonic measurements at ambient pressure and storage capacity as a function of Terzaghi’s effective pressure. Results from Duda (2011) (D11) and Ahrens et al. (2018) (A18) for axial pore-pressure method (axpp), triaxial deformation tests and hydrostatic deformation tests on dry samples for Wilkeson sandstone and Brace (1965) (B65) from linear compressibility measurements for Westerly granite are shown for comparison. The legend applies to all subplots. (axpp: axial pore-pressure method, axVp: axial pore-flow method, rQp: radial pore-flow method)

Full size image
Fig. 7
figure 7

Apparent porosity (see Eq. 3) as a function of Terzaghi’s effective pressure for a Wilkeson sandstone and b Westerly granite in comparison with connected porosity as determined at ambient pressure (dashed horizontal line). Legend applies to (a) and (b). (axpp: axial pore-pressure method, axVp: axial pore-flow method, rQp: radial pore-flow method)

Full size image
Fig. 8
figure 8

Relation between normalized penetration depth (see Sect. 4.3) and pressure perturbation as expressed by amplitude ratio \(\delta _{\mathrm {axpp}}\) observed in axial oscillatory pore-pressure tests. The scaling relation for homogeneous materials between the characteristic duration of a pressure perturbation and its penetration depth is indicated by the red shaded areas, where we considered an amplitude ratio of 0.1 as the limit for a significant perturbation

Full size image

4 Discussion

Frequency dependence of effective permeability is known to result from a switch in dominance of forces affecting the fluid motion in conduits; with increase in flow rate, fluid inertia becomes relevant for the fluid-velocity profile in the individual conduits manifesting itself in the transition from frequency-independent intrinsic or Darcy permeability \(k_{\mathrm {D}}\) at low frequencies to so-called dynamic permeability decreasing with increase in frequency [e.g., (Johnson et al. , 1987; Pazdniakou and Adler , 2013)]. For tubes, the Poiseuille-type profile, typical for the dominance of viscous forces, gives way to rather flat profiles when inertia forces dominate. The spatial scaling parameter, viscous skin depth \(\delta _{\mathrm {vis}} = \sqrt{2\mu /(\rho _{\mathrm {f}\omega )}}\), determines the thickness of the stagnant boundary layer at the tube’s wall, where, \(\rho _{\mathrm {f}}\) denotes fluid density. Biot [see e.g., (Kurzeja and Steeb 2012)] introduced the critical frequency

$$\begin{aligned} \omega _{\rm {\mathrm {c,tube}}} = \frac{\mu }{\rho _{\mathrm {f}R^2}} \end{aligned}$$
(6)

that marks the onset of the decrease in effective permeability with frequency. For permeable media, the challenge arises to relate the characteristic length scale of the microscopic conduits (e.g., pore throat radii) with macroscopic (measurable) hydraulic properties. The frequently used conventional Biot frequency

$$\begin{aligned} \omega _{\mathrm {c}} = \frac{\mu \phi ^2}{\rho _{\mathrm {f}}\alpha _\infty k_{\mathrm {D}}} \end{aligned}$$
(7)

results from introducing an effective hydraulic radius, i.e., \(R \rightarrow R_{\mathrm {hyd}} = \sqrt{\alpha _\infty k_{\mathrm {D}}}/\phi\) in (6), where \(\alpha _\infty\) denotes the geometrical tortuosity (e.g., (Ghanbarian et al. 2013)), i.e., the ratio of the average length of flow paths to the length of the geometric flow path (e.g., sample dimension). See, for example, Malinouskaya et al. (2008) or Bernabé and Bruderer (1998) for alternatives and their discussion. Using viscosity \(10^{-3}\) Pa.s and density \(10^3\) kg/m\(^3\) of water, porosity of \(\sim\)1 and 10 %, upper bounds of intrinsic Darcy permeability of \(10^{-17}\) and \(10^{-16}\) m\(^2\) yields \(\alpha _\infty \omega _{\mathrm {c}} > 10^8\) Hz and \(> 10^7\) Hz for Wilkeson sandstone and Westerly granite, respectively. Thus, the frequencies of \(<0.1\) Hz applied in our experiments are several orders of magnitude smaller than the critical frequency for sensible values in tortuosity [e.g., (Ghanbarian et al. 2013)]. Therefore, we exclude contributions of inertia forces to the observed frequency dependence of the hydraulic parameters deduced in this study.

In the following, we discuss our results from two perspectives. First, we evaluate their order-of-magnitude plausibility and their agreement with previously published data. For permeability, results of stationary Darcy-flow experiments serve as benchmark. To assess the validity of the gained values of specific storage capacity, they are converted to drained moduli according to equation (4) and compared to independent constraints for this parameter. The observed differences between hydraulic properties derived from the various methods exceed experimental uncertainty and the effects of the obvious method-related differences, such as the different stress states for axial and radial symmetry. Thus, in a second step, the systematic aspects of the differences between parameters determined by the various methods are interpreted in the light of anisotropy and heterogeneity of the samples’ pore space.

4.1 Validation of Determined Hydraulic Parameters

Permeability values gained by the axial pore-pressure method exhibit close agreement with the results of the Darcy tests, but those of the other two methods exhibit differences that are significant in the light of the deduced uncertainties (“Appendix” A.4), in particular for Westerly granite (Fig. 2). The direct comparison of our results from measurements on the same sample excludes contributions from natural sample-to-sample variability for the two axial methods, possibly affecting the study by David et al. (2018b). A comparison of our permeability results to literature data, however, has to account for the natural variability of the tested rocks, for which reported porosity may serve as indicator (see “appendix” A.1). The connected porosity values of 8–9 % and around 0.5 % quoted for WS and WG samples in previous studies, respectively, are lower than those of our samples. These differences in porosity provide a plausible explanation why the previously reported permeability values, e.g., for WS from axial pore-pressure method by Duda (2011) and a yet unpublished inter-laboratory test, and for WG from hydraulic pulse tests by Brace et al. (1968) and constant-flow rate tests by David et al. (1994), are lower, up to two orders of magnitude for WS and one order of magnitude for WG, and exhibit a more pronounced dependence on pressure (Fig. 2).

We convert results for specific storage capacity to drained bulk moduli (Eq. 4) to allow for a comparison with our ultrasonic measurements and literature data. The ultrasonic measurements at ambient pressure yield magnitudes of drained bulk moduli not reached by all of the pressure-dependent, hydraulic-measurement-based moduli even at elevated pressure difference (Fig. 6a,b). Of the reasons typically invoked for discrepancies between static and dynamic elastic parameters (e.g., Fjaer 2019), strain dependence seems particularly relevant considering the nonlinearity documented by the pressure dependence of the various determined parameters. We observe that previously reported static bulk moduli of WS (Duda 2011: axial pore-pressure method and triaxial deformation tests; Duda 2011, Ahrens et al. 2018: hydrostatic deformation tests on dry samples) and of WG (Brace 1965: linear compressibility measurements) tend to be slightly higher and show a stronger dependence on Terzaghi’s effective pressure than those deduced in this study from the two axial methods (Fig. 6). As for permeability, this relation between our moduli results and those reported in the literature match with the expected anti-correlation with porosity (see Eq. 4). The order-of-magnitude agreement among moduli does not include those derived from the radial pore-flow method that fall below those of the other methods by up to an order of magnitude for both rocks (Fig. 6a, b). We discuss the relation between this discrepancy and specifics of the rocks’ pore spaces in the next section.

4.2 Relation Between Microstructure and Effective Transport Properties: Method-Related Differences as Evidence for Anisotropy and Heterogeneity

The solutions of the diffusion equation underlying the evaluation strategy of all three oscillatory methods assume isotropy and homogeneity, and in turn, test results should be independent of method and period for a isotropic and homogeneous medium. Homogeneity is required to the extent that a one-dimensional description is valid and that the spatial derivatives inherent in the partial differential equation actually have significance on a spatial scale smaller than the sample size. The studied rock samples are about two orders of magnitude larger than their grains with whose size the length of individual pores probably scales. Yet, network properties also depend on characteristic length scales of clusters of pores and proximity to the percolation limit. Several previous studies showed the association of period dependence of hydraulic properties (Fig. 3) with sample heterogeneity (e.g., Bernabé et al. 2004; Fischer et al. 1992; Rabinovich et al. 2015; Song and Renner 2006a).

4.2.1 Anisotropy

The rock-physical characterization of the blocks and the prepared samples suggest a limited sample-to-sample variability for both rocks. In particular, the connected porosity of samples used for axial and radial flow nominally differs by only 0.1 % for WS and 0.06 % for WG, i.e., less than the uncertainty. The ultrasonic measurements on the block of WS at ambient conditions, however, yield an anisotropy of about 10 %, with the fast propagation direction perpendicular to the drilling direction for sampling (see Tables 1 and 2 and also Ahrens et al. 2018). The difference in permeability between the two flow-methods is independent of period, with radial permeability larger by a factor of four to five than axial permeability (Figs. 2 and 3), qualitatively consistent with the velocity anisotropy assuming that microcracks, preferentially oriented perpendicular to the sample axis, are responsible for either observation. Bedding may contribute to the observed anisotropy. Yet, while hydrostatic loading tends to reduce anisotropy related to microfractures, it does not affect bedding-related anisotropy. The convergence of permeability values determined by the radial pore-flow method to those of the axial pore-flow method observed for increasing pressure (Fig. 2a) thus suggests that microfractures play a significant role for the anisotropy of the physical properties of Wilkeson sandstone.

4.2.2 Heterogeneity

The oscillatory testing provides two principal indicators for heterogeneity, period dependence of the effective hydraulic parameters and apparent porosity. The permeability of Westerly granite gained by the radial pore-flow method exceeds the ones of both axial methods for short periods, while the results for the different methods converge for long periods (Fig. 3). Appreciating that the period affects the hydraulic penetration depth, this period-dependence corresponds to an increased permeability around the central bore, comparable to what is addressed as a negative skin in in situ well testing (e.g., Butler 1997; Matthews and Russel 1967). Drilling of the central borehole may have caused spatially confined damage in the form of microcracks that resulted in a systematic form of heterogeneity of the hydraulic properties.

To test the plausibility of the presence of a distinct hydraulic zone around the borehole, hydraulic properties were determined using a borehole diameter 1 mm larger than the true diameter of 1.5 mm, an extreme model for a zone of some tens of a millimeter thickness with enhanced permeability. Using the larger borehole radius indeed reduces the effective permeability by approximately half an order of magnitude bringing the results of the radial pore-flow method closer to the ones of the two axial methods. Thus, the presence of a more permeable zone around the borehole could be responsible for the discrepancy between the methods. However, the effect may probably not be strong enough to explain all of the difference between the two flow methods.

Equation (3) suggests that the apparent porosity of a homogeneous medium should exceed its connected porosity. For Wilkeson sandstone, this expected relation holds independent of applied method, period, or pressure (Fig. 7a). The moderate difference between specific storage capacity values determined by the three methods suggests a well interconnected pore-space in which all pores contribute to fluid flow independent of the specifics of the source for the excited pore-pressure perturbation. In contrast, for Westerly granite, apparent porosity determined by the axial pore-pressure method falls up to one and a half orders of magnitude below the connected porosity, while the two flow methods yield apparent porosity values that exceed the connected porosity (Fig. 8). This apparent contradiction to Eq. (3) suggests method-related differences regarding the part of a heterogeneous network of pores accessed by the pressure perturbations. The axial pore-pressure method investigates the flow path connecting the two sample ends, that in a heterogeneous medium may constitute a rather limited fraction of the entire network, in percolation theory called the backbone (Stauffer and Aharony 1994). The flow methods investigate all of the pore space connected to the respective upstream surface, irrespective of how far the pore channels reach into the sample (Fig. 9).

Fig. 9
figure 9

Relation between ultrasonic p-wave velocities of dry samples and total porosity (a) and geometrical density (b), between p-wave velocities of saturated samples and connected porosity (c) and correlation of connected porosity to total porosity with identity line (black) (d). The three rock types are labeled by different symbols and different samples of the same material with different colors. Data of Ahrens et al. (2018) (A17), Duda (2011) (D11), Interlab-test (Int, see Acknowledgement), and Brace (1965) (B65) are shown for comparison. The legends in a) and c) apply to all parts of Fig.

Full size image

4.2.3 Scale Effects

The measured initial connected porosity of the Westerly granite is small (\(\sim 0.4~\%\)), and the pore network is therefore perceived to be near the percolation limit, the smallest porosity for which an interconnected cluster of pores exists, the backbone. Above and at the percolation limit, the correlation length controls the size dependence of sample properties. Samples with a size above the correlation length exhibit transport properties that are size independent, controlled solely by the conduits forming the backbone. For sample sizes below the correlation length, clusters with finite size (e.g., dead ends) may significantly contribute to total flow. Within the concept of percolation theory, permeability and its variability are thus expected to increase with decreasing sample size, because the amount of pores not belonging to the backbone but contributing to the flow increases (e.g., Hunt, 2005). Song and Renner (2006a) investigated the influence of sample length on hydraulic properties for heterogeneous samples and found the anti-correlation between permeability and sample length predicted by percolation theory.

In our experiments, the borehole surface, i.e., the penetration area, for the radial method is only about 60 % of that for the axial method, the sample-end face. Reducing the penetration area should lead to a reduction in permeability according to percolation theory. Yet, we observe the opposite, the radial pore-flow method gives higher permeability than the axial pore-pressure method (Fig. 2). The small difference in penetration area of the two sample types may actually be outweighed by the more significant difference in “sample length,” 75 mm for the axial and 13.5 mm (thickness of the hollow cylinder) for the radial flow geometry that we advocate as an explanation for higher permeability in radial than in axial direction for Westerly granite besides drilling induced damage. For Wilkeson sandstone, the higher permeability in radial than in axial direction is, in principle, thus also in agreement with the expectation from percolation theory. Yet, considering that the two axial methods suggest “homogeneity” based on the absence of period dependence—plausible for a rock with 10 % porosity—we favor the explanation by “intrinsic anisotropy.”

4.3 Heterogeneity and the Concept of Penetration Depth

The stationary flow in a Darcy test is solely controlled by permeability. In contrast, transient flow in oscillatory hydraulic tests is controlled by hydraulic diffusivity \(D=k/(\eta s)\), the ratio of transport and storage properties. According to the scaling relation for diffusion in homogeneous materials (e.g., Carslaw and Jaeger 1959), the penetration depth \(r_{\mathrm {p}}\) of pressure diffusion can nominally be estimated by \(r_{\mathrm{{p}}}\sim \sqrt{DT}\) for oscillatory tests, since period T corresponds to the characteristic time scale. This scaling relation between the duration of a pressure perturbation and its penetration depth is based on the spatial–temporal characteristics of the occurrence of “significant” pressure perturbations, typically on the order of a tenth or more of the source strength depending on flow geometry (e.g., Weir 1999). However, the scaling argument does not apply to the axial pore-pressure method that yields information regarding the complete sample length irrespective of the chosen oscillation period because Fourier analysis allows for the detection of rather weak signals. In our study, the lowest amplitude ratios \(\delta _{\mathrm {pp}}\) analyzed are actually on the order of \(5\cdot 10^{-4}\), i.e., much smaller than the significant perturbation underlying the scaling analysis. Thus, the concept of penetration depth has limited applicability to the axial pore-pressure method because its results always represent the pathway(s) connecting the sample ends irrespective of oscillation period, that, however, affects for heterogeneous samples to what extent branches of the pore space, including dead ends, contribute to transport and in particular to storage in addition to the most effective pathway.

To observe a significant fraction of the upstream-pressure perturbation at the downstream end of a homogeneous sample, in theory, requires the penetration depth to exceed sample length L. We thus investigated the relation between “nominal” normalized penetration depth, estimated as \(r_{\mathrm{{p}}}/L\sim \sqrt{DT}/L\), and amplitude ratio, i.e., the fraction of the upstream perturbation observed at the downstream end (Fig. 8). For Wilkeson sandstone, the pairs of \(r_{\mathrm{{p}}}/L\) and amplitude ratio fall in the range consistent with the scaling relation, i.e., either \(r_{\mathrm{{p}}}/L1\) and \(\delta _{\mathrm{{pp}}}10^{-1}\) or \(r_{\mathrm{{p}}}/L1\) and \(\delta _{\mathrm{{pp}}}10^{-1}\), where we considered 1/10 of the excitation at the source the threshold for a significant response (e.g., Weir 1999). For Westerly granite, the pairs do not fall inside these fields. We suppose that such a consistency check with the standard scaling relation, valid for a homogeneous material, constitutes a further tool for the evaluation of sample homogeneity. Here, the observations for Wilkeson sandstone and Westerly granite distinctly differ in their respective relations to the prediction of the standard scaling relation classifying the sandstone as hydraulically homogeneous and the granite as heterogeneous, in accord with the results of the other two explored “heterogeneity indicators,” period dependence and apparent porosity.

5 Conclusions

We applied three different oscillatory methods, the axial pore-pressure and pore-flow methods and the radial flow method, and steady-state flow measurements on samples of two rocks, Wilkeson sandstone and Westerly granite, representing prototypic materials for the hydraulically distinct classes of granular materials with intermediate permeability and microfractured materials with low permeability. All but one method (radial flow) were applied on one sample, an approach that avoided contributions from natural sample-to-sample variability to the differences of the results from the various methods.

The order of magnitude of the determined hydraulic parameters, permeability, specific storage capacity, and hydraulic diffusivity, is in accord with previously reported results and independent constraints, but parameters gained by the three oscillatory methods exhibit differences that are significant in the light of our error analysis. Only for Wilkeson sandstone, a part of these differences seems related to hydraulic anisotropy. We suppose that the various methods investigate different parts of a heterogeneous pore space, and their combination thus permits an extended characterization of porous materials. The oscillatory testing provides three principal indicators for heterogeneity, a) period dependence of the hydraulic parameters, b) the relation between apparent porosity, i.e., the specific storage capacity normalized by fluid compressibility, and connected porosity, and c) the consistency of pairs of nominal penetration depth and associated pressure perturbation with the standard scaling relation for diffusion processes. All three indicators classify the tested samples of Wilkeson sandstone as hydraulically fairly homogeneous and those of Westerly granite as hydraulically heterogeneous despite its petrographical homogeneity, presumably as an expression of its pore structure’s proximity to the percolation limit. Heterogeneity affects the permeability results of the various methods notably less than their results for specific storage capacity, documenting the crucial distinction between conducting and storing pore space.