## Abstract

The seismic ray theory in anisotropic inhomogeneous media is studied based on non-Euclidean geometry called Finsler geometry. For a two-dimensional ray path, the seismic wavefront in anisotropic media can be geometrically expressed by Finslerian parameters. By using elasticity constants of a real rock, the Finslerian parameters are estimated from a wavefront propagating in the rock. As a result, the anisotropic parameters indicate that the shape of wavefront is expressed not by a circle but by a convex curve called a superellipse. This deviation from the circle as an isotropic wavefront can be characterized by a roughness of wavefront. The roughness parameter of the real rock shows that the shape of the wavefront is expressed by a fractal curve. From an orthogonality of the wavefront and the ray, the seismic wavefront in anisotropic media relates to a fractal structure of the ray path.

## 1. Introduction

In seismological studies, the seismic ray theory for a high-frequency wave is based on Fermat's principle that renders the travel time minimum (Červený 2001; Bóna & Slawiński 2003). Generally, a shortest ray path connected to two different points becomes a straight line, i.e. geometrically, the seismic wave along the ray path propagates in Euclidean space. However, a real seismic ray path is not straight but curved in nature. Geometrically, this situation indicates that the Earth's interior can be regarded as a non-Euclidean space. Therefore, the seismic ray theory should be studied based on the non-Euclidean geometry.

The geometrical properties of seismic ray path are given by a metric tensor. When the seismic wave propagates in inhomogeneous media, the ray velocity depends on its position. In this case, the ray velocity leads a Riemannian metric (Červený 2001; Bóna & Slawiński 2002). Moreover, when the ray velocity depends on both position and direction in the anisotropic inhomogeneous media, this ray velocity leads a metric in a non-Riemannian space called Finsler space (Antonelli *et al*. 1993; Bao *et al*. 2000). Then, the seismic ray path and the seismic elementary wavefront can be interpreted as a geodesic curve and an indicatrix in the Finsler space, respectively. The geometrical aspects of the anisotropic ray theory have been studied (Červený 2001; Antonelli *et al*. 2003*a*,*b*; Bucataru & Slawiński 2005; Stavrinos & Arvanitis 2005; Yajima & Nagahama 2006*a*,*b*, 2007) and a Finslerian Lagrangian (*m*th root metric) has been introduced (Antonelli *et al*. 2003*a*,*b*). However, a relationship between this metric and a seismic wavefront has not been discussed. Therefore, in this study, we classify the shape of a seismic elementary wavefront by a geometrical parameter in the Lagrangian. Moreover, the Finslerian parameter in the *m*th root metric is obtained from the anisotropic wavefront of peridotite. This paper is an extended version of the previous papers (Yajima & Nagahama 2006*a*,*b*, 2007).

In §2, a geometrical viewpoint of the seismic ray in anisotropic inhomogeneous media is introduced, based on the Finsler geometry. Then, in §3, we give a geophysical interpretation of the Finsler metric, *m*th root metric. In §4, we estimate a concrete value of Finslerian parameters in the *m*th root metric from a seismic wavefront propagating in a real rock. Finally, by using the geometrical parameters, a roughness of the anisotropic wavefront is discussed.

## 2. Seismic ray path in anisotropic media and Finsler geometry

Let us regard the Earth's surface and interior as a three-dimensional manifold *M*. The coordinate system on *M* is (*x*^{i})=(*x*^{1}, *x*^{2}, *x*^{3})≡(*x*, *y*, *z*), where *i*=1, 2, 3. *x* and *y* are directed to the Earth's surface and *z* is directed to the Earth's interior. The ray path is a parametric curve of time *t*: *x*^{i}=*x*^{i}(*t*). The direction of the tangent to the ray is denoted by , where . The pair is the coordinate system on the tangent bundle *TM* over *M*. Throughout this paper, Einstein's summation convention is used.

### (a) Finsler geometry and seismic ray path in anisotropic media

The travel time of seismic ray from point *A*(*t*=*t*_{1}) to point *B*(*t*=*t*_{2}) is given by(2.1)where d*s* is a length of the seismic ray path. The fundamental function is homogeneous of degree 1 in (2.2)where the function *V* is called the group velocity function or ray velocity function that denotes the energy flow of the seismic wave (Červený 2001). Then, the Lagrangian for the ray path is defined as ** L**=

*L*

^{2}/2. Geometrical properties of the Earth's interior are given by the metric tensor(2.3)This metric tensor is Finslerian because of the direction dependency of the ray velocity function. In order to require the regularity of the metric tensor, we consider only P-wave. The metric tensor of the P-wave is regular (Bucataru & Slawiński 2005). Therefore, the metric tensor of the P-wave has the inverse matrix. With respect to the ray velocity in anisotropic inhomogeneous media, the metric tensor depends on the direction of the ray. Therefore, the Earth's interior where the seismic wave propagates can be regarded as a Finsler space

*F*=(

*M*,

*L*). On the other hand, for isotropic media, the ray velocity depends only on the position

*x*

^{i}and so the metric tensor (2.3) is a function of the position or a constant. In this case, the ray path is a parametric curve in the Riemannian or the Euclidean space.

The ray is obtained by using a variational principle called Fermat's principle that renders the travel time minimum. Three-dimensional Euler–Lagrange equation results from the variation of the energy function *L*^{2} that fixes the parametrization of the geodesics(2.4)where *E*_{i} is the Euler's vector in the Finsler space (Kawaguchi 1962). Geometrically, the Euler–Lagrange equation is equivalent to the geodesic equation on (Bucataru & Miron 2007)(2.5)where(2.6)where *g*^{ik} are components of inverse of the metric tensor *g*_{ik}. In this case, the solution of the second-order differential equations (2.5) can be interpreted as the unified parametric curve on the six-dimensional tangent bundle *TM*, where *I*=1, …, 6, *i*=1, 2, 3 and *i*^{*}=4, 5, 6. Then, the adapted frame on the tangent bundle *TM* is expressed by (Ikeda 1993)(2.7)(2.8)where we put(2.9)(2.10) is called the nonlinear connection that implies an interaction between the (*x*)- and -fields. Geometrically, the nonlinear connection on *TM* prescribes the direct sum structure of the horizontal distribution *H*_{u}*TM* and the vertical distribution *V*_{u}*TM* in the tangent space *T*_{u}*TM* for *u*∈*TM*: *T*_{u}*TM*=*H*_{u}*TM*⊕*V*_{u}*TM*. With respect to the adapted frames (2.7) and (2.8), there are two types of connection coefficients of a connection on *TM*, where and are the horizontal and vertical coefficients of the connection, respectively (Bucataru & Miron 2007). These connections preserve the horizontal and vertical distributions and act similarly on both. Then, the covariant derivative *D* for an arbitrary vector *V*^{i} is defined as (Ikeda 1993)(2.11)Moreover, since we consider the ray path in the framework of Finsler space, the relationship holds (e.g. Antonelli *et al*. 1993). Then, the seismic ray path *c*(*t*) projected from a total curve on the tangent bundle is given by a solution curve of the geodesic equation in the Finsler space (Antonelli & Bucataru 2003; Bucataru & Miron 2007)(2.12)Therefore, from the tangent bundle viewpoint, the seismic ray path *c*(*t*)=(*x*^{i}(*t*)) is a three-dimensional projected curve from the tangent bundle into the base manifold *M*.

## 3. Hamiltonian and Lagrangian for ray path in anisotropic media

In this section, we derive a concrete expression of the Lagrangian for the seismic ray path in anisotropic inhomogeneous media. At first, in the following, we consider the case that the continuum media are regarded as the factorized anisotropic inhomogeneous media (Červený 2001), i.e. the ray velocity function is expressed by a separable form: . Then, the function *L* is rewritten as(3.1)where we put and the Finslerian metric tensor for the function is denoted by .

### (a) Fundamental function of the Hamiltonian for seismic ray path

Let us derive the anisotropic part in the fundamental function of the Lagrangian (3.1). The propagation of the body wave propagating in minerals is governed by the elastodynamic equation (Musgrave 1970; Antonelli *et al*. 2003*b*)(3.2)where *ρ* is the mass density at a point *x*^{j}, *u*^{i} is a vector describing the displacement of the continuum, and *σ*^{ij} is the stress tensor. By the context of linear elasticity, the stress tensor is given by(3.3)where *c*^{ijkl} are the elasticity constants at a point *x*^{j} and is the strain tensor.

Let us consider the solution of the elastodynamics equation (Červený 2001)(3.4)where *A*^{i} is the wavefront amplitude, the function *F* represents a high-frequency signal of wave, and level sets of the function *ψ* imply the wavefront. Then, substituting (3.3) and (3.4) into (3.2), one can obtain the equation(3.5)where *P*_{i}≡∂*ψ*/∂*x*^{i} is called the slowness vector. For the anisotropic media, the gradient of the function *ψ* is perpendicular to the wavefront in the Hamiltonian sense (Bucataru & Slawiński 2005): *g*_{ij}(∇_{g}*ψ*)^{i} d*x*^{j}/d*t*=0, where *g*_{ij} is the Finslerian metric tensor and the gradient is (∇_{g}*ψ*)^{i}=*g*^{ij}∂*ψ*/∂*x*^{j}. In this case, the coordinates on the cotangent bundle *T*^{*}*M* are denoted by (*x*^{i}, *p*_{i})=(*x*, *y*, *z*, *p*_{x}, *p*_{y}, *p*_{z}).

For convenience, let us rewrite equation (3.5) as(3.6)where we put the matrix *Γ*^{ik} called the Christoffel matrix as(3.7)Then, the condition for existence of non-trivial solutions for equation (3.6) is called the Christoffel equation (Musgrave 1970; Červený 2001)(3.8)Solving this eigenvalue equation, one can obtain three velocity functions that give different types of seismic body waves in anisotropic media (Červený 2001). The wave with the highest phase velocity is called the quasi-compressional wave (qP-wave). The remaining two waves are quasi-shear waves (qS1- and qS2-waves). We denote the phase velocity of qP-, qS1- and qS2-waves as *v*_{qP}, *v*_{qS1} and *v*_{qS2}, respectively.

Next, we consider the components of the elasticity tensor. The real seismic wave propagates in a rock aggregated of minerals. The mineral's physical properties are influenced by the crystal structures. Therefore, we study the crystal systems of rock mineral in linear elasticity theory. Generally, there are eight symmetry classes of elasticity tensors for minerals or crystals: triclinic; monoclinic; orthorhombic; tetragonal; trigonal; hexagonal; cubic; and isotropic systems (Böhlke & Brüggemann 2001). Mathematically, the symmetries of structures in three dimensions have discrete symmetric groups. Moreover, in linear elasticity, not all the symmetry classes have discrete symmetry groups, for example transversely isotropic continua have symmetry groups isomorphic with O(2). In the following, we shall denote the components of the elastic tensor of the crystal by capital letters *C*^{IJ}≔*c*^{ijkl}, where *I*, *J*=1, …, 6 and *i*, *j*, *k*, *l*=1, 2, 3 and we put *I*≔*ij* and *J*≔*kl* as follows (Červený 2001):(3.9)Here, we consider the case that the indices 1, 2, 3 of the elastic tensor *c*^{ijkl} are coincident with the indices of coordinates of the crystal axis, (*X*^{1}, *X*^{2}, *X*^{3}).

Especially, the simplest anisotropic case of broad geophysical applicability is called a transverse isotropy or hexagonal system (Thomsen 1986). This has one distinct direction (usually, but not always, vertical), while the other two directions are equivalent to each other. Therefore, we assume the seismic propagation in the (*x*^{1}, *x*^{3})-plane, i.e. we consider two-dimensional manifold *M* that corresponds to the Earth's section and the coordinates are rewritten as(3.10)(3.11)(3.12)where the *y*-axis is equivalent to the *z*-axis because of the transversely isotropic medium and we put *y*≡*z*, and *p*_{y}≡*p*_{z}. With this simplification, the Christoffel equation (3.8) can be solved for the hexagonal crystal system. For the elastic constants of the hexagonal system, three phase velocity functions can be obtained (Daley & Hron 1977; Thomsen 1986)(3.13)(3.14)(3.15)where the phase angle of the phase velocity direction *ϕ* is defined by tan *ϕ*≡*p*_{x}/*p*_{z} and we put(3.16)These phase velocity functions *v*_{qP}, *v*_{qS1} and *v*_{qS2} describe a wave surface called the slowness surface of qP-, qS1- and qS2-waves in the hexagonal system, respectively. In particular, convenient anisotropic parameters *ϵ* and *η* have been introduced (Thomsen 1986)(3.17)(3.18)In a weakly anisotropic media, the anisotropy parameters *ϵ* and *η* satisfy |*ϵ*|≪1 and |*η*|≪1 (Thomsen 1986). Then, the phase velocity of the qP-wave for the hexagonal media (3.13) can be rewritten by a simple expression(3.19)where .

Next, let us consider a fundamental function , which generates the Hamiltonian defined on . The slowness surface of the qP-wave at fixed point *x*_{0} is represented by the end of radius vector *r* with length *v*_{qP}. Then, we put the polar coordinates on the cotangent space (3.20)The case *p*_{x}=*p*_{z}=0 does not appear since the framework is the slashed tangent bundle. Therefore, substituting (3.20) into (3.13), we obtain the following expression:(3.21)The function *D*(*p*_{x}, *p*_{z}) is defined as(3.22)where *κ*=(*C*^{11}−*C*^{44})^{2}, and *μ*=(*C*^{33}−*C*^{44})^{2}. Then, according to the Okubo method (e.g. Antonelli *et al*. 1993; Bao *et al*. 2000), we replace the unit vector *p*_{i} by . As a result, the fundamental function of qP-wave in the hexagonal media can be obtained as(3.23)Consequently, we obtain the result.

*The fundamental function of Hamiltonian for the qP*-*wave propagating in minerals that have a hexagonal crystal system is given by the function* (3.24)*where, for the hexagonal crystal system, the coefficients a*^{ij} *are a*^{11}=*C*^{11}+*C*^{44}, *a*^{13}=*a*^{31}=0, *a*^{33}=*C*^{33}+*C*^{44} and *a*^{i2}=*a*^{2i}=0. *The function D*(*p*_{h}) *is expressed by* (*3.22*). *Then*, *the Hamiltonian of qP*-*wave is given by* .

Moreover, in case of the weakly anisotropic media, we obtain the following.

*The fundamental function of qP*-*wave for the weakly anisotropic media is expressed by*(3.25)

In case of other crystal systems except for triclinic and monoclinic systems, the phase velocity functions of the qP-wave have the similar expressions of hexagonal system (Musgrave 1970). Therefore, the fundamental function of the qP-wave for orthorhombic, tetragonal, trigonal and cubic systems can be also given by the same way as the hexagonal system.

The fundamental functions (3.24) and (3.25) are homogeneous of degree one in *p*_{i}. Therefore, the fundamental function can be rewritten as(3.26)where the non-degenerate metric tensor *h*^{ij} satisfies the relationships(3.27)The metric tensor *h*^{ij} for Hamiltonian (3.24) and (3.25) is Finslerian because the direction is dependent on the metric tensor. Moreover, the fundamental function of the qS1- and qS2-waves can be obtained similarly. However, the Hamiltonian metric of the qS2-wave is singular when . This implies that the Legendre transformation may not be well defined somewhere and we cannot obtain the Lagrangian from the Hamiltonian. Therefore, for convenience, we consider only the qP-wave in the following.

### (b) Fundamental function of the Lagrangian for seismic ray path

The Lagrangian for the seismic ray path can be obtained from the Hamiltonian. Let us consider the smooth map called the Legendre transformation (Hrimiuc & Shimada 1996): , where the map is locally denoted by(3.28)Then, the Lagrangian is given by(3.29)In our case, the qP-wave is only considered and so the Hamiltonian is regular. Therefore, the Lagrangian can be obtained from the Hamiltonian (3.26) by the Legendre transformation. Hence, the anisotropic part of the Lagrangian for the crystal systems is given by , where *h*_{ij} is the inverse matrix of *h*^{ij} and *x*^{i} is a holonomic coordinate. The metric tensor *h*_{ij} is Finslerian because the metric tensor depends on the direction of the ray.

From the fundamental function , the wavefront at a fixed point *x*_{0} is expressed by an indicatrix on the tangent space : , i.e. . Here, the ray direction vector in the holonomic coordinate system represents distances along axes in the tangent space . Then, the axes of this equation of indicatrix can be transformed into different axes using relationships for transformation of axes: , where is the cosine of the angle between the old and new axes. Moreover, a suitable transformation of the axes can remove the cross product terms in , i.e. except for the monoclinic and triclinic symmetries, the matrix (*h*_{ij}) can be diagonalizable (Lovett 1989). Therefore, the fundamental function of the Lagrangian for crystal systems can be rewritten by the canonical form(3.30)where and are principal coefficients of the matrix (*h*_{ij}) and the ray direction vector is directed to the principal axis.

Next, let us consider the fundamental function of the Lagrangian for rocks based on the functions for crystal systems . The texture of rocks as polycrystal is expressed by the internal geometrical structure (e.g. the preferred orientation of crystals in deformed rocks). In order to express the alignment of orientated crystals, we apply an orientation parameter called the scalar order parameter *S* (de Gennes 1974)(3.31)where the bracket 〈 〉 denotes the statistical average of alignment and *ξ* is the angle between the primary axis of each crystals and the direction of average orientation. Especially, in the isotropic case, 〈cos^{2}*ξ*〉=1/3 and the parameter *S* is zero. By using the orientation parameter, the indicatrix of rocks as polycrystals can be expressed by a generalization of indicatrix of single crystal, i.e. the indicatrix of rocks is a convex curve called a superellipse (Kindlmann 2004; Jankun-Kelly & Mehta 2006) and can be given by the following form:(3.32)where the coefficients *α* and *β* express the anisotropy of the preferred orientation of aligned minerals along the Earth's surface and the Earth's interior, respectively. The parameter *m* is connected to the statistical orientation parameter *S* (Jankun-Kelly & Mehta 2006): *m*=2(1−*S*)^{−γ}. The parameter *γ* can express the sharpness of the superellipse, i.e. the shape becomes sharp for *γ*>1 and the shape is smooth for *γ*=1 (Kindlmann 2004; Jankun-Kelly & Mehta 2006). In our case, the seismic wavefront is smooth. Therefore, we put *γ*=1. From this expression of the indicatrix, the fundamental function of the Lagrangian for crystal systems is changed into the generalized function for rocks(3.33)where we put the coefficients (3.34)The coefficients *α* and *β* give the axial anisotropy. Geometrically, the fundamental function (3.33) is called the *m*th root metric (Matsumoto & Shimada 1978; Antonelli & Shimada 1991; Antonelli *et al*. 1993). The metric tensor given by the *m*th root metric is Finslerian because the metric tensor depends on the direction of the ray. Consequently, from (3.1) and (3.33), we have the following result.

*The Lagrangian for the seismic ray in factorized anisotropic inhomogeneous media is given by a generalization of the Lagrangian for minerals into that for rocks*. *In this case*, *the geometric properties of the ray and the wavefront are given by the Finsler metric function L*(3.35)*where the coefficients* *are defined by* (*3.34*). *The parameter m expresses the statistical orientation of wavefront on rock textures*. *According to the paper* (*Antonelli et al. 2003b*), *this metric is called the seismic Finsler metric*.

Similar to the crystals, the typical geometrical symmetry for rock media in the Earth can be classified into five types, i.e. axial, spherical, orthorhombic, monoclinic and triclinic symmetries (Paterson & Weiss 1961). On the other hand, the Lagrangian is symmetric by exchanging for when *α*=*β*. In other words, the polar representation of this Lagrangian is invariant of finite rotational group action, for example, cyclic and dihedral symmetries (Gielis 2003). Therefore, except for the monoclinic and triclinic systems, the symmetry of oriented crystals in rocks (Paterson & Weiss 1961) can hold for this Lagrangian.

## 4. Discussions

### (a) Indicatrix of the seismic Finsler metric

The nature of wave propagation is explained by Huygens' principle, i.e. each point on a wavefront can be treated as a source point of a secondary wavelet (elementary wavefront). The envelope of these elementary wavefronts constructs a new wavefront (Guenther 1990). In our case, the source point is not a point of minerals but a domain of rocks as a collection of minerals.

Geometrically, the shape of the wavefront at a fixed point *x*_{0} is expressed by the indicatrix . For the seismic Finsler metric (3.35), the anisotropy of ray velocity is given by the function (4.1)where *θ* is the ray angle defined by . Then, in the case of *m*=2 and *α*=*β*, the indicatrix is a circle. In this isotropic case, the phase velocity *v* is also a function of the position and so the slowness surface is coincident with the group velocity surface. On the other hand, for the parameter *m*>2, the indicatrix describes a convex curve. Therefore, the shape of the seismic wavefront is anisotropic. Hence, the anisotropy of the seismic wavefront can be expressed by the parameter *m* in seismic Finsler metric (3.35).

### (b) Geometric parameters in the seismic Finsler metric and wavefront of anisotropic rock

Let us discuss a relationship between the anisotropic parameters *m*, *α* and *β* in the seismic Finsler metric (3.35) and the wave surface of anisotropic rocks. The seismic qP-wave propagates globally in the Earth's upper mantle whose rock type is thought to be a peridotite. The peridotite is aggregated of orthorhombic olivine that shows a lattice preferred orientation with strong maxima of a crystallographic *b*-axis close to a pole of the foliation (Barroul & Kern 1996). In this case, the peridotite has a transverse isotropic structure, i.e. when we define a Cartesian coordinate (*X*^{1}, *X*^{2}, *X*^{3}) in the peridotite, the *b*-axis is an axis of rotational symmetry and the other axes (*a*- and *c*-axes) are on the (*X*^{1}, *X*^{2})-plane that is normal to this symmetric axis. Therefore, the continuum exhibits isotropic properties on this plane. As a result, we can consider only the two-dimensional case, i.e. we can put *X*^{1}≡*X*^{2} and the spatial *X*^{3}-axis is directed to the *b*-axis of the peridotite. Therefore, it is possible to study the wavefront of ray velocity in a plane that contains the axis of rotational symmetry.

The elastic constants for the peridotite are obtained based on the elastic constants for the olivine. According to the paper (Bunge 1982), we transform the elastic constants for the olivine as a single crystal (Kumazawa & Anderson 1969) into the elastic constants for the peridotite with a preferred orientation of *b*-axis(4.2)where the density of olivine is given by *ρ*=3.311×10^{3} kg m^{−3}. The phase velocity *v* of the seismic qP-wave for the transverse isotropic media can be obtained by equation (3.13). Then, a relationship between the ray velocity *V* and the phase velocity *v* in the two-dimensional plane is given by (Berrymann 1979; Thomsen 1986)(4.3)where the ray angle *θ* can be computed from the phase angle *ϕ*(4.4)Substituting the elastic constants (4.2) into the function of the phase velocity (3.13), the wavefront of peridotite is obtained from (4.3) and (4.4) (figure 1). In figure 1, the crystal *b*-axis is parallel to direction of depth, *z*-axis. Then, velocity of the slow qP-wave is directed to the olivine *b*-axis. On the other hand, the seismic Finsler metric (3.35) is divided into the position *f*(*x*^{i}) and the direction part independently. Therefore, we consider the indicatrix described by the direction part: . By comparing the wavefront obtained from (4.3) with the indicatrix of the seismic Finsler metric, the *m*-value, *α* and *β* in seismic Finsler metric are *m*=2.1, *α*=8.85 km s^{−1} and *β*=7.60 km s^{−1}, respectively. Thus, the shape of wavefront in peridotite is not elliptic but a superellipse.

### (c) Roughness of the wave surface, ray path and Finsler parameter

Geometrically, the indicatrix corresponds to a Minkowski norm defined over each tangent space (Antonelli *et al*. 1993). For the *m*th root metric, the Minkowski norm is the modified *L*_{p}-norm in the tangent space (Ait-Haddou *et al*. 2000). For *m*=2, the Minkowski norm is given by the standard Euclidean norm. Moreover, when *m* takes real value, the *m*th powered *L*_{p} norm is the Hausdorff measure of tangent space, i.e. the dimension of the indicatrix is measured by the real-valued *m*. In this case, the shape of the indicatrix relates to a fractal structure (Jonot 1991). Therefore, the *m*th root metric with non-integer value *m* means that a measurement in the tangent space is fractal. For the peridotite, the value of parameter *m* takes *m*=2.1 which is slightly different from 2. This non-integer value *m* shows a deviation from an elliptic wavefront whose length is measured by the Euclidean measurement.

Let us discuss the non-integer value *m* as a deviation from the shape of the isotropic wavefront to that of the anisotropic wavefront, i.e. the indicatrix for the anisotropic rocks is smooth only in Finsler space, whereas it is not smooth in the Riemannian or Euclidean space. In order to discuss the deviation, let us introduce the roughness parameter (Miklashevich 2003) that is determined by a ratio of the Riemannian arc length of the indicatrix to the Euclidean arc length of the indicatrix *l*.

In the theory of Finsler geometry, one has (Bao *et al*. 2000) the following.

*The indicatrix length in the two-dimensional Minkowski plane can be computed as*(4.5)*where* *and U is the standard unit circle*(4.6)

In our case, the Riemannian arc length of the indicatrix is determined by the anisotropic part of the Lagrangian for the rocks in (3.35). Then, by definition (2.3), the metric tensor is Finslerian and the determinant is(4.7)Then, substituting (3.33) and (4.7) into (4.5) with the variable transformation , we can obtain the length of indicatrix for the anisotropic rocks as(4.8)The fundamental function of the isotropic media is given by *m*=2 and *α*=*β* in the Lagrangian for rocks (3.35). Then, similar to the anisotropic case, we obtain the length of the indicatrix for the isotropic media(4.9)For the peridotite, the values of *m*, *α* and *β* are *m*=2.1, *α*=8.85 km s^{−1} and *β*=7.60 km s^{−1}, respectively. Therefore, substituting these values into the integral (4.8), the ratio of the length of indicatrix to the length of indicatrix *l* can be obtained. As a result, the roughness of the wave surface is(4.10)This value shows that the wavefront of the peridotite is slightly rough. This deviation from the isotropic wavefront is caused by the small difference between the isotropic value *m*=2 and the anisotropic value *m*=2.1. Therefore, the anisotropy of media where the seismic wave propagates can be represented by the roughness of wave surface.

From a viewpoint of the fractal geometry, a roughness of surfaces is related to the fractal dimension (Mandelbrot 1982). Therefore, the value *R*=1.10 indicates that the shape of wavefront of peridotite is expressed as a fractal curve such as the Koch curve. Moreover, by the Legendre transformation of the seismic Finsler metric, the fundamental function of Hamiltonian *H* in the cotangent space gives a wavefront (slowness surface): *H*=1. In this case, the shape of slowness surface is also expressed by a fractal curve. Generally, a direction of ray path is orthogonal to the slowness surface (Chapman 2004; Bucataru & Slawiński 2005). Therefore, the orthogonality gives a fractal structure of the seismic ray path in anisotropic media. This confirms that a fractal measurement of ray trajectory in anisotropic media such as crystals is expressed by Finslerian distance, and is concordant with the previous remark in the paper (Argyris *et al*. 2000).

## 5. Conclusion

In this paper, the seismic ray theory is studied in terms of differential geometry. When the seismic wave propagates in anisotropic inhomogeneous media, the ray velocity is a function of direction. In this case, the ray path is a geodesic curve in Finsler space. The Finslerian Lagrangian for a two-dimensional ray path is given by a Lagrangian for rocks, which is a generalization of Lagrangian for minerals. Especially, the anisotropic parameters in the Finslerian Lagrangian are estimated from the wavefront of real rocks. With respect to the peridotite, the anisotropic parameters indicate that the shape of the wavefront is expressed by a superellipse. This deviation from isotropic wavefront is characterized by the roughness of the wavefront. The roughness parameter of peridotite shows that the shape of the wavefront is similar to the fractal curve. Moreover, by the orthogonality of the wavefront and the ray path, the wavefront relates to a fractal structure of the seismic ray path.

## Acknowledgements

One of the authors (T.Y.) was financially supported by research fellowships of the Japan Society for the Promotion of Science for Young Scientists and is now supported by the Global Centers of Excellence (COE) Program: ‘Global Education and Research Center for Earth and Planetary Dynamics’ of Tohoku University.

## Footnotes

- Received November 5, 2008.
- Accepted February 10, 2009.

- © 2009 The Royal Society