The compressional wave reflection coefficient R(θ) given by the Zoeppritz equations is simplified to the following: [Formula: see text] The first term gives the amplitude at normal incidence (θ = 0), the second term characterizes R(θ) at intermediate angles, and the third term describes the approach to critical angle. The coefficient of the second term is that combination of elastic properties which can be determined by analyzing the offset dependence of event amplitude in conventional multichannel reflection data. If the event amplitude is normalized to its value for normal incidence, then the quantity determined is [Formula: see text] [Formula: see text] specifies the normal, gradual decrease of amplitude with offset; its value is constrained well enough that the main information conveyed is [Formula: see text] where [Formula: see text] is the contrast in Poisson’s ratio at the reflecting interface and [Formula: see text] is the amplitude at normal incidence. This simplified formula for R(θ) accounts for all of the relations between R(θ) and elastic properties first described by Koefoed in 1955.

I derive an approximate formula for the plane P‐wave reflection coefficient as a function of ray‐parameter. The approximation shows that the behavior of the P‐wave reflection coefficient at nonnormal angles of incidence is mainly controlled by two parameters: (1) [Formula: see text], the fluid‐fluid reflection coefficient (i.e., the reflection coefficient when the S‐wave velocities in both media are set to zero) and (2) Δμ/ρ, the ratio of the contrast in shear moduli to the average bulk density. I also show that the other formulas for the P‐wave reflection coefficient given by R. Bortfeld and R. Shuey can be approximately derived from my formula. I give numerical examples to demonstrate the accuracy of the formula. Using a least‐squares inversion of theoretical values for the reflection coefficients, I demonstrate that, in a linear inversion of amplitude versus offset data, Δμ/ρ is better estimated than the contrast in Poisson’s ratio, Δσ. Finally, comparing the exact reflection coefficient with Shuey’s approximation for typical shale and Class 1 gas‐sand models, I show that Shuey’s approximation under‐estimates the value of |Δσ| for such reflections.

The Geostack technique is a method of analyzing seismic amplitude variation with offset (AVO) information. One of the outputs of the analysis is a set of direct hydrocarbon indicator traces called “fluid factor” traces. The fluid factor trace is designed to be low amplitude for all reflectors in a clastic sedimentary sequence except for rocks that lie off the (mudrock line.) The mudrock line is the line on a crossplot of P‐wave velocity against S‐wave velocity on which water‐saturated sandstones, shales, and siltstones lie. Some of the rock types that lie off the mudrock line are gas‐saturated sandstones, carbonates, and igneous rocks. In the absence of carbonates and igneous rocks, high amplitude reflections on fluid factor traces would be expected to represent gas‐saturated sandstones. Of course, this relationship does not apply exactly in nature, and the extent to which the mudrock line model applies varies from area to area. However, it is a useful model in many basins of the world, including the one studied here. Geostack processing has been done on a 3-D seismic data set over the Mossel Bay gas field on the southern continental shelf of South Africa. We found that anomalously high amplitude fluid factor reflections occurred at the top and base of the gas‐reservoir sandstone. Maps were made of the amplitude of these fluid factor reflections, and it was found that the high amplitude values were restricted mainly to the gas field area as determined by drilling. The highest amplitudes were found to be located roughly in the areas of best reservoir quality (i.e., highest porosity) in areas where the reservoir is relatively thick.

To overcome the weaknesses of conventional prestack amplitude variation with angle inversion based on various linear or quasi-linear approximations, we have conducted a nonlinear inversion method using the exact Zoeppritz matrix (EZAI). However, the inversion using the exact Zoeppritz matrix was highly nonlinear and often unstable, if not properly treated. To tackle these issues, we have used an iteratively regularizing Levenberg-Marquardt scheme (IRLM), which regularizes the inversion problem within an algorithm that minimizes the misfit between the observed and the modeled data at the same time by incorporating the Tikhonov regularization method. As a result, the new EZAI method solved using the IRLM scheme is feasible for seismic data sets with large incidence angles, even up to or beyond the critical angle as well as strong parameter contrasts. Single and multilayered synthetic examples were used to test these features. These tests also showed that EZAI is robust on noisy gathers for parameter extraction and has weak dependence on the initial model. For the influence of inaccurate amplitudes, dominant frequencies, and phase angles, we found that EZAI is less sensitive to the variation in amplitude and phase shifts than to the dominant frequencies. Specifically, the inversion results of EZAI for P- and S-wave velocities and density were reliable if the inaccurate range for the amplitude was within 20% or the angle of the phase shift was no more than 20°. The superiority of EZAI makes it a very promising method for the estimation of subsurface elastic parameters.

Recent advances in seismic data acquisition and processing allow routine extraction of offset-/angle-dependent reflection amplitudes from prestack seismic data for quantifying subsurface lithologic and fluid properties. Amplitude-variation-with-offset (AVO) inversion is the most commonly used practice for such quantification. Although quite successful, AVO has a few shortcomings primarily due to the simplicity in its inherent assumptions, and for any quantitative estimation of reservoir properties, they are generally interpreted in combination with other information. In recent years, waveform-based inversions have gained popularity in reservoir characterization and depth imaging. Going beyond the simple assumptions of AVO and using wave equation solutions, these methods have been effective in accurately predicting the subsurface properties. Developments of these waveform inversions have so far been along two lines: (1) the methods that use a locally 1D model of the subsurface for each common midpoint and use an analytical solution to the wave equation for forward modeling and (2) the methods that do not make any 1D assumption but use an approximate numerical solution to the wave equation in 2D or 3D for forward modeling. Routine applications of these inversions are, however, still computationally demanding. We described a multilevel parallelization of elastic-waveform inversion methodology under a 1D assumption that allowed its application in a reasonable time frame. Applying AVO and waveform inversion on a single data set from the Rock Springs Uplift, Wyoming, USA, and comparing them with one another, we also determined that the waveform-based method was capable of obtaining a much superior description of subsurface properties compared with AVO. We concluded that the waveform inversions should be the method of choice for reservoir property estimation as high-performance computers become commonly available.

The inversion of prestack seismic data using amplitude variation with offset (AVO) has received increased attention in the past few decades because of its key role in estimating reservoir properties. AVO is mainly governed by the Zoeppritz equations, but traditional inversion techniques are based on various linear or quasilinear approximations to these nonlinear equations. We have developed an efficient algorithm for nonlinear AVO inversion of precritical reflections using the exact Zoeppritz equations in multichannel and multi-interface form for simultaneous estimation of the P-wave velocity, S-wave velocity, and density. The total variation constraint is used to overcome the ill-posedness while solving the forward nonlinear model and to preserve the sharpness of the interfaces in the parameter space. The optimization is based on a combination of Levenberg’s algorithm and the split Bregman iterative scheme, in which we have to refine the data and model parameters at each iteration. We refine the data via the original nonlinear equations, but we use the traditional cost-effective linearized AVO inversion to construct the Jacobian matrix and update the model. Numerical experiments show that this new iterative procedure is convergent and converges to a solution of the nonlinear problem. We determine the performance and optimality of our nonlinear inversion algorithm with various simulated and field seismic data sets.

This paper uses a Bayesian approach for inverting seismic amplitude versus offset (AVO) data to provide estimates and uncertainties of the viscoelastic physical parameters at an interface. The inversion is based on Gibbs' sampling approach to determine properties of the posterior probability distribution (PPD), such as the posterior mean, maximum a posteriori (MAP) estimate, marginal probability distributions, and covariances. The Bayesian formulation represents a fully nonlinear inversion; the results are compared to those of standard linearized inversion. The nonlinear and linearized approaches are applied to synthetic test cases which consider AVO inversion for shallow marine environments with both unconsolidated and consolidated seabeds. The result of neglecting attenuation in the seabed is investigated, and the effects of data factors such as independent and systematic errors and the range of incident angles are considered. The Bayesian approach is also applied to estimate the physical parameters and uncertainties from AVO data collected at two sites along a seismic line in the Baltic Sea with differing sediment types; it clearly identifies the distinct seabed compositions. Data uncertainties (independent and systematic) required for this analysis are estimated using a maximum‐likelihood approach.

We have developed a Markov chainMonte Carlo (MCMC) method for joint inversion of seismic data for the prediction of facies and elastic properties. The solution of the inverse problem is defined by the Bayesian posterior distribution of the properties of interest. The prior distribution is a Gaussian mixture model, and each component is associated to a potential configuration of the facies sequence along the seismic trace. The low frequency is incorporated by using facies-dependent depositional trend models for the prior means of the elastic properties in each facies. The posterior distribution is also a Gaussian mixture, for which the Gaussian component can be analytically computed. However, due to the high number of components of the mixture, i.e., the large number of facies configurations, the computation of the full posterior distribution is impractical. Our Gaussian mixture MCMC method allows for the calculation of the full posterior distribution by sampling the facies configurations according to the acceptance/rejection probability. The novelty of the method is the use of an MCMC framework with multimodal distributions for the description of the model properties and the facies along the entire seismic trace. Our method is tested on synthetic seismic data, applied to real seismic data, and validated using a well test.

In geophysical inverse problems, the posterior model can be analytically assessed only in case of linear forward operators, Gaussian, Gaussian mixture, or generalized Gaussian prior models, continuous model properties, and Gaussian-distributed noise contaminating the observed data. For this reason, one of the major challenges of seismic inversion is to derive reliable uncertainty appraisals in cases of complex prior models, non-linear forward operators and mixed discrete-continuous model parameters. We present two amplitude versus angle inversion strategies for the joint estimation of elastic properties and litho-fluid facies from pre-stack seismic data in case of non-parametric mixture prior distributions and non-linear forward modellings. The first strategy is a two-dimensional target-oriented inversion that inverts the amplitude versus angle responses of the target reflections by adopting the single-interface full Zoeppritz equations. The second is an interval-oriented approach that inverts the pre-stack seismic responses along a given time interval using a one-dimensional convolutional forward modelling still based on the Zoeppritz equations. In both approaches, the model vector includes the facies sequence and the elastic properties of P-wave velocity, S-wave velocity and density. The distribution of the elastic properties at each common-mid-point location (for the target-oriented approach) or at each time-sample position (for the time-interval approach) is assumed to be multimodal with as many modes as the number of litho-fluid facies considered. In this context, an analytical expression of the posterior model is no more available. For this reason, we adopt a Markov chain Monte Carlo algorithm to numerically evaluate the posterior uncertainties. With the aim of speeding up the convergence of the probabilistic sampling, we adopt a specific recipe that includes multiple chains, a parallel tempering strategy, a delayed rejection updating scheme and hybridizes the standard Metropolis-Hasting algorithm with the more advanced differential evolution Markov chain method. For the lack of available field seismic data, we validate the two implemented algorithms by inverting synthetic seismic data derived on the basis of realistic subsurface models and actual well log data. The two approaches are also benchmarked against two analytical inversion approaches that assume Gaussian-mixture-distributed elastic parameters. The final predictions and the convergence analysis of the two implemented methods proved that our approaches retrieve reliable estimations and accurate uncertainties quantifications with a reasonable computational effort.

A general method, suitable for fast computing machines, for investigating such properties as equations of state for substances consisting of interacting individual molecules is described. The method consists of a modified Monte Carlo integration over configuration space. Results for the two-dimensional rigid-sphere system have been obtained on the Los Alamos MANIAC and are presented here. These results are compared to the free volume equation of state and to a four-term virial coefficient expansion.