We outline a method of computing synthetic seismograms for stratified, azimuthally anisotropic, viscoelastic earth models. This method is an extended form of the Kennett algorithm that is efficient for multioffset vertical seismic profiling. The model consists of a stack of homogeneous plane layers, and the response is computed iteratively by successive inclusion of deeper layers. In each layer, the 6¿6 system matrix A is diagonalized numerically; this permits treatment of triclinic materials, i.e., those with the lowest possible symmetry. Jacobi iteration is an efficient way to diagonalize A because the entries of A change little from one wavenumber to the next. When the material properties are frequency dependent, the wavenumber loops are inside the frequency loop, and the computation is slow even on a supercomputer. When the material parameters are frequency independent, it is better to make frequency the deepest loop, with diagonalization of A outside the loop, in which case vectorization gives a relatively rapid computation. Temporal wraparound is avoided by making use of complex frequencies, and spatial aliasing is avoided by using a generalized Filon's method to evaluate both the wavenumber integrals. Various methods of generating anisotropic elastic constants from microlayers, cracks, and fractures and joints are discussed. Example computations are given for azimuthally isotropic and azimuthally anisotropic (AA) earth models. Comparison of computations using single and double wavenumber integrations for a realistic AA model shows that single wave number integration often gives incorrect answers especially at near offsets. Errors due to use of a single wavenumber integration are explained heuristically by use of wave fronts diagrams for point and line sources. ¿ American Geophysical Union 1990