Biphasic media with a dynamic interaction between fluid and solid phases must be taken into account to accurately describe seismic wave amplitudes in subsurface and reservoir geophysical applications. Consequently, the modeling of the wave propagation in heteregeneous porous media, which includes the frequency-dependent phenomena of the fluid-solid interaction, is considered for 2D geometries. From the Biot-Gassmann theory, we have deduced the discrete linear system in the frequency domain for a discontinuous finite-element method, known as the nodal discontinuous Galerkin method. Solving this system in the frequency domain allows accurate modeling of the Biot wave in the diffusive/propagative regimes, enhancing the importance of frequency effects. Because we had to consider finite numerical models, we implemented perfectly matched layer techniques. We found that waves are efficiently absorbed at the model boundaries, and that the discretization of the medium should follow the same rules as in the elastodynamic case, that is, 10 grids per minimum wavelength for a P0 interpolation order. The grid spreading of the sources, which could be stresses or forces applied on either the solid phase or the fluid phase, did not show any additional difficulties compared to the elastic problem. For a flat interface separating two media, we compared the numerical solution and a semianalytic solution obtained by a reflectivity method in the three regimes where the Biot wave is propagative, diffusive/propagative, and diffusive. In all cases, fluid-solid interactions were reconstructed accurately, proving that attenuation and dispersion of the waves were correctly accounted for. In addition to this validation in layered media, we have explored the capacities of modeling complex wave propagation in a laterally heterogeneous porous medium related to steam injection in a sand reservoir and the seismic response associated to a fluid substitution.

The quantitative estimation of rock physics properties is of great importance in any reservoir characterization. We have studied the sensitivity of such poroelastic rock physics properties to various seismic viscoelastic attributes (velocities, quality factors, and density). Because we considered a generalized dynamic poroelastic model, our analysis was applicable to most kinds of rocks over a wide range of frequencies. The viscoelastic attributes computed by poroelastic forward modeling were used as input to a semiglobal optimization inversion code to estimate poroelastic properties (porosity, solid frame moduli, fluid phase properties, and saturation). The sensitivity studies that we used showed that it was best to consider an inversion system with enough input data to obtain accurate estimates. However, simultaneous inversion for the whole set of poroelastic parameters was problematic due to the large number of parameters and their trade-off. Consequently, we restricted the sensitivity tests to the estimation of specific poroelastic parameters by making appropriate assumptions on the fluid content and/or solid phases. Realistic a priori assumptions were made by using well data or regional geology knowledge. We found that (1) the estimation of frame properties was accurate as long as sufficient input data were available, (2) the estimation of permeability or fluid saturation depended strongly on the use of attenuation data, and (3) the fluid bulk modulus can be accurately inverted, whereas other fluid properties have a low sensitivity. Introducing errors in a priori rock physics properties linearly shifted the estimations, but not dramatically. Finally, an uncertainty analysis on seismic input data determined that, even if the inversion was reliable, the addition of more input data may be required to obtain accurate estimations if input data were erroneous.

Reliable quantification of carbon dioxide (${\mathrm{CO}}_2$) properties and saturation is crucial in the monitoring of ${\mathrm{CO}}_2$ underground storage projects. We have focused on quantitative seismic characterization of ${\mathrm{CO}}_2$ at the Sleipner storage pilot site. We evaluate a methodology combining high-resolution seismic waveform tomography, with uncertainty quantification and rock physics inversion. We use full-waveform inversion (FWI) to provide high-resolution estimates of P-wave velocity $V_{\mathrm{P}}$ and perform an evaluation of the reliability of the derived model based on posterior covariance matrix analysis. To get realistic estimates of ${\mathrm{CO}}_2$ saturation, we implement advanced rock physics models taking into account effective fluid phase theory and patchy saturation. We determine through sensitivity tests that the estimation of ${\mathrm{CO}}_2$ saturation is possible even when using only the P-wave velocity as input. After a characterization of rock frame properties based on log data prior to the ${\mathrm{CO}}_2$ injection at Sleipner, we apply our two-step methodology. The FWI result provides clear indications of the injected ${\mathrm{CO}}_2$ plume being observed as low-velocity zones corresponding to thin ${\mathrm{CO}}_2$ filled layers. Several tests, varying the rock physics model and ${\mathrm{CO}}_2$ properties, are then performed to estimate ${\mathrm{CO}}_2$ saturation. The results suggest saturations reaching 30%–35% in the thin sand layers and up to 75% when patchy mixing is considered. We have carried out a joint estimation of saturation with distribution type and, even if the inversion is not well-constrained due to limited input data, we conclude that the ${\mathrm{CO}}_2$ has an intermediate pattern between uniform and patchy mixing, which leads to saturation levels of approximately $25\%\pm 15\%$. It is worth noting that the 2D section used in this work is located 533 m east of the injection point. We also conclude that the joint estimation of ${\mathrm{CO}}_2$ properties with saturation is not crucial and consequently that knowing the pressure and temperature state of the reservoir does not prevent reliable estimation of ${\mathrm{CO}}_2$ saturation.

Partially saturated rocks are considered to be major sources of seismic wave velocity dispersion and attenuation in recorded real data. From the physical description of partially saturated gas-water and oil-water reservoirs, we use upscaling theories to compute an equivalent frequency-dependent porous medium. These homogenization methods are associated with mesoscale description of attenuation and dispersion coming from wave-induced flow phenomena. To compute wave propagation, we use numerical codes in the frequency domain that allow us to take into account all the frequency-dependent parameters without approximation or local time steps. In this way, the Biot slow compressional wave is well modeled and its partially diffusive, partially propagative behavior is completely considered. The attenuation and dispersion of the waves in such media are coming partly from the wave mode conversion into diffusive slow waves, not visible on seismograms. But the amplitude of propagative P- and S-waves are mainly affected by these energy losses at interfaces. Using full waveform modeling, we investigate the amplitude versus offset (AVO) attributes with respect to saturation and frequency. For a simple three-layer case, we compute poroelastic wave propagation, extract maximum amplitude with respect to angle, and, through a least-square fitting method, we obtain the AVO attributes for PP- and PS-reflected events. Due to the influence of mesoscale induced-flow phenomena and relatively to the regime of the Biot slow wave, the main results show a strong variability of the AVO attributes with the frequency and a lower variability with the saturation for reflected PP or PS events. We show that gas-water and oil-water systems have similar behaviors. Strong differences between patchy saturation and effective fluid phase theories are highlighted, especially at high frequency, for PP events and for gas-water systems. Then, we conclude that these AVO attributes carry information about the saturation that can be used to estimate the saturation variations in time-lapse studies.

CO₂ saturations are estimated at Sleipner using a two‐step imaging workflow. The workflow combines seismic tomography (full‐waveform inversion) and rock physics inversion and is applied to a two‐dimensional seismic line located near the injection point at Sleipner. We use baseline data (1994 vintage, before CO₂ injection) and monitor data that was acquired after 12 years of CO₂ injection (2008 vintage). P‐wave velocity models are generated using the Full waveform inversion technology and then, we invert selected rock physics parameters using an rock physics inversion methodology. Full waveform inversion provides high‐resolution P‐wave velocity models both for baseline and monitor data. The physical relations between rock physics properties and acoustic wave velocities in the Utsira unconsolidated sandstone (reservoir formation) are defined using a dynamic rock physics model based on well‐known Biot–Gassmann theories. For data prior to injection, rock frame properties (porosity, bulk and shear dry moduli) are estimated using rock physics inversion that allows deriving physically consistent properties with related uncertainty. We show that the uncertainty related to limited input data (only P‐wave velocity) is not an issue because the mean values of parameters are correct. These rock frame properties are then used as a priori constraint in the monitor case. For monitor data, the Full waveform inversion results show nicely resolved thin layers of CO₂–brine saturated sandstones under intra‐reservoir shale layers. The CO₂ saturation estimation is carried out by plugging an effective fluid phase in the rock physics model. Calculating the effective fluid bulk modulus of the brine–CO₂ mixture (using Brie equation in our study) is shown to be the key factor to link P‐wave velocity to CO₂ saturation. The inversion tests are done with several values of Brie/patchiness exponent and show that the CO₂ saturation estimates are varying between 0.30 and 0.90 depending on the rock physics model and the location in the reservoir. The uncertainty in CO₂ saturation estimation is usually lower than 0.20. When the patchiness exponent is considered as unknown, the inversion is less constrained and we end up with values of exponent varying between 5 and 20 and up to 33 in specific reservoir areas. These estimations tend to show that the CO₂–brine mixing is between uniform and patchy mixing and variable throughout the reservoir.

The estimation of quantitative rock physics properties is of great importance for reservoir characterization and monitoring in ${\mathrm{CO}}_2$ storage or enhanced oil recovery as an example. We have combined the high-resolution results of full-waveform inversion (FWI) methods with rock physics inversion. Because we consider a generic and dynamic rock physics model, our method is applicable to most kinds of rocks for a wide range of frequencies. The first step allows determination of viscoelastic effective properties, i.e., quantitative seismic attributes, whereas the rock physics inversion estimates rock physics properties (porosity, solid frame moduli, fluid phase properties, or saturation). This two-step workflow is applied to time-lapse synthetic and field cases. The sensitivity tests that we had previously carried out showed that it can be crucial to use multiparameter inputs to accurately recover fluid saturations and fluid properties. However, due to the limited data availability and difficulties in getting reliable multiparameter FWI results, we are limited to acoustic FWI results. The synthetic tests are conclusive even if they are favorable cases. For the first time-lapse fluid substitution synthetic case, we first characterize the rock frame parameters on the baseline model using P-wave velocity estimations obtained by acoustic FWI. Then, we obtain an accurate estimation of fluid bulk modulus from the time-lapse P-wave velocity. In the Marmousi synthetic case, the rock frame properties are accurately recovered for the baseline model, whereas the gas saturation change in the monitor model is not estimated correctly. On the field data example (time-lapse monitoring of an underground blowout in the North Sea), the estimation of rock frame properties gives results on a relatively narrow range, and we use this estimation as a starting model for the gas saturation inversion. We have found that the estimation of the gas saturation is not accurate enough, and the use of attenuation data is then required. However, the uncertainty on the estimation of baseline rock frame properties is not critical to monitor gas saturation changes.