A seismic scattered wave field numerical simulation method based on a pre-gsor method

Abstract

The present application belongs to the technical field of seismic scattering wave field numerical simulation, and relates to a seismic scattering wave field numerical simulation method based on a PreGSOR method, comprising the following steps: first, introducing a complex wave number background Green function with a virtual component into a generalized successive over-relaxation iteration method for solving an L-S integral equation to obtain a generalized successive over-relaxation iteration algorithm for solving a complex wave number L-S equation, energy loss caused by the virtual component of the Green function is compensated by the virtual component of the scattering site, and the solution of the original L-S equation remains unchanged; second, introducing a diagonal precondition operator to weight the coefficient matrix to reduce the condition number of the coefficient matrix and accelerate the convergence speed of the iteration series to obtain a convergent Born series suitable for a seismic strong velocity disturbance medium. The present method solves the problem that the traditional iteration series method diverges, converges slowly or stagnates when solving the LS equation of the seismic strong velocity disturbance due to the too large condition number of the coefficient matrix of the discrete equation set.

Description

Technical Field

[0001] This invention belongs to the field of numerical simulation technology of seismic scattered wavefield, specifically relating to a method for numerical simulation of seismic scattered wavefield using the complex wave number LS integral equation, and particularly relating to a numerical iterative method for solving the dual equation system of the complex wave number LS equation using the preconditioned generalized successive overrelaxation iteration (PreGSOR) method with preconditioned operator constraints. Background Technology

[0002] In recent years, with the increasing complexity and precision of geological exploration problems, traditional multiple-coverage acquisition techniques have been replaced by wide-azimuth, wide-bandwidth, and high-density (two wide and one high) acquisition techniques. Consequently, the acquired reflection seismic data contains more information related to complex geological structures. Studies have shown that the coexistence of heterogeneous geological bodies of different scales in complex mountainous areas leads to extremely complex reflection seismic signals, specifically manifested as interference between multiple wave groups, weak reflection signal energy, and low signal-to-noise ratio. Even more concerning, in some areas, the reflection data lacks conventional reflection waves altogether, containing only scattering and diffraction—fluctuations typically considered noise in conventional processing. Therefore, it is essential to study reflection seismic data processing and imaging techniques based on scattering theory. As the foundation of seismic migration imaging and full-waveform inversion, forward modeling of seismic scattered wavefields is of great significance for understanding the distribution and properties of subsurface heterogeneous bodies.

[0003] The Helmholtz equations are commonly used to describe electromagnetic scattering, seismic scattering, and other problems. According to perturbation theory, the Helmholtz equations describing scattering problems can be transformed into Fredholm-type integral equations of the second kind, namely the Lippermann–Schwinger (LS) equations. When using the Lippermann–Schwinger (LS) equations to solve seismic wave scattering problems, the LS equations can be adopted... Numerical methods such as the Galerkin method, the method of moments, and iterative methods are used to solve the LS equations. When using iterative methods to solve the LS equations, the widely used Born scattering series can be obtained. However, the Born series only converges when the scatterer size is small or when scattering is weak. How to improve the convergence of the Born series and extend it from weak scattering to strong scattering to meet practical needs has always been a core focus of scattering problem research.

[0004] Several research findings have been published by scholars both domestically and internationally regarding improvements to the convergence properties of the Born scattering series. Osnabrugge et al. obtained a background Green's function with a damping factor by introducing an imaginary component into the background wavenumber and selected an appropriate pre-solution factor to ensure the convergence of the Born series iteration. This method demonstrated high accuracy and fast computational efficiency in optical scattering numerical simulations. Jakobsen and Wu re-derived the renormalized De Wolf scattering series using the De Wolf series and the T-matrix propagation operator, and analyzed its convergence characteristics in detail. Huang et al. reinterpreted the convergent Born series obtained by Osnabrugge from a renormalization perspective and analyzed the convergence properties of this method in forward modeling of low-to-medium frequency seismic wavefields under strongly perturbed scattering media through numerical examples, achieving good application results. Jakobsen and Wu derived the renormalized scattering series solution based on renormalization group theory, showing good convergence in highly inhomogeneous scattering media. Xu Yangyang et al. applied a generalized over-relaxation iterative algorithm for solving strong optical scattering problems to calculate the numerical solution of the Lippmann–Schwinger (LS) problem. By analyzing the applicability and convergence of this algorithm under reflection seismic conditions, they obtained a convergent Born scattering series. However, under conditions of strong velocity perturbation media, the convergence speed of the iterative series remains slow or even stalls at high computation frequencies.

[0005] When using the LS equations to solve seismic scattering problems, current techniques lack a numerical analysis of the mathematical properties of the coefficient matrix and its impact on the convergence of iterative series, specifically addressing the dual equations of the LS equations in seismic scattering. Furthermore, no convergent numerical iterative algorithm or numerical simulation results of the scattered wavefield are provided. In addition, under strongly velocity-perturbed media, the numerical results of high-frequency Born series significantly affect the forward modeling, migration imaging, and inversion results of the seismic scattered wavefield, and improving the convergence of high-frequency Born series remains a research challenge. Current techniques also lack numerical simulation results of the scattered wavefield using high-frequency convergent Born series under strongly velocity-perturbed media. Summary of the Invention

[0006] The purpose of this invention is to provide a numerical simulation method for seismic scattering wavefield based on a generalized successive over-relaxation iterative method with preconditions, in order to solve the problem of divergence and stagnation of convergence speed of high-frequency Born scattering series in strongly seismic scattering media.

[0007] The objective of this invention is achieved through the following technical solution:

[0008] First, the complex wavenumber background Green's function with an imaginary component is introduced into the generalized successive over-relaxation iterative method for solving the LS integral equation, resulting in a generalized successive over-relaxation iterative algorithm for solving the complex wavenumber LS equation. Then, a diagonal preconditioner is introduced to perform row weighting on the coefficient matrix to reduce the condition number of the coefficient matrix and ensure the convergence of the Born scattering series. Finally, the single-frequency scattered wavefield at all frequencies is calculated and inverse FFT is performed to obtain the numerical simulation results of the time-domain seismic scattered wavefield.

[0009] The numerical simulation method for seismic scattered wavefield based on the generalized successive over-relaxation iterative method with pre-conditional constraints, as described in this invention, includes the following steps:

[0010] A numerical simulation method for seismic scattered wavefields based on the PreGSOR method includes the following steps:

[0011] a. Read the velocity model and set the source wavelet, source, receiver point parameters, and background velocity;

[0012] b. Uniformly divided integration region Ω V j Let j be the j-th partition element; take the piecewise constant function space. L2(Ω) is a Hilbert space; basis functions e j (r) satisfies:

[0013]

[0014] c. Based on the interpolation basis function e in step b j (r), construct the interpolation operator:

[0015]

[0016] Where, r j V represents j The center of , v(r) is any continuous function;

[0017] d. Construct the complex wavenumber LS equation and use the interpolation operator P n Using weakly singular kernels for complex wavenumber LS equations The method involves numerical discretization, yielding the matrix form of the linear equation system as follows:

[0018]

[0019] Among them, u=u(r i (i = 1, 2, ..., n) represents the wave field to be determined in the frequency domain. Let s(w) represent the incident wave field, and r represent the source wavelet. s Indicates the location of the earthquake source. It is composed of weight functions The complex wave number Green's operator constitutes It is a complex scattering potential The diagonal operator formed;

[0020] e. Solve the wave field u(r) using a generalized successive over-relaxation iterative algorithm with pre-conditional constraints on the linear equation system in step d. The iterative format is as follows:

[0021]

[0022] or

[0023]

[0024] Where γ is the preconditioner operator, taken as... For the diagonal preconditioner, the coefficient matrix is ​​left-multiplied by the preconditioner γ and then row-weighted; α n As a relaxation factor, it is used to minimize the residual ||r| during the iteration process. n ||To obtain, consider The result of the nth iteration, taking α n Make ||r n || Take the minimum value in the Hilbert space;

[0025]

[0026] f. Substitute the single-frequency source scattered wave field calculated in step e into the real wave number LS equation to obtain the single-frequency scattered wave field at all detector points.

[0027] g. Calculate the scattered wave field at the detector point at all frequencies, and then perform an inverse Fourier transform to obtain the forward modeling record of the scattered wave field of a single shot.

[0028] Further, in step d, the complex wavenumber LS equation is constructed, and the interpolation operator P is used. n Using weakly singular kernels for complex wavenumber LS equations The method for numerical discretization includes the following steps:

[0029] d1. Define the velocity perturbation function v0 is the constant background velocity, and v(r) is the medium velocity. According to the principle of wave field superposition, the Helmholz equation satisfied by the scalar wave equation in the frequency domain is rewritten in the form of an equivalent source:

[0030]

[0031] Where ω is the angular frequency, s(ω)δ(rr) s ) indicates that it is excited by r sThe point source is u(r), which is the wave field to be determined. By adding iεu(r) (ε>0) to both sides of equation (3), we obtain

[0032]

[0033] Pick Where |O| max =max r |O(r)|;

[0034] d2. Define the complex background wavenumber Then there is the complex wave number background Green's function. satisfy

[0035]

[0036] in, The analytical solution is:

[0037]

[0038] d3, corresponding to the complex background wavenumber There is a complex scattering position Based on Green's theorem, the solution to equation (4) can be expressed in integral form, yielding the complex wavenumber LS equation:

[0039]

[0040] Incorporating the linear integral operator with the Dirac function, the matrix chain multiplication form of the complex wavenumber LS equation is:

[0041]

[0042] in,

[0043] d4. Define the integral operator Using interpolation operator P n If the continuous part of the integrand is replaced by its interpolation function, then the integral operator K... Approximate K n It can be represented as:

[0044]

[0045] Wherein, weight function Then the complex wave number LS equation The approximate system of equations can be expressed as:

[0046]

[0047] Using the matrix chain multiplication form (8) of the complex wavenumber LS equation, the matrix (operator) of the above equation is expressed as:

[0048]

[0049] Where I is the identity matrix, It is composed of weight functions The complex wave number Green's operator constitutes It is a complex scattering potential The diagonal operator formed.

[0050] Compared with the prior art, the beneficial effects of the present invention are:

[0051] This invention proposes a pre-conditional, generalized successive over-relaxation iterative algorithm for solving the complex wavenumber LS integral equation for numerical simulation of seismic scattered wavefields. This method addresses the problem of excessively large condition number in the coefficient matrix of the discrete equation system, leading to divergence, slow convergence, or even stagnation in the iterative series, when solving the LS equation for strong seismic velocity perturbations using traditional iterative series methods. By introducing the complex wavenumber Green's function with an imaginary component into the generalized over-relaxation iterative method for solving the LS integral equation, and further employing a pre-conditional operator to reduce the condition number of the coefficient matrix, the convergence speed of the over-relaxation iterative series is accelerated, and the convergence performance of the iterative series is improved. This invention has the following advantages:

[0052] 1. This invention weakens the singularity of the Green's function by selecting an appropriate attenuation factor ε and employing a complex wavenumber Green's function with an imaginary component, thereby localizing the energy of the Green's function and causing it to decay exponentially. Furthermore, the energy loss caused by the imaginary component of the Green's function is compensated by the imaginary component of the scattering potential, preserving the solution to the Helmholtz equation unchanged.

[0053] 2. Introduce the diagonal preconditioner γ to left-multiply the coefficient matrix of the dual equation system of the complex wavenumber LS equation. Perform row weighting, so that Each row has approximately the same ∞-norm, which reduces the condition number of the coefficient matrix;

[0054] 3. This invention employs a generalized successive over-relaxation iterative algorithm to solve the dual equations of the complex wavenumber LS equation after preconditioning constraints, thereby obtaining a generalized successive over-relaxation iterative algorithm with preconditioning operator constraints. This algorithm can obtain convergent series solutions at all calculation frequencies under arbitrarily large strong velocity perturbation media. Attached Figure Description

[0055] To more clearly illustrate the technical solutions of the embodiments of the present invention, the accompanying drawings used in the embodiments will be briefly introduced below. It should be understood that the following drawings only show some embodiments of the present invention and should not be regarded as a limitation on the scope. For those skilled in the art, other related drawings can be obtained based on these drawings without creative effort.

[0056] Figure 1 This is a flowchart of the overall process of the numerical simulation method for seismic scattered wavefield based on the PreGSOR method of the present invention.

[0057] Figure 2 For numerical simulation of velocity model;

[0058] Figures 3a-3b This represents the source-scattered wave field at a frequency of f = 30 Hz.

[0059] Figures 4a-4b This represents the source-scattered wave field at a frequency of f = 60 Hz.

[0060] Figures 5a-5b The curve shows the normalized residual as a function of the number of iterations.

[0061] Figures 6a-6b The scattered wave field at the detector point at a frequency of f = 30 Hz;

[0062] Figures 7a-7b The scattered wave field at the detector point at a frequency of f = 60 Hz;

[0063] Figures 8a-8b The results are numerical simulations of the finite difference and complex wavenumber LS equations using the Pre-GSOR iterative method. Detailed Implementation

[0064] The present invention will be further described below with reference to embodiments:

[0065] The present invention will now be described in further detail with reference to the accompanying drawings and embodiments. It should be understood that the specific embodiments described herein are merely illustrative of the invention and not intended to limit it. Furthermore, it should be noted that, for ease of description, the accompanying drawings show only the parts relevant to the present invention, and not all of the structures.

[0066] It should be noted that similar reference numerals and letters in the following figures indicate similar items; therefore, once an item is defined in one figure, it does not need to be further defined and explained in subsequent figures. Furthermore, in the description of this invention, terms such as "first," "second," etc., are used only to distinguish descriptions and should not be construed as indicating or implying relative importance.

[0067] This invention presents a forward modeling method for seismic scattered wavefields based on a generalized successive over-relaxation iterative algorithm with pre-conditional operator constraints. First, an attenuation factor is introduced and the Helmholz equation is rewritten to include an imaginary component in the background Green's function. Then, the complex wavenumber LS equation is obtained based on the representation theorem. The generalized successive over-relaxation iterative algorithm is then used to solve the forward modeling of the scattered wavefield. The complex wavenumber LS equations are discrete, and a diagonal preconditioner is introduced to row-weight the coefficient matrix to reduce the condition number of the coefficient matrix, thereby accelerating the convergence speed of the iterative series and ensuring convergence.

[0068] This invention provides a numerical simulation method for seismic scattered wave fields based on the PreGSOR method, comprising the following steps:

[0069] a. Read the velocity model, set the source wavelet, source, receiver point parameters, and background velocity; the model size is Nx*Nz;

[0070] b. Uniform subdivision of the integration region (velocity model) Ω V j Let j be the j-th partition element; take the piecewise constant function space. L2(Ω) is a Hilbert space; basis functions e j (r) satisfies:

[0071]

[0072] c. Based on the interpolation basis function e in step b j (r), construct an orthogonal constant interpolation operator:

[0073]

[0074] Where, r j V represents j The center of , v(r) is any continuous function.

[0075] d. Construct the complex wavenumber LS equation and use the interpolation operator P n Using weakly singular kernels for complex wavenumber LS equations Numerical discretization is performed. A velocity perturbation function is defined. v0 is the constant background velocity, and v(r) is the medium velocity. Based on the principle of wave field superposition, the Helmholz equation satisfied by the scalar wave equation in the frequency domain is rewritten in the form of an equivalent source:

[0076]

[0077] Where ω is the angular frequency, s(ω)δ(rr) s ) indicates that it is excited by r s The point source is u(r), and the wave field to be determined is u(r).

[0078] Adding iεu(r) (ε>0) to both sides of formula (3) yields

[0079]

[0080] Pick Where |O| max =max r |O(r)|.

[0081] Define complex background wavenumber Then the complex wave number background Green's function satisfy

[0082]

[0083] in, The analytical solution is:

[0084]

[0085] In actual calculations, a suitable constant 'a' needs to be selected based on the perturbation intensity of the model to ensure... It is not zero within the valid calculation area.

[0086] Corresponding to complex background wavenumber There is a complex scattering position Based on Green's theorem, the solution to equation (4) can be expressed in integral form, yielding the complex wavenumber LS equation:

[0087]

[0088] Incorporating the linear integral operator with the Dirac function, the matrix chain multiplication form of the complex wavenumber LS equation is:

[0089]

[0090] in,

[0091] Define the integral operator Using interpolation operator P n If the continuous part of the integrand is replaced by its interpolation function, then the integral operator K... Approximate K n It can be represented as:

[0092]

[0093] Wherein, weight function Then the complex wave number LS equation The approximate system of equations can be expressed as:

[0094]

[0095] Using the matrix chain multiplication form (8) of the complex wavenumber LS equation, the matrix (operator) of the above equation is expressed as:

[0096]

[0097] Where I is the identity matrix, It is composed of weight functions The complex wave number Green's operator constitutes It is a complex scattering potential The matrix form of the system of linear equations formed by the diagonal operators is:

[0098]

[0099] Among them, u=u(r i (i = 1, 2, ..., n) represents the wave field to be determined in the frequency domain. Let s(w) represent the incident wave field, and r represent the source wavelet. s Indicates the location of the earthquake source.

[0100] e. Solve the wave field u to be determined using the generalized successive over-relaxation iterative algorithm (PreGSOR) with preconditional constraints on the linear equation system in step d. The iterative format is as follows:

[0101]

[0102] or

[0103]

[0104] Where γ is the preconditioner operator, taken as... For the diagonal preconditioner, the coefficient matrix is ​​left-multiplied by the preconditioner γ and then row-weighted. α n As a relaxation factor, it is used to minimize the residual ||r| during the iteration process. n ||To obtain, consider The result of the nth iteration, taking α n Make ||r n || Take the minimum value in the Hilbert space;

[0105]

[0106] f. Substitute the single-frequency source scattered wave field calculated in step e into the real wave number LS equation to obtain the single-frequency scattered wave field at all detector points.

[0107] g. Calculate the scattered wave field at the detector point at all frequencies, and then perform an inverse Fourier transform to obtain the forward modeling record of the scattered wave field of a single shot.

[0108] To better illustrate the implementation effect of the above specific algorithm process, a specific example is given below:

[0109] Example 1

[0110] The overall flow of the forward modeling method for seismic scattered wavefields based on the PreGSOR method of this invention is as follows: Figure 1 As shown.

[0111] a. Reading speed model ( Figure 2 The model size is 1850m × 1050m, the source wavelet is a Ricker wavelet with a dominant frequency of 15Hz, and the source location is r. s = (900m, 0m), the detector is placed on the interface of z = 0m, and the background velocity is v0 = 1500m / s.

[0112] b. Divide the integration region Ω into a uniform partition and construct constant interpolation basis functions.

[0113] c. Construct constant interpolation operators based on the basis functions from step b.

[0114] d. Construct the complex wavenumber LS equations. Based on the condition number of the coefficient matrix in the dual system of the real wavenumber LS equations, and to ensure the convergence of the algorithm while obtaining an effective numerical solution, the attenuation factor ε is taken as 0.03|O|. max Using the interpolation operator P n Using weakly singular kernels for complex wavenumber LS equations Numerical discretization was performed to obtain the dual equations of the complex wavenumber LS equations, and preconditioning operators were then used.

[0115] e. Solve the wave field u(r) using the PreGSOR iterative algorithm on the linear equations in step d. Figure 3a and Figure 3b Let represent the real and imaginary parts of the source-scattered wave field at a frequency of 30 Hz. Figure 4a and Figure 4b Let represent the real and imaginary parts of the source-scattered wave field at a frequency of 60 Hz. Figure 5a and Figure 5b The figures show the curves of the normalized residuals of the iterative series as a function of the number of iterations at frequencies of 30Hz and 60Hz, respectively. Compared with the generalized over-relaxation iterative method of the real wavenumber LS equation, the PreGSOR iterative method of the complex wavenumber LS equation has better convergence and also has good convergence at high frequencies.

[0116] f. Substitute the single-frequency source scattered wave field calculated in step e into the real wave number LS equation to obtain the single-frequency scattered wave field at all detector points. Figure 6a and Figure 6bLet be the real and imaginary parts of the scattered wave field at the detector point at a frequency of 30Hz. Figure 7a and Figure 7b Let represent the real and imaginary parts of the scattered wave field at the detector point at a frequency of 60 Hz. It can be seen that the PreGSOR iterative numerical solution of the complex wavenumber LS equation has a high degree of matching with the numerical solution obtained by the direct method (Singular Value Decomposition (SVD)) of the real wavenumber LS equation.

[0117] g. Calculate the scattered wavefield at the receiver point at all frequencies, and then perform an inverse Fourier transform to obtain the forward modeling record of the scattered wavefield of a single shot, such as... Figure 8b As shown. The results of this example are compared with those of the finite difference forward modeling simulation, as follows. Figure 8a As shown, the numerical simulation results obtained using the PreGSOR method with complex wavenumber LS equations do not suffer from a loss of computational accuracy, such as... Figure 8b As shown.

[0118] The seismic scattered wave field forward modeling method based on the generalized successive over-relaxation iterative algorithm with preconditioning operator constraints proposed in this invention can be further applied to seismic scattered wave migration imaging and full waveform inversion to calculate the migration Green's function and forward and inverse Green's function of the full wave field.

[0119] Note that the above description is merely a preferred embodiment of the present invention and the technical principles employed. Those skilled in the art will understand that the present invention is not limited to the specific embodiments described herein, and various obvious changes, readjustments, and substitutions can be made without departing from the scope of protection of the present invention. Therefore, although the present invention has been described in detail through the above embodiments, the present invention is not limited to the above embodiments, and may include many other equivalent embodiments without departing from the concept of the present invention, the scope of which is determined by the scope of the appended claims.

Claims

1. A numerical simulation method for seismic scattered wavefields based on the PreGSOR method, characterized in that, Includes the following steps: a. Read the velocity model and set the source wavelet, source, receiver point parameters, and background velocity; b. Uniformly divided integration region Ω V j Let j be the j-th partition element; take the piecewise constant function space. L2(Ω) is a Hilbert space; basis functions e j (r) satisfies: c. Based on the interpolation basis function e in step b j (r), construct the interpolation operator: Where, r j V represents j The center of , v(r) is any continuous function; d. Construct the complex wavenumber LS equation and use the interpolation operator P n Using weakly singular kernels for complex wavenumber LS equations The method involves numerical discretization, yielding the matrix form of the linear equation system as follows: Among them, u=u(r i (i = 1, 2, ..., n) represents the wave field to be determined in the frequency domain. Let s(w) represent the incident wave field, and r represent the source wavelet. s Indicates the location of the earthquake source. It is composed of weight functions The complex wave number Green's operator constitutes It is a complex scattering potential The diagonal operator formed; e. Solve the wave field u(r) using a generalized successive over-relaxation iterative algorithm with pre-conditional constraints on the linear equation system in step d. The iterative format is as follows: or Where γ is the preconditioner operator, taken as... For the diagonal preconditioner, the coefficient matrix is ​​left-multiplied by the preconditioner γ and then row-weighted; α n As a relaxation factor, it is used to minimize the residual ||r| during the iteration process. n ||To obtain, consider The result of the nth iteration, taking α n Make ||r n || Take the minimum value in the Hilbert space; f. Substitute the single-frequency source scattered wave field calculated in step e into the real wave number LS equation to obtain the single-frequency scattered wave field at all detector points. g. Calculate the scattered wave field at the detector point at all frequencies, and then perform an inverse Fourier transform to obtain the forward modeling record of the scattered wave field of a single shot.

2. The numerical simulation method for seismic scattered wavefield based on the PreGSOR method according to claim 1, characterized in that: Step d: Construct the complex wavenumber LS equation and use the interpolation operator P n Using weakly singular kernels for complex wavenumber LS equations The method for numerical discretization includes the following steps: d1. Define the velocity perturbation function v0 is the constant background velocity, and v(r) is the medium velocity. According to the principle of wave field superposition, the Helmholz equation satisfied by the scalar wave equation in the frequency domain is rewritten in the form of an equivalent source: Where ω is the angular frequency, s(ω)δ(rr) s ) indicates that it is excited by r s The point source is u(r), which is the wave field to be determined. By adding iεu(r) (ε>0) to both sides of equation (3), we obtain Pick Where |O| max =max r |O(r)|; d2. Define the complex background wavenumber Then there is the complex wave number background Green's function. satisfy in, The analytical solution is: d3, corresponding to the complex background wavenumber There is a complex scattering position Based on Green's theorem, the solution to equation (4) can be expressed in integral form, yielding the complex wavenumber LS equation: Incorporating the linear integral operator with the Dirac function, the matrix chain multiplication form of the complex wavenumber LS equation is: in, d4. Define the integral operator Using interpolation operator P n If the continuous part of the integrand is replaced by its interpolation function, then the integral operator K... Approximate K n It can be represented as: Wherein, weight function Then the complex wave number LS equation The approximate system of equations can be expressed as: Using the matrix chain multiplication form (8) of the complex wavenumber LS equation, the matrix (operator) of the above equation is expressed as: Where I is the identity matrix, It is composed of weight functions The complex wave number Green's operator constitutes It is a complex scattering potential The diagonal operator formed.

Images