Wave-induced Fluid Flow Research Articles

Seismic wave attenuation and modulus dispersion in sandstones

We have conducted laboratory experiments over the 1–200 Hz band to examine the effects of viscosity and permeability on modulus dispersion and attenuation in sandstones and also to examine the effects of partial gas or oil saturation on velocities and attenuations. Our results have indicated that bulk modulus values with low-viscosity fluids are close to the values predicted using Gassmann’s first equation, but, with increasing frequency and viscosity, the bulk and shear moduli progressively deviate from the values predicted by Gassmann’s equations. The shear moduli increase up to 1 GPa (or approximately 10%) with high-viscosity fluids. The P- and S-wave attenuations ([Formula: see text] and [Formula: see text]) and modulus dispersion with different fluids are indicative of stress relaxations that to the first order are scaling with frequency times viscosity. By fitting Cole-Cole distributions to the scaled modulus and attenuation data, we have found that there are similar P-wave, shear and bulk relaxations, and attenuation peaks in each of the five sandstones studied. The modulus defects range from 11% to 15% in Berea sandstone to 16% to 26% in the other sandstones, but these would be reduced at higher confining pressures. The relaxations shift to lower frequencies as the viscosity increased, but they do not show the dependence on permeability predicted by mesoscopic wave-induced fluid flow (WIFF) theories. Results from other experiments having patchy saturation with liquid [Formula: see text] and high-modulus fluids are consistent with mesoscopic WIFF theories. We have concluded that the modulus dispersion and attenuations ([Formula: see text] and [Formula: see text]) in saturated sandstones are caused by a pore-scale, local-flow mechanism operating near grain contacts.

Seismic attenuation in partially saturated Berea sandstone submitted to a range of confining pressures

AbstractUsing the forced oscillation method, we measure the extensional‐mode attenuation and Young's modulus of a Berea sandstone sample at seismic frequencies (0.5–50 Hz) for varying levels of water saturation (~0–100%) and confining pressures (2–25 MPa). Attenuation is negligible for dry conditions and saturation levels <80%. For saturation levels between ~91% and ~100%, attenuation is significant and frequency dependent in the form of distinct bell‐shaped curves having their maxima between 1 and 20 Hz. Increasing saturation causes an increase of the overall attenuation magnitude and a shift of its peak to lower frequencies. On the other hand, increasing the confining pressure causes a reduction in the attenuation magnitude and a shift of its peak to higher frequencies. For saturation levels above ~98%, the fluid pressure increases with increasing confining pressure. When the fluid pressure is high enough to ensure full water saturation of the sample, attenuation becomes negligible. A second series of comparable experiments reproduces these results satisfactorily. Based on a qualitative analysis of the data, the frequency‐dependent attenuation meets the theoretical predictions of mesoscopic wave‐induced fluid flow (WIFF) in response to a heterogeneous water distribution in the pore space, so‐called patchy saturation. These results show that mesoscopic WIFF can be an important source of seismic attenuation at reservoir conditions.

Shear properties of heterogeneous fluid-filled porous media with spherical inclusions

An exact analytical solution is presented for the effective dynamic shear modulus in a heterogeneous fluid-filled poroelastic medium containing spherical inclusions. The complex and frequency-dependent properties of the derived shear modulus are solely caused by the physical mechanism of mesoscopic-scale wave-induced fluid flow whose scale is assumed to be smaller than wavelength but larger than the size of pores. Our model consists of three phases: a spherical inclusion, a shell of matrix material with different mechanical and/or hydraulic properties and an outer region of effective homogeneous medium of infinite extent. This three-phase model represents a self-consistent model or an approximate model of composite having periodically distributed inclusions. The behaviors of both the inclusion and the matrix are described by Biot's equations (1941) with standard conditions of Deresiewicz and Skalak (1963) at the inclusion-matrix interface. The effective medium is regarded as an equivalent elastic or viscoelastic material with complex and frequency-dependent moduli to be determined. The derived effective shear modulus is used to quantify the shear-wave attenuation and velocity dispersion. For the problem of fluid patchy saturation (inclusions with the same solid frame as the matrix but with a different pore fluid from that in the matrix), the gas pocket does not affect the shear attenuation and dispersion characteristic of the water-filled matrix medium at all. For the problem of double porosity structure (inclusions having a different solid frame than the matrix but the same pore fluid as the matrix), with the increase of frequency the heterogeneous medium transitions from a low-frequency state having drained inclusions and drained matrix with no pore pressure difference to a higher-frequency state having undrained inclusions and undrained matrix with no fluid communication at the inclusion's surface. The relaxation frequency at which the maximum value of inverse quality factor occurs moves to frequencies by two orders of magnitude lower if the size of a unit cell increases by one order of magnitude. Stiff inclusions imbedded in a relatively soft matrix can cause significant and observable attenuation at seismic frequency bands, but softer inclusions imbedded in a relatively stiff matrix cause very weak attenuation. The mixed heterogeneity in both the solid frame and pore fluid also has important influences on the frequency-dependent shear wave attenuation.

Scattering of SH-wave for a non-welded contact fluid saturated natural fracture

Closed-form time-domain formulas have been obtained for reflected and transmitted coefficients of SH-wave across a saturated fracture. SH-wave incident upon a fluid filling saturated fracture at arbitrary angle containing the scattering attenuation and wave-induced fluid flow attenuation mechanism. The fracture linear slip theory is introduced under the saturated condition and the dynamic shear compliance is expressed by the viscous fluid mechanics. A finite jump displacement and velocity across a fracture is expressed as a function of wave frequency and geometrical parameters of fracture. The frequency-dependent wave reflection and transmission in wave equations are given by the discontinuity boundary conditions. Coupling relationship between scattering and intrinsic attenuation can be explained partially based on our model.

Seismic Attenuation and Velocity Dispersion to Discriminate Gas Hydrates and Free Gas Zone, Makran Offshore, Pakistan

Gas hydrates gained a remarkable attention as an unconventional energy resource recently. In order to interpret gas hydrates (part of fluid) and free gas saturated zone accurately, it is essential to implement new technique related to seismic attenuation and velocity dispersion. P wave attenuation and velocity dispersion in porous media made promising imprints for exploration of gas hydrates. The most prominent phenomenon for attenuation and velocity dispersion in porous media is wave induced fluid flow in which wave inhomogeneities are larger than pore size but smaller than wavelength. Numerical simulation technique is applied to analyze frequency dependent velocity dispersion and attenuation in gas hydrates and free gas layer in Makran offshore of Pakistan. Homogeneous and patchy distribution patterns of gas hydrates and free gas within pore spaces of host sediments at lower and higher frequency regime are considered. It is noted that the attenuation and velocity dispersion increase with the increase in gas hydrates saturation. The maximum attenuation is observed at 66% saturation of gas hydrates in the area under investigation. However, in case of water and gas mixture the maximum attenuation and velocity dispersion occur at low gas saturation (~15%). Therefore, based on our numerical simulation, velocity dispersion and attenuation can be used as seismic attributes to differentiate various gas saturations and gas hydrates saturation for Makran offshore area of Pakistan.

Energy dissipation of P- and S-waves in fluid-saturated rocks: An overview focusing on hydraulically connected fractures

An important characteristic of fractured rocks is their high seismic attenuation, which so far has been mainly attributed to wave-induced fluid flow (WIFF) between the fractures and the embedding matrix. The influence of fracture connectivity on seismic attenuation has, however, recently, only been investigated. Numerical compressibility and shear tests based on Biot’s quasi-static poro-elastic equations illustrate that an important manifestation of WIFF arises in the presence of fracture connectivity. The corresponding energy loss, which can be significant for both P- and S-waves, is mainly due to fluid flow within the connected fractures and is sensitive to the permeabilities as well as the lengths and intersection angles of the fractures. Correspondingly, this phenomenon contains valuable information on the governing hydraulic properties of fractured rocks and hence should be accounted for whenever realistic seismic models of such media are needed.

Modulus defect, velocity dispersion and attenuation in partially-saturated reservoirs of Jurassic sandstone, Indus Basin, Pakistan

Partially saturated reservoirs are one of the major sources of seismic wave attenuation, modulus defect and velocity dispersion in real seismic data. The main attenuation and dispersion phenomenon is wave induced fluid flow due to the heterogeneity in pore fluids or porous rock. The identification of pore fluid type, saturation and distribution pattern within the pore space is of great significance as several seismic and petrophysical properties of porous rocks are largely affected by fluid type, saturation and fluid distribution pattern. Based on Gassmann-Wood and Gassmann- Hill rock physics models modulus defect, velocity dispersion and attenuation in Jurassic siliclastic partially-saturated rocks are studied. For this purpose two saturation patterns - uniform and patchy - are considered within the pore spaces in two frequency regimes i.e., lower frequency and higher frequency. The results reveal that at low enough frequency where saturation of liquid and gas is uniform, the seismic velocity and bulk modulus are lower than at higher frequency where saturation of fluid mixture is in the form of patches. The velocity dispersion and attenuation is also modeled at different levels of gas saturation. It is found that the maximum attenuation and velocity dispersion is at low gas saturation. Therefore, the dispersion and attenuation can provide a potential way to predict gas saturation and can be used as a property to differentiate low from high gas saturation.

Attenuation of sonic waves in water-saturated alluvial sediments due to wave-induced fluid flow at microscopic, mesoscopic and macroscopic scales

Attenuation of sonic waves in water-saturated alluvial sediments due to wave-induced fluid flow at microscopic, mesoscopic and macroscopic scales

Seismic attenuation due to heterogeneities of rock fabric and fluid distribution

AbstractThe study of the influence of rock fabric and porefluid distribution on the seismic wavefield is important for the prediction and detection of reservoir properties such as lithology and fluid saturation. Wave-induced local fluid flow (WILFF), which is affected by local heterogeneities of the pore structure and fluid saturation, is believed to be the main mechanism to explain the measured attenuation levels at different frequency bands. These two types of heterogeneity affect seismic waves as a combined WILFF process. In this work, we consider a double-porosity system, each part with a different compressibility and patchy saturation, and derive the wave equations from Hamilton's principle. A plane-wave analysis yields the properties of the classical P-wave and those of the four slow waves. The examples show that patchy fluid saturation dominates the peak frequency of the relaxation mechanism. The relation between seismic anelasticity (velocity dispersion and attenuation) and saturation depends on frequency and on the geometrical features of the two heterogeneities. The proposed theory constitutes the comprehensive description for wave propagation process through reservoirs rocks of shallow Earth and porous media in general, to estimate fluid content and distribution.

Water saturation effects onP-wave anisotropy in synthetic sandstone with aligned fractures

The seismic properties of rocks are known to be sensitive to partial liquid or gas saturation, and to aligned fractures. P-wave anisotropy is widely used for fracture characterization and is known to be sensitive to the saturating fluid. However, studies combining the effect of multiphase saturation and aligned fractures are limited even though such conditions are common in the subsurface. An understanding of the effects of partial liquid or gas saturation on P-wave anisotropy could help improve seismic characterization of fractured, gas bearing reservoirs. Using octagonal-shaped synthetic sandstone samples, one containing aligned penny-shaped fractures and the other without fractures, we examined the influence of water saturation on P-wave anisotropy in fractured rocks. In the fractured rock, the saturation related stiffening effect at higher water saturation values is larger in the direction across the fractures than along the fractures. Consequently, the anisotropy parameter ‘?’ decreases as a result of this fluid stiffening effect. These effects are frequency dependent as a result of wave-induced fluid flow mechanisms. Our observations can be explained by combining a frequency-dependent fractured rock model and a frequency-dependent partial saturation model.

Seismic reflection dispersion due to wave-induced fluid flow in heterogeneous reservoir rocks

We have investigated the impact of wave-induced fluid flow, including Biot flow and mesoscopic flow, on the signatures of seismic reflectivity in heterogeneous reservoir rocks. We have incorporated the dynamic poroelastic responses of mesoscopic flow into the classical Biot theory. The resulting effective Biot media could capture the characteristics of velocity dispersion and wave attenuation in heterogeneous poroelastic media. On the basis of this effective Biot media, an approach was developed to compute the poroelastic reflection at arbitrary angles and frequencies from the boundary of two heterogeneous porous media. The computed poroelastic reflections not only depended on the elastic properties’ contrast and incident angle, but also relied on the fluid mobility and observational frequency. For a typical sand-shale reflector, with the given rock and fluid properties, we found that the effect of mesoscopic flow causes P-wave reflection amplitude variations with the frequency being as high as 40% and a maximum phase shift as high as 16° at the seismic exploration frequency band. In addition, it was found that the amplitude variation with offset intercept and the gradient at the poroelastic interface were impacted by the mesoscopic flow and had a decreasing trend with frequency. Therefore, ignoring the impact of mesoscopic flow could possibly lead to uncertainty in seismic imaging as well as quantitative interpretation of reservoir properties. In comparison, the Biot flow-induced seismic dispersion effect, which occured at a very high-frequency range, was almost negligible.

Modeling attenuation and dispersion in porous heterogeneous rocks with dynamic fluid modulus

Based on poroelasticity analysis, we developed a new concept of frequency-dependent dynamic fluid modulus (DFM) to understand the velocity dispersion and wave attenuation due to wave-induced fluid flow. Conventional applications of Gassmann’s equation require a complete homogeneity inside the porous media (or equivalently, at zero frequency) and closed boundary condition. We first analyzed the fluid effect on the bulk modulus in homogeneous porous media with a nonclosed condition. The partial drainage of pore fluid causes additional pore volume change under applied stress. An incoming fluid flow stiffens the porous system, and an outgoing flow softens it. Such a phenomenon can still be effectively formulated as a closed system with Gassmann’s equation, by introducing a DFM, which adds a flow term into the original fluid modulus. We further proved that in heterogeneous porous media, the wave-induced internal fluid flow caused additional bulk volume deformation. It equals the amount of the fluid flow from a soft phase to a stiff phase times the difference of Skempton coefficients between the two phases. Again, with a frequency-dependent complex number DFM, the use of Gassmann’s equation can be extended into heterogeneous porous media at nonzero frequencies to fully characterize its viscoelastic behavior on rock’s bulk modulus. We evaluated three examples to show how to model the P-wave attenuation and velocity dispersion in porous media due to microscopic and mesoscopic scale heterogeneities, as well as in the presence of multiple sets of heterogeneities. Finally, we demonstrated that the DFM can be conveniently and deterministically inverted from the measured lab data. The inverted results can be used to identify the possible mechanisms behind the observed wave attenuation and velocity dispersion, and they were good indicators of the degree and distribution of heterogeneities inside the rocks.

Including poroelastic effects in the linear slip theory

Numerical simulations of seismic wave propagation in fractured media are often performed in the framework of the linear slip theory (LST). Therein, fractures are represented as interfaces and their mechanical properties are characterized through a compliance matrix. This theory has been extended to account for energy dissipation due to viscous friction within fluid-filled fractures by using complex-valued frequency-dependent compliances. This is, however, not fully adequate for fractured porous rocks in which wave-induced fluid flow (WIFF) between fractures and host rock constitutes a predominant seismic attenuation mechanism. In this letter, we develop an approach to incorporate WIFF effects directly into the LST for a 1D system via a complex-valued, frequency-dependent fracture compliance. The methodology is validated for a medium permeated by regularly distributed planar fractures, for which an analytical expression for the complex-valued normal compliance is determined in the framework of quasistatic poroelasticity. There is good agreement between synthetic seismograms generated using the proposed recipe and those obtained from comprehensive, but computationally demanding, poroelastic simulations.

An analytical study of seismoelectric signals produced by 1-D mesoscopic heterogeneities

The presence of mesoscopic heterogeneities in fluid-saturated porous rocks can produce measurable seismoelectric signals due to wave-induced fluid flow between regions of differing compressibility. The dependence of these signals on the petrophysical and structural characteristics of the probed rock mass remains largely unexplored. In this work, we derive an analytical solution to describe the seismoelectric response of a rock sample, containing a horizontal layer at its center, that is subjected to an oscillatory compressibility test. We then adapt this general solution to compute the seismoelectric signature of a particular case related to a sample that is permeated by a horizontal fracture located at its center. Analyses of the general and particular solutions are performed to study the impact of different petrophysical and structural parameters on the seismoelectric response. We find that the amplitude of the seismoelectric signal is directly proportional to the applied stress, to the Skempton coefficient contrast between the host rock and the layer, and to a weighted average of the effective excess charge of the two materials. Our results also demonstrate that the frequency at which the maximum electrical potential amplitude prevails does not depend on the applied stress or the Skempton coefficient contrast. In presence of strong permeability variations, this frequency is rather controlled by the permeability and thickness of the less permeable material. The results of this study thus indicate that seismoelectric measurements can potentially be used to estimate key mechanical and hydraulic rock properties of mesoscopic heterogeneities, such as compressibility, permeability, and fracture compliance.

Effect of local fluid flow on the propagation of elastic waves in a transversely isotropic double-porosity medium

This study considers the propagation of elastic waves in a fluid-saturated double-porosity solid. Constitutive relations are derived for the double porosity medium having the transverse isotropy with vertical axis of symmetry. The equations of motion are solved for the propagation of harmonic waves. Then, the amount of wave-induced fluid flow at any point is expressed in terms of normal strain components in the composite medium. For propagation in a plane, the equations of motion decouple into two independent systems. One of these represents the particle motion normal to the plane and identifies the SH-type motion for solid particles. The other system represents the in-plane motion of particles and governs the propagation of four coupled waves. Among these four waves, three are quasi-longitudinal (qP1, qP2, qP3) waves and one is a quasi-transverse (qSV) wave. A numerical example is solved to calculate the velocity and attenuation anisotropy of the five waves. Effect of local fluid flow is observed on the wave velocities as well as the polarization of the constituent particles in double porosity medium. Effects of wave-frequency, propagation direction, pore-fluid viscosity, permeability and radius of spherical inclusions are analysed on the velocities and attenuation.

Seismic modelling study of P-wave attenuation and velocity dispersion in patchy-saturated porous media

Seismic wave propagation in patchy-saturated porous media is studied by numerical simulation in time domain at the seismic frequency band (1–1000 Hz). The models consist of hundreds of representative elementary volumes (REVs), where the REV is partially saturated with water and gas pockets. Seismic modelling experiments are implemented in a traditional way, with 'periodic' boundary conditions applied to get rid of undrained boundary conditions at the outer edges of the REVs. The characteristics of confining pressure, induced pore pressure, solid particle velocities and Darcy filtration velocities are analysed. The snapshots show that strong pore pressure gradients are generated across the interface between gas and water phases, and significant fluid flow occurs. The conversion of a fast P-wave into a dissipating slow P-wave takes place at seismic frequencies, and the converted slow P-wave diffuses strongly in both gas- and water-saturated phases. These numerical results can help us to understand the loss mechanism at seismic frequencies. Then, P-wave attenuation and velocity dispersion of a heterogeneous REV are calculated during traditional seismic modelling at seismic frequencies. The numerical results show good agreement with theoretical predictions obtained from patchy saturation theory. Furthermore, the effects of different fluid distributions on P-wave attenuation and velocity dispersion are analysed numerically. A series of experiments are implemented by considering large, small and random gas-patchy inclusions. The decrease of gas pocket size makes the peak frequency move towards high frequencies. Random distribution of gas patches may affect both the peak attenuation and peak frequencies. Seismic attenuation caused by Biot global flow, elastic scattering and wave-induced fluid flow (WIFF) associated with patchy saturation are computed numerically. The results show that the contribution of Biot's global flow and scattering to the overall attenuation are negligible. Mesoscopic WIFF may dominate the loss mechanism in patchy-saturated media at seismic frequencies (1–1000 Hz).

Experimental observation of water saturation effects on shear wave splitting in synthetic rock with fractures aligned at oblique angles

Fractured rocks are known to exhibit seismic anisotropy and shear wave splitting (SWS). SWS is commonly used for fractured rock characterization and has been shown to be sensitive to fluid type. The presence of partial liquid/gas saturation is also known to affect the elastic properties of rocks. The combined effect of both fractures and partial liquid/gas saturation is still unknown. Using synthetic, silica-cemented sandstones with aligned penny-shaped voids, we conducted laboratory ultrasonic experiments to investigate the effect fractures aligned at an oblique angle to wave propagation would have on SWS under partial liquid/gas saturation conditions. The result for the fractured rock shows a saturation dependence which can be explained by combining a fractured rock model and a partial saturation model. At high to full water saturation values, SWS decreases as a result of the fluid bulk modulus effect on the quasi-shear wave. This bulk modulus effect is frequency dependent as a result of wave-induced fluid flow mechanisms, which would in turn lead to frequency dependent SWS. This result suggests the possible use of SWS for discriminating between full liquid saturation and partial liquid/gas saturation.

Seismic attenuation and velocity dispersion in fractured rocks: The role played by fracture contact areas

ABSTRACTThe presence of fractures in fluid‐saturated porous rocks is usually associated with strong seismic P‐wave attenuation and velocity dispersion. This energy dissipation can be caused by oscillatory wave‐induced fluid pressure diffusion between the fractures and the host rock, an intrinsic attenuation mechanism generally referred to as wave‐induced fluid flow. Geological observations suggest that fracture surfaces are highly irregular at the millimetre and sub‐millimetre scale, which finds its expression in geometrical and mechanical complexities of the contact area between the fracture faces. It is well known that contact areas strongly affect the overall mechanical fracture properties. However, existing models for seismic attenuation and velocity dispersion in fractured rocks neglect this complexity. In this work, we explore the effects of fracture contact areas on seismic P‐wave attenuation and velocity dispersion using oscillatory relaxation simulations based on quasi‐static poroelastic equations. We verify that the geometrical and mechanical details of fracture contact areas have a strong impact on seismic signatures. In addition, our numerical approach allows us to quantify the vertical solid displacement jump across fractures, the key quantity in the linear slip theory. We find that the displacement jump is strongly affected by the geometrical details of the fracture contact area and, due to the oscillatory fluid pressure diffusion process, is complex‐valued and frequency‐dependent. By using laboratory measurements of stress‐induced changes in the fracture contact area, we relate seismic attenuation and dispersion to the effective stress. The corresponding results do indeed indicate that seismic attenuation and phase velocity may constitute useful attributes to constrain the effective stress. Alternatively, knowledge of the effective stress may help to identify the regions in which wave induced fluid flow is expected to be the dominant attenuation mechanism.

Sensitivity of S-wave attenuation to the connectivity of fractures in fluid-saturated rocks

Biot’s equations of poroelasticity were solved to study the effects of fracture connectivity on S-wave attenuation caused by wave-induced fluid flow at the mesoscopic scale. The methodology was based on numerical quasistatic pure-shear experiments performed on models of water-saturated rocks containing pairs of either connected or unconnected fractures of variable inclination. Each model corresponded to a representative elementary volume of a periodic medium. Inertial terms were neglected, and hence, the observed attenuation was entirely due to wave-induced fluid flow at the mesoscopic scale. We found that when fractures are not connected, fluid flow in the embedding matrix governs S-wave attenuation, whereas fluid flow through highly permeable fractures, from one fracture into the other one, may dominate when fractures are connected. Each of these energy-dissipation phenomena has a distinct characteristic frequency, with the S-wave attenuation peak associated with flow through connected fractures occurring at higher frequencies than that associated with flow in the embedding matrix. Exploring a range of geometric arrangements of either connected or unconnected fractures at different inclinations, we also observed that the magnitude of S-wave attenuation at both characteristic frequencies shows a strong dependence on fracture inclination. For comparison, we performed quasistatic uniaxial compressibility tests to compute P-wave attenuation in the same models. We found that the attenuation patterns of S-waves tend to differ fundamentally from those of P-waves with respect to fracture inclination. The attenuation characteristics of P- and S-waves in fractured media are thus, largely complementary. With respect to fracture connectivity, we observed that S-wave attenuation tends to follow a specific pattern, indeed, more consistently than that of the P-waves. Our results point to the promising perspective of combining estimates of attenuation of P- and S-waves to infer information on fracture connectivity as well as on the effective permeability of fractured media.

Finite-element harmonic experiments to model fractured induced anisotropy in poroelastic media

Fractures in a fluid-saturated poroelastic–Biot–medium can be modeled as very thin highly permeable and compliant layers within a porous background. A Biot medium containing a dense set of aligned fractures behaves as an effective transversely isotropic and viscoelastic (TIV) medium at the macroscale when the predominant wavelength is much larger than the average distance between fractures. One important mechanism in Biot media at seismic frequencies is wave-induced fluid flow generated by fast compressional waves at mesoscopic-scale heterogeneities, generating slow diffusion-type Biot waves. In this work, we present and analyze a collection of time-harmonic finite element experiments that take into account the effects of the presence of aligned fractures and interlayer fluid flow occurring at the mesoscale, allowing us to determine the complex and frequency dependent stiffnesses of the effective TIV medium at the macroscale.These numerical upscaling experiments are defined as boundary value problems on representative samples of the fractured material, with boundary conditions associated with compressibility and shear tests, which are solved using the finite element (FE) method. The FE space chosen to discretize each component of the solid displacement vector is that of globally continuous piecewise bilinear functions, while for the fluid phase the vector part of the Raviart–Thomas–Nedelec space of zero order is employed. We present results on the uniqueness of the solution of the continuous and discrete problems, and derive optimal a priori energy error estimates. First, the numerical results are validated with those of a theory valid for fluid flow perpendicular to the fracture layering and independent of the loading direction, so that the attenuation mechanism can be represented by a single relaxation function. Then, the methodology is applied to cases for which no analytical solutions are available, such as a fractured Biot medium saturated with brine and patches of CO2 and a brine saturated sample of uniform background and fractures with fractal variations in their petrophysical properties.