Abstract
Computation of a normal incidence synthetic seismogram for a system of vertically stacked homogeneous plane layers can be approached by methods familiar in communication theory (Robinson and Treitel, 1966). In this approach, waves in adjacent layers are related by a 2 X 2 matrix. Waves in nonadjacent layers can be related by applying these elementary matrices recursively--a computational task equivalent to computing the ordered product of the elementary matrices for the intervening layers. In the absence of total internal reflection, the relation of waves in the first and last layers obtained in this way characterizes the propagation of seismic energy in the system. The characterization is complete in the Fourier transform domain when matrix elements are evaluated at N/2 frequencies, where N is the number of interfaces. (If the elements were represented individually, N/2 + 1 frequencies would be required. N/2 is a prorated value based on a more efficient representation.)Analysis of the computational requirements of the standard algorithm (described in the preceding) indicates computer run time scales as N 2 for N large. This paper describes a new algorithm for which the scaling is N(1n 2 N) 2 and, therefore, is faster for large systems.Examples of seismograms computed with the new algorithm are given. Comparative timings on a VAX 11/780 computer indicate the new algorithm has about a one order of magnitude speed advantage when N = 512.
