An Efficient Mesh for Wavefront Construction Methods and its Numerical Properties for Anisotropic Media

The wavefront construction method is an efficient way to model wave propagation in anisotropic media. This method is based on asymptotic ray theory, and it explicitly tracks the propagation of wavefronts through the earth model. The first step is to trace a fan of rays through the earth model, initializing the wavefront by constructing a mesh from the set of all ray points at a specified travel time. The wavefront is thus a mesh composed of quadrilateral cells defined by four neighboring rays. Rays are computed by kinematic ray tracing using Runge-Kutta methods. The basic geometry, or coordinate system, used to select the rays comprising the initial mesh has significant affects on the precision and performance of the numerical calculation. A common approach is to trace rays with regular increments in the azimuthal and declination takeoff angles. These two angles, along with travel time, define a ray coordinate system that is commonly used for implementations of conventional ray tracing for wavefront construction schemes. While the ray coordinate system is straightforward to implement and has been used extensively, it has some drawbacks (Gibson et al., 2002). The most important is that the derivatives of Cartesian coordinates on the ray with respect to the azimuthal takeoff angle vanish near the poles of the coordinate system, leading to potential numerical errors. Our implementation applies paraxial methods, and the poor estimates of derivatives can seriously degrade the performance of the algorithm. Also, specifying an initial ray field with even increments in takeoff angles leads to large concentrations of rays near the poles, which is numerically inefficient. To overcome these important limitations of the ray coordinate system, we apply a new mesh generation algorithm that utilizes a cubic gnomonic mesh. A cubic gnomonic mesh maps points chosen at regular intervals on the surface of a cube surrounding the source point to the focal sphere. In essence, the initial directions of the rays are defined by drawing a line segment from the source point through each initial mesh point on the cubic surface. Then all of the calculations of derivatives for the paraxial method are performed on the geometry defined by this new set of ray coordinates. The ray density is uniform in all directions, and the derivatives never vanish since there are no poles in the new ray coordinate system. We will compare the numerical results of ray coordinate and cubic gnomonic coordinate systems to illustrate the advantages of the new system. We will also describe the implementation of the new scheme and show examples of applications to anisotropic earth models.