Full-waveform Modeling Research Articles

Retrieval of the P wave reflectivity response from autocorrelation of seismic noise: Jakarta Basin, Indonesia

AbstractWe autocorrelate the continuously recorded seismic wavefield across a dense network of seismometers to map the P wave reflectivity response of the Jakarta Basin, Indonesia. The proximity of this mega city to known active faults and the subduction of the Australian plate, especially when the predominance of masonry construction and thick sedimentary basin fill are considered, suggests that it is a hot spot for seismic risk. In order to understand the type of ground motion that earthquakes might cause in the basin, it is essential to obtain reliable information on its seismic velocity structure. The body wave reflections are sensitive to the sharp velocity contrasts, which makes them useful in seismic imaging. Results show autocorrelograms at different seismic stations with reflected‐wave travel time variations, which reflect the variation in basement depth across the thick sedimentary basin. We also confirm the validity of the observed autocorrelation waveforms by conducting a 2‐D full waveform modeling.

Estimating geomechanical parameters from microseismic plane focal mechanisms recorded during multistage hydraulic fracturing

We extend a full-waveform modeling method to invert source focal-plane mechanisms for microseismic data recorded with dual-borehole seismic arrays. Combining inverted focal-plane mechanisms with geomechanics knowledge, we map the pore pressure distribution in the reservoir. Determining focal mechanisms for microseismic events is challenging due to poor geometry coverage. We use the P-wave polarities, the P- and S-wave similarities, the SV/P amplitude ratio, and the SH/P amplitude ratio to invert the focal-plane mechanisms. A synthetic study proves that this method can effectively resolve focal mechanisms with dual-array geometry. We apply this method to 47 relatively large events recorded during a hydraulic fracturing operation in the Barnett Shale. The focal mechanisms are used to invert for the orientation and relative magnitudes of the principal stress axes, the orientation of the planes slipping in shear, and the approximate pore pressure perturbation that caused the slip. The analysis of the focal mechanisms consistently shows a normal faulting stress state with the maximum principal stress near vertical, the maximum horizontal stress near horizontal at an azimuth of N60°E, and the minimum horizontal stress near horizontal at an azimuth of S30°E. We propose a general method that can be used to obtain microseismic focal-plane mechanisms and use them to improve the geomechanical understanding of the stimulation process during multistage hydraulic fracturing.

Full-waveform modeling and inversion of physical model data

Because full elastic waveform inversion requires considerable computation time for forward modeling and inversion, acoustic waveform inversion is often applied to marine data for reducing the computational time. To understand the validity of the acoustic approximation, we study data collected from an ultrasonic laboratory with a known physical model by applying elastic and acoustic waveform modeling and acoustic waveform inversion. This study enables us to evaluate waveform differences quantitatively between synthetics and real data from the same physical model and to understand the effects of different objective functions in addressing the waveform differences for full-waveform inversion. Because the materials used in the physical experiment are viscoelastic, we find that both elastic and acoustic synthetics differ substantially from the physical data over offset in true amplitude. If attenuation is taken into consideration, the amplitude versus offset (AVO) of viscoelastic synthetics more closely approximates the physical data. To mitigate the effect of amplitude differences, we apply trace normalization to both synthetics and physical data in acoustic full-waveform inversion. The objective function is equivalent to minimizing the phase differences with indirect contributions from the amplitudes. We observe that trace normalization helps to stabilize the inversion and obtain more accurate model solutions for both synthetics and physical data.

2D full-waveform modeling of seismic waves in layered karstic media

We have developed a new propagator-matrix scheme to simulate seismic-wave propagation and scattering in a multilayered medium containing karstic voids. The propagator matrices can be found using the boundary element method. The model can have irregular boundaries, including arbitrary free-surface topography. Any number of karsts can be included in the model, and each karst can be of arbitrary geometric shape. We have used the Burton-Miller formulation to tackle the numerical instability caused by the fictitious resonance due to the finite size of a karstic void. Our method was implemented in the frequency-space domain, so frequency-dependent [Formula: see text] can be readily incorporated. We have validated our calculation by comparing it with the analytical solution for a cylindrical void and to the spectral element method for a more complex model. This new modeling capability is useful in many important applications in seismic inverse theory, such as imaging karsts, caves, sinkholes, and clandestine tunnels.

Fast computation of eikonal and transport equations on graphics processing units computer architectures

Eikonal models have been widely used for traveltime computations in the field of seismic imaging, but they are often criticized for having low accuracy and poor resolution of the output image. Including amplitude information can provide higher model resolution and accuracy of the images. We have developed a new approach for computing eikonal traveltimes and amplitudes, and we implemented it for multicore computer processing unit and graphics processing unit computer architectures. Traveltimes and amplitudes are computed simultaneously in iterations of the 3D velocity model. This is achieved by extending a fast sweeping method by computing amplitudes directly after the traveltimes with upwind numerical difference schemes. By performing the extra computations simultaneously with the traveltimes, the additional cost for the amplitude and raypaths is low. We tested our method on synthetic 3D data sets to compute traveltimes, amplitudes, and raypaths, from which the Helmholtz Green’s function was assembled. Using a grid of [Formula: see text] nodes, the computations were performed in less than 1 s. The proposed method could work as a feasible alternative to full waveform modeling in seismic applications, which suffers from demanding computations because it requires several order of magnitudes shorter computing times.

Interferometric salt-flank estimation

Accurately determining the position of salt-sediment interfaces is important in oil and gas exploration, salt mining, and [Formula: see text] or oil and gas storage operations. We used seismic waves from a vertical seismic profile conducted above a salt dome in Texas to image a 2D vertical slice of the salt flank that was longer than can be obtained with traditional depth migration methods alone. This was possible by combining interferometric imaging methods with data from 100 3C receiver positions over a 1524 m (5000 ft) spread in a borehole close to the flank. Interferometry was used to construct virtual source gathers at the receiver positions within the salt body. This velocity-independent redatuming process effectively moved the sources from the surface to the borehole, removed complications of the overlying geology, and provided favorable illumination of the suspected flank position. Full-waveform elastic forward modeling helped us to understand and confirm the flank position. From the processed virtual source data, we created a 2D image in the source-receiver plane. In effect, we have created an imaging problem geometrically equivalent to that of a simple 2D surface seismic line rotated by 90°, with the requirement that we knew the velocity model only within the salt body. Interferometric imaging of the salt flank reduced the picking ambiguity and provided a direct depth image that was unaffected by overburden effects where the salt-flank location could be interpreted more confidently. Ultimately, we were able to image the salt flank over 1371 m (4500 ft) of depth with a lateral resolution in the order of 45–60 m (150–200 ft). The original data were obtained with more traditional imaging methods and cost considerations in mind, but the emergence of interferometric imaging methods suggested significant improvements to the acquisition strategy that would optimize the illumination of the subsurface.

Imaging high-pressure rock exhumation in eastern Taiwan

Increasingly detailed studies of out crops of high-pressure rock terranes in combination with rapidly evolving numerical modeling studies have given rise to a number of possible explanations for the processes by which these rocks are exhumed. Imaging actively exhuming high-pressure terranes remains one of the fundamental, but elusive, tasks that could advance the understanding of how these important rocks reach Earth’s surface. Seismic tomography along the active arc-continent collision in eastern Taiwan images a high P- and S-wave velocity zone that extends from the shallow subsurface beneath a high-pressure metamorphic terrane to ∼50 km depth. We present a petrophysical analysis of this high-velocity zone that indicates the presence of rock types common to high-pressure terranes. The high-velocity zone is seismically active throughout. We determine focal mechanisms for 57 earthquakes, and carry out full waveform modeling on 10; these have double-couple focal mechanisms with a compensated linear vector dipole component up to 20.6%. We suggest that the high-velocity zone comprises an exhuming high-pressure terrane. Focal mechanisms for earthquakes within it indicate that shear faulting dominates in the deformation, but high fluid pressure may also play a role.

Full Waveform Modelling for Converted Waves Seismic Reflections in Mountainous and Marine Environment

The application of seismic waves allows us to achieve adequate results by compressional wave (P-wave) surveys alone. However, in the presence of gas P-wave transmission disrupts and obscures underlying targets. Many reservoirs don’t present sufficient impedance contrast to the overburden and not reflect P-wave strongly to produce an impedance image. High impedance rock such as basalt or hard volcanic rocks are difficult to image with P-wave. To overcome these challenges shear-wave (S-wave) or converted wave (P-S) surveys are usedfor last 20 years by making the use of down going P waves converting to upcoming S waves at the mode conversion boundaries.The processing of converted waves requires studying asymmetric reflection at the conversion point, difference in geometries and conditions of source and receiver, and the partitioning of energy into orthogonally polarized components. Interpretation of P-S sections incorporates the identification of P-S waves, full waveform modelling, correlation with P-wave sections and depth migration. The objectives of this study is to model P-S wave reflections in onshore and offshore environment and to examine the major differences in processing of P and P-S wave surveys together with the identifying converted mode reflections by P-wave sources in anisotropic media. To achieve these, realistic mountainous and marine environment models have been developed and synthetic seismograms are generated by full waveform modelling technique. First a mountain foothill model was studied. A Kirchhoff-based technique that includes anisotropic velocities is used for depth migration of P-S waves. The results from depth imaging show that P-S section help in distinguishing amplitude associated with hydrocarbons from those caused by localized stratigraphic changes. Marine model shows a good correlation with identified converted waves. In addition, the full waveform elastic modellingproves useful in finding an appropriate balance between capturing high-quality P-wave data as well as P-S data challenges in a survey. Key Words: Converted-waves (P-S); P-S Wave; Kirchhoff migration; Depth migration; Gas clouds; Shale diapers

Genetic waveform modeling for the crustal structure in Northeast Japan

We propose a genetic algorithm (GA) search procedure for waveform modeling of local crustal earthquakes for optimal one-dimensional (1-D) crustal velocity model. Both waveforms and travel-time data are used for the structure determination. The use of travel times in model evaluation improves the waveform modeling performance in the sense of computation speed and accuracy. We applied this method to broadband waveforms of a local crustal earthquake (M 4.2) in Northeast Japan. P-wave velocities of the crustal model are found to be 4.95±0.30, 5.9±0.02, and 6.51±0.20km/s for a surface layer, upper crust and lower crust, respectively. The surface layer thickness and the Conrad and Moho depths are found to be 3.01±0.8, 17.77±0.4 and 34.59±1.0km, respectively. For epicentral distances <200km, our synthetic waveforms match the observed ones generally well. Early arrivals are mainly observed at stations near the Pacific coast in the forearc area having a thinner crust. In contrast, delayed arrivals appear at stations near the volcanic front and back-arc areas where low-velocity anomalies exist due to the effect of the Pacific slab dehydration and the hot upwelling flows in the mantle wedge. In general, our results agree well with the main tectonic setting of the study area, which confirms the reliability of the proposed approach. Despite a 1-D velocity model is too simple to represent the complex crustal structure, it is still required for the conventional routine analysis of seismology, such as earthquake location and source parameter studies. The current approach is considered as a step toward the genetic full waveform modeling for the 3-D velocity model estimation.

Order of magnitude increase in subducted H2O due to hydrated normal faults within the Wadati-Benioff zone

It is widely proposed that the oceanic mantle is hydrated by outer rise normal faults, and carries large amounts of water to the deep mantle. However, the extent of oceanic mantle hydration is poorly constrained by existing observations, and is a major source of uncertainty in determining the total water delivered to the mantle. Full waveform modeling of dispersed P-wave arrivals from events deep within the Wadati-Benioff zone of northern Japan shows that hydrated fault zone structures are present at intermediate depths. Analysis of the P-wave coda associated with events 5-35 km below the top of the slab gives an overall indication of the bulk hydration of the subducting oceanic mantle, and can be explained by a 40-km-thick layer that is 17%-31% serpentinized. This suggests that the top of the oceanic mantle is 2.0- 3.5 wt% hydrated, subducting 170-318 Tg/m.y. of water per meter of arc beneath northern Japan. This order-of-magnitude increase in the estimated H 2 O fl ux in this arc implies that over the age of the Earth, the equivalent of as many as 3.5 present-day oceans of water could be subducted along the Kuril and Izu-Bonin arcs alone. These results offer the fi rst direct measure of the lower lithosphere hydration at intermediate depths, and suggest that regassing of the mantle is more vigorous than has previously been proposed.

The 29 September 1969, Ceres, South Africa, Earthquake: Full Waveform Moment Tensor Inversion for Point Source and Kinematic Source Parameters

The Ceres earthquake of 29 September 1969 is the largest known earthquake in southern Africa. Digitized analog recordings from Worldwide Standardized Seismographic Network stations (Powell and Fries, 1964) are used to retrieve the point source moment tensor and the most likely centroid depth of the event using full waveform modeling. A scalar seismic moment of 2.2–2.4×1018 N·m corresponding to a moment magnitude of 6.2–6.3 is found. The analysis confirms the pure strike‐slip mechanism previously determined from onset polarities by Green and Bloch (1971). Overall good agreement with the fault orientation previously estimated from local aftershock recordings is found. The centroid depth can be constrained to be less than 15 km. In a second analysis step, we use a higher order moment tensor based inversion scheme for simple extended rupture models to constrain the lateral fault dimensions. We find rupture propagated unilaterally for 4.7 s from east‐southwest to west‐northwest for about 17 km (average rupture velocity of about 3.1 km/s). Online Material: Waveform and spectral fits.

Examining the processing differences between P and P-S waves in a Rocky Mountain Foothills model

For the past two decades, most of the shear-wave (S-wave) or converted wave (P-S) acquisitions were performed with P-wave source by making the use of downgoing P-waves converting to upgoing S-waves at the mode conversion boundaries. The processing of converted waves requires studying asymmetric reflection at the conversion point, difference in geometries and conditions of source and receiver, and the partitioning of energy into orthogonally polarised components. Interpretation of P-S sections incorporates the identification of P-S waves, full waveform modeling, correlation with P-wave sections and depth migration. The main applications of P-S wave imaging are to obtain a measure of subsurface S-wave properties relating to rock type and fluid saturation (in addition to the P-wave values), imaging through gas clouds and shale diapers, and imaging interfaces with low P-wave contrast but significant S-wave changes. This study examines the major differences in processing of P and P-S wave surveys and the feasibility of identifying converted mode reflections by P-wave sources in anisotropic media. Two-dimensional synthetic seismograms for a realistic rocky mountain foothills model were studied. A Kirchhoff-based technique that includes anisotropic velocities is used for depth migration of converted waves. The results from depth imaging show that P-S section help in distinguishing amplitude associated with hydrocarbons from those caused by localised stratigraphic changes. In addition, the full waveform elastic modeling is useful in finding an appropriate balance between capturing high-quality P-wave data and P-S data challenges in a survey.

Influence of frequency and saturation on AVO attributes for patchy saturated rocks

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.

Examining the processing differences between P and P-S waves in a Rocky Mountain Foothills model

For the past two decades, most of the shear-wave (S-wave) or converted wave (P-S) acquisitions were performed with P-wave source by making the use of downgoing P-waves converting to upgoing S-waves at the mode conversion boundaries. The processing of converted waves requires studying asymmetric reflection at the conversion point, difference in geometries and conditions of source and receiver, and the partitioning of energy into orthogonally polarised components. Interpretation of P-S sections incorporates the identification of P-S waves, full waveform modeling, correlation with P-wave sections and depth migration. The main applications of P-S wave imaging are to obtain a measure of subsurface S-wave properties relating to rock type and fluid saturation (in addition to the P-wave values), imaging through gas clouds and shale diapers, and imaging interfaces with low P-wave contrast but significant S-wave changes. This study examines the major differences in processing of P and P-S wave surveys and the feasibility of identifying converted mode reflections by P-wave sources in anisotropic media. Two-dimensional synthetic seismograms for a realistic rocky mountain foothills model were studied. A Kirchhoff-based technique that includes anisotropic velocities is used for depth migration of converted waves. The results from depth imaging show that P-S section help in distinguishing amplitude associated with hydrocarbons from those caused by localised stratigraphic changes. In addition, the full waveform elastic modeling is useful in finding an appropriate balance between capturing high-quality P-wave data and P-S data challenges in a survey.

Assessing the seismic AVA detectability of CO 2 storage sites using novel time-lapse attributes

Monitoring of stored carbon dioxide (CO 2 ) in subsurface reservoirs is fundamental to operation and management of the storage site, and is a requirement of some national and international legislation. As a consequence, effectiveness of monitorability (the ability to observe the evolving location of subsurface CO 2 ) for any given level of investment in monitoring technology is a significant investment uncertainty that must be assessed along other components of the storage-site selection criteria (e.g. capacity, injectivity and storage economies). We develop a workflow to assess the time-lapse seismic detectability of changes in subsurface aquifer reservoirs by analysing expected changes in seismic amplitude variation with angle (AVA) in the field. Laboratory measurements are used to calculate the seismic response of the reservoir at different saturations and pressures. We include the scattering effect of material above and below the reservoir by using a finite-difference, full-waveform modelling approach AVA analysis then assimilates local site effects into the detectability assessment. We show that performing waveform modelling which includes local geological heterogeneities above and below the reservoir interval is essential to assess the storage site monitroability. In order to quantify expected time-lapse changes in the seismic response, we introduce a new set of robust time-lapse attributes based on time–frequency decomposition. The attributes effectively separate amplitude and phase changes (time-shifts) of time-lapse seismic records, and allow us to quantify their repeatability against the background noise. Furthermore, the frequency-dependent nature of the attributes provides a quantification of the frequency–domain effects of time-lapse changes. The approach is employed to assess the detectability of supercritical CO 2 in two analogue storage sites in the near-shore UK North Sea. Analysis of laboratory measurements and AVA responses indicate the contrasting monitorability of the two sites, which helps decision making about further site investigation and development. Application of the approach to hydrocarbon reservoir monitoring is straightforward.

Sensitivity of compaction-induced multicomponent seismic time shifts to variations in reservoir properties

Pore-pressure variations inside producing reservoirs result in excess stress and strain that influence the arrival times of reflected waves. Inversion of seismic data for pressure changes requires better understanding of the dependence of compaction-induced time shifts on reservoir pressure reduction. Using geomechanical and full-waveform seismic modeling, we investigate pressure-dependent behavior of P-, S-, and PS-wave time shifts from reflectors located above and below a rectangular reservoir embedded in a homogeneous half-space. Our geomechanical modeling algorithm generates the excess stress/strain field and the stress-induced stiffness tensor as linear functions of reservoir pressure. Analysis of time shifts obtained from full-waveform synthetic data shows that they vary almost linearly with pressure for reflectors above the reservoir, but become nonlinear for reflections from the reservoir or deeper interfaces. Time-shift misfit curves computed with respect to noise-contaminated data from a reference reservoir for a wide range of pressure reductions display well-defined global minima corresponding to the actual pressure. In addition, we evaluate the influence of the reservoir width on time shifts and the possibility of constraining the width using time-lapse data. We also discuss the impact of moderate perturbations in the strain-sensitivity coefficients (i.e., third-order stiffnesses) on time shifts and on the accuracy of pressure inversion. Our feasibility analysis indicates that the most stable pressure estimation from noisy data is provided by multicomponent time shifts from reflectors below the reservoir. For multicompartment reservoirs, time shifts can be accurately modeled by linear superposition of the excess stress/strains computed for the individual compartments.

Multi-offset ground-penetrating radar imaging of a lab-scale infiltration test

Abstract. A lab scale infiltration experiment was conducted in a sand tank to evaluate the use of time-lapse multi-offset ground-penetrating radar (GPR) data for monitoring dynamic hydrologic events in the vadose zone. Sets of 21 GPR traces at offsets between 0.44–0.9 m were recorded every 30 s during a 3 h infiltration experiment to produce a data cube that can be viewed as multi-offset gathers at unique times or common offset images, tracking changes in arrivals through time. Specifically, we investigated whether this data can be used to estimate changes in average soil water content during wetting and drying and to track the migration of the wetting front during an infiltration event. For the first problem we found that normal-moveout (NMO) analysis of the GPR reflection from the bottom of the sand layer provided water content estimates ranging between 0.10–0.30 volumetric water content, which underestimated the value determined by depth averaging a vertical array of six moisture probes by 0.03–0.05 volumetric water content. Relative errors in the estimated depth to the bottom of the 0.6 m thick sand layer were typically on the order of 2%, though increased as high as 25% as the wetting front approached the bottom of the tank. NMO analysis of the wetting front reflection during the infiltration event generally underestimated the depth of the front with discrepancies between GPR and moisture probe estimates approaching 0.15 m. The analysis also resulted in underestimates of water content in the wetted zone on the order of 0.06 volumetric water content and a wetting front velocity equal to about half the rate inferred from the probe measurements. In a parallel modeling effort we found that HYDRUS-1D also underestimates the observed average tank water content determined from the probes by approximately 0.01–0.03 volumetric water content, despite the fact that the model was calibrated to the probe data. This error suggests that the assumed conceptual model of laterally uniform, one-dimensional vertical flow in a homogenous material may not be fully appropriate for the experiment. Full-waveform modeling and subsequent NMO analysis of the simulated GPR response resulted in water content errors on the order of 0.01–0.03 volumetric water content, which are roughly 30–50% of the discrepancy between GPR and probe results observed in the experiment. The model shows that interference between wave arrivals affects data interpretation and the estimation of traveltimes. This is an important source of error in the NMO analysis, but it does not fully account for the discrepancies between GPR and the moisture probes observed in the experiment. The remaining discrepancy may be related to conceptual errors underlying the GPR analysis, such as the assumption of uniform one-dimensional flow, a lack of a sharply defined wetting front in the experiment, and errors in the petrophysical model used to convert dielectric constant to water content.

A Stochastic Analysis of Cross‐Hole Ground‐Penetrating Radar Zero‐Offset Profiles for Subsurface Characterization

Cross‐hole ground‐penetrating radar (GPR) zero‐offset profiling (ZOP) is a well‐established technique for the measurement of the one‐dimensional (1‐D) vertical distribution of soil dielectric permittivity and has often been used in time‐lapse mode for the monitoring of natural or man‐made infiltration and soil moisture redistribution processes in the vadose zone. However, in spite of its widespread use, the quantitative interpretation of ZOPs in terms of dielectric permittivity profiles is known to be fraught with difficulties. Often a simplified approach is adopted that translates directly ZOP travel times into electromagnetic (EM) velocities and these, in turn, into dielectric permittivities and soil volumetric moisture content. This approach is known to lead to over‐smoothed moisture content profiles, which are a consequence of the ZOP measurement scale, controlled by critical refractions along fast layers and averaging effects within the first Fresnel zone. Such smooth profiles are not necessarily compatible with the true soil moisture content that should be reproduced by water flow models. In this paper we present a stochastic inversion approach that aims at solving these issues. The approach is based on a forward stochastic simulation that generates the expected travel times by means of a two‐dimensional (2‐D) full waveform modeling, thus reproducing all physical processes that contribute to profile smoothing. We first assess the robustness of this approach on synthetic data, where the true dielectric permittivity profile is known. Then we apply the technique to two different case studies, where the results of the proposed technique show supporting evidence from independent information on the sites' stratigraphy. The proposed approach proves to be capable of reconstructing sharp dielectric profiles, in addition to assigning relevant uncertainty bounds derived from the expected errors in travel time picking. These results may prove very useful in a correct assessment of the hydrological implications of the measured dielectric permittivity, and thus moisture content, profile.

Modeling and analysis of compaction-induced traveltime shifts for multicomponent seismic data

Modeling of time shifts associated with time-lapse (4D) seismic surveys is helpful in evaluating reservoir depressurization and inverting for subsurface stress. Using coupled geomechanical and full-waveform seismic modeling, we study the influence of compaction-induced stress and strain around a simplified reservoir on compressional (P), shear (S), and mode-converted (PS) waves. We estimate compaction-induced time shifts and analyze their dependence on reflector depth and pressure drop inside the reservoir. Time shifts between synthetic baseline and monitor surveys are obtained by processing techniques that are potentially applicable to field data. Although P-wave time-shift lags for reflectors in the overburden are indicative of induced anisotropy, they are two to three times smaller than S-wave time-shift leads for reflectors beneath the reservoir. We also investigate the contributions of the deviatoric and volumetric stains to the time shifts for all three modes. Time shifts for S- and PS-waves are strongly influenced by elevated volumetric and deviatoric strains inside the reservoir. Almost constant S-wave time shifts for a range of offsets and source locations indicate that the contribution of stress-induced velocity anisotropy to shear-wave signatures is weak because the symmetry is close to elliptical. Our modeling also shows that mild tilt of a rectangular reservoir, or its replacement with an elliptically shaped reservoir of the same aspect ratio, has little influence on time shifts. Potentially, the developed methodology can be applied to estimate compaction-induced stress fields using simple compartmentalized reservoir models.

Efficient and accurate modeling of ocean bottom seismometer data using reciprocity

Seismic experiments in which the number of sources is considerably larger than the number of receivers occur regularly. An important example is the collection of crustal scale seismic data using ocean bottom seismometers and marine sources. We describe a method to accurately and efficiently compute synthetic seismograms for such experiments by using finite differences and reciprocity. We show numerically how to decompose an explosive source into its equivalent body force components using the staggered-grid finite-difference technique with a fourth-order approximation for the spatial derivative and a second-order approximation for the temporal derivative. This decomposition results in a source configuration where the equivalent body forces are defined in 12 points around the point where the ex-plosive source is applied. We then use the derived equivalent body forces for the explosive source and seismic reciprocity theorems to convert the common receiver gather to a common shot gather. The method is tested on a structurally complex elastic model of the crust and the results show that it is accurate within floating point precision. The synthetic data are compared to data from a real ocean bottom seismometer experiment conducted across a continent-ocean transition zone. A good fit in terms of traveltime is observed for many of the prominent seismic phases. The amplitude fit of these arrivals is not always as good as the traveltime fit. This indicates that full-waveform modeling of such data can provide useful information about the subsurface that cannot be obtained from traveltime modeling. If enough data are available, the modeling method can be used in full-waveform inversion.