Tal-Ezer Method

When the velocity structure varies arbitrarily in the lateral direction, the eigenvalues and eigenvectors of matrix A can no longer be obtained by inspection. It would therefore seem that a matrix diagonalization would have to be performed before each propagation to equation (11). However, recent work by Tal-Ezer(1986) indicates how the calculation of matrix exponential can be done without having to resort to expensive matrix diagonalizations. The solution is based on a Chebychev expansion for the function e^x given by

$$e^x=\sum_{k=0}^\infty {C_k J_k(R) T_k\left ({x\over R}\right )}, \vert x\vert<R$$ (18)

where x is imaginary and C₀=1 and C_k=2 for k>0. The polynomials T_k(x) satisfy the recurrence relations

T₀(x) = 1,

T₁(x) = x,

and

T_k+1(x)=T_k-1(x)+2xT_k(x)

(Tal-Ezer, 1986). By analogy with equation (18), the exponent in equation (4) is expanded according to

$$\pmatrix {P \cr {\partial P \over \partial z}}_{z+dz}=\exp[A... ...z\over R})}{\pmatrix {P \cr {\partial P \over \partial z}}_{z}}$$ (19)

This expansion is valid when the eigenvalues of [Adz] are purely imaginary and R must be chosen large enough to span the range of the eigenvalues of [Adz]. It was shown in Tal-Ezer(1986) that for k > R, the series expansion converges exponentially. The number of terms required in the sum in equation (19) will therefore always be finite. Equation (19) serves as the basis for implementing the generalized phase-shift migration. First the range R of the eigenvalues of [Adz] need to be estimated. Based on the case of laterally uniform velocity, Kosloff showed the estimate r=wdz/c_min with c_min denoting the lowest velocity in the strip (z,z+dz), is sufficient for stable results. The Bessel functions J_k(R) are computed next. The solution

$$\pmatrix { P \cr {\partial P\over \partial z}}_{z+dz}$$

is then calculated recursively according to following: The first two values of T₀, T₁ needed to initialize the recursion are given by

$$T_0\pmatrix { P \cr {\partial P\over \partial z}}_z= \pmatrix { P \cr {\partial P\over \partial z}}_z$$

and

$$T_1\pmatrix { P \cr {\partial P\over \partial z}}_z= [{Adz \over R}]\pmatrix { P \cr {\partial P\over \partial z}}_z$$

The next polynomial is generated by the recursion formula

$$T_{k+1}\pmatrix { P \cr {\partial P\over \partial z}}_z= T_{... ...dz \over R}]T_k\pmatrix { P \cr {\partial P\over \partial z}}_z$$

These steps are repeated until a sufficient number of terms have been calculated in the sum(19). Then the solution is carried out in the next lower level. When $P(x,z,\omega)$ has been calculated completely, the final migrated section is cumulated by

$$P_{mig}(x,z)=\sum_\omega P(x,z,\omega)$$

(Kosloff and Baysal, 1983).