Empirical rock physics relationships on carbonate dry-frame elastic properties

Laboratory tests were conducted to analyze the ultrasonic velocity response to the pressure change in dry carbonate rocks from the Weyburn oilfield, Canada. Twenty-four samples are from seven wells with helium porosities ranging from 1% to 29%. Thin-section images, SEM and mercury intrusion porosimetry were performed to show their inner structures and pore throat size distributions. P- and S-wave velocities (Vp and Vs) measurements were first done under hydrostatic loading and then while unloading, with confining pressures varying between 3 and 35 MPa. The results indicate that Vp and Vs in these samples follow a linear relation independent of the pressure change. The ratio Vp/Vs is more responsive to pressure change irrespective of the pore volume. One-third of the carbonate samples show abnormal Vp/Vs reduction with the increase in the effective pressure. The pressure dependence of velocities (PDV) of Weyburn carbonate rocks varies widely even for samples from the same formation with similar sedimentary history. Samples with loosely packed crystals and/or relatively large dominant pore diameter have higher PDV. The exponential empirical model V = A-CeDPe\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$A-C\text{e}^{DP_{\text{e}}}$$\end{document} was tested; therein, V is the elastic wave velocity, Pe is the effective confining pressure, and A, C and D are the best fitting coefficients determined by curve fitting. The model gives good fits for most of the Weyburn carbonate samples. From a statistical point of view, there is no difference between the Vp- and Vs-derived exponential coefficient D.


Introduction
The change in stress condition influences the elastic frame moduli of naturally occurring rocks and accordingly results in the change in their elastic wave velocities. The pressure dependence of velocities (PDV) of upper crust rocks is of vital importance for a variety of applications, such as pore pressure investigation of hydrocarbon reservoirs (Xu et al. 2006) and time-lapse monitoring of underground fluid production or injection (Pervukhina et al. 2010;Perozzi et al. 2017a, b;Zeng et al. 2019). Nonlinear response of elastic wave velocities to the effective pressure change has been observed for a long time in a wide variety of rocks, such as mudstones and sandstones (e.g., Han et al. 1986;Freund 1992;Khaksar et al. 1999;Sayers 1999Sayers , 2002Ong et al. 2016;Gao et al. 2019;Li et al. 2019a). From these studies, it is clear that with increasing effective pressure, the elastic wave velocities of rocks first quickly increase and then get slower before reaching a plateau (Fig. 1). Phenomenologically, this is explained by the progressive closure with confining pressure of small aspect ratio and compressible crack-like pores (Walsh 1965;Batzle et al. 1980). For hypothetical elliptical cracks, those with the smallest aspect ratios will close under the lowest pressures. Upon each closure, the material becomes incrementally stiffer and more rigid and its elastic wave speeds correspondingly rise (e.g., the numerical modeling work by Wang et al. 2016). The rate of change depends on the distribution of the crack-like porosity.

Edited by Jie Hao and Chun-Yan Tang
Ideally, all these cracks are closed once a certain confining pressure is obtained, whereupon further compression results in the almost linear stress-strain relationships for rocks.
For operating convenience, a number of simplified mathematical expressions have been used to provide a descriptive relationship between the wave speeds and the confining pressure. These usually assume a trial equation whose coefficients are calculated in linear or nonlinear least squares regression of the observed velocities versus pressure. The most popular one is a combination of exponential and linear items based on velocity measurements (Eberhart-Phillips et al. 1989): where V is the elastic wave velocity, P e is the effective confining pressure and A, B, C and D are the best fitting coefficients determined by curve fitting. From Eq. (1), the linear item can be ignored at relative low pressure (P e < 100 MPa) according to the analysis of Shapiro (2003), resulting in: This empirical model gives satisfactory approximations to the measurements of sandstones. Another type of empirical model describing the nonlinearity of rocks comes from the strain-stress measurements (Zimmerman et al. 1986), such as (Melendez-Martinez and Schmitt 2013): where a, b and c are determined by curve fitting. This empirical model provided a superior fit of linear strains with pressure in low-permeability unconventional reservoir shales. For the analysis of PDV, Eq. (3) can be converted to the velocity form with similar variation trend: Figure 1 shows the fittings of Eqs. (1), (2) and (4) based on the velocities of a sandstone sample, with all giving good performances.
Relative to siliciclastics, carbonate rocks have a much more heterogeneous pore network (Xia et al. 2019), which results in larger variances of the acoustic properties for a given porosity (e.g., Assefa et al. 2003;Kenter et al. 2007;Verwer et al. 2008;Weger et al. 2009;Li et al. 2019b). First, in carbonate rocks the velocity-porosity relation shows wide scatter (Anselmetti and Eberli 1993;Njiekak and Schmitt 2019;Li et al. 2020). This is caused among others by the composition of the mineral grains, the texture of the rock frame and the pore characteristics. Second, the PDV of carbonate rocks still lack consistent understanding . Generally, the velocity-effective pressure curve of carbonate rocks is very flat; the observed hysteresis effect following unloading in carbonates underground indicates low PDV (Huffman 2002). Gomez et al. (2007) analyzed the pore structure and the PDV of 50 carbonate samples. The authors found that the PDV of samples dominated by interparticle/intercrystalline pores is stronger than those of samples displaying mostly stiff pores, such as moldic pores, intraparticle pores and vugs. The observations of Guo et al. (2012) on carbonate outcrops from Sichuan Basin, China, indicated that the PDV of dolomites with crack structures is stronger than that of limestones.
As for the Weyburn reservoir formations investigated in our study, Njiekak and Schmitt (2019) examined the influences of pore characteristics and pore fluids on the stress dependence of velocities on a small set of carbonate samples (two limestones and one dolomite). The study was developed from an earlier extensive laboratory testing of 24 Weyburn samples (Schmitt et al. 2012). On the studied subset of samples, the authors also tested the empirical model of Todd and Simmons (1972) for effective stress coefficient for seismic velocity. Moreover, Melendez-Martinez and  investigated the elastic anisotropy change with confining pressure on five samples from the same Weyburn reservoir.
Here, we focus on PDV under dry conditions on the larger set of the Weyburn reservoir samples including 22 carbonate samples and 2 sandstone samples. We will first provide an overview of the geological background and textural features  (2) and (4). The measurements shown here are from a sandstone sample (sample P47, porosity = 0.1437, clay content = 0.13) in Han et al (1986) of the 24 studied Weyburn samples. Results of P-and S-wave velocities acquired at different confining pressures will then be given. Using these measurements, the PDV of the dry carbonate rocks will be analyzed in three different ways in order to seek useful empirical trends that may be applied in estimating rock properties under differing conditions: the V p -V s linear relation, the change in V p /V s ratio and the normalized velocity. Empirical velocity-pressure relations (Eqs. (2) and (4)) will then be tested on the studied carbonate rocks. Finally, the potential controlling factors on PDV for carbonate rocks will be discussed. Overall, this study aims to provide a new set of experimentally measured P-and S-wave speeds in a suite of carbonates, to assess the PDV of these carbonates and the factors controlling this and finally to test some existing empirical velocity-pressure relations that are often employed in the literature.

Geological background
The carbonate samples used in our study are from the Weyburn oilfield in western Canada as shown in Fig. 2a (Fig. 2b). Above the Weyburn carbonates and evaporite formation is the relatively impermeable Triassic Lower Watrous formation, which is a regional sub-Mesozoic unconformity (e.g., Whittaker et al. 2004).

Sample characterization
In our study, we are using 24 geologically representative core samples ( For example, sample 13-O was cut obliquely to the rock beds; i.e., the main axis of the tested cylindrical plug is oblique to the bedding plane. When no letter is attached to a sample index, it means that the axis of the core sample is parallel to the vertical, i.e., normal to the bedding plane. The mineral compositions of the samples were determined through XRD and XRF whole rock analysis. Thinsection images and scanning electron microscopy (SEM) of the samples were performed to show their inner structures. Pore throat size distributions were assessed by mercury porosimetry. A He pycnometer was used to determine the grain densities and the rock porosity. The air permeability as measured at room temperature of most of the samples is lower than 2mD. Their porosity from He pycnometer ranges from 1% to 29% with a mean value of about 11%. Further details of these measurements are available in Njiekak et al (2018). In the following, sample characteristics will be described for each formation unit. A summary of the sample characteristics is given in Table 1. The data and characterizations provided here all originate from the earlier reporting by Schmitt et al. (2012).

Watrous red bed unit
The Watrous Red Bed unit is represented by samples 39 and 43. They have a porosity of 7% to 11% and a grain density varying from 2.74 g/cm 3 to 2.76 g/cm 3 , respectively (see Table 4 for the supporting data in Appendix 1). Porosity is mainly intercrystalline and the grain size is up to 100 µm as depicted in Fig. 4a. Quartz and feldspar are the main minerals. Anhydrite, dolomite, calcite and clay-size crystals usually fill the space between the quartz and feldspar grains (Fig. 5a). Low content in anhydritic and carbonate cement is usually indicative of high pore volume and permeability (up to 26 mD) as in sample 43. The dominant pore throat diameter is around 4 µm (Fig. 6a).

Top midale evaporite unit
Three samples from the Top Midale Evaporite unit were investigated. Overall, these samples have a tightly packed and dolomite-rich matrix, which is crosscut by discontinuous thin anhydritic layers (Figs. 4b and 5b). Their grain density ranges between 2.73 g/cm 3 and 2.89 g/cm 3 (Table 4), and their porosity, which is intercrystalline, ranges from 1 to 3%. The main pore throat diameter ranges between 0.015 µm and 0.7 µm (Fig. 6b, c). Permeability is lower than 2 mD in the three samples.

Marly unit
There are ten samples from the Marly units, including eight dolostone and two limestone samples. The grain density varies between 2.64 g/cm 3 and 2.83 g/cm 3 , and the permeability is also smaller than 2 mD (Table 4). Besides some locally distributed vugs or channels (e.g., Fig. 4c), there are mainly intercrystalline pores (Fig. 4d). The main pore throat diameter in those samples varies from 0.02 µm to 1.5 µm (Fig. 6d, e). Based on their porosity, they can be subdivided into two groups: 1 The tight Marly samples ( Fig. 5d) with porosity up to 8% and composed of samples 9, 16 and 31. Their tight matrixes occasionally have sparse channels and vugs, which are often filled with anhydrite ( Fig. 4d). 2 The loose Marly samples with porosity ranging from 17% to 27%, such as samples 20, 26 and 46. Their frames commonly display alternating loose domains (Fig. 5c).

Vuggy unit
The Vuggy unit is represented by eight limestone samples (Table 4). They are mainly wackstones and packstones, except sample 42 which is a grainstone (Njiekak and Schmitt 2019;their Fig. 1e, f and g). Porosity in those samples ranges from 4% to 16% and may be vuggy ( Fig. 5f, g), moldic, intergranular and intragranular. Their permeability is generally good, ranging from 16 to 50 mD (Table 4). Recrystallization usually leads to the reduction of porosity by filling the pore space (Fig. 5e). The dominant pore throat diameters in the Vuggy samples vary from 0.01 µm to 8 µm ( Fig. 6f, g), and some can reach 35 µm (Njiekak and Schmitt 2019; see sample OL in their Fig. 2, sample OL = sample 42 in this paper).

Frobisher Evaporite unit
Sample 38 is from a lime-rich interval within the Frobisher Evaporite unit. It has a tight matrix mostly made of calcite and dolomite crystals (Figs. 4f and 5h). Other minerals include quartz and feldspar.
The grain density of sample 38 is 2.73 g/cm 3 . Its porosity is 9% and is intercrystalline, with the dominant pore throat size around 0.1 µm (Fig. 6h). Its permeability is lower than 2 mD (Table 4).

Velocity measurements
Experimental setup (Fig. 7) and procedures for the P-and S-wave velocity (V p and V s ) determinations were given in Njiekak et al (2013). Briefly, the tested samples were all cylindrical plugs of 38.1 or 25.4 mm in diameter and of about 20.1 to 58.0 mm in length (Table 4). Once placed inside the pressure vessel, the sample's pore space was subject to vacuum for about 8 h (pore fluid pressure, P p = 0 MPa). Then, the hydraulic oil pressure within the vessel (confining pressure, P c ) increased from P c = 3 MPa to P c = 35 MPa, stepped by 2 or 3 MPa. During this hydrostatic loading process, the V p and V s of the dry sample at each pressure point were measured with the standard pulse transmission at a central frequency equal to 1 MHz. The main source of errors in the velocity measurements is associated with the picking of travel times from the acquired raw waveforms (e.g., Njiekak et al 2013, their Sect. 4.4). Overall, when considering all the sources of uncertainties in our measurements, the percentage errors are up to 0.7% and 0.4% for P-and S-wave velocities, respectively. The velocity measurements were repeated at the same pressure points while unloading. For these measurements done under vacuum dry and hydrostatic conditions, the Terzaghi effective pressure is the same as the confining pressure.
The measurements show that the velocities nonlinearly decrease with the increase in porosity, but it is hard to use one simple relation to fit the wide scatterings ( Fig. 8). Such scatter in the velocity-porosity cross-plot for carbonate rocks is consistent with observations made in previous studies (e.g., Anselmetti and Eberli 1993;Assefa et al. 2003;Rogen et al. 2005;Verwer et al. 2008;Weger et al. 2009). Also, a broadly accepted rationalization is that the complex textures and related pore structures of carbonate rocks are responsible for this wide dispersion (Weger et al. 2009;. For comparative purposes, we are using a differential effective medium (DEM) model (Berryman 1992; see detail in Appendix 2) to calculate the velocities of hypothetical rock frames composed entirely of a single carbonate mineral but with different pore aspect ratios (ARs). For our calculations, the rock matrix is assumed to be either pure calcite or dolomite. The values of the input parameters are shown in Table 2. From the DEM modeling, the velocity-porosity relation depends on the pore aspect ratio. The latter defines the stiffness of the rock skeleton. The velocity-porosity relation is nonlinear, especially at low ARs (AR < 0.1). For the Fig. 2 The location of Weyburn oilfield and its stratigraphic section. a The location of Weyburn oilfield, from Underground Gas Storage: Worldwide Experiences and Future Development in the UK and Europe, Volume 313, Riding and Rochelle (2009), Subsurface characterization and geological monitoring of the CO 2 injection operation at Weyburn, Saskatchewan, Canada, pp. 229, with permission granted according to the Geological Society of London fair use policy; b the stratigraphic section of the Weyburn oilfield, modified from Proceedings of the 6th International Conference on Greenhouse Gas Control Technologies, Volume I, Whittaker and Rostron (2003), Geologic storage of CO 2 in a carbonate reservoir within the Williston Basin, Canada: an update, pp. 387, Copyright (2003), with permission from Elsevier ◂ same porosity, samples with lower effective AR have lower velocities. Irrespective of the dominant mineral composition, the measured velocities overall fall within the velocity envelopes corresponding to AR = 0.05 and AR = 0.2 (Fig. 8). This means the effective ARs representing the stiffness of the dry skeleton of these samples vary between 0.05 and 0.2, pointing out significant differences in their pressure sensitivity. Sample 45 is dolomite-dominated but, unlike samples 13 and 29 from the same unit and also dolomite-rich, it falls out of the velocity envelopes representing AR = 0.05 and AR = 0.2. A possible explanation for this lies in the presence of clay-like crystals in sample 45 (Fig. 9). Those clay-like crystals, which likely originated from feldspars breakdown, led to lower velocities.
The measured V p during the hydrostatic loading process on the 24 samples is shown in Fig. 10 (see detailed results in Table 5). Their formation units and the samples' dominant minerals are also reported in the same figure. The measurements during the unloading process are not shown here. In fact, the velocities measured during the unloading process were usually higher than those recorded during the loading process. This is typical for such experiments when there is not enough time left between measurements during the unloading cycle. The observed discrepancy results from the Fig. 3 Core pictures of the Weyburn field from seven wells, a from Well 10-11-6-14W2, b from Well 2-16-6-13W2, c from Well 6-1-6-13W2, d from Well 2-10-6-14W2, e from Well 8-36-6-14W2, f from Well 4-23-7-14W2 and g from Well 4-16-7-13W2 fact that the pace at which pores close during the loading is faster than the pace at which they open during the unloading cycle. Indeed, with enough time left between measurements during the unloading cycle, all the measured velocities will match the "loading" velocities unless the tested sample at some point experiences mechanical damage. We can confirm here that no mechanical alteration of the samples occurred. In fact, the same loading process was applied to a subset of the tested samples at a later stage and similar velocities were obtained.
As a whole, the measured velocities increase with the effective confining pressure, but the PDV are quite different from each other. Wide variations in PDV also exist among samples of the same formation, i.e., samples originating from a similar depositional environment (Fig. 10). All this makes it difficult to analyze the effects of the mineral composition or the rock texture on the PDV by solely relying on the velocity-effective pressure plot.

Pressure dependence of velocities
For the pressure dependence of velocities (PDV) analysis, we choose as maximum and minimum effective pressures, 35 and 6 MPa, respectively. This is because, for most of the studied carbonate rocks, high-quality velocity data were obtained from P c = 6 MPa upward and the P-and S-wave velocities barely changed when the effective pressure was higher than 35 MPa. For consistency, in the following sections only the "vertical" cylindrical plugs are considered in our analyses, unless otherwise stated. For example, sample 13 and sample 13-V are all vertical samples collected from the same core material. However, the former was considered together with vertical samples from other cores in the statistical treatment of our data. Moreover, velocities recorded on samples 13-V, 13-O and 13-H are quite different from those measured on sample 13 (Fig. 10). The discrepancy is likely related to textural differences between samples 13-V, 13-O and 13-H (all cut from a single big rock fragment) and sample 13. In fact, the core from which the four samples were obtained displays strong local textural variations.

V p -V s linear relation
Well-known V p -V s relations for carbonate rocks were developed in Castagna et al. (1993), wherein the empirical relation for dolomites is linear while that for limestones requires a second polynomial. The measurements of the studied 22 carbonate samples from the Weyburn oilfield are slightly different from those reported in Castagna et al. (1993). Our newly obtained carbonates dataset, regardless of the rock mineralogy, follows a linear V p -V s relation (Fig. 11): Also, linear fittings of the V p -V s relation corresponding to different pressure conditions (e.g., at P e = 6 and 35 MPa;  (2003), for a given rock, the fitting coefficient D in Eqs. (1) and (2) should be in the first approximation the same for both V p and V s , that is: which can be further rewritten as follows: It should be noted that fitting parameters in Eq. (7) are independent of the pressure, which means that this linear relation between P-and S-wave velocities is suitable for different pressure conditions. As a result, the same empirical model for S-wave velocity prediction can be used for certain formation even though the formation depths are quite different.

The change in V p -V s ratio
The V p /V s is a commonly used attribute for fluid and lithology identification by geophysicists. V p /V s variation with effective pressure is also suggested as an indicator of overpressure and shallow water flow (SWF) zones (Dvorkin et al.    . For sandstones, V p /V s decreases with the increase in effective pressure under liquid-saturated conditions, and the opposite is true under dry or gas-saturated conditions (Dvorkin et al. 1999;Dutta 2002). Testing the suitability of this rule in carbonate formations needs to be done with caution because of the complexity of their textures.
In order to quantitatively describe the V p /V s variation with effective pressure, the change in V p -V s ratio (Δ(V p /V s )) in this paper is defined as: For dry carbonate samples in our measurements, Δ(V p / V s ) < 0 means the V p /V s normally increases, while Δ(V p / V s ) > 0 means the V p /V s abnormally decreases. As reported at the beginning of Sect. 4.2, the low and high effective pressures considered for the PDV analysis of these carbonate samples are 6 and 35 MPa, respectively. The results are shown in Fig. 12. Unlike the V p -V s linear relation discussed in the above section, the V p /V s is more responsive to pressure change. According to our measurements, when the main axis of the sample is normal to the bedding plane, the V p /V s of 18 carbonate samples normally increases with the pressure increasing from 6 to 35 MPa, while the V p /V s of the other three samples slightly decreases. However, when the main axis of the sample is oblique or parallel to the bedding plane, the V p /V s of only three carbonate samples normally increases with the pressure, while the V p /V s of the other seven subcores abnormally decreases. Thus, sample orientation likely played a role in the abnormal behavior of V p /V s recorded on one-third of the studied Weyburn carbonate samples. This is quite a high number of samples, and the observed abnormal trends cannot be explained solely by measurement errors. To put our results into perspective, the measurements on dry carbonate samples published by Fournier et al. (2014) and Hairabian et al. (2014) are also included in Fig. 12. The additional data include 251 carbonate samples that are variably porous and which display a large variety of pore types and textures. It should be emphasized that both Fournier et al. (2014) and Hairabian et al. (2014) described that their samples were "vertically oriented, cylindrical samples." The uncertainties in their velocity measurements were around 0.5% to 1%.
Overall, the Δ(V p /V s ) scatters widely and this is independent of the pore volume (Fig. 12). Among the dry verticaloriented carbonate samples, 229 samples (about 83.9% percent of samples analyzed here) show a Δ(V p /V s ) smaller than 5% (i.e., |Δ(V p /V s )|< 0.05), whereas 101 samples (about 37.0% percent of samples analyzed here) show an abnormal decrease in V p /V s with the increase in the effective pressure. These different trends clearly imply that using V p /V s as a tool to predict abnormal (pore) pressure in carbonate gas reservoirs requires careful calibration.

Normalized velocity
Since the amount of velocity change is dependent on the pressure range, the quantitative comparison of the PDV of different carbonate samples should be constrained within the same pressure window. For our measurements, the chosen lower and upper effective pressure bounds for the PDV analysis are 6 and 35 MPa, respectively. Indeed, within this pressure range, the P-and S-wave velocities of all the studied rock specimens were reliably collected. To quantitatively analyze the PDV, the normalized velocities (V pn and V sn ) are defined in the following: and As shown in Fig. 1, rock velocity usually increases with pressure; thus, the values of the normalized velocities defined in Eqs. (9) and (10) are higher than 1.0. The  Fig. 7 Schematics of the pressure vessel and the data acquisition system for the P-and S-wave velocities measurements. Two different pump systems controlled the confining and the pore pressures. For the measurements in this paper, the pore pressure is zero normalized velocities quantitatively represent the PDV of these dry carbonate samples because they are calculated at the same pressure window. As shown in Fig. 13, the normalized velocities do not show clear dependence to the pore volume. V pn is slightly higher than V sn . Compared to the linear V p -V s relation in Sect. 4.2.1 (Fig. 11), the linear fit between V pn and V sn shows stronger scattering, but it is still a significant linear relation. Indeed, the critical linear correlation coefficient (R 0.01, 20 ) for a fitting with 20 degrees of freedom (i.e., 22 data points as in Fig. 13) and with a 0.99 significance level is 0.537. The latter is much smaller than the correlation coefficient R obtained here, which is 0.825. The V pn and V sn have similar distributions: (i) V pn ranges between 1.009 and 1.142 with the mean value of 1.043, while V sn ranges between 1.004 and 1.095 with the mean value of 1.032; (ii) the V pn and V sn of more than 72% samples are smaller than 1.05, which means the increase in velocity is less than 5%. According to the measurements of Freund (1992) on 88 dry clastic rocks, the mean values of V pn and V sn calculated by the velocities at the effective pressure of 8 MPa and 40 MPa are 1.196 and 1.127, respectively. Thus, the normalized velocity for carbonate rocks is obviously smaller than that of sandstones at similar pressure range. This means that the PDV of carbonate rocks is much weaker than that of sandstones, and velocities of carbonate rocks are less sensitive to abnormal pore pressure. Pore pressure prediction of carbonate reservoirs based on velocity-effective stress relation may not be as reliable as that in sandstone formations.

Empirical velocity-pressure models
To test the suitability of the two types of empirical models reported in Sect. 1 for carbonate rocks, we use the elastic wave velocities measured on the 22 Weyburn carbonates during the loading process. For the exponential model (Eqs. (1) and (2)), the measured velocities at different pressures are fitted by Eq. (2), i.e., without the linear item. This is because our measurements were carried out under relatively low effective pressures. For the two empirical models (Eqs. 2 and 4), we use the same fitting method and the same database to do the fittings. The average value of the squared correlation coefficients (R 2 ) for the fittings of the carbonate samples is higher than 0.97 (Table 3). Predicted velocities on all the carbonate samples based on the two empirical models are compared to our measurements in Fig. 14.
As shown in Table 3 and Fig. 14, the fitting results indicate that the two models are both effective in describing the PDV for carbonate rocks at relative low pressure (P e < 60 MPa). Only few samples did not follow the variation trends obtained through these two empirical models. The low quality of the waveforms recorded on those few samples likely explains such discrepancy. Also, velocities measured at different effective pressures on this subset of samples sometimes show very irregular fluctuations. As shown in Fig. 14, the exponential model (Eq. (2)) slightly performs better than the power function model (Eq. (4)). Also, the R 2 of fittings from the exponential model has slightly higher means and lower standard deviations (Table 3). However, for many applications, the differences between the fitting results of these two models are too small and can be ignored safely.
The coefficient of the exponent (D) in Eq. (2) represents the level of PDV. The fitted values of D derived from V p and V s are statistically the same (Fig. 15), regardless of the

Fig. 8
The V p (circle points) and V s (triangle points) at P e = 35 MPa of the dry Weyburn samples whose main axis is normal to the bedding plane. The blue color represents samples whose dominant mineral is dolomite, while red color is for calcite-rich samples. The solid and dashed curves are obtained using differential effective medium concepts (Berryman 1992) with pore aspect ratios of 0.2 and 0.05, respectively. For these calculations, the rock matrix is assumed to be pure calcite (red curves) and dolomite (blue curves), respectively. The arrows indicate the velocities of three samples from the Top Midale Evaporite (samples 13, 29, 45). Sample 45 that contains noticeable amount of quartz and anhydrite according to the XRD data dominant mineral composition or the orientation of the sample main axis. This is in line with sandstones data listed in Table 2 of Khaksar et al. (1999, see also Fig. 15). Thus, both our observations and those of Khaksar et al. (1999) support the suggestion of Shapiro (2003) that all coefficients D are identical in the first approximation. The D values for our set of carbonate samples are obviously more scattered than those of the sandstones reported in Khaksar et al. (1999), which confirms the complex behavior of the PDV of carbonate rocks.

Mineral composition and PDV
While seismic velocity of natural rocks, especially clastic rocks, strongly depends on their pore volumes, it is quite different for carbonate rocks. Neither the V p -V s ratio change with pressure ( Fig. 12), nor the normalized velocity (Fig. 13) shows apparent correlations with the pore volume. This is because pore volume reduction induced by pressure increase is quite small for most carbonate rocks. The stiffness of the   Fig. 10 The measured P-wave velocity for the 24 dry Weyburn samples during the loading process. The samples are grouped according to their formation units in the graph horizontal axis. For each column, it shows the velocity increases for each sample during the loading process (confining pressure goes up from 3 to 35 MPa). The dominant mineral in each sample is indicated by the color bar above the sample index rock skeleton, more specifically some factors controlling it such as the mineral composition and the pore structure, would play more influential roles on PDV than the pore volume.
In order to assess the influence of the mineralogy on the PDV of rocks, the V pn at P e = 35 MPa is plotted against the pore volume in Fig. 16; therein, rocks with different dominant minerals are distinguished. For the 24 studied samples, 11 samples are dolomite-rich, 11 samples are calcite-rich and 2 samples are quartz-rich (sandstone). As discussed above, the normalized velocity change seems independent of the pore volume (Fig. 13). Thus, it is not surprising that the V pn -φ cross-plot is quite scattered (Fig. 16). The linear fitting coefficient is 0.57 and 0.43 for dolomite-rich and calcite-rich samples, respectively, which are smaller than the critical linear correlation coefficient (R 0.01, 9 = 0.735) for a fitting with nine degrees of freedom (i.e., 11 data points) and 0.99 significance level. There is no noticeable difference in the V pn -φ between the dolomite-and calcite-rich samples.
The V pn -φ cross-plot of the sandstones from our dataset cannot be interpreted with enough certainty due to less data points. Therefore, for more reliability, the measurements of dry clastic sandstones by Freund (1992) are added here (Fig. 16). On the one hand, as a whole the V pn of the sandstones is much higher than that of carbonate rocks, which shows that the sandstones have stronger PDV than the carbonate rocks. On the other hand, the V pn -φ cross-plot of sandstones displays a wide scattering. The carbonate rocks data points show some scattering as well but to a smaller degree. In light of experimental studies conducted on sandstones, the presence of clay minerals can be a factor causing the scattering of their PDV (Eberhart-Phillips et al. 1989;Freund 1992;Sams and Andrea 2001). As for carbonate rocks, a similar observation has not been made so far. In our study, overall there is no evidence to suggest that the observed V pn -φ scattering for carbonate samples is due to non-carbonate minerals (quartz, feldspar, anhydrite) typically present in their solid frames.

Pore structure and PDV
The V pn data were grouped by formation unit to assess how pore geometries may influence the PDV under similar sedimentary conditions (Fig. 17). The black dashed lines in the middle subplots of Fig. 17 indicate the average curves for different formations. The left and right subplots are, respectively, the pore throat distributions of samples showing maximum and minimum V pn in the same formation. It is easy to find that the dominant pore throat diameters for maximum and minimum V pn have significant differences. For the same formation, the variation of PDV is probably caused by the matrix texture and the pore structure.
For each unit, the sample with minimum V pn has a small dominant pore throat diameter. The latter is usually not larger than 0.1 µm (the right subplots of Fig. 17a, b and  d), except for the loose Marly unit. A small dominant pore throat diameter is often typical for tightly packed matrixes.
Moreover, the average value of the normalized velocity change in the loose Marly unit is 1.09 at P e = 35 MPa (Fig. 17c), while it is smaller than 1.06 for other formations (Fig. 17a, b and d). The PDV of the loose Marly unit is significantly larger than those of the other formations due to its loosely packed matrix. Interestingly too, for the loose Marly unit, the samples with the maximum and minimum V pn , samples 20 and 25, respectively, almost have the same main pore throat diameter (1 to 2 µm; Fig. 17c). Yet, the difference in their pore volume, 27% and 20% for samples 20 and 25, respectively, is possibly the main reason behind their different PDV. Thus, within a carbonate formation made of loosely packed material, pore volume can well be considered as a PDV main controlling factor.

Fig. 11
The linear relations of P-and S-wave velocities recorded on the dry carbonate specimens at P e = 6 MPa (red color) and P e = 35 MPa (blue color), respectively. The circle points represent the samples whose dominant mineral is dolomite, while the triangle points represent calcite-rich samples. The two fitting lines nearly overlap, which indicates that the linear relation of V p and V s is independent of the pressure condition V p V s ∆ Fig. 12 Changes in V p /V s on dry carbonate rocks following the increase in the effective stress. The histogram subplot represents the number of samples with different levels of V p /V s changes. The red points are based on our measurements. Measurements carried out under similar conditions by Fournier et al. (2014) and Hairabian et al. (2014) are also included. Here, pore volume change caused by increasing effective stress is neglected

Conclusions
Elastic compressional and shear velocities were measured on 24 dry carbonate specimens following loading and unloading processes under hydrostatic condition. The rock samples were selected from seven wells of the Weyburn oilfield, and they display a wide range of pore volume, pore type and depositional texture. Making use of the data collected on the carbonate rocks forming the reservoir, the pressure dependence of velocities (PDV) under dry conditions were investigated. The results indicate that: (1) P-and S-wave velocities in the dry Weyburn carbonate specimens follow a linear relation that is independent of the pressure change. The ratio V p /V s is more responsive to pressure change irrespective of the pore volume of the carbonate rocks. Log differential intrusion, mL/g S12 Pore throat diameter, µm P e , MPa Pore throat diameter, µm Pore throat diameter, µm P e , MPa Pore throat diameter, µm Table 4 The location, density, porosity and permeability of all samples* *Porosity is measured by the helium porosimeter, Perm = the air permeability, Rom = grain density. It should be noted that some samples were cut to perform measurements following different angles (see anisotropy study by Melendez-Martinez and Schmitt 2013). Their sizes are smaller (diameter = 25.4 mm, length = 20.1-37.0 mm). For these samples, their indexes include a letter, such as V, O and H meaning vertical, oblique and horizontal, respectively. For example, sample 13-O was cut obliquely to the rock bedding plane; i.e., the main axis of the tested cylindrical plug is oblique to the bedding plane. When no letter is attached to a sample index, it means that it was cut vertically, i.e., normal to the bedding plane Appendix 2: The differential effective medium (DEM) model The differential effective medium (DEM) model is widely used to predict the effective bulk (K * ) and shear (μ*) moduli of two-phase composites. It is a coupled system of two ordinary differential equations (Berryman 1992): with the initial conditions (11) K * (0) = K m * (0) = m where K and μ are the bulk and shear modulus, respectively, m and f stand for the solid matrix and pore fluid, respectively, and y is the pore volume fraction. The parameters P *f and Q *f for elliptical pores of arbitrary aspect ratio (AR) are given by and T 1 and T 2 can be calculated as follows (Berryman 1980): and wherein with A, B and R given by: where v * is the effective Poisson's ratio of the composite mineral. The functions β and f are given by: