Seismic Inversion Chaiwoot Boonyasiriwat October 21, 2020 - - PowerPoint PPT Presentation

Seismic Inversion

Chaiwoot Boonyasiriwat

October 21, 2020

▪ Petroleum was generated in organic-rich sedimentary rocks and then migrated upward and was trapped in various structures such as anticline and fault. ▪ Images of subsurface structures can be used to interpret where petroleum reservoirs are located in the earth.

2

Petroleum Exploration

Image from http://media-3.web.britannica.com/eb-media/69/669-004-7B2A0E51.jpg

▪ Mechanical waves generated by several sources propagate into the earth and reflect back to the surface. ▪ The wavefields recorded by receivers at the surface are the observed data used to estimate the earth structures. ▪ Seismic exploration can be performed on both land and marine environments.

3

Seismic Exploration

Subsurface Structure Seismic Data

4

Marine Seismic Acquisition

https://www.rigzone.com/training/insight.asp?insight_id=303&c_id=

Air gun source

https://9ca7a3fb47-custmedia.vresp.com/ 303b3a0f74/Bolt_e_source_airgun.jpg http://filemaker2-server.cbl.umces.edu/ sensorimages/9663-TC4032.gif

Hydrophone

5

Land Seismic Acquisition

http://static.panoramio.com/photos/large/27525469.jpg

Geophone

http://img.diytrade.com/cdimg/1296480/15832570/0/1286613918/LAND _GEOPHONE_STRING.jpg

Vibroseis

http://www.oilinuganda.org/wp-content/media/2013/06/Seismic-Illustration1.bmp

A seismic velocity model v can be estimated by first- arrival travel time T by solving the constrained minimization problem subject to where To is the observed data, xr and xs are the receiver and source locations, respectively, and s = 1/v is slowness.

6

Traveltime Inversion

Using the method of Lagrange multiplier, we first form the Lagrangian as where When ,

7

Traveltime Inversion

We have Setting yields the state problem Setting yields the adjoint problem

8

Traveltime Inversion

Since , The gradient of the objective functional is the product between the slowness s and the Lagrange multiplier . The state and adjoint problems can be solved using ▪ Ray tracing ▪ Fast marching method (Sethian and Popovici, 1999) ▪ Fast sweeping method (Zhao, 2004) ▪ Fast iterative method (Jeong and Whitaker, 2008)

9

Traveltime Inversion

Let t be unit vector along ray path. We have where is slowness So or

10

Ray Tracing

Rawlinson et al. (2008, p. 207)

Taking gradient of the last equation and using the definition of tangent vector t yields ray equation The ray equation is the Euler-Lagrange equation of the Lagrangian which is called the Fermat integral. Fermat’s principle of stationary time: ray paths correspond to extremal curves of Fermat integral.

11

Ray Equation

▪ Ray paths are characteristic curves of the Hamiltonian which can be written in various form. ▪ In the Hamiltonian formalism of ray tracing, the momentum is corresponding to the slowness vector defined as . ▪ The eikonal equation can then be written as with a Hamiltonian

12

Hamiltonian Formalism

▪ The 6-vector is called the canonical coordinate vector in a 6D phase space. ▪ is a hypersurface in the phase space. ▪ On the hypersurface, we have ▪ The Einstein summation convention is always used. ▪ Along a characteristic, we have

13

Characteristic System

Cerveny (2001, p. 104)

We then obtain the first 6 ODEs of the system The last equation for T can be obtained by which can be written in vector form as

14

Characteristic System

Cerveny (2001, p. 104)

▪ This system of 7 first-order ODEs are the Hamilton equations and can be solved for the characteristic curve r, the slowness vector p, and the traveltime T. ▪ Note that the first 6 equations are independent of the last equation. ▪ The parameter u depends on the form of H.

15

Characteristic System

Cerveny (2001, p. 104)

▪ Starting from a general Hamiltonian ▪ The corresponding characteristic system is

16

General Ray Tracing Systems

Cerveny (2001, p. 105)

▪ When , ▪ When , ▪ When ,

17

Examples of Ray Tracing Systems

Cerveny (2001, p. 105)

Because , the 3 systems become ▪ When , ▪ When , ▪ When ,

18

Examples of Ray Tracing Systems

Cerveny (2001, p. 106)

Initial conditions Initial point is source point Initial momentum is slowness vector where is the shooting angle Let’s use the system with In 2D, the system becomes

19

Shooting Method: Initial-Value Problem

General form of Hamilton-Jacobi equation (HJE) where is the generalized coordinates is the Hamilton’s principal function is the Hamiltonian is the generalized momentum Hamilton-Jacobi equation is first-order nonlinear PDE while Euler-Lagrange equations are second-order PDE. It can be shown that and, therefore, is the classical action.

20

Hamilton-Jacobi Equation

https://en.wikipedia.org/wiki/Hamilton-Jacobi_equation

Static Hamilton-Jacobi equation (HJE) can be written as So the eikonal equation is a special case of HJE when

21

Static Hamilton-Jacobi Equation

Gremaud and Kuster (2006)

On a 2D Cartesian grid, the Hamilton-Jacobi equation can be discretized as Let Godunov-flux Hamiltonian (first-order upwind): Lax-Friedrichs-flux Hamiltonian

22

Numerical Hamiltonian

Gremaud and Kuster (2006, p. 2)

We consider the 2D boundary-value problem “The slowness field can be unbounded, corresponding to the presence of obstacles.” (Gremaud and Kuster, 2006)

23

2D Problem

In 2D, the upwind scheme can be written as where In practice, we use a square grid with

24

Upwind FD Scheme

Sethian and Popovici (1999)

Godunov upwind scheme can be written as where

25

Godunov Scheme

Zhao (2004)

Let Traveltime at a node can be updated in one of these cases: Case 1: Case 2: Case 3: Case 4:

26

Update Cases

▪ Apply the Godunov upwind scheme to all grid points except the boundary points and replace the

• ld value if the new value is smaller.

▪ 4 sweeping schemes are used:

27

Fast Sweeping Method

Zhao (2004)

▪ Inversion of first-arrival traveltime data can only provide a smoothed velocity model. ▪ Inversion of full-waveform data can provide a higher- resolution velocity model. ▪ Suppose the wave propagation is governed by ▪ Inserting a planewave into the wave equation yields the Helmholtz equation where s is slowness.

28

Frequency-domain FWI

In 2D, wavefield is governed by The direct scattering problem can be solved using the finite-difference frequency-domain (FDFD) method:

2D Helmholtz Problem

Absorbing boundary conditions are used at the boundary For example, the absorbing boundary condition at the bottom is Using the finite difference method, we obtain

2D Helmholtz Problem

x z

Example Result

Ajo-Franklin (2005, p. 13)

▪ Consider the constrained minimization problem subject to ▪ The corresponding Lagrangian is where

32

Frequency-domain FWI

The gradient of objective functional with respect to slowness is when and in this case

33

Frequency-domain FWI

Setting yields the state equation Since setting yields the adjoint equation

34

Frequency-domain FWI

Example Result: Wave Path

Ajo-Franklin (2005, p. 14)

▪ Ajo-Franklin, J. B., 2005, Frequency-domain modeling techniques for the scalar wave equation: An introduction, Earth Resources Laboratory Technical Report, MIT. ▪ Cerveny, V., 2001, Seismic Ray Theory, Cambridge University Press. ▪ Gremaud, P. A., and C. M. Kuster, 2006, Computational study of fast methods for the eikonal equation: SIAM Journal on Scientific Computing, 27, 6, 1803-1816. ▪ Jeong, W., and R. T. Whitaker, 2008, A fast iterative method for eikonal equations: SIAM Journal of Scientific Computing, 30(5), 2512–2534. ▪ Rawlinson, N., J. Hauser, and M. Sambridge, 2008, Seismic ray tracing and wavefront tracking in laterally heterogeneous media: Advances in Geophysics, 49, 203-273. ▪ Sethian, J. A., and A. M. Popovici, 1999, 3-D traveltime computation using the fast marching method: Geophysics 64(2). ▪ Zhao, H., 2004, A fast sweeping method for eikonal equations: Mathematics

• f Computation, 74, 603-627.

36

References