Comparison of Fluid Pressure Wave between Biot Theory and Storativity Equation

Compressibilities of pore fluid and rock skeleton affect pressure profile and flow velocity of fluid in aquifers. Storativity equation is often used to characterize such effects. The equation suffers from a disadvantage that at infinite large frequency, the predicted velocity of fluid pressure wave is infinitely large, which is unrealistic because any physical processes need certain amounts of time. In this paper, Biot theory is employed to investigate the problem. It is shown that the key equations of Biot theory can be simplified to storativity equation, based on low-frequency assumption. Using Berea sandstone as an example, we compare phase velocity and the quality factor between Biot theory and storativity equation. The results reveal that Biot theory is more accurate in yielding a bounded wave velocity. At frequency lower than 100 kHz, Biot theory yields a wave velocity 8 percent higher than storativity equation does. Apparent permeability measured by fluid pressure wave (such as Oscillatory Hydraulic Tomography) may be 14 percent higher than real permeability measured by steady flow experiments. If skeleton is rigid, Biot theory at very high frequencies or with very high permeabilities will yield the same velocity as sound wave in pure water. The findings help us for better understanding of the physical processes of pore fluid and the limitations of storativity equation.


Introduction
Fluid in the subsurface is very important for hydrogeologists and petroleum engineers, as significant portions of fresh water and hydrocarbon are stored in rock pores. Porous rocks consist of two phases, i.e., solid and fluid. Voided solid is referred to as skeleton, also called dry rock. Mechanically, fluid differs from skeleton in two aspects: (1) it is often more compressible than skeleton; (2) fluid has zero shear modulus and thus can seep between pores and reach far away but skeleton cannot. For static fluid, zero shear modulus also causes pressure in different directions to equal.
In hydrogeology, fluid is often treated to be incompressible for simplicity, which yields the Laplacian equation as the integration of Darcy law and fluid mass conservation [1]. As a reminder, Darcy law may not be accurate if seepage is very slow [2]. Another fact is that Darcy law ignored fluid acceleration which consists of the unsteady and inertial terms (associated with the temporal and spatial derivatives of fluid Lagrange velocity, respectively). The latter term may be important near wellbores such that the effect of non-Darcy flow can be considerable [3].
With the consideration of compressibilities of fluid and skeleton, a storativity equation [4,5] is more advanced than the aforementioned Laplacian equation. Mathematically, the storativity equation is a diffusion equation. At low frequencies, the equation yields a finite velocity of fluid pressure wave. Recall that any physical processes need certain amounts of time. As such, the equation is physically reasonable at low frequencies. However, at infinite high frequency, the equation yields infinite large velocity of fluid pressure wave, which is unrealistic. In this regard, the storativity equation is not accurate enough.
Pure fluid allows a single type of compressional (P) wave (fast P wave or sound wave). In contrast, there are two types of P wave in fluid-saturated rocks (fast P wave and slow P wave). Fast P wave travels with high speed and small attenuation which have been observed in seismology [6,7], while slow P wave is diffusive with low velocity and large attenuation. Plona [8] observed slow P wave in water-saturated sintered glass beads. Slow P wave is fluid pressure driving groundwater flow, and the associated skeleton deformation is often negligible. Hydraulic fracturing [9,10] is an extreme case that slow P wave has high amplitude to enlarge fissures in tight shale.
Biot [11,12] proposed an advanced theory to characterize waves in fluid saturated rocks, based on single porosity, elastic skeleton, and viscous fluid. In the theory, the dilatational wave of the first kind was fast P wave, while the dilatational wave of the second kind was slow P wave. In this regard, Biot theory unified sound wave in acoustics and groundwater flow in hydrogeology.
For fast P wave, skeleton is active while fluid is passive. Solid and fluid tend to have synchronous motions [26], and friction between them is relatively small. In contrast, for slow P wave, fluid is active while skeleton is passive. Fluid and solid have out-of-phase motions, and friction between them is very large. Li et al. [26] showed that for fast P wave with frequency at 1-100 Hz, skeleton pressure is one order of magnitude higher than fluid pressure and the two pressures have a very slight phase difference from each other.
Historically, transient pumping tests (water is extracted or injected very quickly via pumping wells) have been used to estimate aquifer properties [27]. In the recent years, tomographic analysis based on the storativity equation was used to invert for hydraulic conductivity/permeability. Cardiff et al. [28] developed a method named Oscillatory Hydraulic Tomography (OHT), in which periodic pumping with different frequencies served as the stimulation and fluid pressure responses were recorded at different locations. What essentially in OHT is stationary wave due to boundary conditions at two ends. Zhou et al. [29] conducted sandbox experiments (fiber-optic pressure transducers were used to record fluidpressure change due to oscillatory pumping) and showed that OHT can be used to estimate hydraulic conductivity. Overall, OHT is based on the assumption that the storativity equation is accurate.
Similar to OHT, Becker et al. [30] measured fluidpressure change and fracture displacement using transducers and fiber optic distributed acoustic sensors. In their field experiment, receivers were in a borehole while water-head oscillations were stimulated in a companion borehole 30 m away.
In the production life of hydraulically fractured reservoirs, gas is very compressible, and therefore, storativity is an important parameter for flow velocity and pressure profile of gas. Al-Rbeawi [31] analyzed pressure behaviors and flow regimes of gas in hydraulically fractured gas reservoirs, in which the gas pressure profile with time was analytically modeled and the effect on gas pressure by fracture storativity was well considered.
Fluid pressure drives flow in porous media, and there are two approaches to solving the problems of groundwater flow. One is directly solving storativity equation subjected to specified boundary condition and initial condition. The other approach is solving the eigenvalue problem, i.e., using plane/monochromatic waves, and the final solution is the superpo-sition of the individual plane waves. The first approach is advantageous in getting flow velocity and pressure profile, but neither the velocity nor attenuation of fluid pressure propagation is well known. The second approach has the assets that both velocity and attenuation of fluid pressure wave are got as functions of frequency but has the liabilities that flow velocity and pressure profile are dependent on the specific boundary condition and initial condition. This paper uses the second approach (the eigenvalue approach) because the solution of an eigenvalue problem has a general meaning, i.e., independent of 1D, 2D, or 3D. Figure 1 is a schematic propagation of fluid pressure in pores (denoted as P p ) which is caused by the change of well pressure (P w ) by periodic pumping.
Hydrogeologists and petroleum engineers often use a storativity equation with the first approach to simulate flow in the subsurface. The novelty of this paper is using Biot theory to study the propagation of fluid pressure in pores. As Biot [11,12] theory involves a set of partial differential equations (that are difficult to directly solve), the second approach is simpler and significantly facilities comparing fluid pressure wave between Biot theory and storativity equation. By the comparison, the error of the storativity equation, e.g., in measuring permeability [28,29], can be obtained. In addition, we are interested with a scientific question on the condition under which fluid pressure wave will have the same velocity as sound wave in pure water.
At first, the key equations of Biot theory are presented in Section 2. Then in Section 3, the equations are simplified to the storativity equation. In Section 4. Berea sandstone is used as an example to compare Biot theory with storativity equation. In Discussion, we estimate the error of the OHT method and analyze an extreme case. The last section is summary and conclusions.

Biot Theory
2.1. Constitutional Relations. Equations (2.11-2.12) in [11] stated the constitutional relations of solid and fluid, respectively. Solid stress and fluid pressure are linear combinations of volumetric strains of solid and fluid. Reversely, volumetric strains of solid and fluid are linear combinations of solid stress and fluid pressure. In a 1D scenario, the constitutional relation of solid is that skeleton strain on the right hand side of equation (1) is determined by the mean normal stress (the confining pressure), P C , and fluid pressure (P p ) on the left hand side; please refer to equation (7.41) in [32].
where α B = ðβ eff − β s Þ/β eff is the Biot-Willis coefficient [32], β eff and β s are the compressibility coefficients of skeleton and solid material, respectively, and u is the skeleton displacement in the x direction. The left hand side of equation (1) represents effective stress. Note that α B = 1 is often used in hydrogeology or soil mechanics [5].

Geofluids
Using the normal stress in the x direction (σ 11 ), the above relation can expressed as follows [32]: where K eff = 1/β eff is bulk modulus of skeleton and G is shear modulus of skeleton.
In the lowest part of equation (2.11) in [11], volumetric strain of fluid is linearly related to the confining pressure (P C ) and fluid pressure (P p ), which is identical to fluid mass conservation as follows [25], except that the latter is the derivative of the former with respect to time: where q is the Darcy flux rate in the x direction, η is the rock porosity, and β f is the compressibility coefficient of fluid. Equation (3) is more accurate than Li et al.'s [26] because the compressibility coefficients have been revised according to [33,34].

Momentum Equations of Fluid and
System. For slow P wave, fluid is active and skeleton is passive, such that squirt is minor. At low frequencies, the term of fluid acceleration associated with fluid density (ρ f ) and the term of relative acceleration associated with the coupling density (ρ 12 ) are both small in equation (6.7) in [11]. For these reasons, the fluid momentum equation, i.e., the lower part of equation (6.7) in [11] is simplified to where μ and k D are the fluid viscosity and Darcy permeability, respectively, and v = ∂u/∂t is the Lagrange velocity of skeleton but approximated by the Euler derivative. The solid momentum equation is the upper part of equation (6.7) in [11]. Adding the upper part with the lower part of equation (6.7) in [11] cancels the internal force between fluid and solid, yielding the momentum equation of rock system as follows: where ρ s and ρ f are the skeleton density and fluid density, respectively. Evidently, equation (5) is Newton's second law for the rock system. Equations (2)-(5) involve four quantities, namely, σ 11 , P p , u, and q. Solving the four unknowns in terms of plane waves yields a wavenumber equation (the eigenvalue equation) which has two branches of solution. The first branch is fast P wave while the second branch is slow P wave. In each branch, phase velocity (v p ) and attenuation (1/Q p ) are functions of frequency.
The quality factor (denoted as Q p ) is a dimensionless quantity describing the ratio of energy loss to the total (mechanical) energy in one cycle/period of wave propagation [35]. For a fixed frequency (a monochromatic wave), the factor is low when the attenuation of the wave by friction is large and vice versa.
In order, equations (6)-(9) approximately represent the constitutional relation of rock skeleton, the fluid mass conservation, the fluid momentum equation, and the momentum equation of rock system. Eventually, all the four quantities will vanish and be trivial. In short, the condition of zero rock stress results in a trivial solution, because  Figure 1: A schematic propagation of fluid pressure in pores (P p ). A pumping well changes the well pressure (P w ) as sinusoidal function of time (t), which propagates (with attenuation) into the aquifer. Note that this paper solves the eigenvalue problem and the resulting functions of wave velocity and attenuation with argument frequency are independent of 1D, 2D, or 3D. 3 Geofluids equations (6)-(9) are overdetermined for three unknown quantities (u * , q * , and P p ). In addition, once the condition is implemented, shear modulus (G) will not be a parameter any longer.
To seek a nontrivial solution, we notice that at low frequencies, the right hand side of equation (9) associated with solid and fluid accelerations will vanish. Consequently, the condition of zero rock stress, i.e., equation (9), is approximately satisfied at low frequencies. Eventually, equations (6)-(8) result in the following storativity equation: where S = ηβ f + β eff − β s − ηβ s is storativity (specific storage). Equation (10) is similar to a diffusion equation in heat conductance [36].
The plane wave solution of equation (10) yields phase velocity (v p ) and the quality factor (Q p ) as functions of angular frequency (ω).
Biot [11,12] theory used the Euler derivative (∂q/∂t) to approximate the Lagrange acceleration (dq/dt). The latter, however, consists of the former and the inertial terms (involving the spatial derivate). The inertial term is precisely the non-Darcy flow effect in [3]. In a 1D scenario, the inertial term vanishes such that Biot theory is accurate. Nevertheless, in the case where the spatial derivative is not small any longer [3], Biot theory may not be accurate; such a complicated case is beyond the scope of this paper.

Illustrative Example
Berea sandstone is a classic sandstone [37]. Its parameters are listed in Table 1. Velocities of P and S waves in the dry sandstone were measured as 3300 and 2100 m/s, respectively [38]. Hence, the bulk and shear modules are calculated to be 10.58 and 9.30 GPa, respectively. With the parameters, slow P wave [11] has v p and Q p as depicted in Figures 2-4. The two quantities in storativity equation (10) are calculated via equation (11) and also plotted in Figures 2-4 for comparison. Figure 3 is the amplified version of Figure 2 at frequency lower than 1 Hz, in order to resolve the small difference of v p between the two models.
As shown in Figure 2, Biot [11] theory and storativity equation (10) differ significantly at frequency higher than 100 kHz. The latter has phase velocity increasing to infinite large, while the former has phase velocity approaching to a bounded value. At frequency lower than 1 Hz, the latter has a velocity lower than the former by 7-8 percent.
In Figure 4, Q p is constantly at 0.5. At frequency lower than 100 kHz, Biot [11] theory and storativity equation (10) yield the same Q p . For frequency above 100 kHz, Q p   [11] theory and storativity equation (10), for water-saturated Berea sandstone and with f denoting frequency. 4 Geofluids differs increasingly between the two models. At very high frequencies, the former has a vanishing attenuation, well consistent with constant v p in Figure 2. This consistency is reasonable in accordance with Kramers-Kronig relations [7].

Discussion
As shown in Figure 2, phase velocity yielding from storativity equation (10) will become infinite large with the increase of frequency, which is unrealistic because any physical processes need certain accounts of time. The disadvantage inherent in storativity equation arises from the simplification of Biot [11] theory. Ignoring the compressibility coefficient of solid material (β s = 0) and setting the Biot-Willis coefficient α B = 1, Dome-nico and Schwartz [5] derived the storativity equation as follows: The derivation in Appendix A in [5] differs from our above approach in that the former explicitly used mass conservation of solid in Euler approach. In contrast, our derivation of equation (3) used Lagrangian approach to analyze poroelasticity such that the mass conservation of solid material is implicitly included in the equation. Comparing our equation (10) with in equation (12) from [5], ours is slightly more accurate in considering β s .
The error of storativity equation in flow velocity can be calculated using plane wave, P p = exp ½iωðt − ðx/v p ÞÞ − ðωx/2v p Q p Þ. With the consideration that the skeleton Lagrange velocity is typically much smaller than the fluid Lagrange velocity, substituting this plane wave into equation (4) yields: In Figure 3, at frequency lower than 100 kHz, Biot [11] theory has v p 7-8 percent faster than the storativity equation (10) has. According to equation (13), flow velocity in the former is 7-8 percent lower than that in the latter. In short, Biot theory has a higher velocity of pressure wave, a longer wavelength, a smaller pressure gradient, and a lower flow (Lagrange) velocity than storativity equation has.
At frequency lower than 100 kHz, with the same parameters, Biot theory has wave velocity 7-8 percent  [11] theory, for water saturated Berea sandstone and with f denoting frequency.

Geofluids
higher than storativity equation has. According to equation (11), in order for the latter to yield the same velocity of fluid pressure wave as the former does, either its apparent storativity (S) has to be 13 percent lower than the real one, or its apparent permeability (k D ) has to be higher than the real permeability by 14 percent. The real permeability is measured (via Darcy law) by flux rate and pressure gradient in steady flow experiments. As such, the above result suggests that the apparent permeability indirectly measured via storativity equation [28,29] may be 14 percent higher than the real permeability.
For water-saturated Berea sandstone, Biot theory has an asymptotical velocity of 1294 m/s at infinite large frequency. We set the skeleton and solid material to be rigid, to investigate phase velocity and the quality factor of fluid pressure wave in [11] (Figures 5 and 6, respectively). At very high frequencies, fluid pressure wave has the same velocity as sound wave (fast P wave) in pure water (1474 m/s) [42]. Therefore, it is skeleton compressibility that decreases phase velocity of fluid pressure wave.
Dimensional analysis of Biot [11] theory reveals that phase velocity and the quality factor of fluid pressure wave are controlled by dimensionless angular frequency (Ω) as follows: where fluid density rather than skeleton density is used because fluid pressure wave depends largely on fluid property. According to equation (14), a very high ω with a bounded k D is equivalent to a very high k D with a bounded ω; both increase Ω to a very large number. In other words, the high frequency limit in Figures 5 and 6 is equivalent to the high permeability limit. Therefore, for rigid rock with a very high permeability, fluid pressure wave propagates with the same velocity as sound wave in pure water does. Similarly, a very low ω with a bounded k D is equivalent to a very low k D with a bounded ω; both decreases Ω to a very small number. In other words, the low frequency limit is equivalent to the low permeability limit in which groundwater pressure transmits very slowly.
Steady flow is driven by static pressure gradient which has a vanishing wave velocity at zero frequency as depicted in Figure 3. Besides hydrogeology, storativity equation is also used in hydraulic fracturing [43]. As Biot [11] theory is more accurate than storativity equation, this study provides a useful reference for people focusing on the equation.
In reservoir engineering, capillary force in a matrix (or very fine fissures) is not ignorable. Actually, the matrix capillary pressure is used purposely to suck water into the matrix, thus forcing petroleum out of the matrix [44,45]. The capillary force represents the pressure difference between two fluid phases (otherwise the net force on the interface between the two phases would not vanish and would yield an infinite large acceleration of the interface, which is unrealistic). This paper is focused on a single phase of fluid (pressure is continuous everywhere in the fluid) and therefore is limited in not considering the matrix capillary pressure between two fluid phases, which is an interesting topic of future study.

Summary and Conclusions
(1) Under the assumption of low frequency, the constitutional relation of rock skeleton, the fluid mass conservation, and the fluid momentum equation in Biot theory can be simplified to storativity equation. For water-saturated Berea sandstone, at frequency lower than 100 kHz, both phase velocity and the quality factor of fluid pressure wave are close between the equation and Biot theory (2) At frequency lower than 100 kHz, the storativity equation has a phase velocity 7-8 percent lower than Biot theory has. The apparent permeability indirectly measured via storativity equation such as OHT may be 14 percent higher than the real permeability directly measured by steady flow experiments (3) At a frequency of 1 Hz, Biot theory has a phase velocity of approximately 3 m/s ( Figure 3) which is very slow compared with the sound velocity in pure water (1474 m/s), but is far over the velocity of groundwater seepage. If skeleton is rigid, Biot theory at very high frequencies or with very high permeabilities will yield the same velocity as sound wave in pure water Nomenclature G: Shear modulus of skeleton k D : Darcy permeability K ef f : Bulk modulus of skeleton (K eff = 1/β ef f ) P C : The confining pressure of rock P p : Fluid pressure in pores P w : Fluid pressure in a wellbore  Figure 6: Same as Figure 5 except for the quality factor (Q p ). 6 Geofluids q: Darcy flux rate in the x direction q * : Darcy flux rate in a vector Q p : The quality factor of a plane P wave S: Storativity (specific storage, S = ηβ f + β ef f − β s − ηβ s ) u: Skeleton displacement in the x direction u * : Skeleton displacement in a vector v: Lagrange velocity of skeleton (v = ∂u/∂t is the approximation by the Euler derivative) v p : Phase velocity of a plane P wave α B : The Biot-Willis coefficient (α B = ðβ eff − β s Þ/β eff ) β ef f : Compressibility coefficient of skeleton β f : Compressibility coefficient of fluid β s : Compressibility coefficient of solid material η: Rock porosity μ: Fluid viscosity ω: Angular frequency Ω: Dimensionless angular frequency (Ω = ρ f ðk D /μÞω) ρ f : Fluid density ρ s : Skeleton density ρ 12 : The coupling density in Biot theory σ 11 : The normal stress in the x direction.

Additional Points
Highlights. Biot theory at low frequencies can simplify to storativity equation. At frequency lower than 100 kHz, Biot theory yields wave velocity higher than storativity equation does. Apparent permeability measured by OHT may be 14 percent higher than real permeability

Conflicts of Interest
The authors declare that they have no conflicts of interest.