## Abstract

Using analytical and numerical methods, we analyse the Raj–Ashby bicrystal model of diffusionally accommodated grain-boundary sliding for finite interface slopes. Two perfectly elastic layers of finite thickness are separated by a given fixed spatially periodic interface. Dissipation occurs by time-periodic shearing of the viscous interfacial region, and by time-periodic grain-boundary diffusion. Although two time scales govern these processes, of particular interest is the characteristic time *t*_{D} for grain-boundary diffusion to occur over distances of order of the grain size. For seismic frequencies *ωt*_{D}≫1, we find that the spectrum of mechanical loss *Q*^{−1} is controlled by the local stress field near corners. For a simple piecewise linear interface having identical corners, this localization leads to a simple asymptotic form for the loss spectrum: for *ωt*_{D}≫1, *Q*^{−1}∼const.*ω*^{−α}. The positive exponent *α* is determined by the structure of the stress field near the corners, but depends both on the angle subtended by the corner and on the orientation of the interface; the value of *α* for a sawtooth interface having 120^{°} angles differs from that for a truncated sawtooth interface whose corners subtend the same 120^{°} angle. When corners on an interface are not all identical, the behaviour is even more complex. Our analysis suggests that the loss spectrum of a finely grained solid results from volume averaging of the dissipation occurring in the neighbourhood of a randomly oriented three-dimensional network of grain boundaries and edges.

## 1. Introduction

Motivated by the problem of seismic attenuation, the mechanical loss spectrum of fine-grained mantle minerals has been measured at high temperatures in forced torsional–oscillation experiments (e.g. Gribb & Cooper 1998; Jackson *et al.* 2004; Sundberg & Cooper 2010). According to these experiments, within the seismic frequency range, the mechanical loss ℒ (inverse of the quality factor *Q*) varies with angular frequency *ω* according to a power law: ℒ∝*ω*^{−α} with 0.2<*α*<0.35. Gribb & Cooper (1998, §4) summarize the experimental evidence supporting diffusionally accommodated grain-boundary sliding as the attenuation mechanism in these experiments. To their considerations, we can add the following argument of dynamical similarity (Morris & Jackson 2009*a*): values of ℒ measured as a function of frequency for different grain sizes and temperatures define a single curve when graphed against the dimensionless frequency *ωη*′/*μ*. Here *μ* is the grain rigidity, and *η*′ is the steady-state viscosity for Coble creep measured in the independent uniaxial compression tests of Faul & Jackson (2007). Because *η*′ is controlled by grain-boundary diffusion, it follows that ℒ is also. (As discussed by Gribb & Cooper (1998), the experiments are designed to eliminate dislocation damping; grain sizes are kept sufficiently small that, within individual grains, dislocation numbers are negligible for the experimental levels of shear stress.)

Although the experimental results can be fitted by spring–dashpot models containing a sufficient number of elements, Gribb & Cooper (1998) and Cooper (2002) argue that the power-law form of the spectrum can be explained more physically by accounting for the spatial variation of stress within grains. To test that explanation, they use the bicrystal model of grain-boundary sliding described by Raj (1975); two isotropic Hookean layers are separated by a fixed prescribed non-plane interface upon which the shear stress vanishes, and across which the normal velocity is discontinuous owing to (grain boundary) diffusion along the interface. To determine the loss spectrum for that model, Gribb & Cooper (1998) solve the initial-value problem determining the response for a step change in applied stress, and then find the loss spectrum by Laplace transformation. As shown in their fig. 10, the agreement between theory and experiment supports their explanation of the observed power-law spectrum.

That agreement is called into question by recent work. Morris & Jackson (2009*b*) repeat that calculation using the same assumption of infinitesimal grain-boundary topography. The new solution differs from the earlier one in two essentials: the loss spectrum is now obtained directly by imposing a sinusoidally varying boundary stress; and an explicit asymptotic form is obtained giving ℒ at high frequencies, i.e. for *ωη*′/*μ*≫1. According to the new solution, ℒ decreases much more slowly than Gribb & Cooper predicted; it does not even follow a power law, but instead decays inversely with the logarithm of frequency. Thus, although experimental evidence points to grain-boundary sliding as the explanation for the high-temperature attenuation background, detailed analysis of the simplest (bicrystal) model predicts a spectrum that is qualitatively different from that observed.

Despite this result, we argue here that useful lessons can still be drawn from the bicrystal model, provided the effect of finite interface slopes is included. Using analytical and numerical methods, we show that for the sawtooth and truncated sawtooth interfaces (modes 1 and 2 sliding surfaces of Raj & Ashby (1971)), the bicrystal model does indeed predict a power-law spectrum when the slope is finite: for *ωη*′/*μ*≫1, the mechanical loss ℒ∝(*ωη*′/*μ*)^{−α} at high frequencies. Although the exponent *α* is uniquely determined by the angle subtended by the corner on these piecewise linear interfaces, the constant of proportionality in the loss relation depends on the orientation of the interface. In the limit of vanishing slope, the new result is consistent with the scaling found by Morris & Jackson (2009*b*). Further, using a model problem, we argue that, at high frequencies, the total dissipation rate within the sample is determined as the sum of contributions from each corner on the interface. The magnitude of individual contributions depends on the angle subtended by the corner, and on the stress amplitude at the corner. Because the latter proves to depend on the orientation of the interface, so too does the dissipation.

Although we establish this result for the bicrystal model, we expect a similar result to apply in a three-dimensional sample. This result is, of course, a refinement of the Gribb & Cooper explanation. It provides a definite picture of where dissipation is occurring, and implies that it is fruitless to seek a simple theory making quantitative predictions of ℒ for a three-dimensional sample. For, although the dissipation is localized, its magnitude and scaling with frequency depends on the geometry of grains and corners. Cross sections of experimental samples (e.g. Barnhoorn *et al.* 2007, fig. 1*e*) do not resemble that of a regular hexagonal array, and the three-dimensional geometry is likely to be even more complex. Consequently, to be physically instructive, any future studies of diffusionally accommodated grain-boundary sliding in regular hexagonal arrays would need to be carefully motivated. In addition to providing this (negative) guide to model building, our analysis also provides quantitative results suitable for testing numerical studies.

Following the statement of the boundary-value problem (b.v.p.) in §2, and the outline of our numerical method in §3, we use scaling in §4 to show that, for a sawtooth interface, ℒ∝(*ωη*′/*μ*)^{−α} for *ωη*′/*μ*≫1. In particular, equation (4.4) gives the formula relating the power-law exponent *α* to the angle subtended by the corner on the sawtooth interface. In §5, we show that our numerical results agree quantitatively with that power law. Because cross sections of experimental samples show a range of corner angles, some corresponding to triple junctions, with others apparently corresponding to kinks in the grain boundary, we then consider an interface having two different corner angles. According to our numerical solutions, the slope of the loss spectrum then decreases gradually with increasing frequency; at high frequency, the behaviour of the mechanical loss appears to be controlled by the corner having the strongest singular stress behaviour. Consequently, one should not expect a single power law to fit the entire range of experimental frequencies. This result might account for the range of *α*-values found in experimental studies. In §6, we summarize our chief results and conclusions.

Throughout this work, dimensional variables are denoted by asterisks.

## 2. Boundary-value problem

Figure 1 shows the geometry of the bicrystal model. Two perfectly elastic grains with rigidity *μ* and Poisson ratio *ν* are separated by an interface . The interface is periodic with a wavelength 2*π*/*ξ*, where *ξ* is the wavenumber. Because samples in the attenuation experiments are subjected to small strains *O*(10^{−6}) (Jackson *et al.* 2004), we assume the interface position to be time-independent, given by a function *f**(*x**). Unit vectors in the coordinate directions are denoted by and ; unit tangent and unit normal vectors of the interface are denoted by and , respectively. Along the upper and lower boundaries at *y**=±*a*/*ξ*, the imposed displacement varies sinusoidally in time with angular frequency *ω** and amplitude *U*_{0}, i.e. . The grains are assumed to be undergoing plane deformation and the *x* and *y* components of the displacement vector **u*** are denoted by *u**(*x**,*y**) and *v**(*x**,*y**), respectively. Similarly, the Cartesian components of the stress and strain tensors are denoted by *σ**_{ij}(*x**,*y**) and *e**_{ij}(*x**,*y**), respectively.

On the grain interface **𝒮**_{I}, we impose the following constitutive equations:
2.1a
and
2.1b
The parameters ℓ, *η*, *v*, *D*, *k* and *T* denote boundary thickness, boundary viscosity, molecular volume, grain-boundary diffusivity, Boltzmann constant and temperature, respectively. Equation (2.1a) states that the shear stress exerted across 𝒮_{I} is proportional to the discontinuity in tangential velocity across 𝒮_{I}. As described by Raj & Ashby (1971), the thin disordered boundary phase acts as if it is a liquid film of uniform viscosity *η* and constant thickness ℓ. We may note that the steady-state creep viscosity *η*′ is the manifestation of grain-scale diffusion, whereas the boundary viscosity *η* is, at least for a high-angle boundary, a manifestation of diffusion at the scale ℓ (Ashby 1972, p. 511). Equation (2.1b) is obtained by combining Fick’s Law with interfacial mass balance. The volumetric flow rate *j** (per unit *z*-length) along the interface due to grain-boundary diffusion is related to the normal stress by *j**=(*v*ℓ*D*/*kT*)(d*σ**_{nn}/d*s**) in a form analogous to Fick’s Law (Lifshitz 1963; Raj & Ashby 1971). Using that definition of the volumetric flow rate and invoking interfacial mass balance
2.2
leads to the second constitutive equation. According to that equation (2.1b), mass flows along the interface from regions under compression to regions in tension.

Following Mosher & Raj (1974) and Raj (1975), we define the sliding time scale *t*_{η} and the diffusive time scale *t*_{D} by
2.3a,b

Physically, *t*_{η} and *t*_{D} are, respectively, the time scales on which the two sides of equations (2.1a) and (2.1b) balance, if derivatives along the interface scale with its wavelength. We note that, if we identify 2*π*/*ξ* with the grain dimension *d*, the time scale *t*_{D} is within a factor of 2 of the Maxwell time *η*′/*μ* based on the Coble creep viscosity. According to Morris & Jackson (2009*a*, fig. 3), experiments conducted in the seismic frequency range lie within 0.1<*ω*_{*}*η*′/*μ*<10^{8}, with most cases occurring at *ω*_{*}*η*′/*μ*≫1. Consequently, though our numerical results will cover the whole range of dimensionless frequencies, the limiting behaviour at large dimensionless frequencies is of particular interest. At those high frequencies, matter can diffuse along the grain boundary only over a short distance compared with the grain size, before the time-oscillatory stress reverses. Balancing terms in equation (2.1b), we find that matter diffuses over a distance of order of the diffusion length defined as follows:
2.4
From the identity ℓ_{D}*ξ*=1/(*ω*_{*}*t*_{D})^{1/3}, it follows that, for *ω*_{*}*t*_{D}≫1, ℓ_{D}≪*d*, as claimed.

Dimensionless variables (without asterisks) are defined as follows:
2.5a
2.5b
2.5c
2.5d
2.5e
and
2.5f
In equation (2.5c), *ε* is the characteristic slope of the interface.

The dimensionless b.v.p. is as follows: in grains 1 and 2,
2.6a
on *y*=±*a*,
2.6b
2.6c
on *y*=*εf*(*x*),
2.6d
2.6e
2.6f,g
on *x*=2*π* and *x*=0,
2.6h
2.6i
In equation (2.6d), we define the viscosity parameter
2.7
When (fixed frequency), the interface becomes effectively inviscid, i.e. *σ*_{ns}=0. See Morris & Jackson (2009*b*) for further discussion of the b.v.p.

Problem (2.6) is linear and time-separable because the interface is fixed. Consequently, the solution of problem (2.6) for a time-periodic boundary displacement is also time-periodic with the same angular frequency *ω*.

Letting
2.8
and *γ*(*t*)=e^{iωt}/*a*, we define the sample shear modulus by
2.9
We see that *G* is independent of *t*, because, on the right side of (2.9), both the numerator and denominator vary as e^{iωt}. Indeed, to within a factor of *a*, *G* is equal to the Fourier transform of *τ*(*t*).

The mechanical loss ℒ is defined, as usual, by the equation
2.10
If the material can be modelled as a network of springs and dampers, the quantity defined in equation (2.10) is equal to the ratio of the loss per cycle to 4*π* times the mean strain energy stored within the grains (Bland 1960; O’Connell & Budiansky 1978).

We initially consider the two types of interface illustrated in figure 2. These interfaces can be represented using piecewise linear functions defined by
2.11
where the specific values *m*=1/2 and *m*=1/4 correspond to types *S* and *TS* interface, respectively. To relate the characteristic slope *ε* to the interface slope angle *φ*, we use equations (2.12*a*,*b*) for types *S* and *TS* interface, respectively.
2.12
These interfaces are found in a regular array of hexagonal grains. In that array, the slope angles for types *S* and *TS* interfaces have values *φ*=30^{°} and *φ*=60^{°}, respectively.

## 3. Numerical method

Solving b.v.p. (2.6) directly using the conventional finite-element method is challenging because boundary condition (2.6e) requires approximation of the second derivative of normal stress *d*^{2}*σ*_{nn}/d*s*^{2}. As a result of the stress concentration given in equation (4.2), numerical approximation of the term *d*^{2}*σ*_{nn}/d*s*^{2} will incur a large numerical error and requires an excessively fine mesh near the corners.

To circumvent that difficulty, we use the following method, based on that of Sethian & Wilkening (2003). Using the principle of superposition, we decompose problem (2.6) into two separate b.v.p.’s. By doing so, we can recast the original two-dimensional b.v.p. into a one-dimensional partial differential equation (PDE) defined along the interface . That PDE is defined by a composite operator embedded with a spatial differential operator originating from equation (2.6e). B.v.p. (2.6) is solved if the eigenvalues and the eigenfunctions of the composite operator are found. To avoid calculating stress derivatives, the eigenvalues and the eigenfunctions are found indirectly using a constructed ‘pseudo-inverse’ of the composite operator. The solution procedure is described in appendix A and the details are given in Sethian & Wilkening (2003) and Lee (2010).

## 4. Asymptotes to the loss spectrum

To derive the form for these asymptotes, we need the mechanical energy balance. According to Morris & Jackson (2009*b*) and Lee & Morris (2010), for the bicrystal system shown in figure 1, the external power supplied at the sample boundaries is either dissipated at the grain interface or stored as strain energy within the perfectly elastic grains, i.e.
4.1a
4.1b
and
4.1c
define the strain energy function *W*(*t*) and the dissipation rate . Here, 𝒱 is the combined volume of grains 1 and 2, and *τ* is the *x*-averaged shear stress defined in equation (2.8). As noted in §2, we are taking the grain interface to be time-independent throughout this work.

Before considering the power-law behaviour that is the main topic of this work, we note two results from previous papers. First, according to Morris & Jackson (2009*b*, eqn 53), for *ω*≪1, ℒ∝*ω*^{−1}. This result can be interpreted as follows: for , the quality factor *Q*=ℒ^{−1} is proportional to *ω*, as one might expect from Taylor’s theorem. Second, owing to the slip viscosity in equation (2.1a), ℒ may have a local maximum describing the loss allowed by elastically accommodated grain-boundary sliding. As discussed by Morris & Jackson (2009*b*), for , that local maximum occurs at a large frequency, *ω*=*O*(*ℳ*^{−1}). At these very large frequencies, the background loss caused by diffusion becomes negligibly small, so that the structure of the resulting loss maximum is as described by Lee & Morris (2010).

The power-law spectrum discussed in §1 occurs for (*ω* fixed and large). Let us consider how b.v.p. (2.6) now simplifies. According to equation (2.6d), the shear stress vanishes on the interface: *σ*_{ns}=0. The mass balance expressed by equation (2.6e) also simplifies. According to equation (2.4), the terms on the left side of equation (2.6e) balance on the dimensionless length scale given by ℓ_{D}*ξ*=*ω*^{−1/3}. Because this scale vanishes with increasing *ω*, at any fixed distance from a corner, diffusion along the interface becomes negligibly small, and equation (2.6e) simplifies to [*u*_{n}]=0. According to this discussion, for *ω* fixed and large, and at distance *r* from the corner that is fixed (possibly small), the interfacial conditions simplify to [*u*_{n}]=0=*σ*_{ns}. These are the boundary conditions imposed by Picu & Gupta (1996) in their local analysis of the stress state near a triple junction. According to their analysis, the interfacial normal stress *σ*_{nn} varies with distance *r* measured from the corner according to
4.2
The stress exponent *λ* is independent of material properties, and depends only on corner angle; it satisfies the condition 1>*λ*>0. The first inequality ensures that the strain energy is finite, and the second inequality follows because stress is singular at a corner.

We use this stress field to estimate the dissipation and strain energy. Because diffusion acts to smooth the stress singularity at dimensionless distance *r*_{ℓ}∼*ω*^{−1/3}, we estimate the corresponding integrals by excluding a small neighbourhood of radius *r*_{ℓ} centred on the corner. This cut-off length *r*_{ℓ} determines the form of the loss spectrum. Using equation (4.2) to evaluate equations (4.1b) and (4.1c), we find that
4.3a
and the dissipation per cycle
4.3b
We note that, in equation (4.3a), the integration is carried out over an annular region, so that the area element scales as *r* d*r*. Because *λ*<1, we see that *W* approaches a limit as *r*_{ℓ}→0; the strain energy *W* is not concentrated near the corner. By contrast, the dissipation is focused into the corner region, and its magnitude is controlled by the cut-off scale. Substituting for *r*_{ℓ}, we find that *Υ*∼*ω*^{2(λ−1)/3}. Using the energetic interpretation of mechanical loss ℒ given below (2.10), we obtain
4.4a
where
4.4b

Because *λ* depends on corner angle, so too does *α*. Equation (4.4) holds for both interfaces shown in figure 2, with one exception. A type *S* interface with slope angle *φ*=±45^{°} coincides with the principal axes of stress for simple shear (Lee & Morris 2010). As a result, at the high frequencies at which the simplified boundary conditions apply, grains can deform under simple shear. The entire stress field is then independent of *r*, and the stress exponent *λ*=0. Substituting that value into (4.4), we find that ℒ∼*ω*^{−2/3}. We note that although, in this special case, the Picu & Gupta analysis still predicts a non-zero value for the stress exponent, the boundary conditions ensure that the amplitude of the corresponding eigenfunction is zero. We return to this point in §5.

Figure 3 summarizes the results given above. If the frequencies defining each region of the spectrum are widely separated (i.e. ), the mechanical loss ℒ should scale accordingly as defined in the figure.

In addition to predicting the high-frequency asymptote to the loss spectrum for an inviscid interface, the scaling argument above also implies that the stress near a corner should be self-similar. For, within the corner region, both terms in the interfacial mass balance (2.6e) must be of comparable magnitude; moreover, the stress within that inner region must match the outer stress field given by Picu & Gupta. Using equation (4.2), and the cut-off scale ℓ_{d}, we see that values of the interfacial normal stress *σ*_{nn} computed without approximation as a function of distance *r* from the corner should define a single curve when graphed using the similarity variables *σ*_{nn}*ω*^{−λ/3} and *rω*^{1/3}. This prediction of self-similarity allows another test of the arguments underlying the power-law spectrum; it is verified in §5.

## 5. Comparison with numerical solutions

We show results for 0.1<*ω*<10^{8}, corresponding roughly to the range of dimensionless frequencies encountered in the experiments (e.g. Morris & Jackson 2009*a*, fig. 3). Results are given for the Poisson ratio *ν*=0.3, comparable to that measured in olivine (Christensen 1996); our conclusions are insensitive to this choice.

Figure 4 shows that the mechanical loss spectrum is sensitive to the slope angle. For this figure, we have set *ℳ*=0, so that the interface is inviscid. First, consider the top curve (slope angle *φ*=0.36^{°}); for *ω*>0.2, that curve agrees closely with the small-slope, high-frequency asymptote given by Morris & Jackson (2009*b*, eqn 39). Because that portion of the curve has been obtained by two independent methods, without use of adjustable constants, the agreement provides a test of our numerical method; it also confirms the analysis of Morris & Jackson. The remaining curves (*φ*≥18^{°}) show that, for the range of *ω* shown, ℒ decreases strongly with increasing slope angle; specifically, ℒ decreases by about a factor of 10 when *φ* increases from 0.36^{°} to 30^{°}.

To verify the power-law scaling given by equation (4.4), we note that, for the larger values of *φ*≥18^{°} shown in the figure, ℒ varies as *ω*^{−α} for *ω*≫1. In table 1, we give the values of *α* obtained by fitting equation (4.4a) to the computed spectrum. The *λ* values shown in column 3 of that table are calculated using equation (4.4). Because the normal stress distribution for a type *S* interface is an odd function with respect to the corner, these stress exponents *λ* can be compared with the eigenvalues *λ*_{PG} associated with an anti-symmetric eigenfunction given by Picu & Gupta (1996, fig. 5). Comparing columns 3 and 4 of the table, we see that the computed stress exponents *λ* agree closely with those obtained from the Picu & Gupta analysis, except when *φ*=45^{°}. As explained below equation (4.4), for that special case, *λ*=0 and ℒ∼*ω*^{−2/3}. That prediction is verified in column 2 of the table. We do not display the corresponding values of *λ* and *λ*_{PG} because, as discussed in §4 above, they correspond to different eigenfunctions in this case.

Figure 5 shows the interfacial normal stress *σ*_{nn} near a corner as a function of distance *r* along the interface, with *ω* as a parameter. The figure verifies the self-similarity of the stress field. We also note that, for the type *S* interface, *σ*_{nn} is an odd function of distance along the interface; for this reason, *λ*_{PG} values cited in table 1 were obtained using the curve given in fig. 5 of Picu & Gupta for an anti-symmetric stress field. (In their figure, the curve labels are interchanged; the solid line should correspond to the anti-symmetric eigenfunction.) Figures 6 and 7 show the relation between the loss-maximum occurring when *ℳ*≠0, and the background spectrum discussed above.

Figure 6 shows the rigidity *G* computed as a function of angular frequency *ω* with as a parameter for a type *S* interface with *φ*=30^{°}.

Figure 6*a* shows ℒ as a function of *ω* for a viscous interface. In the curve for *ℳ*=10^{−8}, all the features summarized in figure 3 are present: for *ω*≪1, the mechanical loss ℒ varies as *ω*^{−1}; for 1≪*ω*≪10^{5}, ℒ follows the power-law asymptote discussed above; the local maximum caused by elastically accommodated grain-boundary sliding is found at *ω*∼10^{−8}; thereafter, ℒ varies as *ω*^{−1}, as shown in figure 3. At the local maximum ℒ≃0.05, approximately equal to the value found in Lee & Morris (2010, fig. 9) for the same values of the control parameters. Although we do not show the loss spectra for other values of *ℳ* we note that, once the maximum is clearly visible, its height is independent of *ℳ*, because the loss due to diffusionally accommodated sliding is then small at the peak frequency.

The curve for *ℳ*=10^{−3} is included to show that, when the sliding time scale and the diffusion time scale are not widely separated, the loss decreases rapidly with increasing frequency, except for a short plateau covering a couple of decades in frequency.

Figure 6*b* shows the sample rigidity |*G*|. From the curve for *ℳ*=10^{−8}, we see that the response consists of two regions of constant |*G*| separated by transition regions. The first plateau covers the range 10^{2}<*ω*<10^{7}. Within this frequency range, ℒ follows the power-law asymptote and the shear stress vanishes over most of the interface; because only normal stresses act on the interface, |*G*| is less than the unit rigidity of the grains. The second plateau occurs for *ω*>10^{8}. At these high frequencies, the grains behave as if they are welded at the interface, i.e. [*u*_{n}]=0 and [*u*_{s}]=0, and |*G*|→1. Similar behaviour is predicted by the small-slope analysis (Morris & Jackson 2009*b*).

Figure 7 shows the corresponding results for the type *TS* interface. They are included to show that the slowly varying region in the mechanical loss spectrum depends on corner orientation, as well as on the angle subtended by the corner. For this type *TS* interface with *φ*=60^{°}, the subtended angle is identical to that of the type *S* interface discussed in figure 6. The orientation is different, however. Figure 7 shows that, in the power-law regime, ℒ decreases more rapidly in the present case. This more rapid decay reflects the parity of the most singular allowable stress eigenfunction. According to figure 5*a*, for the type *S* interface, *σ*_{nn} is an *odd* function of distance along the interface, whereas for a type *TS* interface, *σ*_{nn} is nearly an *even* function. Using Picu & Gupta (1996, fig. 5), we find that, for a symmetric stress eigenfunction *σ*_{nn}, the stress exponent *λ*=0.45. The same value is obtained by fitting the values of ℒ shown in figure 7 to equation (4.4). We conclude that although, at high frequencies, dissipation is concentrated near grain corners, we cannot predict the loss spectrum without accounting for the orientation of grain boundaries.

In cross sections of experimental samples, corner angles of differing sizes occur. It is interesting to see how two corners subtending different angles affect the loss spectrum. Because the strain energy *W* is insensitive to local stress behaviour, the mechanical loss ℒ can be found by summing the contribution of the dissipation *Υ* from each region surrounding a corner. Consequently, the mechanical loss ℒ behaviour in the slowly varying region is a summation of the power-law scaling associated with each corner. The constants of each scaling are determined by the respective constants of proportionality found in the Picu & Gupta local stress description. Our scaling analysis suggests that the mechanical loss behaviour in polycrystals at sufficiently high frequencies, i.e. , will be controlled by the corner having the largest stress exponent *λ*.

To test this prediction, we consider an interface illustrated in the inset of figure 8 by the solid line. Along the interface, there are two different corners *C*_{1} and *C*_{2} having angles *ϕ*_{1}=175^{°} and *ϕ*_{2}=107^{°}, respectively. For these two corners *C*_{1} and *C*_{2}, the local analysis by Picu & Gupta (1996) predicts the strongest stress exponents *λ* to be 1 and 0.5, respectively. The behaviour of the mechanical loss ℒ at sufficiently high frequencies is therefore expected to be controlled by *C*_{1}.

Figure 8 shows the mechanical loss spectrum obtained for the interface shown in the inset. There are two main features in the figure. First, the behaviour of the mechanical loss is consistent with the above prediction and appears to approach the logarithmic scaling, i.e. corresponds to a stress exponent *λ*=1 at *C*_{1}. The graph is truncated at *ω*=5×10^{8} owing to a lack of numerical resolution at higher frequencies. Second, the slope decreases gradually with frequency in the slowly varying region due to the diminishing effect on the loss spectrum from the other corner *C*_{2}. To show that the effect of *C*_{2} indeed diminishes with increasing frequency *ω*, we also graph the scaling ℒ∼*ω*^{−0.33} produced by *C*_{2}.

This result is also consistent with the behaviour of the mechanical loss ℒ found in experiments. Because corner angles in triple junctions vary spatially within polycrystals, a gradual decrease in the slope of the mechanical loss spectrum caused by the diminishing effect from corners having smaller stress exponents *λ* is also expected to be observed in experiments. This may explain the behaviour seen in Morris & Jackson (2009*a*, fig. 3), in which the measured quality factor *Q*=ℒ^{−1} becomes decreasingly sensitive to *ω* at higher frequencies.

## 6. Conclusion

Using a bicrystal model, we have made an analytical and numerical study of grain-boundary sliding along the prescribed spatially periodic interface separating two Hookean layers. This work differs from the analysis of Morris & Jackson (2009*b*) because the interface is now allowed to be finite. Like the model analysed by Morris & Jackson (2009*b*), our model contains two timescales; as defined by equation (2.3), these are the timescale *t*_{D} defined by grain-boundary diffusion, and the timescale *t*_{η} characterizing a deformation for which the stress due to slip along the viscous interface is comparable with the stress within the bulk of the crystal. For reasons explained in Morris & Jackson (2009*a*,*b*), we emphasize the behaviour occurring when *t*_{D}≫*t*_{η}. The time scales are then widely separated, and the shear stress on the interface is effectively zero at the angular forcing frequencies *ω* of seismic interest.

Although the qualitative features predicted by the small-slope analysis of Morris & Jackson (2009*b*) are also present when the interface has finite slope, the magnitude of the mechanical loss ℒ now depends strongly on slope. We find that for *ωt*_{D}≫1, the loss factor is still controlled by stress concentrations near corners on the interface. However, the contribution of each corner now depends on the angle subtended by the corner, and on the local orientation of the interface relative to the principal axes of applied stress.

To demonstrate these new effects made possible by the finite interface slope, we have used the example of a piecewise linear interface whose corners subtend identi- cal angles. For *ωt*_{D}≪1, the mechanical loss now follows the power-law asymptote given by ℒ∼const.*ω*^{−α}. The exponent *α* depends on the orientation of the sliding surfaces, and on the angle subtended by the corner (compare figures 6*a* and 7*a*); it is, however, bounded by the condition 0<*α*≤2/3. Although Raj & Ashby (1971) viewed the bicrystal model as being representative of orthogonal sliding surfaces within a crystal, the dependence of ℒ on the orientation of the sliding surface means that no single interface can be taken as representative of an experimental sample. For example, the mechanical loss for a (highly idealized) experimental sample consisting of a plane array of regular hexagons would still depend strongly on the orientation of that array relative to the principal axes of applied stress.

From the example of the simplest piecewise linear interface, we learned that, with increasing *ω*, the loss factor for a small-slope interface decays more slowly than that for an interface having larger slope. That observation suggested our second example: a primary sawtooth interface having on each rising face a slight kink. In this case, the interface has two types of corner, each type subtending a different angle, and we expect the slope of loss spectrum in a loglog plot to decrease with increasing frequency, as the effect of the corners having weaker stress concentration diminishes. That idea is consistent with experimental observations; we speculate that the differing values of power-law exponent *α* reported experimentally might be a consequence of differing ranges of *ωt*_{D} being studied in the various experiments. Moreover, for a two-corner interface, we expect that, with increasing *ωt*_{D}, the mechanical loss should be controlled ultimately by the corner having the largest stress exponent. Our numerical results are consistent with that prediction.

In the above discussion, we have assumed that, at seismic frequencies, the shear stress acting on the grain boundary is small. However, our model does include the boundary viscosity, and at frequencies the mechanical consequently has a local maximum. Although the precise magnitude of the background loss due to dif- fusion depends on the details described above, our assumption *t*_{D}≫*t*_{η} means that this local maximum occurs at such large frequencies that the diffusive loss is then negligibly small in any case. Consequently, as discussed in the context of figure 6, the results of Lee & Morris (2010) concerning the local maximum apply directly.

When the bicrystal model is scaled to the conditions of the experiments discussed in §1, the loss spectrum at seismic frequencies is controlled by local stress behaviour near corners. This result should also apply in polycrystals. Because our analysis of the bicrystal model shows that the quantitative behaviour of the loss spectrum depends on the corner angle and interface orientation, we also expect that to be true in polycrystals. Consequently, it seems unlikely that the computed loss spectrum for an regular array of hexagons in plane strain would quantitatively match experiments. We believe that, for quantitative comparison, the next model must include a random distribution of crystal orientations and corner angles at triple junctions. To address these complications, one may have to resort to homogenization techniques.

## Acknowledgements

We are grateful to Prof. Ian Jackson and Prof. Tarek Zohdi for their helpful comments and discussions. We also thank the reviewers for their valuable comments that have helped us to improve the presentation. L.C.L. was supported in part by a Committee on Research Faculty Research Grant to S.J.S.M. from the University of California. J.K. was supported by the National Science Foundation through grant DMS-0955078 and by the Director, Office of Science, Computational and Technology Research, US Department of Energy under contract DE-AC02-05CH11231.

## Appendix A. Solution procedure

Because interface is time-independent and the b.v.p. given in equation (2.6) is linear, the principle of superposition applies. We decompose that b.v.p. into two separate b.v.p.’s, which we denote here as b.v.p.^{(1)} and b.v.p.^{(2)}. These two b.v.p.’s share the same geometry shown in figure 1. Using superscripts (1) and (2) to denote, respectively, variables associated with b.v.p.^{(1)} and b.v.p.^{(2)}, the stress field *σ*_{ij}, strain field *e*_{ij} and the displacement fields *u*, *v* of problem (2.6) can be obtained by superposing the solution of the two b.v.p.’s, i.e.
A 1a
A 1b
and
A 1c

To simplify the notation, we use *g*_{n} and *g*_{s} here to denote, respectively, the normal displacement jump [*u*_{n}] and the tangential displacement jump [*u*_{s}] across the interface . Interfacial stresses and displacement jumps are also denoted using 2×1 vector of functions *σ*_{n}=[*σ*_{nn},*σ*_{ns}]^{T} and *g*=[*g*_{n},*g*_{s}]^{T}, respectively.

The plane elastostatic equation in equation (2.6a), the periodic boundary conditions in equations (2.6h) and (2.6i) and the requirement that the normal and tangential stresses across the grain boundary are continuous in equations (2.6*f*,*g*) all apply in b.v.p.^{(1)} and b.v.p.^{(2)}. The other boundary conditions are now stated. In b.v.p.^{(1)}, the boundary conditions at *y*=±*a* are
A 2a,b
and the boundary conditions along the interface are
A 3a,b
Conversely in b.v.p.^{(2)}, boundary conditions at *y*=±*a* are
A 4a,b
whereas boundary conditions along the interface are
A 5a,b

By inspection of b.v.p.^{(1)}, the two grains do not interact with one another through the interface . Hence, the two grains move rigidly across one another and the displacement field ** u** of the upper grain and the lower grain are

*ie*

^{iωt}and −

*ie*

^{iωt}, respectively. These displacement fields satisfy all equations given in b.v.p.

^{(1)}, and the resulting normal displacement jump and tangential displacement jump across the interface 𝒮

_{I}are, respectively, A 6 We also note that, in b.v.p.

^{(1)}, the stress field .

To solve b.v.p.^{(2)}, we use eigenfunction expansion. In essence, we reduce a two-dimensional problem given in b.v.p.^{(2)} to a one-dimensional problem defined along interface . We define a linear operator S that maps the given displacement jumps *g*^{(2)} onto the interfacial stresses . Note that S solves for when *g*^{(2)} is prescribed along the interface . Because stresses in b.v.p.^{(1)} are zero, the interfacial stresses in b.v.p.^{(2)} are equivalent to that in the original b.v.p., i.e. . The operator S is defined as follows:
A 7
We also define the differential operator L as
A 8
In equation (A 8), L operates separately on functions *σ*_{ns} and *σ*_{nn}; multiplying *σ*_{ns} by and taking the second derivative of *σ*_{nn} with respect to *s*. Using the definitions given in equations (A 7) and (A 8), and noting that , we find, from the constitutive equations (2.6d) and (2.6e) of the original b.v.p., that . Applying the principle of superposition to that equation, the two-dimensional elasticity problem is absorbed into the operators leaving a single equation governing the time-evolution of the interfacial displacement jumps,
A 9
The r.h.s. term in equation (A 9) can be calculated using equation (A 6).

Time evolution of the interfacial gap *g*^{(2)} defined in equation (A 9) can be obtained by eigenfunction expansion once the eigenvalues *γ*_{k} and the eigenfunctions *Z*_{k}(*s*) associated with the composite operator LS are known, i.e.
A 10
Using *N*_{z} eigenfunctions, the solution to the homogeneous part of equation (A 9) (i.e. with ) is given by a separable form
A 11
where the subscript *h* refers to the homogeneous solution and *β*_{k} are coefficients determined by the initial condition . The coefficients *β*_{k} can be found by requiring them to satisfy
A 12
Letting *Φ*_{Z} be a 1×*N*_{z} vector containing these eigenfunctions,
A 13
and *Φ**_{Z} be the adjoint operator of *Φ*_{Z} so that is a *N*_{z}×1 vector of scalars defined as
A 14
the coefficients ** β**=[

*β*

_{1},

*β*

_{2},…,

*β*

_{NZ}]

^{T}, upon solving equation (A 11) for

*β*

_{k}, can be written as A 15 Substituting equation (A 15) into equation (A 11), the latter equation can be written compactly as A 16 where E(

*t*) is defined as the evolution operator, or propagator A 17 and

**is a diagonal matrix defined as A 18 The solution to the inhomogeneous PDE given in equation (A 9) can then be obtained using Duhamel’s principle, A 19 Hence, b.v.p.**

*Λ*^{(2)}is solved, if the eigenvalues

*γ*

_{k}and the eigenfunctions

*Z*

_{k}of LS defined in equation (A 10) are found.

The steady-state response of *g*^{(2)} can be obtained by setting the first r.h.s. term in equation (A 19) to zero (because it vanishes as ), and setting the lower integration limits in the second r.h.s. term from 0 to i.e.
A 20
The subscript ss is used here to denote steady-state solution. Substituting equations (A 17) and (A 6) into equation (A 20) and then evaluate the resulting integral, the steady-state response of the displacement jump in b.v.p.^{(2)} becomes
A 21a
where its frequency response is given as
A 21b
and *D* is a *N*_{z}×*N*_{z} diagonal matrix with its *k*th component given by
A 21c
Noting that , as explained above in equation (A 7), the steady-state response of the interfacial stresses can be calculated using the operator S, i.e.
A 22
Integrating the *x*-projection of along the interface then leads to the *x*-averaged shear stress *τ* defined in equation (2.8). The mechanical loss ℒ can thereafter be calculated using *τ* as described in the main text. Thus, for any given interface, the mechanical loss spectrum can be obtained by computing the eigenvalues *γ*_{k} and the eigenfunctions *Z*_{k}(*s*) of the operator LS. We thus reduce problem (2.6) into an eigenvalue problem.

To avoid computing the second derivative of *σ*_{nn}, a pseudo-inverse of LS is used, instead, to find the eigenvalues and eigenfunctions. The pseudo-inverse A has the same eigenfunctions *Z*_{k}(*s*) as LS, and its eigenvalues *ζ*_{k} are related to those of LS by
A 23
The pseudo-inverse A is constructed using finite-element methods; details are given in Lee (2010) and Sethian & Wilkening (2003).

- Received August 24, 2010.
- Accepted November 24, 2010.

- This journal is © 2010 The Royal Society