Fast Computation of Fourier Integral Operators

We introduce a general purpose algorithm for rapidly computing certain types of oscillatory integrals which frequently arise in problems connected to wave propagation and general hyperbolic equations. The problem is to evaluate numerically a so-called Fourier integral operator (FIO) of the form $\int e^{2\pi i \Phi(x,\xi)} a(x,\xi) \hat{f}(\xi) \mathrm{d}\xi$ at points given on a Cartesian grid. Here, $\xi$ is a frequency variable, $\hat f(\xi)$ is the Fourier transform of the input $f$, $a(x,\xi)$ is an amplitude and $\Phi(x,\xi)$ is a phase function, which is typically as large as $|\xi|$; hence the integral is highly oscillatory at high frequencies. Because an FIO is a dense matrix, a naive matrix vector product with an input given on a Cartesian grid of size $N$ by $N$ would require $O(N^4)$ operations. This paper develops a new numerical algorithm which requires $O(N^{2.5} \log N)$ operations, and as low as $O(\sqrt{N})$ in storage space. It operates by localizing the integral over polar wedges with small angular aperture in the frequency plane. On each wedge, the algorithm factorizes the kernel $e^{2 \pi i \Phi(x,\xi)} a(x,\xi)$ into two components: 1) a diffeomorphism which is handled by means of a nonuniform FFT and 2) a residual factor which is handled by numerical separation of the spatial and frequency variables. The key to the complexity and accuracy estimates is that the separation rank of the residual kernel is \emph{provably independent of the problem size}. Several numerical examples demonstrate the efficiency and accuracy of the proposed methodology. We also discuss the potential of our ideas for various applications such as reflection seismology.Comment: 31 pages, 3 figure

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

Elastic wave propagation in anisotropic media : source theory, traveltime computations and migration

Thesis (Ph. D.)--Massachusetts Institute of Technology, Dept. of Earth, Atmospheric, and Planetary Sciences, 1992.Title as it appears in the M.I.T. Graduate List, June 1992: Advances in the theory of elastic wave propagation in anisotropic media, source theory, traveltime computations and migration.Includes bibliographical references (p. 214-222).by Arcangelo Gabriele Sena.Ph.D

Offset-continuation stacking: Theory and proof of concept

CAPES - COORDENAÇÃO DE APERFEIÇOAMENTO DE PESSOAL DE NÍVEL SUPERIORCNPQ - CONSELHO NACIONAL DE DESENVOLVIMENTO CIENTÍFICO E TECNOLÓGICOThe offset-continuation operation (OCO) is a seismic configuration transform designed to simulate a seismic section, as if obtained with a certain source-receiver offset using the data measured with another offset. Based on this operation, we have introduced the OCO stack, which is a multiparameter stacking technique that transforms 2D/2.5D prestack multicoverage data into a stacked common-offset (CO) section. Similarly to common-midpoint and common-reflection-surface stacks, the OCO stack does not rely on an a priori velocity model but provided velocity information itself. Because OCO is dependent on the velocity model used in the process, the method can be combined with trial-stacking techniques for a set of models, thus allowing for the extraction of velocity information. The algorithm consists of data stacking along so-called OCO trajectories, which approximate the common-reflection-point trajectory, i.e., the position of a reflection event in the multicoverage data as a function of source-receiver offset in dependence on the medium velocity and the local event slope. These trajectories are the ray-theoretical solutions to the OCO image-wave equation, which describes the continuous transformation of a CO reflection event from one offset to another. Stacking along trial OCO trajectories for different values of average velocity and local event slope allows us to determine horizon-based optimal parameter pairs and a final stacked section at arbitrary offset. Synthetic examples demonstrate that the OCO stack works as predicted, almost completely removing random noise added to the data and successfully recovering the reflection events.The offset-continuation operation (OCO) is a seismic configuration transform designed to simulate a seismic section, as if obtained with a certain source-receiver offset using the data measured with another offset. Based on this operation, we have introduced the OCO stack, which is a multiparameter stacking technique that transforms 2D/2.5D prestack multicoverage data into a stacked common-offset (CO) section. Similarly to common-midpoint and common-reflection-surface stacks, the OCO stack does not rely on an a priori velocity model but provided velocity information itself. Because OCO is dependent on the velocity model used in the process, the method can be combined with trial-stacking techniques for a set of models, thus allowing for the extraction of velocity information. The algorithm consists of data stacking along so-called OCO trajectories, which approximate the common-reflection-point trajectory, i.e., the position of a reflection event in the multicoverage data as a function of source-receiver offset in dependence on the medium velocity and the local event slope. These trajectories are the ray-theoretical solutions to the OCO image-wave equation, which describes the continuous transformation of a CO reflection event from one offset to another. Stacking along trial OCO trajectories for different values of average velocity and local event slope allows us to determine horizon-based optimal parameter pairs and a final stacked section at arbitrary offset. Synthetic examples demonstrate that the OCO stack works as predicted, almost completely removing random noise added to the data and successfully recovering the reflection events.815V387V401CAPES - COORDENAÇÃO DE APERFEIÇOAMENTO DE PESSOAL DE NÍVEL SUPERIORCNPQ - CONSELHO NACIONAL DE DESENVOLVIMENTO CIENTÍFICO E TECNOLÓGICOCAPES - COORDENAÇÃO DE APERFEIÇOAMENTO DE PESSOAL DE NÍVEL SUPERIORCNPQ - CONSELHO NACIONAL DE DESENVOLVIMENTO CIENTÍFICO E TECNOLÓGICOsem informaçãosem informaçã

Efficient seismic imaging with the double plane wave data

Seismic imaging is critical in providing the image of the Earth’s subsurface, and it plays an important role in hydrocarbon explorations. Obtaining high resolution images with accurate reflectivities and accurate positions of subsurface structures is the goal for exploration geophysicists. Reverse time migration (RTM), which solves the two-way wave equation, can resolve all wavefield propagation phenomena. In geologically complex regions, RTM has been proven to outperform other imaging methods in correctly revealing the subsurface structures. However, implementing the traditional pre-stack shot profile RTM is computationally expensive. Time consuming wavefield propagation processes need to be performed for each shot gather to obtain high resolution images. The traditional RTM can become extremely expensive with increasing shot numbers. In this dissertation, I focus on improving the migration efficiency of the RTM using the double plane wave (DPW) data, which are the fully decomposed plane wave data. Three RTM methods are developed to migrate the DPW data, all of which can improve the migration efficiency comparing to the traditional shot profile RTM. Two of the methods utilize the adjoint state method, and they are known as the time domain DPW-based RTM and the frequency domain DPW-based RTM. A third migration method using the DPW data is derived under the Born approximation. This method employs the frequency domain plane wave Green’s functions for imaging, and it is named as frequency domain DPW RTM. Among the three proposed RTM methods, the frequency domain DPW RTM is the most efficient. Comparing to the traditional shot profile pre-stack RTM, the frequency domain DPW RTM can increase migration efficiency of RTM by an order of magnitude, making the frequency domain DPW RTM a preferable option for migrating large seismic datasets. All of the three proposed migration methods can image subsurface structures with given dips, which makes them target-oriented imaging methods. The proposed methods are beneficial to migration velocity analysis. To improve the resolution of migration results, a least squares RTM method using the DPW data is proposed. A Born modeling operator that predict the DPW data at the surface and its adjoint operator, which is a migration operator, are derived to implement the least squares RTM. Both of the operators require only a limited number of plane wave Green’s functions for the modeling and the migration processes. The proposed least squares RTM substantially increases the efficiency of the least squares migration. In the DPW domain, the applicability of the reciprocity principle is also investigated. The reciprocity principle can be applied to the seismic data that are processed with proper seismic processing flow. Utilizing the reciprocity principle, a DPW dataset transformed from one-sided shot gathers can approximate a DPW dataset transformed from split-spread shot gathers. Therefore, I suggest that one-sided acquisition geometries should be extended to the largest possible offsets, and the reciprocity principle should be invoked to improve subsurface illumination. Migration efficiency can be further improved with the help of the reciprocity principle.Geological Science

Imaging of vertical seismic profiling data using the common-reflection-surface stack. Abbildungsverfahren für seismische Daten aus Bohrlochmessungen mit der Common-Reflection-Surface Stapelung

Diese Dissertation beschäftigt sich mit der Entwicklung eines automatisierten, datenorientierten Abbildungsverfahrens, das auf der sogenannten Common-Reflection-Surface (CRS) Stapelung basiert. Durch die Miteinbeziehung von benachbarten Experimenten bei der Rekonstruktion eines Einzelexperiments ergibt sich ein verbessertes Signal-zu-Rauschen Verhältnis und eine starke Bereinigung von Mehrdeutigkeiten. Hauptaugenmerk liegt auf der Adaption der Methode für Bohrloch- und Mehrkomponentendaten

Approximate inversion of generalized Radon transforms

Generalized Radon transforms (GRT) serve, for instance, as linear models for seismic imaging in the acoustic regime. They occur when the corresponding inverse problem is linearized about a known background compression wave speed (Born approximation). The resulting GRT is completely determined by this background velocity. In this work, we present an implementation of approximate inversion formulas for this class of GRTs proposed and analyzed in [Inverse Problems, 34 (2018), 014002, 114001], where we restrict ourselves to layered background velocities in 2D. In a series of numerical experiments, we intensively test our implementation, reproducing theoretical predictions. Further, we drive the validity of the linearization to its limits