Modeling and inversion of seismic data using multiple scattering, renormalization and homotopy methods

Seismic scattering theory plays an important role in seismic forward modeling and is the theoretical foundation for various seismic imaging methods. Full waveform inversion is a powerful technique for obtaining a high-resolution model of the subsurface. One objective of this thesis is to develop convergent scattering series solutions of the Lippmann-Schwinger equation in strongly scattering media using renormalization and homotopy methods. Other objectives of this thesis are to develop efficient full waveform inversion methods of time-lapse seismic data and, to investigate uncertainty quantification in full waveform inversion for anisotropic elastic media based on integral equation approaches and the iterated extended Kalman filter. The conventional Born scattering series is obtained by expanding the Lippmann-Schwinger equation in terms of an iterative solution based on perturbation theory. Such an expansion assumes weak scattering and may have the problems of convergence in strongly scattering media. This thesis presents two scattering series, referred to as convergent Born series (CBS) and homotopy analysis method (HAM) scattering series for frequency-domain seismic wave modeling. For the convergent Born series, a physical interpretation from the renormalization prospective is given. The homotopy scattering series is derived by using homotopy analysis method, which is based on a convergence control parameter $h$ and a convergence control operator $H$ that one can use to ensure convergence for strongly scattering media. The homotopy scattering scattering series solutions of the Lippmann-Schwinger equation, which is convergent in strongly scattering media. The homotopy scattering series is a kind of unified scattering series theory that includes the conventional and convergent Born series as special cases. The Fast Fourier Transform (FFT) is employed for efficient implementation of matrix-vector multiplication for the convergent Born series and the homotopy scattering series. This thesis presents homotopy methods for ray based seismic modeling in strongly anisotropic media. To overcome several limitations of small perturbations and weak anisotropy in obtaining the traveltime approximations in anisotropic media by expanding the anisotropic eikonal equation in terms of the anisotropic parameters and the elliptically anisotropic eikonal equation based on perturbation theory, this study applies the homotopy analysis method to the eikonal equation. Then this thesis presents a retrieved zero-order deformation equation that creates a map from the anisotropic eikonal equation to a linearized partial differential equation system. The new traveltime approximations are derived by using the linear and nonlinear operators in the retrieved zero-order deformation equation. Flexibility on variable anisotropy parameters is naturally incorporated into the linear differential equations, allowing a medium of arbitrarily anisotropy. This thesis investigates efficient target-oriented inversion strategies for improving full waveform inversion of time-lapse seismic data based on extending the distorted Born iterative T-matrix inverse scattering to a local inversion of a small region of interest (e. g. reservoir under production). The target-oriented approach is more efficient for inverting the monitor data. The target-oriented inversion strategy requires properly specifying the wavefield extrapolation operators in the integral equation formulation. By employing the T-matrix and the Gaussian beam based Green’s function, the wavefield extrapolation for the time-lapse inversion is performed in the baseline model from the survey surface to the target region. I demonstrate the method by presenting numerical examples illustrating the sequential and double difference strategies. To quantify the uncertainty and multiparameter trade-off in the full waveform inversion for anisotropic elastic media, this study applies the iterated extended Kalman filter to anisotropic elastic full waveform inversion based on the integral equation method. The sensitivity matrix is an explicit representation with Green’s functions based on the nonlinear inverse scattering theory. Taking the similarity of sequential strategy between the multi-scale frequency domain full waveform inversion and data assimilation with an iterated extended Kalman filter, this study applies the explicit representation of sensitivity matrix to the the framework of Bayesian inference and then estimate the uncertainties in the full waveform inversion. This thesis gives results of numerical tests with examples for anisotropic elastic media. They show that the proposed Bayesian inversion method can provide reasonable reconstruction results for the elastic coefficients of the stiffness tensor and the framework is suitable for accessing the uncertainties and analysis of parameter trade-offs

Computational Inverse Problems for Partial Differential Equations

The problem of determining unknown quantities in a PDE from measurements of (part of) the solution to this PDE arises in a wide range of applications in science, technology, medicine, and finance. The unknown quantity may e.g. be a coefficient, an initial or a boundary condition, a source term, or the shape of a boundary. The identification of such quantities is often computationally challenging and requires profound knowledge of the analytical properties of the underlying PDE as well as numerical techniques. The focus of this workshop was on applications in phase retrieval, imaging with waves in random media, and seismology of the Earth and the Sun, a further emphasis was put on stochastic aspects in the context of uncertainty quantification and parameter identification in stochastic differential equations. Many open problems and mathematical challenges in application fields were addressed, and intensive discussions provided an insight into the high potential of joining deep knowledge in numerical analysis, partial differential equations, and regularization, but also in mathematical statistics, homogenization, optimization, differential geometry, numerical linear algebra, and variational analysis to tackle these challenges

A computational framework for infinite-dimensional Bayesian inverse problems: Part II. Stochastic Newton MCMC with application to ice sheet flow inverse problems

We address the numerical solution of infinite-dimensional inverse problems in the framework of Bayesian inference. In the Part I companion to this paper (arXiv.org:1308.1313), we considered the linearized infinite-dimensional inverse problem. Here in Part II, we relax the linearization assumption and consider the fully nonlinear infinite-dimensional inverse problem using a Markov chain Monte Carlo (MCMC) sampling method. To address the challenges of sampling high-dimensional pdfs arising from Bayesian inverse problems governed by PDEs, we build on the stochastic Newton MCMC method. This method exploits problem structure by taking as a proposal density a local Gaussian approximation of the posterior pdf, whose construction is made tractable by invoking a low-rank approximation of its data misfit component of the Hessian. Here we introduce an approximation of the stochastic Newton proposal in which we compute the low-rank-based Hessian at just the MAP point, and then reuse this Hessian at each MCMC step. We compare the performance of the proposed method to the original stochastic Newton MCMC method and to an independence sampler. The comparison of the three methods is conducted on a synthetic ice sheet inverse problem. For this problem, the stochastic Newton MCMC method with a MAP-based Hessian converges at least as rapidly as the original stochastic Newton MCMC method, but is far cheaper since it avoids recomputing the Hessian at each step. On the other hand, it is more expensive per sample than the independence sampler; however, its convergence is significantly more rapid, and thus overall it is much cheaper. Finally, we present extensive analysis and interpretation of the posterior distribution, and classify directions in parameter space based on the extent to which they are informed by the prior or the observations.Comment: 31 page

Probabilistic estimation of reservoir properties by means of wide-angle AVA inversion and a petrophysical reformulation of the Zoeppritz equations

We apply a target-oriented amplitude versus angle (AVA) inversion to estimate the petrophysical properties of a gas-saturated reservoir in offshore Nile Delta. A linear empirical rock-physics model derived from well log data provides the link between the petrophysical properties (porosity, shaliness and saturation) and the P-wave, S-wave velocities and density. This rock-physics model, properly calibrated for the investigated reservoir, is used to re-parameterize the exact Zoeppritz equations. The so derived equations are the forward model engine of a linearized Bayesian AVA-petrophysical inversion that, for each data gather, inverts the AVA of the target reflections to estimate the petrophysical properties of the reservoir layer, keeping fixed the cap-rock properties. We make use of the iterative Gauss-Newton method to solve the inversion problem. For each petrophysical property of interest, we discuss the benefits introduced by wide-angle reflections in constraining the inversion and we compare the posterior probability distributions (PPDs) analytically obtained via a local linearization of the inversion with the PPDs numerically computed with a Markov Chain Monte Carlo (MCMC) method. It results that the porosity is the best resolved parameter and that wide-angle reflections effectively constrain the shaliness estimates but do not guarantee reliable saturation estimates. It also results that the local linearization returns accurate PPDs in good agreement with the MCMC estimates

HMCLab: a framework for solving diverse geophysical inverse problems using the Hamiltonian Monte Carlo method

The use of the probabilistic approach to solve inverse problems is becoming more popular in the geophysical community, thanks to its ability to address nonlinear forward problems and to provide uncertainty quantification. However, such strategy is often tailored to specific applications and therefore there is a lack of a common platform for solving a range of different geophysical inverse problems and showing potential and pitfalls. We demonstrate a common framework to solve such inverse problems ranging from, e.g, earthquake source location to potential field data inversion and seismic tomography. Within this approach, we can provide probabilities related to certain properties or structures of the subsurface. Thanks to its ability to address high-dimensional problems, the Hamiltonian Monte Carlo (HMC) algorithm has emerged as the state-of-the-art tool for solving geophysical inverse problems within the probabilistic framework. HMC requires the computation of gradients, which can be obtained by adjoint methods, making the solution of tomographic problems ultimately feasible. These results can be obtained with "HMCLab", a tool for solving a range of different geophysical inverse problems using sampling methods, focusing in particular on the HMC algorithm. HMCLab consists of a set of samplers and a set of geophysical forward problems. For each problem its misfit function and gradient computation are provided and, in addition, a set of prior models can be combined to inject additional information into the inverse problem. This allows users to experiment with probabilistic inverse problems and also address real-world studies. We show how to solve a selected set of problems within this framework using variants of the HMC algorithm and analyze the results. HMCLab is provided as an open source package written both in Python and Julia, welcoming contributions from the community.Comment: 21 pages, 4 figure

Quantitative seismic interpretation in thin-bedded geology using full-wavefield elastic modelling

Refleksjonsseismikk brukes til å lage seismiske «bilder» av den øverste delen av jordskorpen, blant annet med tanke på leting etter reservoarer for olje, gass, karbonlagring og geotermisk energi. I tillegg til å gi grunnlag for en strukturell tolkning, kan de seismiske dataene brukes til å kvantifisere egenskapene til det faste materialet og væskeinnholdet i bergartene. Et viktig verktøy i slik kvantitativ seismisk tolkning er analyse av såkalt AVO: amplitudenes variasjon med avstanden mellom kilde og mottaker (offset). Tynne geologiske lag gir utfordringer for AVO-modellering og tolkning, fordi lagtykkelsen vil kunne være mindre enn oppløsningen i de seismiske dataene. En problemstilling som tas opp i denne avhandlingen er nettopp hvordan man kan gjøre nøyaktig seismisk (forover) modellering i medier med tynne lag. En konvensjonell tilnærming innen AVO- modellering og inversjon er å bruke såkalt konvolusjonsmodellering. Denne metoden tar imidlertid bare hensyn til de primære seismiske refleksjonene og er derfor unøyaktig når modellene har tynne lag. To bedre alternativer er endelig-differanse-modellering og reflektivitetsmetoden. Reflektivitetsmetoden er en delvis analytisk modelleringsmetode for horisontalt lagdelte medier og er beregningsmessig billigere enn endelig-differansemodellering, der beregningene er basert på et tett samplet rutenett (grid). Jeg viser i avhandlingen at reflektivitetsmetoden er godt egnet for AVO-modellering i lagdelte medier. Seismiske data har en båndbegrenset karakter. En konsekvens er at beregning av reservoaregenskaper fra seismiske data generelt ikke er entydig, noe som særlig kommer til uttrykk for lagdelt geologi med tynne lag. Probabilistiske inversjonsmetoder, som for eksempel bayesianske metoder, tar hensyn til denne flertydigheten ved å forutsi sannsynligheter, noe som gjør det mulig a kvantisere usikkerheten. I avhandlingen kombinerer jeg seismisk modellering med bayesiansk klassifisering og inversjon. Modelleringen er utført med reflektivitetsmetoden og er basert på det komplette elastiske bølgefeltet. Formålet er å adressere to konkrete kvantitative seismiske tolkningsproblemer: 1) kvantifisering av usikkerhet i bayesiansk porevæske-klassifisering i nærvær av tynne lag med høy impedans, forårsaket av kalsittsementering i sandstein, og 2) estimering av reservoaregenskapene til turbiditt-reservoarer karakterisert ved alternerende lag av sandstein og skifer. I den første anvendelsen viser jeg i en modelleringsstudie at kalsitt-sementerte lag kan gi en detekterbar refleksjonsrespons, noe som kan påvirke amplituden målt ved reservoartoppen og dermed forstyrre AVO-målingen. Den observerte effekten øker usikkerheten ved porevæske-klassifisering basert på AVO-attributter, som jeg har demonstrert i en case-studie. Følgelig øker sannsynligheten for en falsk hydrokarbon-indikasjon betydelig i nærvær av kalsittsementerte lag. I den andre anvendelsen presenterer jeg en bayesiansk inversjon som tar AVO-skjæringspunktet og gradienten målt på toppen av et reservoar som inngangsdata og estimerer sannsynlighetstetthetsfunksjonen til forholdstallene «net-to-gross» og «net-pay-to-net». Metoden ble anvendt på syntetiske data og AVO-attributtkart fra Jotunfeltet på norsk kontinentalsokkel. Det ble funnet at AVO-gradienten korrelerer med reservoarets net-togross forhold, mens AVO-skjæringspunktet er mest følsomt for typen porevæske. Etter inversjon genererte jeg kart over de mest sannsynlige verdiene av forholdene net-to-gross og net-pay-to-net, samt kart over net pay og usikkerhetene. Disse kartene kan bidra til å identifisere potensielle soner med høy reservoarkvalitet og hydrokarbonmetning.Reflection seismics is used to image the subsurface for the exploration of oil and gas, geothermal or carbon storage reservoirs, among others. In addition to the structural interpretation of the resulting seismic images, the seismic data can be interpreted quantitatively with the goal to obtain rock and fluid properties. An essential tool in quantitative seismic interpretation is the analysis of the amplitude variation with offset (AVO). Thin-bedded geology below the seismic resolution poses challenges for AVO modelling and interpretation. One problem addressed in this thesis is accurate seismic forward modelling in thin-bedded media. Primaries-only convolutional modelling, commonly used in conventional AVO modelling and inversion, is prone to failure in the presence of thin beds. Better alternatives are finite-difference modelling or the reflectivity method. The reflectivity method is a semi-analytic modelling method for horizontally layered media and is computationally cheaper than finite-difference modelling on densely sampled grids. I show in this thesis that the reflectivity method is well-suited for the AVO modelling of layered media. The band-limited nature of seismic data is one reason for the non-unique estimation of reservoir properties from seismic data, especially in thin-bedded geology. Probabilistic inversion methods, such as Bayesian methods, honour this non-uniqueness by predicting probabilities that allow the uncertainty to be quantified. In this thesis, I integrate full-wavefield elastic seismic modelling by the reflectivity method with Bayesian classification and inversion. The objective is to address two concrete quantitative seismic interpretation problems: 1) the uncertainty quantification of Bayesian pore-fluid classification in the presence of thin high-impedance layers caused by calcite cementation in sandstone, and 2) the estimation of reservoir properties of turbidite reservoirs characterised by sand-shale interbedding. In the first application, I show through a modelling study that calcite-cemented beds lead to detectable reflection responses that can interfere with the target reflection at the reservoir top and thereby perturb the AVO behaviour. The observed effect increases the uncertainty of pore-fluid classification based on AVO attributes, as demonstrated by a case study. Consequently, the probability of a false hydrocarbon indication is significantly increased in the presence of calcite-cemented beds. In the second application, I present a Bayesian inversion that takes the AVO intercept and gradient measured at the top of a reservoir as input and estimates the probability density function of the net-to-gross ratio and the net-pay-to-net ratio. The method was applied to synthetic data and AVO attribute maps from the Jotun field on the Norwegian Continental Shelf. It was found that the AVO gradient correlates with the net-to-gross ratio of the reservoir, while the AVO intercept is most sensitive to the type of pore fluid. After inversion, maps of the most-likely values of the net-to-gross ratio, net-pay-to-net ratio, net pay and the uncertainty could be generated. These maps help to identify potential zones of high reservoir quality and hydrocarbon saturation.Doktorgradsavhandlin

Seismic Inversion and Uncertainty Analysis using Transdimensional Markov Chain Monte Carlo Method

We use a transdimensional inversion algorithm, reversible jump MCMC (rjMCMC), in the seismic waveform inversion of post-stack and prestack data to characterize reservoir properties such as seismic wave velocity, density as well as impedance and then estimate uncertainty. Each seismic trace is inverted independently based on a layered earth model. The model dimensionality is defined as the number of the layers multiplied with the number of model parameters per layer. The rjMCMC is able to infer the number of model parameters from data itself by allowing it to vary in the iterative inversion process, converge to proper parameterization and prevent underparameterization and overparameterization. We also use rjMCMC to enhance uncertainty estimation since it can transdimensionally sample different model spaces of different dimensionalities and can prevent a biased sampling in only one space which may have a different dimensionality than that of the true model space. An ensemble of solutions from difference spaces can statistically reduce the bias for parameter estimation and uncertainty quantification. Inversion uncertainty is comprised of property uncertainty and location uncertainty. Our study revealed that the inversion uncertainty is correlated with the discontinuity of property in such a way that 1) a smaller discontinuity will induce a lower uncertainty in property at the discontinuity but also a higher uncertainty of the location of that discontinuity and 2) a larger discontinuity will induce a higher uncertainty in property at the discontinuity but also a higher ``certainty'' of the location of that discontinuity. Therefore, there is a trade-off between the property uncertainty and the location uncertainty. To our surprise, there is a lot of hidden information in the uncertainty result that we can actually take advantage of due to this trade-off effect. On the basis of our study using rjMCMC, we propose to use the inversion uncertainty as a novel attribute in an optimistic way to characterize the magnitude and the location of subsurface discontinuities and reflectors

Seismic imaging: a practical approach

In the geophysics of oil exploration and reservoir studies, the surface seismic method is the most commonly used method to obtain a subsurface model in 2 or 3 dimensions. This method plays an increasingly important role in soil investigations for geotechnical, hydrogeological and site characterization studies regarding seismic hazard issues. The goal of this book is to provide a practical guide, using examples from the field, to the application of seismic methods to surface imaging. After reviewing the current state of knowledge in seismic wave propagation, refraction and reflection seismic methods, the book aims to describe how seismic tomography and fullwave form inversion methods can be used to obtain seismic images of the subsurface. Through various synthetic and field examples, the book highlights the benefit of combining different sets of data: refracted waves with reflected waves, and body waves with surface waves. With field data targeting shallow structures, it shows how more accurate geophysical models can be obtained by using the proposed hybrid methods. Finally, it shows how the integration of seismic data (3D survey and VSP), logging data (acoustic logging) and core measurements, combined with a succession of specific and advanced processing techniques, enables the development of a 3D high resolution geological model in depth. In addition to these examples, the authors provide readers with guidelines to carry out these operations, in terms of acquisition, as well as processing and interpretation. In each chapter, the reader will find theoretical concepts, practical rules and, above all, actual application examples. For this reason, the book can be used as a text to accompany course lectures or continuing education seminars. This book aims to promote the exchange of information among geologists, geophysicists, and engineers in geotechnical fields