## Abstract

We present direct finite strain and rheological estimates of naturally deformed mantle materials using field observations in the Twin Sisters ultramafic body of Washington State. Folded and elongated (boudinaged or ductilely thinned) orthopyroxenite dikes within the host dunite provide strain markers that allow us to characterize the finite strain over an ∼15,000 m^{2} area. Using dynamic instability analysis on the folded orthopyroxenite dikes, orthopyroxene is calculated to have ∼25 times the effective viscosity of olivine-rich host rocks (based on a power-law exponent of 3.0 for this dislocation creep process). Detailed mapping also indicates that inclusion of up to 15% orthopyroxene in the olivine-rich host rocks does not affect their viscosity, constraining the rheological behavior of two-phase, mantle material.

## INTRODUCTION

Knowledge of the rheology of Earth's lithosphere is critical to quantitative understanding of tectonic processes. Increasingly sophisticated models are able to characterize the three-dimensional development of Earth systems, but these models require accurate rheological estimates. To date, the rheology of the lithosphere has been quantified through the use of experimental rock deformation studies (e.g., Brace and Kohlstedt, 1980; Kohlstedt et al., 1995). While the experimental results form the basis for our quantitative understanding of rheology, the experimental approach is limited by extrapolation of laboratory-based results to geologically relevant spatial and temporal scales, natural external conditions (pressure, temperature), and compositional heterogeneity of naturally deformed materials (Paterson, 1987, 2001).

Further, models typically describe the bulk behavior of distinct lithospheric polyphase layers by the mechanical behavior of a single mineral. Olivine, for example, as the volumetrically dominant mineral in the lithospheric mantle (e.g., Hirth, 2002), is thought to control its bulk behavior.

Field-based constraints on rheology have lagged behind data from rock deformation experiments due to the inherent complexity (composition, water content, grain size, preexisting fabrics, etc.) of geological materials. This complexity makes it difficult to constrain the influence of these variables on the rheology of naturally deformed materials. Field-based—or “natural laboratory”—approaches are particularly relevant, however, because they analyze naturally deformed materials at geological conditions.

The Twin Sisters ultramafic body in Washington State is an excellent natural laboratory to estimate the relative rheologies of important constituents of the lithospheric mantle (orthopyroxenite, dunite, harzburgite). Using the arc length (wavelength)/thickness ratio of folded orthopyroxene dikes, we determined the relative effective viscosity of the orthopyroxenite versus the olivine-rich rocks (dunite and harzburgite). Second, we determined the effect of the second phase (orthopyroxene) on the rheology of a typical, two-phase mantle composition (harzburgite).

## GEOLOGICAL BACKGROUND

The Twin Sisters Range is part of the North Cascade Mountains in the Cordillera of western North America and extends from northern Washington State to southwest British Columbia. Consisting of a faulted assemblage of accreted oceanic terranes (e.g., Davis et al., 1978), the range exposes several large peridotite bodies, including the Ingalls, Fidalgo, and Twin Sisters bodies. Unlike the Ingalls and Fidalgo ophiolites, the undated Twin Sisters body consists solely of ultramafic rocks. Located 6 km southwest of Mount Baker and ∼30 km east of Bellingham (Fig. 1), the Twin Sisters ultramafic body is 16 km long parallel to its NNW strike, is 5–6 km wide, and is estimated from gravity modeling to extend as a sheet to 2 km depth (Thompson and Robinson, 1975). Temperature estimates for the Twin Sisters range from ∼800 °C (Toy et al., 2010; based on two-pyroxene and olivine + spinel geothermometry) to ∼1000 °C (Onyeagocha, 1978; based on olivine + spinel geothermometry).

We completed detailed structural mapping (Fig. 1) in a 100 × 150 m area located west of the crest of the Twin Sisters Range and approximately midway between the northern and southern ends. This area was chosen for the presence of numerous dikes, continuity of glaciated exposure, and lack of serpentinization and vegetation. The distribution, shape, and three-dimensional geometry of orthopyroxene dikes, faults, and chromite veins were mapped using a laser surveying station at centimeter-scale accuracy. The area is dominated by meter-scale, alternating bands of dunite (>95% olivine) and harzburgite (∼15% orthopyroxene; ∼85% olivine) (Fig. 1). Foliations (∼163, 75°E) and lineations (pitching ∼40°S in the foliation plane) are consistent throughout the field area. Shape-preferred orientations of elongate chromite grains in the dunite and orthopyroxene grains in the harzburgite define the foliation. Preferred alignment of individual and aggregates of chromite grains defines the lineation. Compositional bands in the harzburgite and dunite are parallel to foliation (Fig. 1).

There is a bimodal grain-size distribution in olivine grains in dunites (>90% olivine) and harzburgites (>80% olivine), with grain size ranging from 2 to 10 mm. Olivine grains exhibit triple junctions and subgrain development. Locally, particularly adjacent to a boundary with an orthopyroxene dike, a finer grain size (∼0.25 mm) population is observed. Lattice preferred orientation (LPO) analysis of olivine grains of different grain sizes in host dunites and harzburgites shows that the [100] axes of that mineral consistently form a point maximum plunging shallowly to the SE, parallel to the lineation observed in the field (Fig. 1). The [010] axes define a partial girdle that plunges moderately to the ESE. The [001] axes define a maximum that also lies within a partial girdle in the NE quadrant of the stereogram. These patterns, typical for the Twin Sisters Range (Christensen, 1971, 2002), suggest operation of the (0kl)[100] slip systems in dislocation creep and deformation temperatures of ∼1000 °C (e.g., Carter and Ave Lallemant, 1970).

Two sets of deformed dikes are observed in the field area: (1) older, dominantly orthopyroxene-bearing, folded dikes; and (2) younger, clinopyroxene-bearing, planar dikes. The folded orthopyroxene dikes (Fig. 2) are the critical marker for this study and show a distinctive pattern with respect to foliation. The axial planes of the folds generally strike NW-SE, parallel to the foliation, whereas fold hinges are strongly girdled within this axial plane (Fig. 1). The variation of fold hinge orientation suggests: (1) that all the dikes were not originally parallel; and (2) the finite strain ellipsoid characterizing bulk deformation is probably oblate, since a prolate finite strain ellipsoid would result in parallel fold hinges. Fold limbs that strike in a NW-SE direction (parallel to foliation) are boudinaged or ductilely thinned with respect to the folded portions of the dike. Importantly, the compositional variation within the host rock (dunite versus harzburgite) has no discernible effect on the folding or boudinage/thinning behavior of the orthopyroxene dikes (Fig. 1). Alternatively stated, when crossing a dunite/harzburgite contact, there is no change at any scale in the geometry of the orthopyroxene dikes (e.g., shape of folded dike; thickness and orientation of straight dike).

Orthopyroxene grains within the dikes are coarse (∼10 mm) and display undulose extinction, subgrain development, and a strong, often bent, cleavage. Locally, finer grains (0.05 mm) are observed along orthopyroxene grain boundaries, especially close to the boundary of orthopyroxene dikes with adjacent dunite or harzburgite. The undulose extinction, subgrains, and new grains suggest deformation by dislocation creep.

## FINITE STRAIN ANALYSIS

Finite strain was calculated using the methodology of Talbot (1970), which characterizes the stretches and orientations of all deformed markers. The stretch (change in length or *S* = *l*_{f}/*l*_{o}, where *l*_{f} is the final length, and *l*_{o} is the original length) was calculated for a variety of orthopyroxene dikes. In cases of shortening (folding), the fold shapes were projected into the profile plane that is perpendicular to the fold hinge. The final length of the dike (*l*_{f}) was simply measured in the profile plane. The original length of the dike (*l*_{o}) was determined by numerically “unfolding” the dike. We first calculated the hinge-to-hinge arc length and amplitude of the folds in the profile plane, from which we could calculate the length along the sinusoidal form. This length was the original length of the dike (*l*_{o}). Because we did not know the original thickness of the dike prior to deformation, we could not independently determine the amount of layer-parallel shortening, and thus the shortening values represent minimum estimates. The relative stiffness of the orthopyroxene dikes (see following) suggests that this assumption is not a major factor.

When elongated, the orthopyroxene dikes occasionally were boudinaged, but more generally underwent ductile thinning (a similar result to that obtained by Talbot (1970) for different materials). This thickness was compared to thicknesses in the folded areas of the same dike in order to calculate the stretch. If layer-parallel shortening occurred during folding, the original thickness would be less than assumed, and thus the stretch values are maximum estimates.

The method of Talbot (1970) works on a radial coordinate system, shown graphically in Figure 3A. For any individual dike, the present azimuthal orientation of the dike is given by the angle in real space, and the stretch is plotted relative to a unit length. Each dike is plotted twice in this graph, on both ends of the azimuthal line (e.g., 160 and 340, in Fig. 3A). If the dikes are undeformed, they plot in a circle. Deformed dikes trace out an ellipse in two dimensions, which corresponds to the finite strain ellipse. As shown in Figure 3B, the dikes in the field area describe the finite strain ellipse with an elliptical ratio 1.596, oriented at 151° (within the variation of foliation measurements in the field area).

The three-dimensional finite strain ellipsoid can also be calculated, if: (1) the maximum finite stretch direction (*S*_{1}) is assumed to lie parallel to the lineation direction; and (2) no volume loss has occurred (*S*_{1}*S*_{2}*S*_{3} = 1). Because the minimum finite stretch direction (*S*_{3}) lies approximately in the horizontal plane, it can be directly measured from the deformed dikes (*S*_{3} = 0.286). Thus, a relationship between the two other principal stretch directions is obtained (*S*_{1}*S*_{2} 0.286 = 1). Using a numerical simulation, the relative magnitudes of the finite strain axes were iteratively arrived at, until a horizontal ellipse of the correct aspect ratio of 1.596 was obtained. The best-fitting answer, which defines an oblate ellipsoid, was: *S*_{1} = 3.15, *S*_{2} = 1.11, and *S*_{3} = 0.286. Thus, the Twin Sisters has an oblate-shaped finite strain ellipsoid, the long axis of which plunges 40°SE. This result is consistent with thinning of the planar dikes in both the horizontal and vertical directions. We emphasize that this calculation is a first-order estimate, but it nonetheless characterizes the deformation of the ultramafic body after injection of the orthopyroxenite dikes.

## RHEOLOGICAL ANALYSIS

Like most single-layer buckle folds observed in nature (e.g., Chapple, 1978), the folded orthopyroxenite dikes display an approximately linear relationship between arc length and layer (dike) thickness. This relationship has been widely explained using the elementary theory of dynamic instabilities for a folded layer that is embedded in a less viscous medium (e.g., Fletcher, 1974; Smith, 1977; Ramberg, 1981). By fitting a line to our fold arc length (λ) – thickness (*H*) data for the orthopyroxene dikes, we get the equation *y* = 8.48(±2.59)*x* + 31.56 (±77.66) (Fig. 4). The parentheses designate the 95% confidence interval (technically prediction intervals) given for *x* = 50. The slope thus is constrained between 5.89 and 11.07. The intercept, given by the second term in the line equation, should go through zero, which is consistent with the data.

*n*= 3 (e.g., Jin et al., 1994). Using the method outlined in Smith (1977), assuming a thin-plate approximation, the effective viscosity ratio for non-Newtonian folded material can be determined using the following relationship:

*/*

_{m}*H*is the ratio of arc length to thickness;

*n*

_{1}(matrix) and

*n*

_{2}(folded layer) are the stress exponents for the two materials (

*n*= 3 for both dunite and orthopyroxene, respectively, deformed by dislocation creep); and

*m*is the effective viscosity ratio of the two materials.

The *m* value represents the relative viscosities of the orthopyroxene dikes and olivine-rich host rock. Thus, we estimate that orthopyroxene dikes in this area are ∼25 times more viscous than the surrounding olivine-rich rock. Both the folded layer and the matrix are essentially monophase, so that the viscosity of the folded layer relative to the matrix host suggests the relative viscosities of orthopyroxene and olivine, respectively. Thus, orthopyroxene is estimated to be ∼25 times more viscous than olivine.

Our methodology is similar to one used by Trepmann and Stockhert (2009), who determined the relative viscosity of a folded quartz vein relative to its host matrix, a polyphase graywacke. Fletcher and Sherwin (1978) calculated the fold arc length in a different way, by using frequency distribution graphs to determine a mean arc length (also see Hudleston and Holst, 1984). Their approach is generally preferable because it does not bias, for example, the arc length–thickness relations by including folds within a single thick, folded layer. Our data were not collected in a way that would allow us to complete the calculation following the methodology of Fletcher and Sherwin (1978). Although 95% confidence intervals are given for the line, these values cannot be used to evaluate uncertainties in viscosities estimates.

One last difficulty relates to the form of the thin-plate approximation equation (e.g., Smith, 1977). Essentially, a small uncertainty in the ratio of arc length to thickness causes a large uncertainty in the relative viscosity measurements because of the cubic root dependence. This uncertainty is unavoidable in this type of analysis.

## IMPLICATIONS

This study yields three significant results related to the rheology of the lithospheric mantle. First, orthopyroxenite dikes have an effective viscosity ratio ∼25 times the olivine-rich host rock (assuming dislocation creep and a power-law exponent of 3). Because of the composition of the dikes and the dunite host, we interpret this result in terms of the relative viscosities of orthopyroxene and olivine, respectively, at the conditions (pressure, temperature, strain rate) of deformation. This result can be compared with results from two experimental deformation studies of olivine-orthopyroxene aggregates. Ji et al. (2001) reported exactly the opposite relation than that observed in the Twin Sisters rocks; orthopyroxene is less viscous than olivine. However, the deformation microstructures in their experimentally deformed aggregates suggest that dislocation creep operated in the olivine, while diffusion creep operated in the orthopyroxene. Thus, the results of Ji et al. (2001) are not directly comparable to our observation of naturally deformed peridotites deformed by dislocation creep. The experimental results of Hitchings et al. (1989) are consistent with the findings of this study in two ways. Their experiments indicate a power-law exponent of 3 for experiments on orthopyroxene-only samples, and the orthopyroxene-only samples are somewhat stronger than their olivine-only samples, although not to the degree that we observe in the Twin Sisters field area. In addition, experimentally derived flow laws for orthopyroxene and olivine (Lawlis, 1998; Hirth and Kohlstedt, 2003), at the conditions of deformation at the Twin Sisters (900 °C; strain rate = 10^{−15} s^{−1}), do not suggest as significant a viscosity contrast as estimated by this analysis of naturally deformed rocks.

Second, we observe no indication of strain localization in this area, except for late-stage faulting. We note only minor localization on the edges of the dikes and no significant localization across the compositional banding (harzburgite/dunite) within the field area (Fig. 1). This result suggests that the lithospheric mantle deforms in a relatively homogeneous manner in some settings (e.g., Little et al., 2002). Further, the finite strain analysis provides a minimum strain estimate responsible for the observed fabrics; the estimate is a minimum because some deformation may have occurred prior to orthopyroxene dike emplacement.

Third, the homogeneity of deformation allows us to evaluate the extrapolation of an experimentally determined rheology, based on the behavior of a single mineral phase, to polyphase rocks that contain at least two minerals. The field observations and mapping (Fig. 1) provide clear qualitative evidence that compositional differences between dunites to harzburgites (up to 15% orthopyroxene) have no discernible effect on the folding (or thinning/boudinage) behavior of the orthopyroxene dikes. Consequently, we infer that the weaker olivine controls the bulk behavior of the rock, even in olivine-rich rocks that contain up to 15% of the stronger orthopyroxene phase. This result is consistent with theoretical analyses of two-phase behavior (e.g., Handy, 1990, 1994; Treagus, 2002).

One caveat to our study is that orthopyroxene dikes often contain relatively coarse-grain-size (∼10 mm) orthopyroxene grains. The grain size of olivine within the dunite and harzburgite host rocks varies significantly (2–10 mm), but is generally smaller than orthopyroxene within the orthopyroxene dikes. Consequently, there may also be a grain-size effect in our study area, which acts to increase the effective viscosity of the dikes.

## CONCLUSIONS

This study provides quantitative, field-based constraints on large-scale finite strain and the relative viscosity of the two dominant minerals in the lithospheric mantle. Using dynamic instability analysis, orthopyroxenite dikes have an effective viscosity ratio ∼25 times the olivine-rich host rock (assuming dislocation creep and a power-law exponent of 3). Given the relatively monophase composition of both the dikes and the dunite host, we interpret this result in terms of the relative rheology of orthopyroxene and olivine under these conditions. This large difference in mechanical properties could result in significant variations in the rheology of mantle rocks, especially if the orthopyroxene content is relatively high. However, our study also indicates that less than 15% of a second phase (orthopyroxene) does not significantly affect the rheology of olivine-rich material. Thus, for relatively minor amounts of a stronger second mineral (e.g., orthopyroxene), direct field observations verify the commonly applied assumption that the mineral olivine controls naturally deformed olivine-rich lithospheric mantle behavior.

The research was funded by a Packard Foundation grant to Tikoff and National Science Foundation grant EAR-0409522 to Tikoff and Newman. The authors express their gratitude to Paul Riley for mathematical consultation. Thanks are due to Ray Fletcher, Greg Hirth, Bryn Benford, and an anonymous reviewer for constructive reviews of the manuscript. Corky Smith, of Olivine Corporation, is gratefully acknowledged for logistical help and encouragement.