Search

The 3D acoustic wave equation for a heterogeneous medium is used for the seismic modeling of compressional (P-) wave propagation in complex subsurface structures. A combination of finite difference and finite integral transform methods is employed to obtain a “2.5D” solution to the 3D equation. Such 2.5D approaches are attractive because they result in computational run times that are substantially smaller than those for the 3D finite difference method. The acoustic parameters of the medium are assumed to be constant in one of the three Cartesian spatial dimensions. This assumption is made to reduce the complexity of the problem, but still retain the salient features of the approach. Simple models are used to address the computational issues that arise in the modeling. The conclusions drawn can also be applied to the more general fully inhomogeneous problem. Although similar studies have been carried out by others, the work presented here is new in the sense that (i) it applies to subsurface models that are both vertically and laterally heterogeneous, and (ii) the computational issues that need to be addressed for efficient computations, which are not trivial, are examined in detail, unlike previous works. We find that it is feasible to generate true-amplitude synthetic seismograms using the 2.5D approach, with computational run times, storage requirements, and other factors, being at reduced and acceptable levels.

SUMMARY Consider a plane homogeneous wave incident upon an interface between two anelastic halfspaces. Computing the plane-wave displacement reflection and transmission coefficients requires determining the proper signs of the vertical slownesses of all the reflected and transmitted waves. In certain cases, this is not straightforward. Previous work has shown that choosing the signs by applying the elastic radiation condition results in certain vertical slownesses, and hence coefficients, varying discontinuously and unphysically with the angle of incidence, and that choosing the signs so that the vertical slownesses vary continuously can also produce errors. We suggest three approaches for treating these cases. In the first approach, the signs are chosen so that the vertical slownesses vary continuously up to a certain angle of incidence (close to the elastic critical angle), with the elastic radiation condition being applied beyond that. We show that this is actually just an extension of the elastic radiation condition to complex-valued squared slownesses. In the spherical wave case in which a point source lies in the upper half-space, the approach also agrees with the results obtained by applying the saddle point method to approximate the integral for the reflected wavefield. This approximation is just geometrical ray theory. The first approach also produces coefficients with unphysical aspects, but these are confined to the critical zone, where geometrical ray theory is not valid anyway. Outside of this zone, the coefficients compare quite well with the elastic ones. In the second approach, real values of the horizontal slowness are used to compute the coefficients. This results in coefficients which are continuous and have no unphysical aspects, but do not compare as well with the elastic ones. The third approach involves using other paths in the complex plane (of horizontal slowness) to evaluate the coefficients. We show two examples. The first example is an approximation which involves replacing certain vertical slownesses with their complex conjugates before the coefficients are computed. This is equivalent to evaluating the horizontal slowness along a curve just above the real axis. The second example involves using a path which yields coefficients that are almost identical to the elastic ones. The third approach also results in continuous coefficients with no unphysical aspects, which compare well with the elastic ones. The errors mentioned above suggest to us that the problem of how to correctly compute plane-wave anelastic reflection and transmission coefficients, even in the apparently simple case of a homogeneous incident plane wave, still requires further study.

When comparing solutions for the propagation of SH waves in plane parallel layered elastic and viscoelastic (anelastic) media, one of the first things that becomes apparent is that in the elastic case the location of the saddle points required to obtain a high frequency approximation are located on the real p axis. This is true of the branch points also. In a viscoelastic medium this is not typical. The saddle point corresponding to an arrival lies in the first quadrant of the complex p-plane as do the branch points. Additionally, in the elastic case the saddle point and branch points lie on a straight line drawn through the origin (the positive real axis in the complex p-plane), while in the viscoelastic case this is generally not the case and the saddle point and branch points lie in such a manner as to indicate the degree of their complex values. In this paper simple SH reflected and transmitted particle displacement arrivals due to a point torque source at the surface in a viscoelastic medium composed of a layer over a half space will be considered. The path of steepest descent defining the saddle point in the first quadrant will be parameterized in terms of a real variable and the high frequency solutions and intermediate analytic results obtained will be used to formulate more specific constraints and observations regarding saddle point location relative to branch point locations in the complex p-plane. As saddle point determination for an arrival is, in general, the solution of a non-linear equation in two unknowns (the real and imaginary parts of the complex saddle point p 0), which must be solved numerically, the use of analytical methods for investigating this problem type is somewhat limited. Numerical experimentation using well documented solution methods, such as Newton's method, was undertaken and some observations were made. Although fairly basic, they did provide for the design of algorithms for the computation of synthetic traces that displayed more efficient convergence and accuracy than those previously employed. This was the primary motivation for this work and the results from the SH problem may be used with minimal modifications to address the more complicated subject of coupled P-SV wave propagation in viscoelastic media. Another reason for revisiting a problem that has received some attention in the literature was to approach it in a fairly comprehensive manner so that a number of specific observations may be made regarding the location of the saddle point in the complex p-plane and to incorporate these into computer software. These have been found to result in more efficient algorithms for the SH wave propagation and a significant enhancement of the comparable software in the P-SV problem.

The Alekseev–Mikhailenko method (AMM) is the name given to a series of algorithms that use one or more finite spatial transforms to reduce the dimensionality of a wave-propagation problem to that of one space dimension and time. This reduced equation is then solved using finite-difference techniques, and the space–time solution is recovered by applying inverse finite spatial transform(s). In this paper the elastodynamic wave equation that governs the coupled P–Sv motion in an isotropic, vertically inhomogeneous elastic half space is investigated using the AMM. Two types of impulsive body forces that may be used to excite the medium are examined, as is the problem of obtaining accurate transformed finite-difference analogues at the free surface. The second of these is accomplished by introducing the boundary conditions that the shear and normal stress must vanish here and by incorporating their transforms into the transformed elastodynamic equations. The stability criterion for the explicit finite-difference method is given cursory treatment, as detailed discussion of this aspect may be found in many texts that deal with the subject of finite differences.A coal-seam model (two thin, low-velocity layers embedded in a half space) illustrates the method. Both horizontal and vertical seismic traces are computed for this model and the results examined in relation to other seismic-modelling techniques.

A high-frequency approximation for nongeometrical arrivals due to highly concentrated sources adjacent to the interfaces between two different elastic media has been developed following their detection in the real field data by Gutowski et al. (1982). In our approach, the “★” waves are associated with secondary saddle-point contributions of the contour of integration circumventing branch points in the integral representation of transmitted wave field. Our formulas are presented in a simple form suitable for an easy incorporation into any ray synthetic seismogram computation. The range of validity of our approximations is discussed and their accuracy tested by comparing the pertinent ray seismograms with their equivalents produced by the Alekseev-Mikhailenko method. A good match between both sets of seismograms suggests that our high-frequency approximations for “★” waves are well justified as are our formulas derived for the head wave-like arrivals driven by the “★” waves, and seen in the full wave field seismograms and detected in real data as well.

The existence of‘*-waves’has, in recent years, prompted a renewed interest in these non-geometrical arrivals, which are generated by point sources located adajcent to plane interfaces. It has led to the re-evaluation of seismic data aquisition techniques and to the question of how to use this real phenomena in enhancing existing seismic interpretation methods. This paper considers a non-geometrical SH-arrival which is generated by a point torque source unrealistically buried within a half-space. The method of solution is essentially the same as presented in an earlier paper, with the modification that the limitation placed on the distance of the source from the interface has been removed in the saddle point method used to obtain a high-frequency approximate solution. In the earlier paper, a preliminary assumption forced the saddle point, which corresponded to the *-wave arrival, to be real when it is generally complex. However, for offsets removed from the distinct ray, the imaginary part of this complex quantity is negligible. A problem which arose when comparing exact synthetic traces with those obtained using zero-order saddle point methods, was the inability to match either the amplitude or phase of the geometrical arrival in the range of offsets when the *-wave and this corresponding geometrical ray were well separated. For this range of offsets the geometrical arrival was approaching grazing incidence and another term in the saddle point expansion of the integral was necessary to rectify this error. This method is also being used to validate the results for higher order terms obtained using asymptotic ray theory. Analytical formulae are given for both the *-wave and the higher order expansion of the geometrical event, together with a comparison of synthetic seismograms using the method developed here and a numerical integration algorithm.

When an explosive point source is buried less than one wavelength from the free surface in a perfectly elastic, homogeneous, isotropic half-space, a very strong non-geometric arrival, tentatively denoted as S* by Hron and Mikhailenko (1981), exists in the reflected shear wave field. It has all basic characteristics of an ordinary shear body wave (linear polarization, transverse particle motion, and shear wave velocity) and, as such, it can be also reflected or transmitted upon incidence at an interface according to the rules governing ordinary body waves. Because of its very strong amplitude which at most epicentral distances exceeds even the amplitude of a direct P wave, the S* arrival is undoubtedly responsible for a large part of the shear wave energy content seen in the seismic records in oil exploration, which deals almost exclusively with shallow explosive sources. In our paper, we develop a high-frequency approximation to the S* arrival applying a saddle-point method to the integral representing a shear wave reflected from the free surface due to the incidence of a spherical longitudinal wave radiated by the shallow point source. Numerical results presented in the paper show that our high-frequency approximation preserves all basic features of the S* arrival, thereby making it suitable for easy incorporation into any synthetic seismogram computation based on the ray approach.

It has become necessary in seismology to consider more complicated models of the Earth's structure in order to obtain synthetic seismograms that are more consistent with actual field data. Gassmann (1964) and Postma (1955) have presented results dealing with travel-time methods in anisotropic media—in particular, transversely isotropic media. Kinematic properties alone, however, are not enough to conclusively interpret seismic records. Consequently, dynamic properties must be considered producing a need for synthetic seismograms. One of the most efficient methods for obtaining synthetic seismograms is through the use of asymptotic ray theory (Hron and Kanasewich, 1971; Hron, 1973; Hron, Kanasewich and Alpaslan, 1974). A necessary step in the implementation for layered media displaying transverse isotropy is the computation of reflection and transmission coefficients at the interface between two such layers. Reflection coefficients for a free interface and the corresponding surface conversion coefficients must be computed, as well. Theoretical formulas for reflection, transmission, and surface conversion coefficients corresponding to the zero-order approximation of asymptotic theory are presented for the above-mentioned media.