A new Mg-in-plagioclase geospeedometer for the determination of cooling rates of mafic rocks

Available online at www.sciencedirect.com

ScienceDirect Geochimica et Cosmochimica Acta 140 (2014) 691–707 www.elsevier.com/locate/gca

A new Mg-in-plagioclase geospeedometer for the determination of cooling rates of mafic rocks Kathrin Faak a,b,⇑, Laurence A. Coogan b, Sumit Chakraborty a a

Institut fuer Geologie, Mineralogie und Geophysik, Ruhr-Universitaet Bochum, Universitaetsstr. 150, D-44801, Germany School of Earth & Ocean Science, University of Victoria, PO BOX 1700 STN CSC, Victoria, BC V8W 2Y2, Canada

b

Received 31 October 2013; accepted in revised form 6 June 2014; Available online 18 June 2014

Abstract We present a new Mg-in-plagioclase geospeedometer based on the diffusive exchange of Mg between plagioclase and clinopyroxene that allows cooling rates of a wide range of mafic rocks to be determined. The major element composition of plagioclase is shown to play an important role in driving the flux of Mg as well as the partitioning of Mg between plagioclase and clinopyroxene. Therefore, commonly used analytical solutions to the diffusion equation are not applicable in this system and an alternate approach is developed. Practical aspects of the method such as the role of anorthite zoning patterns and influence of shapes of cooling paths and grain geometries in controlling the Mg concentration profile shapes are discussed. Propagation of uncertainties in the knowledge of diffusion and partition coefficients shows that the method is capable of resolving small differences in cooling rates (e.g., log[dT/dt °C yr1] = 3 ± 0.3 in one example). The approach is illustrated by application to two samples from the lower oceanic crust. Ó 2014 Elsevier Ltd. All rights reserved.

1. INTRODUCTION Quantifying the rate of change of temperature of rocks provides important insights into energy transfer processes on Earth and other planetary bodies. This can be achieved by petrological approaches and provides an important link between petrological and geodynamic studies. There are two general approaches used to determine cooling (or less commonly heating) rates: (i) radiogenic dating of minerals with different closure temperatures (Tc) to define a temperature–time path; and (ii) geospeedometry. The former is predicated on the fact that the temperature at which the daughter isotope of a radiogenic parent becomes trapped in a mineral is dependent on the diffusivity, or more

⇑ Corresponding author at: Institut fuer Geologie, Mineralogie und Geophysik, Ruhr-Universitaet Bochum, Universitaetsstr. 150, D-44801, Germany. Tel.: +49 234 32 23228. E-mail address: [email protected] (K. Faak).

http://dx.doi.org/10.1016/j.gca.2014.06.005 0016-7037/Ó 2014 Elsevier Ltd. All rights reserved.

generally, the kinetics of exchange with the surroundings, of the daughter in this phase. The temperature at which the daughter isotope is ‘locked’ into the phase is termed the closure temperature (Dodson, 1973) and by dating various minerals, using different radiogenic systems, a series of temperature–time points can be defined. The second approach, commonly termed geospeedometry (Lasaga, 1983), is based on the extent to which temperature dependent ion exchange reactions (geothermometers) proceed during changes in temperature. For example, if a rock is initially in equilibrium at a high temperature, and then cools rapidly, there will be little time at high temperatures for diffusive modification of the distribution of the ions – this will freeze in a high closure temperature of this thermometer; i.e., measurements of the mineral compositions and their use in the thermometer calibration will give a high temperature. In contrast, slower cooling will allow greater diffusive modification of compositional zoning and lead to a lower closure temperature for the thermometer. Quantifying the cooling rate in such situations requires knowledge of the

692

K. Faak et al. / Geochimica et Cosmochimica Acta 140 (2014) 691–707

temperature dependence of both the exchange reaction and the diffusion coefficient of the ion(s) involved (Lasaga, 1983). Other similar approaches have been applied to structural transitions in glasses and minerals that are kinetically controlled (e.g., Seifert, 1977; Ganguly, 1982; Wilding et al., 1995; Redfern et al., 1996). A number of geospeedometers have been developed and widely applied such as garnet-biotite Mg/Fe exchange (Spear, 1991; Chakraborty and Ganguly, 1991; Lindstrom et al., 1991), olivine-spinel Mg/Fe exchange (Ozawa, 1984), Ca-in-olivine (Coogan et al., 2007), and Zr in oxides (Onorato et al., 1981; Blackburn et al., 2012). A very common mineral association in mafic rocks is plagioclaseclinopyroxene. A geospeedometer applicable to this mineral pair could be applied in a wide variety of settings ranging from mafic and intermediate igneous rocks, to metamorphic rocks (granulites) and even lunar and meteorite samples. Faak et al. (2013) recently calibrated the temperature dependence of Mg exchange between plagioclase and clinopyroxene as a thermometer and determined the diffusion coefficient for Mg in plagioclase building on previous work by LaTourrette and Wasserburg (1998). These experimental data now allow a Mg-in-plagioclase geospeedometer, applicable to plagioclase-clinopyroxene bearing rocks, to be developed – this is the purpose of this contribution. The general approach is first outlined (Section 2) and followed by a detailed description of Mg-in-plagioclase geospeedometry (Section 3). Section 4 explores the applicability of this approach including a discussion of the robustness of the approach. Section 5 presents two case studies using gabbroic rocks from the lower oceanic crust. 2. DIFFUSION PROFILES AND THE EXTRACTION OF COOLING RATES – GENERAL CONSIDERATIONS The equilibrium distribution of a chemical element, i, between two phases, a and b, is governed by the partition coefficient, K a=b i , which is related to the equilibrium constant, K, of a suitably formulated exchange reaction that accounts for stoichiometric balance. K a=b is given by i ¼ C ai =C bi , where C ai etc. are concentrations of i in the K a=b i phase a in weight units (wt% or ppm). For practical convenience, the necessary conversion factors from molar units are absorbed within the constant itself. In cases where K a=b is a function of temperature T, K a=b i i ðT Þ, the equilibrium concentrations of i in a and b will change as a function of temperature. In this study we are interested in the quantity K Pl=Cpx ðT Þ. Mg If the kinetics of element exchange are sluggish relative to the timescale of thermal evolution of the system (e.g., cooling rate), then the concentrations in the crystals lag behind the concentration predicted by K a=b i ðT Þ. This kinetic lag comprises two parts: (i) the transfer of element i between a and b (inter-phase kinetics which may be governed by interfacial reaction rates, or, if the phases are physically separated from each other, by grain boundary transport rates), and (ii) the transfer of i within the phases a and b (volume or lattice diffusion within a crystal). The two processes occur in series and therefore, the slower

one governs the rate of overall exchange (e.g., Lasaga, 1998). How these two aspects of overall kinetics are linked to each other has been quantified and evaluated by Dohmen and Chakraborty (2003). They describe observations that may be used to identify which kinetic process is dominant in any given case. In most applications, a and b are chosen to be in contact with each other, so that the problem is reduced to an evaluation of interfacial element exchange kinetics versus volume diffusion rates. Two commonly made simplifying assumptions allow the distribution of i to be determined analytically based on the grain size and diffusion parameters (D0 and E, the ‘diffusion coefficient at infinite temperature’ and activation energy for diffusion, respectively) for i in a for a simple cooling history and simple grain geometries (Dodson, 1973, 1986). Firstly, it is assumed that the surfaces of the two minerals attain partitioning equilibrium instantaneously relative to the timescales of diffusion across the radius of a crystal; i.e., inter-phase exchange is assumed to be much more efficient than volume diffusion. Secondly, phase b is treated as an infinitive reservoir (Dodson, 1973) for i; this means that diffusion of i in phase a is the rate-limiting process and the concentration of i in phase b can be assumed to be constant. In nature this approximates the case where the concentration of i in b is much greater than in a and diffusion of i in b is sufficiently rapid that the small flux of i into or out of b as temperature changes can be homogenized across the phase. These assumptions have been relaxed to varying degrees in subsequent works. For example, when the surfaces of the two minerals are not in contact, but transport through the grain boundary is infinitely fast and efficient, the fast grain boundary model (FGB) can be used numerically to calculate cooling rates (e.g., Eiler et al., 1992; Jenkin et al., 1994). When the phases are in contact with each other but transport in the adjacent phase is not infinitely fast, it can be explicitly accounted for in the model of Lasaga (1983). However, a system with the two simplifying assumptions of Dodson is more easily amenable to an illustration of the general principles and this is the path we follow below. Based on the two assumptions given above, Fig. 1 illustrates how the content of species i in crystal a will evolve during cooling. For a given cooling history, a temperature (T) – time (t) path can be defined, which here is assumed to be linear (Fig. 1a). If the temperature dependence of the partition coefficient K a=b is known, K a=b i i ðT Þ along this T-t path is also known (Fig. 1b). Here we have implicitly assumed the common form of temperature dependence of partition coefficients dictated by thermodynamics: K(T) = K(1) exp(DH/RT), where K(1) is the partition coefficient at infinite temperature, DH is the standard enthalpy change of the partitioning reaction, and R is the gas constant. Hence, the evolution of the concentration of a species i in phase a in equilibrium with phase b along the chosen T–t path is determined (Fig. 1c), if the assumption of b as an infinite reservoir holds – this is the temperature dependent interface concentration for phase a. At high temperatures (e.g., T1 in Fig. 1), where diffusion coefficients are high, a remains almost homogeneous (Fig. 1d). However, as the diffusion coefficient decreases

K. Faak et al. / Geochimica et Cosmochimica Acta 140 (2014) 691–707

Temperature T

(a)

T1 T2 T3 T4 T5

Time t

(b)

T3

T2

T4

T5

l n Ki

/

T1

1/T Tccore

(c)

T2

T3

Concentrationi

T1

Tcbulk Tcrim T4

T5

core bulk rim equi

Time t

Concentrationi

(d)

rim

core

rim

Distance Fig. 1. Schematic illustration of the evolution of diffusive concentration profiles during cooling on a given time (t) – temperature (T) – path. Panel (a) shows a t–T-plot with an assumed linear decrease of temperature with time, i.e., a constant cooling rate. Panel (b) shows a plot of the partition coefficient K a=b against 1/T, for the i a=b case of decreasing K i with decreasing T. As the different temperatures (T1 to T5) were equally spaced, they are unevenly spaced when plotted as 1/T. Panel (c) shows the evolution of the concentration of element i in mineral a with time on the given t–Tpath. The figure is shown for a given cooling rate and a given grain size of the crystal a, but it is noted that Tccore and Tcrim are a function of cooling rate and grain size and will be different as these parameters are changed. Panel (d) shows the successive evolution of concentration profiles of element i in a grain of mineral a during cooling. (For interpretation of the references to color in this figure, the reader is referred to the online version of this article.).

693

exponentially with decreasing temperature, there will be a temperature at which the core composition can no longer maintain equilibrium with the changing interface composition (e.g., T2 and T3) and phase a becomes zoned in i (Fig. 1d). With progressive cooling, at some temperature (e.g., T4) diffusion ceases to significantly change the concentration of i in a and the concentration becomes frozen (defining a closure profile). Dodson (1973) formally introduced the concept of the closure temperature, Tc, of a cooling geochronological system and provided a mathematical formulation for calculating Tc based on geometry and size of the mineral grains, cooling rate of the host rock, and diffusion kinetic properties of the geochronological system in the mineral of interest. Because a becomes zoned in i, but b is an infinite reservoir and thus remains homogeneous, a Tc can be determined for each point across the crystal of interest Tc(x), decreasing towards the rim (Onorato et al., 1981; Dodson, 1986; Fig. 1c). Thus, instead of a discrete closure temperature, there is a closure temperature profile with lower Tc at the rim (Tcrim) than at the core (Tccore) of the crystal (Fig. 1c) leading to the development of a concentration profile of species i in mineral a that will be zoned (Fig. 1d). Even under the assumption that the interface between a and b always maintains equilibrium, eventually volume diffusion becomes so slow that, within the resolution of the analytical tool of choice, no further changes in the distribution of i within a can occur (e.g., the concentration profiles are identical at T4 and T5). The term ‘rim’ is used here to refer to a small distance into the crystal from the interface where the exact distance will depend on the measurement technique used. The shape of the closure profile records information about the cooling history of the sample provided the assumptions described above were valid. The bulk closure temperature Tcbulk is determined from the volume weighted average of Tc(x) (where Tc(x) can be calculated for each measured concentration of i in a). In simple systems, where only one element is zoned, Tcbulk is equivalent to the mean closure temperature Tc as defined by Dodson (1973), which is calculated for the whole crystal without considering the concentration profile (Ganguly and Tirone, 1999). In this paper we use the notation Tcmean to refer to Dodson’s Tc in order to distinguish it from the other quantities related to closure temperature that are considered in this work. 3. A MATHEMATICAL MODEL FOR MG-INPLAGIOCLASE GEOSPEEDOMETRY The general case of an element undergoing temperature sensitive partitioning between two phases during cooling described in the previous section is expanded here to the specific case of temperature dependent exchange of Mg between plagioclase and clinopyroxene. Firstly, a diffusion model that accounts for gradients in both the concentration of Mg, and of major elements in plagioclase, in determining the diffusive flux of Mg in plagioclase is described. The boundary conditions for this diffusion model are then explored followed by a discussion of the initial conditions. With these three aspects quantified, the distribution of Mg in plagioclase can be determined for any cooling history

694

K. Faak et al. / Geochimica et Cosmochimica Acta 140 (2014) 691–707

given the assumptions outlined in Section 2 of instantaneous surface equilibrium with clinopyroxene, and clinopyroxene acting as an infinite reservoir for Mg. Subsequently, we consider the effects of relaxing these assumptions. 3.1. A diffusion model for Mg in plagioclase – coupling with major components The partitioning of Mg between plagioclase and clinopyroxene has been demonstrated to depend on temperature and the major element composition of the plagioclase (quantified as anorthite content or XAn), which influences the activity of the Mg-bearing component in plagioclase (Blundy and Wood, 1994; Bindeman et al., 1998; Faak et al., 2013). As a result, plagioclase crystals with heterogeneous anorthite contents, that have equilibrated their Mg contents with the surroundings, are zoned in Mg with such zoning correlated to the major element chemistry of the plagioclase (i.e., zoning of XAn). Diffusion rates of major elements (Na, Si, Ca and Al) are so sluggish that limited changes in XAn occur in most geological settings (Grove et al., 1984; Morse, 1984; Liu and Yund, 1992). For rocks held at high enough temperatures for a long enough duration such that XAn could change via diffusion, the more than 5 orders of magnitude faster diffusion of Mg (LaTourrette and Wasserburg, 1998; Faak et al., 2013) implies that Mg would always redistribute itself in response to changes in XAn much more rapidly than diffusion can modify XAn. Thus, when modeling the distribution of Mg in plagioclase the XAn distribution can be considered fixed. As a result, the distribution of Mg in plagioclase is controlled by the changing boundary condition (see below) as well as the fixed (observed) distribution of XAn in the crystal. In terms of modeling the diffusion of Mg out of plagioclase during monotonous cooling, this means that the total flux has two separate components that need to be considered. The first is driven by the concentration gradient generated by the decreasing Mg content at the plagioclase rim (i.e., changing boundary condition) and the second is the diffusive flux associated with the chemical potential gradients generated by the variation of XAn within the plagioclase crystal. Costa et al. (2003) describe the total flux of Mg resulting from these two components as: J Mg ¼ DMg

@C Mg @X An þ LMg A @x @x

where t = time. This equation may be solved numerically to track the evolution of Mg concentration due to diffusion in a zoned plagioclase crystal (Electronic Annex A). Solving Eq. (1) requires knowledge of the diffusion coefficient at any temperature which, in the simplest case, can be determined from the parameters Do and E. Faak et al. (2013) report these parameters determined from experiments that were specifically designed to investigate the diffusive exchange of Mg between plagioclase and clinopyroxene. Their results are consistent with the earlier studies of LaTourrette and Wasserburg (1998). However, Faak et al. (2013) also found that the Mg diffusion coefficient in plagioclase depends on the silica activity of the system, aSiO2 , and provide a relationship of the form m E DPl where m is a constant. They Mg ¼ D0 exp ð RT Þ ðaSiO2 Þ did not observe a dependence of the diffusion rate of Mg in plagioclase on the XAn content of plagioclase. Such a XAn-dependence of DPl Mg was, however, postulated by Costa et al. (2003), based on analogy with the compositional dependence of DPl Sr (Giletti and Casserly, 1994; Cherniak and Watson, 1994), and has been reported by Van Orman et al. (2014). The solution to Eq. (1) depends on whether DPl Mg is considered to be dependent on XAn, aSiO2 , and/or other parameters, as discussed in Electronic Annex A. 3.2. Boundary conditions In order to model the distribution of Mg in plagioclase for any given cooling history using Eq. (1) the boundary conditions defining the interface concentration of Mg in plagioclase as a function of time is required. If there is evidence that plagioclase and clinopyroxene maintained exchange equilibrium at their interface throughout cooling (see Section 2), then Eq. (1) can be solved using the equilibrium partition coefficient to define the boundary conditions at each rim as a function of temperature:  C Pl  Mg K Pl=Cpx ¼ ¼ f ðT Þ ð2Þ  Mg interface C Cpx Mg Pl=Cpx The partition coefficient, K Mg , has been determined experimentally by Faak et al. (2013) and is given by the relationship:

ln K Pl=Cpx ¼ 9219 ½K Mg

where DMg = diffusion coefficient of Mg in plagioclase, CMg = concentration of Mg in plagioclase, x = distance a=b and A = factor to describe the dependence of K i on XAn. Their factor LMg is a phenomenological coefficient, equivalent to (DMgCMg)/RT (see Costa et al., 2003 for details). The resulting diffusion equation needed to describe the flux of Mg to account for both of these driving forces was presented by Costa et al. (2003) as their Eq. (7):   @C Mg @ 2 C Mg @C Mg @DMg ¼ DMg þ 2 @t @x @x @x   A @C Mg @X An @DMg @X An @ 2 X An  þ C Mg þ DMg C Mg DMg 2 @x @x @x @x @x RT

þ ln aSiO2

ð3Þ

If the assumption that clinopyroxene acted as an infinite reservoir is justifiable (see Section 2) then this can be rearranged to give the interface plagioclase composition as a function of the measured clinopyroxene composition: Cpx C Pl Mg ¼ C Mg exp½9219 ½K

þ

ð1Þ

1 16913 ½Jmol1  X An þ 1:6 þ T RT

1 þ 1:6 T

16913 ½Jmol1  X An þ ln aSiO2  RT

ð4Þ

Solving Eqs. (3) and (4) requires knowledge of the silica activity as a function of temperature. In many systems of interest this is either buffered by a given mineral assemblage

K. Faak et al. / Geochimica et Cosmochimica Acta 140 (2014) 691–707

or can at least be constrained to lie between the values given by two buffer reactions (e.g., Carmichael et al., 1970). Faak et al. (2013) show the effect of silica activity on K Pl=Cpx for Mg the case of a system buffered either by pure SiO2 or by coexisting olivine–orthopyroxene (Ol + Opx). For example, at a temperature of 1150 °C this corresponds to ln aSiO2 = 0 and from ln aSiO2 = 0.6, respectively, which changes ln K Pl=Cpx Mg 4.1 to 4.6 (for XAn = 0.6). A more detailed example of how the silica activity can be estimated for a given system is discussed in Section 5, and the assumptions of instantaneous surface equilibration and clinopyroxene acting as an infinite reservoir are scrutinized below. 3.3. Initial conditions Modeling the redistribution of Mg between plagioclase and clinopyroxene during cooling requires an estimate of the initial conditions for the system. How this is done will depend on the system being considered. If the combination of peak temperature, grain size and cooling rate implies that the system would lose its memory of initial conditions entirely by diffusion during the cooling process, then the final results are independent of the choice of an initial concentration distribution. In these cases, any arbitrary initial distribution may be taken; a reasonable approach is to calculate an equilibrium concentration profile of Mg at a sufficiently high temperature using Eq. (4). Such high enough temperatures that can be used as initial temperatures to start the modeling (Tstart) may be identified using several approaches. Dodson (1973, 1976, 1986) derived his equations for closure temperatures assuming that the conditions of loss of memory of initial conditions are fulfilled. He found that such conditions may be identified on the basis of the cooling time constant, s (defined as the time taken for the diffusion coefficient, D, to decrease to e1D, or equivalently, E/RT to increase by 1). Dodson suggested that concentrations at temperatures obtained at 5s before closure are completely erased by diffusion and this may be used to determine a Tstart, from Dodson’s equations, given by 1 E ð  5Þ E RT c

T start ¼ R

ð5Þ

Using the activation energy of 321 kJmol1 from Faak et al. (2013), we obtain Tstart = 1020 °c for a Tc of 900 °C and a Tstart = 650 °C for a Tc of 600 °C. Tc may be obtained, for example, from the calculated Tcbulk (see below). An alternative approach is to use the extended Dodson formulation of Ganguly and Tirone (1999), who specifically evaluated conditions at which the memory of initial conditions is lost in a diffusing system. They identified a dimensionless parameter, M ¼ DðT 0 Þ s=a2 , where D(T0) is the diffusion coefficient at peak temperature, T0, a is the half width of the crystal, and s is the quantity defined above. For a linear change of 1/T with time (t), s can be expressed as s ¼ R=ðE sÞ, with R being the gas constant, E being the activation energy for diffusion, and s being the cooling rate at peak temperature T0. Ganguly and Tirone (1999) found

695

that for values of M > 1.1 the system loses its memory of initial conditions (Fig. 2 in Ganguly and Tirone, 1999) for a slab geometry of the diffusing system, which is considered appropriate for plagioclase crystals. Their equation for M can be recast to obtain the minimum value of peak temperature (T0) for which cooling would erase the memory of initial conditions in any given cooling system: sffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi 1:1 a2 E s T0 ¼ ð6Þ DðT 0 Þ R Using the data from Faak et al. (2013) for the diffusion coefficient of Mg in plagioclase (for a silica activity of 1), a crystal half width of a = 1 mm and a cooling rate of 0.1°Cyr1, it is found that for systems with peak temperatures above 1100 °C (i.e., most mafic igneous systems) the memory of the initial conditions would be erased by diffusion. For metamorphic systems, cooling rates are likely to be much slower, implying that the memory of initial conditions would be erased at even lower temperatures. For instance, for a crystal half width a = 1 mm and a cooling rate s = 1000 °CMyr1 (103 °Cyr1) T0 is 915 °C; for the same grain size, but even slower cooling of 0.1 CMyr1 (107 °Cyr1) T0 is 665 °C. Setting T0, or any temperature higher than that as Tstart would make the choice of initial conditions insignificant; an equilibrium partitioning of Mg, obtained using Eq (4), would be appropriate in these cases. Section 5.1 gives examples of such calculations. Estimating the peak temperature of the system in order to understand the appropriate initial conditions can be achieved in different ways for different systems. In igneous rocks the peak temperature can be estimated based on the temperature of co-saturation of plagioclase and clinopyroxene in a representative parental melt composition. In a metamorphic rock, peak temperature may be estimated from phase equilibria constraints or geothermometry using systems with closure temperatures above the peak temperature. For cases where the memory of initial conditions are not erased and initial equilibrium conditions cannot be assumed, some idea of the initial profile shapes of Mg in plagioclase may still be obtained by studying the profile shapes of other, slower diffusing trace elements (e.g., Costa et al., 2003). 4. EVOLUTION OF MG-CONCENTRATION PROFILES IN PLAGIOCLASE – DISCUSSION OF MODEL SYSTEMS 4.1. Mg-profiles for different zoning patterns of anorthite in a crystal Based on the model described in Section 3 some examples of the evolution of Mg concentration profiles in plagioclase during cooling are given in Fig. 2. These models assume instantaneous surface equilibrium between plagioclase and clinopyroxene and that clinopyroxene acts as an infinite reservoir. The models shown are for different cooling rates (dT/dt = 0.1 and 0.001 °C yr1), and different shapes of XAn profiles in the plagioclase (P1: flat, P2: stepped normal zoning, and P3: stepped reverse zoning).

696

K. Faak et al. / Geochimica et Cosmochimica Acta 140 (2014) 691–707

(g)

(d)

(a) 0.8

XAn

0.7 0.6 0.5

MgO [wt%]

0.4

rim

core

dT/dt 0.1°Cyr-1

(e)

dT/dt 0.1°Cyr-1

(h)

dT/dt 0.1°Cyr-1

dT/dt 0.001°Cyr-1

(f)

dT/dt 0.001°Cyr-1

(i)

dT/dt 0.001°Cyr-1

rim

core

(b)

(c)

0.20

rim

rim

rim

core

rim

0.15 0.10 0.05

MgO [wt%]

0.20 0.15 0.10 0.05

0

500

1000

1500 0

Distance [µm] T0=1200°C

T1=1100°C

500

1000

1500 0

Distance [µm] T2=1000°C

T3=800°C

T4=600°C

500

1000

1500

Distance [µm] diffusion profile

partitioning profile

Fig. 2. Schematic evolution of Mg-profiles in plagioclase during cooling in contact with clinopyroxene (MgO = 14 wt%) for three different assumed plagioclase crystals (length = 1500 lm, plane sheet geometry). The silica activity is assumed to be buffered by the assemblage of Ol + Opx in the system (see Section 5.1 for further discussion). P1 (a–c) has a flat XAn profile, P2 (d–f) and P3 (g–h) have stepped XAn profiles that are normally and inversely zoned. Panels (b), (e) and (h) show the calculated Mg-profiles at different temperatures (red = 1200 °C, pink = 1100 °C, orange = 1000 °C, green = 800 °C and blue = 600 °C) for a linear cooling rate of 0.001 °Cyr1. Panels (c), (f) and (i) shows the calculated Mg-profiles for the same temperatures as above, but for faster linear cooling of 0.1 °Cyr1. Dashed lines show the calculated partitioning profiles using Eq. (4), assuming equilibration between plagioclase and clinopyroxene at the respective temperatures. Some solid lines are thicker only for visibility reasons of overlapping lines. (For interpretation of the references to color in this figure, the reader is referred to the online version of this article.).

The models assume a constant cooling rate starting from an equilibrium distribution of Mg at 1200 °C, which is above Tstart and T0 as defined above (Eq. (5) and (6)), and are run down to 600 °C which is below the closure temperature. Plagioclase P1 (Fig. 2a–c), which is unzoned in XAn, will also be unzoned in Mg if the Mg-content of the plagioclase is in equilibrium with the clinopyroxene (dashed lines in Pl=Cpx Fig. 2b and c). During cooling, K Mg decreases (Eq. (3)) meaning that the plagioclase Mg concentration at the interface decreases; if the diffusion coefficient for Mg in plagioclase is sufficiently high (i.e., at high temperatures) Mg can diffuse down the concentration gradient to the interface and the plagioclase will remain effectively unzoned in Mg (i.e., the zoning will be negligible relative to analytical precision). This is the case down to 1000 °C for a cooling rate of 0.001 °C yr1 for the example given in Fig. 2c but even by 1100 °C this is no longer the case for the faster cooling example (0.1 °C yr1) given in 2b. As temperature

decreases, the flux of Mg from the core of the plagioclase towards its rim can no longer keep pace with the changing plagioclase interface composition and the Mg-content starts to become zoned. With further cooling the diffusive fluxes within the interior of the plagioclase become negligible such that the concentration of Mg in plagioclase does not change (again, within analytical precision). The final concentration profile is convex upwards (i.e., higher concentrations of Mg in the core than at the rims; green and blue lines in Fig. 2b and c) and, importantly, the concentration of Mg in the faster cooled example (Fig. 2b) is much higher than in the slower cooled sample (Fig. 2c). This difference is the basis of being able to extract cooling rates from measuring compositional profiles. In this case of a plagioclase crystal that is not zoned in XAn, Tcmean (mean closure temperature calculated for the whole crystal according to Dodson, 1973) is equivalent to a Tcbulk (volume weighted average of the closure profile Tc(x)). Tc(x) can be determined

K. Faak et al. / Geochimica et Cosmochimica Acta 140 (2014) 691–707

(a) 0.8

XAn

0.6 0.4 0.2 rim

(b) 0.20

MgO [wt%]

from the Mg-concentrations at each point along the final concentration profile using the thermometer of Faak et al. (2013) and the respective XAn at each point. Application of the two methods yields Tcmean = Tcbulk = 978 °C for a cooling rate of 0.1°Cy1 and Tcmean = Tcbulk = 822 °C for a cooling rate of 0.001°Cy1. The situation is different for a crystal that is zoned in XAn. The evolution of the concentration profile of Mg in plagioclase crystals zoned in XAn is a little more complex than that for plagioclase P1 because, as discussed above, the chemical potential of Mg in plagioclase depends on XAn (Blundy and Wood, 1994; Bindeman et al., 1998; Costa et al., 2003; Faak et al., 2013; Eq. (3)). Using Eq. (4) to define the plagioclase starting composition at 1200 °C the Mg concentration profile has the same form as the XAn profile (Fig. 2d and 2 g). During cooling the plagioclase interface Mg contents decrease according to Eq. (4) and the Mg content of the interior of the plagioclase follows. However, unlike in the case of plagioclase P1 the Mg concentration profile is not smoothly curved due to the variability in XAn. The absolute difference in concentration between the plagioclase rim and the core becomes smaller for lower temfor the different peratures, because the difference in K Pl=Cpx Mg XAn at the core and the rim becomes smaller with decreasing temperature (see Eqs. (3) and (4)). In the cases of plagioclase P2 and P3, which are zoned in XAn, diffusion during cooling still leads to Mg-profile shapes that are convex upwards (Fig. 2e, f, h and i). As seen from Eq (4), the distribution of Mg is positively coupled with the XAn content and therefore, plagioclase that is normally zoned in XAn (P2, Fig. 2e and f) will also be normally zoned in Mg when the Mg distribution within the plagioclase is at equilibrium (Fig. 2e and f). In this example the low XAn in the outer portion of the normally zoned plagioclase leads to relatively lower MgO in this portion of the crystal than in the reverse zoned plagioclase. The gradients in XMg set up by changing boundary conditions are thus, relatively smaller leading to a smaller diffusion flux of Mg out of the plagioclase over a given cooling interval. A more direct comparison of the resulting Mg-profile shapes for plagioclase crystals with different XAn-zoning patterns is shown in Fig. 3. Even though all three examples have the same XAn content at the core (XAn = 0.5), the Mg-concentration at the core resulting from diffusion along the same cooling history is significantly different (up to 0.03 wt% MgO; Fig. 3) because each of the three crystals has a different XAn-concentration gradient. This comparison illustrates the importance of analyzing and modeling full zoning profiles rather than just fitting for the Mgcomposition in the core or for the bulk crystal. Clearly the latter approach would lead to spurious results. Also, the calculation of Tcmean following Dodson (1973) is not straightforward in the case of zoned XAn, as the Dodson equation does not account for the additional diffusive flux of Mg due to a gradient in XAn (Eq. (1)). However, a Tcbulk can be determined from the Mg-concentrations of the closure profile and the respective XAn content at each point. In summary, for a linear cooling history, the model described above generally produces Mg-closure-profiles that are convex upwards, and slower cooling rates result

697

core

rim

dT/dt = 0.1°Cyr-1; T = 600°C

0.15 0.10 0.05 0

500 1000 Distance [µm]

1500

Fig. 3. Direct comparison of the Mg-profiles resulting from cooling at the same rate (0.1°Cyr1 from 1200 to 600 °C) for three plagioclase crystals that have different zoning patterns in XAn as shown in part (a). Calculated Mg-profiles for these three plagioclase crystals (b) demonstrate that even though they have equal XAn content at the core, the core Mg-concentrations are significantly different (by up to 0.03 wt% MgO). Thus, spot analyses cannot be used to extract cooling rates.

in lower Mg-concentrations in the interior of the plagioclase crystal. 4.2. Mg-profile shapes as a function of geometry of cooling path All Mg-profiles shown in Fig. 2 and 3 were based on the assumption of a constant cooling rate. The following section illustrates the effect of some alternative cooling paths on the resulting shape of the Mg-profile in plagioclase assuming a homogeneous distribution of XAn in the plagioclase. Fig. 4 shows six different cooling paths (a, d, g, j, m, and p) and the resulting Mg-profiles after cooling along these paths from 1200 to 600 °C along a slower cooling history (Fig. 4b, e, h, k, n and q) and a faster cooling history (Fig. 4c, f, i, l, o, and r). Previous work by Dodson (1973,1976) argued that the exact cooling history is not critical except at temperatures around Tc. Consistent with this, for a generally slower cooling history (middle panel in Fig. 4), and thus, lower Tc, the system is more sensitive to changes in the low T range of the cooling path (Fig. 4h, k, n and q). For faster cooling (lower panel in Fig. 4) and thus, higher Tc, the resulting Mg-profile shape is not very sensitive to changes in the cooling history at lower temperatures (Fig. 4i, l and o), but is sensitive to changes in the cooling path at high temperatures (Fig. 4f and Fig. 4r). However, even for slow cooling, the memory of the high temperature cooling history (but below Tstart) may not be entirely erased by diffusion along the low temperature cooling history (e.g., Fig. 4e).

Temperature [°C]

1200

MgO [wt%]

K. Faak et al. / Geochimica et Cosmochimica Acta 140 (2014) 691–707

0.10 0.08 0.06 0.04 0.02

MgO [wt%]

698

0.10 0.08 0.06 0.04 0.02

(a)

(d)

(g)

(j)

(p)

(m)

1000 Tcbulk=980°C

800 Tc =820°C bulk 600 200.000 400.000 600.000 2000 4000 6000

Time [years] dT/dt 0.001°Cyr-1

200.000 400.000 600.000 2000 4000 6000

Time [years]

200.000 400.000 600.000 2000 4000 6000

Time [years]

200.000 400.000 600.000 2000 4000 6000

Time [years]

200.000 400.000 600.000 2000 4000 6000

Time [years]

200.000 400.000 600.000 6000 2000 4000

Time [years]

(b)

(e)

(h)

(k)

(n)

(q)

(c)

(f)

(i)

(l)

(o)

(r)

dT/dt 0.1°Cyr-1 500 1000 1500 Distance [µm]

500 1000 1500 Distance [µm]

500 1000 1500 Distance [µm]

500 1000 1500 Distance [µm]

500 1000 1500 Distance [µm]

500 1000 1500 Distance [µm]

Fig. 4. Illustration of the effect of different cooling histories on the resulting shape of the Mg-profile. For simplicity, a homogeneous anorthite-content of XAn = 0.4 is assumed for all examples. All plagioclases are 1500 lm long and assumed to maintain interface equilibrium with a Cpx with 14 wt% MgO during cooling. All models were run from 1200 °C to 600 °C and silica activity was assumed to be buffered by the coexistence of Ol + Opx. The upper panel (a, d, g, j, m, p) shows the cooling paths as the evolution of temperature with time. Note that there are two time scales – the grey numbers refer to the cooling history used to produce the Mg-profiles in the middle panel (b, e, h, l, n and q) and the black numbers refer to the lower panel (c, f, i, l, o, and r). Black solid lines refer to the cooling path used to produce the Mg-profiles shown below each cooling path. For easy comparison, results from the two linear examples used in Section 4.1 and Fig. 2 (dT/ dt = 0.001 °Cyr1 and 0.1 °Cyr1, respectively) are repeated in Fig. 4a–c and plotted as grey dashed lines in all other examples. Grey dotted lines in (m), (n) and (o) additionally provide comparison to the cooling path and respective Mg-profiles in (g), (h) and (i). Grey dotted lines in (p), (q) and (r) additionally provide comparison to the cooling path and respective Mg-profiles in (m), (n) and (o). Tcbulk was calculated after Dodson (1973) for a plane sheet geometry, using the activation energy from Faak et al. (2013) and the linear cooling rates of 0.1 °Cyr1 (Tcbulk = 980 °C; black) and 0.001 °Cyr1 (Tcbulk = 820 °C; grey).

As illustrated in Section 4.1, a constant cooling rate always produces convex Mg-closure profiles. This general shape is modified significantly by diffusion along different cooling paths. For a cooling path with progressively slower cooling with time (Fig. 4d), the Mg-profiles become more convex (Fig. 4e and f). Increasing cooling rates with time (Fig. 4g and j) lead to a significant flattening of the profile shape (Fig. 4h and k) unless this change in cooling rate happens below the closure temperature (Fig. 4i and l). Even flatter profiles are produced along cooling paths where the temperature is held constant at a certain temperature above the closure temperature for an extended period and then cooled rapidly to below the closure temperature (Fig. 4m–r). It should be clear from the discussion above, and the variation of profile shape shown in Fig. 3, that measuring and modeling full zoning profiles (for Mg and XAn) is essential if the cooling history is being investigated. Additionally, the detailed shape of the concentration profiles records a memory of the shape of the cooling path, and in favorable circumstances it may be possible to extract more than a single cooling rate from a frozen concentration profile. 4.3. Evaluating the sensitivity and robustness of the Mg-inplagioclase geospeedometer There is a substantial literature on inverting a measured concentration profile into a cooling history (e.g., Dodson, 1973, 1986; Onorato et al., 1981; Lasaga, 1983; Jenkin et al., 1994; Dohmen and Chakraborty, 2003). The discussion here focuses on the Mg-in-plagioclase geospeedometer

and the reader is referred to these sources for a more general treatment of geospeedometry. 4.3.1. Effect of using a ‘rim’ composition instead of the real interface Using an electron microprobe (discussion of the appropriate electron microprobe analytical setup for determination of MgO in plagioclase is given in Electronic Annex B), it is not possible to analyze the composition at the actual interface of the plagioclase crystal and the outermost position where the composition may be measured is several lm away from the interface. We denote this outermost position that is practically accessible to chemical analysis using a given tool as the ‘rim’. Fig. 5 illustrates how different distances between the analyzed rim and the interface affect the retrieved cooling history. A concentration profile (bold grey line in Fig. 5) is calculated for cooling at dT/ dt = 0.1°Cyr1 to room temperature with equilibrium maintained throughout at the interface, i.e., there is effectively no MgO at the interface at the end of cooling. However, a few lm away from this interface the concentration is already significantly higher (e.g., 0.01wt% MgO at a distance of 5 lm), as the diffusive flux within the plagioclase is very small at low temperatures. The concentration profiles (taken to be 500 lm from rim to core) that would be measured with the rim at 2, 5, 10, or 15 lm away from the interface, respectively, may be fitted to obtain cooling rates. It is found that a constant cooling rate of dT/dt = 0.1°Cyr1, i.e., the value used to generate the synthetic profile, is retrieved in each of these cases. The rim composition freezes at the temperature Tcrim (see inset in Fig. 5)

MgO [wt%]

K. Faak et al. / Geochimica et Cosmochimica Acta 140 (2014) 691–707

0.06 Tcrim=880°C Tcrim=850°C

0.04

real profile 2 µm from interface 5 µm from interface

Tcrim=750°C

0.02 Tcrim=660°C

10 µm from interface 15 µm from interface fitted profile

0 5 10 15

0

200

400 600 Distance [µm]

800

1000

Fig. 5. Illustration of the effect of analyzing a rim composition that is not at the exact interface of the plagioclase crystal. The bold grey line represents the real concentration profile, if the plagioclase crystal would have cooled with dT/dt = 0.1 °Cy1 down to room temperature. For simplicity, a homogeneous XAn content of 0.4, MgO in the Cpx = 14 wt%, and a silica activity of the system buffered by the assemblage Ol + Opx was assumed. Tstart was 1200 °C. Profiles that would be measured starting at different distances from the actual interface are shown by different symbols (enlarged in inset). These profiles were modeled using methods described in the text and the parameters noted above. Best fits, shown as black lines, were obtained for dT/dt = 0.1 °Cy1 in all cases. The temperatures at which the respective rim compositions freeze, Tcrim, are shown beside the respective symbols in the inset.

and the profile shape cannot change below this temperature. Therefore, quantitative information on the cooling history below Tcrim cannot be retrieved from such a profile. Pl=Cpx 4.3.2. Propagation of uncertainties in DPl into Mg and K Mg cooling rate estimates The choice of the diffusion coefficient DPl Mg and the partiPl=Cpx tion coefficient K Mg are crucial for the parameters and input conditions of the diffusion model, and uncertainties in these parameters will affect the accuracy of the cooling rates obtained from the model. The propagation of Pl=Cpx uncertainties in DPl into Eq. (1) is not straightforMg and K Mg ward, as their effects will depend, for example, on the fluxes introduced by different XAn distributions in plagioclase. However, an estimate of the effect of the uncertainties in Pl=Cpx on the cooling rate may be obtained from DPl Mg and K Mg assuming a simple scenario of a plagioclase with unzoned XAn, in a system in which the silica activity is well constrained (i.e., the uncertainty in aSiO2 can be neglected), and a cooling history that is linear in 1/T. In this case, Pl=Cpx can be propagated the uncertainties on D0, E and K Mg into the Dodson equation (Dodson, 1973) to obtain the uncertainty on the cooling rate at Tcbulk (e.g., determined from the volume weighted average concentration of Mg along a profile). Errors were propagated using a Monte Carlo approach. The uncertainty on K Pl=Cpx are simulated Mg by drawing randomly from a probability distribution with a standard deviation of ±20 °C around Tcbulk based on the uncertainty on the thermometer (Faak et al., 2013). The uncertainty in DPl Mg comprises the uncertainties in the activation energy E and in the pre-exponential factor D0, but these two parameters cannot vary independently. To

699

propagate the uncertainty on DPl Mg we draw randomly from the probability distribution for E (321,000 ± 40,000 J mol1) determined by Faak et al. (2013) and then determine the value of D0 for this E that best fits the experimental data set of Faak et al. (2013; all of their data points for log DPl Mg were recalculated for a silica activity of 1 for the purpose of these simulations). The inverse procedure, of drawing D0 from its probability distribution and then calculating E was also used. For a cooling rate of 0.1 °Cyr1, or log(dT/dt) = 1, the uncertainties propagate to give log(dT/dt) = 1 ± 0.25 (i.e., the cooling rate is constrained to be within the range 0.05 to 0.18 °Cyr1). For a slower cooling rate of 0.001 °Cyr1 (i.e., log(dT/dt) = 3) the propagated uncertainties yield log(dT/dt) = 3 ± 0.3 (i.e., the cooling rate is constrained to be within the range 0.002 to 0.0004 °Cyr1). These uncertainties do not significantly affect relative cooling rates determined for different samples if the same parameters for the calculation of DPl Mg are used for all samples. If the plagioclase crystal is zoned in XAn it will be crucial, whether the chosen DPl Mg has a significant dependence on XAn (Costa et al., 2003; Faak et al., 2013; Van Orman et al., 2014). To illustrate the effect of a XAn dependence of DPl Mg on the obtained cooling rate, we re-modeled plagioclase crystals P2 and P3 (see Section 3.1), assuming that the diffusion coefficient had a negative correlation with XAn of E m the form DPl  10n ðX An 0:6Þ , Mg ðX An Þ ¼ D0 exp ð RT Þ ðaSiO2 Þ where D0, E and m from Faak et al. (2013) are used, to keep the data comparable to the previous calculations. The factor n describes the XAn dependence of DPl Mg and is assumed to be 3, following Van Orman et al. (2014). As the data from Faak et al. (2013) are mainly determined for XAn = 0.6, this is used as a reference state. For the same cooling rate dT/dt = 0.1°Cyr1, the closure profile for plagioclase P2 (normally zoned in XAn) shows lower Mg-concentrations when it is modeled with a DPl Mg ðX An Þ that is negatively correlated to XAn (solid black line in Fig. 6a). This is expected, because the rims of this plagioclase have XAn < 0.6 and therefore, diffusion of Mg is relatively faster in this area. In turn, plagioclase P3, which is inversely zoned with XAn > 0.6 at the rims, has higher Mg-concentration when modeled with DPl Mg ðX An Þ (solid black line in Fig. 6b). The difference in the cooling rate, dT/dt, that would be retrieved from models with different diffusion coefficients was investigated by varying cooling rates in a model with Pl DPl Mg ðX An Þ until the Mg-profile from the model with DMg was reproduced best. With DPl ðX Þ, the best possible An Mg description of the profile shape calculated using DPl Mg and dT/dt = 0.1 °Cyr1 (dashed black lines in Fig. 6a and b) is obtained for a cooling rate that is a factor of 3 higher (Plagioclase P2, solid grey line in Fig. 6a) or 4 lower (Plagioclase P3, solid grey line in Fig. 6b), although the profile shapes cannot be completely reproduced (compare dashed black lines and solid grey lines in Fig. 6a and b). 4.3.3. How robust is the assumption that clinopyroxene acts as an infinite reservoir? The modeling and discussion above assume that clinopyroxene acts as an infinite reservoir (Dodson, 1973) and

700

K. Faak et al. / Geochimica et Cosmochimica Acta 140 (2014) 691–707

(a)

(b)

0.8 0.6 DMg independent on XAn

XAn

0.6

DMg dependent on XAn

0.20

MgO [wt%]

0.8

XAn

0.4

0.4

0.15 0.10

dT/dt 0.1°Cyr-1

0.05

0

500

dT/dt 0.3°Cyr-1 1000

Distance [µm]

dT/dt 0.1°Cyr-1 1500 0

500

dT/dt 0.025°Cyr-1 1000

1500

Distance [µm]

Pl Fig. 6. Comparison of Mg-concentration profiles obtained from a model where DPl Mg does not depend on XAn with one where DMg ðX An Þ is negatively correlated with XAn (see text for discussion). (a) For the same cooling rate dT/dt = 0.1 °Cy1 the Mg-profile in plagioclase P2 (normally zoned in XAn) shows higher Mg-concentration when diffusion was modeled with DPl Mg ðX An Þ (solid black line), compared to the Pl profile obtained from a model where DPl Mg does not depend on XAn (dashed black line). For DMg ðX An Þ the cooling rate dT/dt has to be adjusted to 0.3 °Cy1 to obtain the best match between the two models (solid grey line). (b) As plagioclase P3 has an inversely zoned XAn profile, the Pl Mg-profile obtained with DPl Mg ðX An Þ has higher Mg-concentration than the Mg-profile obtained with DMg . Therefore, dT/dt has to be 1 Pl decreased to 0.025 °Cy to match the Mg-profile from the model where DMg does not depend on XAn.

hence, it is assumed that the clinopyroxene composition remains fixed during the entire duration of exchange. Because the concentration of Mg in clinopyroxene is much higher than in plagioclase, this will be a good approximation if diffusion allows the redistribution of Mg through the clinopyroxene crystals. For example, assuming equal masses of plagioclase and clinopyroxene, if the plagioclase initially contains 0.2 wt% MgO, and after cooling contains 0.05 wt% MgO, then the bulk MgO content of the clinopyroxene only increases by 0.15 wt% which is a negligible increase under most circumstances. However, Mg diffusion in clinopyroxene is several orders of magnitude slower than in plagioclase (e.g., Zhang et al., 2010; Mueller et al., 2013) meaning that it is likely that the increase in Mg content of clinopyroxene is commonly concentrated in the rims of crystals leading to a larger increase in clinopyroxene Mg content than if the Mg was evenly distributed through the crystals. Additionally, if the mass of plagioclase is much larger than that of clinopyroxene the increase in average Mg content of clinopyroxene would be larger. Whether the clinopyroxene approximated an infinite reservoir for Mg should be investigated when applying the Mg-in-plagioclase geospeedometer by measuring detailed zoning profiles at the rims of the clinopyroxene crystals – if the crystals are unzoned in Mg then the assumption of an infinite reservoir is probably reasonable. 4.3.4. How robust is the assumption of instantaneous surface equilibrium? The discussion above has been predicated on the (commonly used) assumption that plagioclase and clinopyroxene attain surface equilibrium effectively instantaneously; i.e., relative to the rates of diffusion within the phases both the surface reactions, and any inter-granular transport, are very efficient. Dohmen and Chakraborty (2003) provide a detailed analysis of when this assumption is valid. For the

plagioclase–clinopyroxene system, diffusion of Mg is considerably slower in clinopyroxene than in plagioclase (Zhang et al., 2010; Mueller et al., 2013; Faak et al., 2013). In this case, following the analysis of Dohmen and Chakraborty (2003), the kinetics of diffusive uptake of Mg by clinopyroxene would govern the kinetics of the overall exchange process that involves the liberation or incorporation of Mg atoms at the surfaces of the crystals, as well as transport of these atoms through the intergranular medium (even if this is only a few atomic layers wide, as in a grain boundary – exceptions may be where clinopyroxene grain sizes are substantially smaller than plagioclase grain sizes). In a general system where the two minerals are not in direct contact with each other and communicate via a grain boundary network, evaluation of exchange kinetics requires knowledge of: (i) the rate law constants for surface exchange of Mg between the minerals and the grain boundary region; (ii) the partition coefficient of Mg between clinopyroxene and the grain boundary region; and (iii) the diffusion coefficient for Mg on the grain boundaries. However, for the case where the phases are in direct contact with each other (i.e., the distance between the grains is infinitesimally small) the situation becomes considerably simpler and whether equilibrium partitioning is attained at the interface depends entirely on the interfacial exchange kinetic rate (Dohmen and Chakraborty, 2003). In our case this is the kinetics of attachment or detachment of Mg atoms to clinopyroxene; this may occur by reactions involving transport of multiple components rather than a simple cation exchange such as: CaMgSi3 O8 þ CaAl2 SiO6 ¼ CaAl2 Si2 O8 þ CaMgSi2 O6

ð7Þ

Note that such a reaction may be valid irrespective of where Mg sits dominantly in the plagioclase structure, as long as some fraction of the total Mg occupies the T-site and it is likely that Mg dissolves in plagioclase in both the M- and

K. Faak et al. / Geochimica et Cosmochimica Acta 140 (2014) 691–707

the T-sites (Miller et al., 2006; Faak et al., 2013). If the temperature dependence of the exchange kinetics is stronger than the temperature dependence of diffusion rates, then in a cooling system interface kinetics may become rate determining below a certain temperature (i.e., there may be a crossover of rates) and prevent the interface from equilibrating. This would shut down the diffusive exchange process and therefore, adjacency of phases is not necessarily a guarantee for attaining equilibrium partitioning (see also Dohmen and Chakraborty, 2003). Unfortunately, the temperature dependence of attachment/detachment kinetics in clinopyroxene is not known directly and therefore, empirical means have to be used to evaluate whether this was a factor in specific settings. One empirical test of whether surface equilibrium was achieved near-instantaneously is distribution of the element(s) of interest within the exchanging phases (e.g., Dohmen and Chakraborty, 2003). If all rims of different sized grains of a given mineral in a thin section, or domain thereof, have the same concentration (e.g., as revealed by element mapping of multiple grains with different sizes), this is usually a good indicator that interfacial equilibrium was achieved. In the special case of Mg in plagioclase, adjustments would have to be made for the role of XAn content in governing equilibrium Mg concentrations. As shown in Fig. 2, if surface equilibrium is attained throughout then after slow, uniform cooling the Mg profile in plagioclase is predicted to be zoned towards a lower chemical potential of the Mg-component at the rims (lower XMg assuming constant XAn or equivalently lower Tc). If surface equilibrium is not achieved near instantaneously then the phases will tend to be weakly or negligibly zoned in Tc (e.g., Fig. 3 of Dohmen and Chakraborty, 2003). This may also be obtained, however, by a long anneal at a given temperature after a period of cooling or an increase in cooling rate over time (see Fig. 4). Differentiating whether such flat Tc profiles originate due to lack of surface equilibrium or sudden changes in cooling rate (Fig. 3) may require independent information. Detailed observation of whether the curvature of the Tc profile depends on what the adjacent phase is may help in some circumstances. In particular, if interfacial kinetics freezes the element exchange process, it may be expected that all rims of adjacent plagioclase– clinopyroxene pairs in a thin section would freeze at similar temperatures (the temperature of crossover of the diffusionand interfacial- kinetic rates) However, freezing of interface exchange between plagioclase and clinopyroxene does not automatically indicate that plagioclase-plagioclase exchange could not continue and lead to subsequent modification of the Mg distribution within plagioclase. It is also worth noting here that lack of adjacency in a two-dimensional cut through the rock (e.g., a thin section) may not reveal adjacent clinopyroxene in the third dimension and this needs to be borne in mind when interpreting Mg-concentration profiles in plagioclase. 4.3.5. How large are the uncertainties introduced by variable plagioclase geometry? Extracting time scales from 1D (or 2D) models of a process that occurred in 3D can introduce errors, because

701

the diffusive fluxes from dimensions, that are not being modeled, are neglected (e.g., Costa et al., 2003, 2008; Costa and Chakraborty, 2004). In general this leads to an overestimate of the time required to obtain a given extent of diffusive modification of a concentration distribution, i.e., this would lead to an underestimation of the cooling rate. To illustrate the effects of 2D vs. 1D modeling, three plagioclase crystals with aspect ratios of 1:1, 1:2 and 1:3 are considered (Fig. 7), where the short dimension x is 750 lm in all cases.. The Mg-concentration of the 2D model was traced along two perpendicular transects through the middle of the crystal, one parallel to the short dimension x of the crystal (Fig. 7b, e and h) and the other parallel to the long dimension y of the crystal (Fig. 7c, f and i). For comparison, the Mg-profiles along these transects were also modeled in a 1D diffusion model using a plane sheet geometry (grey dashed lines in Fig. 7). In the case of high aspect ratios (here 1:2 and 1:3, Fig. 7d and g), the Mg-profiles along the short dimension x are very similar for a 1D and a 2D model (Fig. 7e and h), suggesting that at aspect ratios above 1:2, diffusion along the longer dimension y does not significantly affect the concentration profile developed along the shorter dimension x. However, for the same cooling rate 1D-modeling of the crystals with high aspect ratio reveals significantly higher Mg-concentrations along the long dimension y, than a 2D model (grey dashed lines in Fig. 7f and i). As the dimensions of the crystal with an aspect ratio of 1:1 are equal, the Mg-profiles along these dimensions are equal, too. In this case, the Mg-concentrations from a 1D plane sheet model are slightly higher than those from a 2D model, as in the 2D model each profile is affected from diffusive fluxes in both dimensions. The uncertainty in the cooling rates resulting from using a 1D plane sheet model to extract cooling rates from the long axis of a plagioclase were determined by changing the cooling rate in the 1D model (grey dotted lines in Fig. 7c, f and i) until the developed Mg-profile along the long dimension matches with the respective Mg-concentration at the core from the 2D model. For the crystal with an aspect ratio of 1:1, the obtained cooling rate decreased from the ’true’ cooling rates of 0.1°Cyr1 to 0.06°Cyr1, for an aspect ratio of 1:2, the best fit cooling rate is 0.025°Cyr1 and for an aspect ratio of 1:3 the cooling rate had to be lowered to 0.01°Cyr1. As discussed above, modeling the diffusive flux only in one dimension and neglecting fluxes from the other dimensions will always underestimate the actual cooling rate. However, for the short dimension of crystals with aspect ratios P 1:2, the difference between the actual cooling rate and the cooling rate obtained from a 1D model becomes insignificantly small (Fig. 7e and h). A caveat to the above discussion is that 2D observations of thin sections are always prone to uncertainty due to cutting effects (e.g., Wallace and Bergantz, 2003). Given the preferred shape of plagioclase crystals, grains that are relatively large in both dimensions seen in a thin section compared to other grains in the same section, are likely to have a shorter 3rd dimension and should generally be

702

K. Faak et al. / Geochimica et Cosmochimica Acta 140 (2014) 691–707 ]

e

c an

m

[µ

e nc

ta Dis300

0.10

MgO [wt%]

MgO [wt%]

st 600 Di 300

0.05 0.00

x[

0.10 0.05

300

tan

ce

Dist

anc

600

MgO [wt%]

ey

]

(b) cross section along x

]

m x [µ ce tan 600 s i D 300

0.10 0.05 0.00 600

600

y[

µm

Dista

1200

nce

[µm1200 ]

y [µm

]

1800

(e) cross section along x

(h) cross section along x

0.10 0.08 0.06 0.04

2D model, dT/dt=0.1°Cyr

2D model, dT/dt=0.1°Cyr-1

2D model, dT/dt=0.1°Cyr

-1

-1

1D model, dT/dt=0.1°Cyr

1D model, dT/dt=0.1°Cyr

1D model, dT/dt=0.1°Cyr

-1

0.02

-1

150

0.12

MgO [wt%]

(g)

600

0.00

Dis

0.12

]

µm

(d)

MgO [wt%]

(a)

x

300

450

600

750

Distance x [µm]

(c) cross section along y

150

(f)

300

450

-1

600

750

Distance x [µm]

150

300

450

600

750

Distance x [µm]

(i) cross section along y

cross section along y

0.10 0.08 0.06 0.04 0.02

1D model, dT/dt=0.06°Cyr 150

300

450

600

1D model, dT/dt=0.025°Cyr

-1

750

300

Distance y [µm]

600

900

1200

Distance y [µm]

1D model, dT/dt=0.01°Cyr

-1

1500

600

1200

-1

1800

Distance y [µm]

Fig. 7. Illustration of the effects of a 2D-diffusion model versus a 1D diffusion model on the resulting Mg-profile shapes on the example of three hypothetical plagioclase crystals with different aspect ratios: (a–c) 1:1, (d–f) 1:2, and (g–i) 1:3. For simplicity, a homogeneous anorthite content of XAn = 0.6 is assumed for all plagioclases. Diffusion was modeled for cooling from 1200 to 600 °C with a constant cooling rate. Silica activity was assumed to be buffered by the coexistence of Ol + Opx and all plagioclase interfaces were assumed to be in equilibrium with Pl=Cpx clinopyroxene throughout cooling (14 wt% MgO); K Mg and DPl Mg were taken from Faak et al. (2013). Black solid lines show Mg-profiles along cross sections through the 2D model for dT/dt = 0.1°Cyr1. Grey dashed lines show the comparison to the respective Mg-profiles produced by a 1D model under the same cooling rate. Grey dotted lines were calculated using 1D models with different cooling rates to match the core compositions of the respective Mg-profiles obtained from the 2D models.

avoided when doing Mg-in-plagioclase geospeedometry. Even if such equi-dimensional crystals can be shown to be elongate in the third dimension, measured profiles would have to be fit using a 2D diffusion model to extract the correct cooling rates. Summarizing, if a 1D diffusion model is being used to fit measured plagioclase MgO contents and to extract cooling rates then the best approach is to analyze and model profiles across the short axis of plagioclase with aspect ratios of P2 (higher will be better especially if there is strong XAn zoning that will affect the diffusive fluxes of Mg). If 2D diffusion models are being used to fit the data then perpendicular profiles through high aspect ratio plagioclase should be measured and modeled but maps of the XAn distribution through the entire crystal will be needed if the crystal has complex Mg zoning. Ideally multiple crystals of different grain sizes in each sample should be analyzed to ensure consistency of the retrieved cooling rates.

4.3.6. Effects arising from dissolution-reprecipitation reactions and structural phase transitions The extraction of thermal histories from diffusion modeling is predicated on diffusion being the mechanism by which the element(s) of interest are transported The occurrence of dissolution-precipitation reactions that can disturb compositions (e.g., Ruprecht and Cooper, 2012) are typically recognized (e.g., irregular boundaries) using tools such as Nomarsky interference microscopy (NDIC, e.g., Singer et al., 1993), cathodoluminescence or backscatter electron imaging (e.g., see a concise discussion in Ginibre et al., 2007). Crystals, or at least the zonal boundaries, that are demonstrably irregular due to such dissolution-precipitation events should be avoided for diffusion modeling. It is sometimes found that plagioclase crystals become depleted in Mg by low-temperature, fluid mediated, dissolution-reprecipitation reactions as well as by diffusion during cooling (e.g., Vanko and Laverne, 1998). Even if there is no step in anorthite content, plagioclase crystals

K. Faak et al. / Geochimica et Cosmochimica Acta 140 (2014) 691–707

that have lost Mg by dissolution-precipitation can often be identified by a substantial depletion in Fe that can be identified through cathodoluminescence, electron microprobe analysis or via Normarsky interference microscopy. The modeling approach developed here relies on the diffusing medium remaining unmodified over relevant length scales during the duration of diffusion. Low temperature unmixing related structural transitions (e.g., the Bogild, Huttenlocher or Voll gaps in anorthitic plagioclase at temperatures below about 800 °C; the peristerite gap at lower temperatures of about 500 °C in albitic plagioclase) are known to result in the production of a variety of mosaic textures in plagioclase (e.g., Carpenter, 1994; Brown and Parsons, 1994) that are typically visible only on the TEM scale. The development of these textures themselves depends on cooling rates, among other factors. There are two ways in which the development of such textures can affect the results of diffusion modeling as described here. Firstly, the partitioning behavior of Mg may be affected. However, the fact that feldspar thermometry and phase

relations have been used reliably in metamorphic rocks (e.g., see Nekvasil, 1994 for a review) suggests that the influence of these phase transitions on the overall Gibbs Free Energy is not so large as to invalidate the partition coefficients. Secondly, and more importantly, the mosaic textures may offer fast diffusion pathways that could short circuit volume diffusion. Measurement of multiple concentration profiles, or if possible, element mapping, and the identification of small scale, high frequency oscillations in compositional profiles may help to identify such fast diffusion pathways where they exist. If evidence of such high diffusion pathways are found, such crystals should be avoided for diffusion modeling. 5. A WORKED EXAMPLE FROM THE LOWER OCEANIC CRUST The new Mg-in-plagioclase geospeedometer was applied to Mg-profiles measured in plagioclase crystals formed in the lower oceanic crust, where plagioclase typically occurs

Concentration MgO [wt%]

(a) 0.14

(b)

P_20

P_20

0.12 0.10 0.08

dT/dt = 3.6 10 °Cyr

0.06

Tcrim = 960°C

-1

0.04 0.02 rim

core 500

-1

dT/dt = 3.6 10 °Cyr -1

Xan

0.7

0.7

0.6

0.6

0.5

0.5

0.4

rim

1000

0.14

rim

1500

(c) Concentration MgO[ wt%]

703

core 500

Xan

0.4

rim

1000

1500

(d)

P_500

P_500

0.7

0.7

0.6

0.12

0.6

Xan

0.5

0.5

0.4

0.10

-1

0.08

Tcrim = 860°C

0.06

dT/dt = 1.1 10 °Cyr -3

Xan

0.4

-1

dT/dt = 1.1 10 °Cyr -3

-1

0.04 0.02 rim

core 2004

rim 00

Distance [µm]

600

rim

core 200

rim 400

600

Distance [µm]

Fig. 8. Measured (black circles) and fitted (black line) Mg-concentration profiles in two different plagioclase crystals P_20 and P_500 (samples 2212–1338 (17 mbsd), and 2218–1132 (520 mbsd), respectively, from the Hess Deep sample suite, Karson et al., 2002). The insets in all panels show the respective XAn profiles, measured along the same traverse as the Mg-profiles. Panels (a) and (b) show the same measured profile in P_20. The fit shown in (a) results from stopping the model at Tcrim = 960 °C, to obtain the best fit for the concentration of Mg at the rim (5 lm into the crystal). The cooling rate obtained from this fit is 0.36 °Cyr1. The fit shown in (b) results, if the same cooling rate (0.36 °Cyr1) is applied, but the model was continued down to lower temperatures (here shown for 527 °C, see text for discussion). Panels (c) and (d) show the same measured profile in P_500. The fit shown in (c) results from stopping the model at Tcrim = 860 °C, to obtain the best fit for the concentration of Mg at the rim. The cooling rate obtained from this fit is 0.0011 °Cyr1. The fit shown in (d) results, if the same cooling rate (0.0011 °Cyr1) is applied, but the model was continued up to lower temperatures (600 °C) at the same cooling rate.

704

K. Faak et al. / Geochimica et Cosmochimica Acta 140 (2014) 691–707

in contact with clinopyroxene as a test of the methodology described above. These samples come from Hess Deep (Karson et al., 2002) and were collected from 20 m and 500 m below the sheeted dike (mbsd) complex meaning that they are expected to have cooled at different rates (e.g., Coogan et al., 2007). Plagioclase crystals with an aspect ratio P1:2 were chosen and the plagioclase composition was measured along the short dimension of the crystal (from one rim to the opposite rim through the core). The analytical conditions used to measure the element concentrations in the samples are outlined in Electronic Annex B. The distance between analyzed spots along the profile was 10 lm. The first and last measurements at the rims of the plagioclase were approximately 4 to 5 lm away from the interface, to avoid contamination of Mg from the adjacent clinopyroxene (although secondary fluorescence is minimal as discussed in Electronic Annex B, contamination by mixed analyses has to be avoided – such mixed analyses can arise, for example, if the contact between the plagioclase and the clinopyroxene is inclined to the plane of the thin section). With these measurement conditions it was possible to determine compositional zoning profiles of high accuracy even for low Mg-concentrations in plagioclase (Fig. 8). The plagioclase in the shallow sample (P_20) has high concentration of Mg in the core (0.13 wt% MgO) and lower concentrations of Mg at the rims (0.04 wt% MgO) (Fig. 8a and b). The XAn profile is also zoned, with higher XAn at the core than at the rims (Fig. 8a and b). In the sample from greater depth below the sheeted dike complex the plagioclase (P_500) has an almost homogeneous XAn content and a flat Mg-profile of around 0.04 wt% MgO (Fig. 8c and d). The different shapes of profiles found in plagioclase crystals from two different depths is systematic in that all crystals in which profiles were measured from a given sample have similar shapes, independent of grain size. 5.1. Calculating cooling rates for rocks from the lower oceanic crust The initial and boundary conditions that are applied need to be chosen to match the inferred geological situation. In the examples presented here, these are the conditions of rocks formed in the lower oceanic crust, i.e., solidification of a basaltic melt around 1200 °C, which produces plagioclase and clinopyroxene in contact with each other, that can exchange Mg via diffusion during continuous cooling. The initial profile was obtained from calculating the equilibrium distribution using Eq. (4) and the observed XAn-zoning profile at an initial temperature P Tstart (Eq. (5)). Using the measured MgOconcentrations in the natural plagioclase crystals together with Eq. (3) and Eq. (5) we obtain Tstart = 1225 °C (plagioclase P_20) and Tstart = 1000 °C (plagioclase P_500), respectively. The final profiles that are obtained following the procedure outlined below are found to be identical if modeling is carried out with initial conditions calculated at temperatures that are higher than these values. Initial conditions set at lower temperatures of 1200 °C (plagioclase P_20) and 980 °C (plagioclase P_500) also yield

the same results. This shows that our method of determining Tstart (a temperature that is high enough so that all memory of the initial condition is erased by diffusion during cooling) is reliable. The boundary conditions are also taken from Eq. (4) using the measured plagioclase Mg-and XAn contents at the rim, and the measured clinopyroxene Mg-content. Detailed measurements of the rim compositions of the clinopyroxene in both samples showed no zoning in Mg and thus, we assume that the clinopyroxene acted as an infinite reservoir. As discussed in Section 3.2, Eq. (4) requires constraints on the silica activity of the system. Although the assemblage olivine–orthopyroxene (Ol + Opx) is not present in the specific thin sections modeled here, it is found in several other samples from the same location. Therefore, aSiO2 in these rocks may be constrained by the reaction Mg2SiO4 + SiO2 = 2MgSiO3 and is given by the relationship proposed by Carmichael et al. (1970): log aSiO2 ¼

DGr 2:303 R T

ð8Þ

where DGr is the Gibbs Free energy of the reaction. The silica activity was calculated for different temperatures using the data set of Ghiorso and Sack (1995) to determine DGr, and the resulting dataset fitted by a 2nd-order polynomial to obtain: aSiO2 ¼ 4:86904  107 T 2 þ 1:51570  103 T  0:618707

ð9Þ

The diffusion coefficient DPl Mg of Faak et al. (2013) was used, since their experiments specifically determined the diffusive exchange of Mg between plagioclase and clinopyroxene in the compositional range of the lower oceanic crust. As a first approximation, a linear cooling history was assumed. As diffusion is a thermally activated process, there should be a closure temperature Tc below which the shape of the developed Mg-concentration profile does not change significantly due to diffusion. It has been shown in Section 2, that Tc is a function of distance from the rim (i.e., Tc(x)) and that Tcrim 6 Tccore. Considering the analytical resolution of profile measurement, we define the ‘rim’ as a region 5 lm away from the interface of the plagioclase. Subsequently, Tcrim can be calculated from the measured Mgand XAn-concentrations at the rim using the thermometer from Faak et al. (2013) and we obtain Tcrim of 960 °C for plagioclase P_20 and 860 °C for plagioclase P_500. This approach does not make any assumptions about the cooling rate and accounts for the effect of zoning in XAn on the diffusion of Mg. The diffusion modeling was stopped at Tcrim. The calculated Tcbulk values are, as expected, a little higher than Tcrim (Tcbulk = 1060 °C for P_20 and Tcbulk = 885 °C for P_500). Within the temperature range spanned by Tstart and Tcrim the cooling rate was adjusted iteratively until the visual misfit between the modeled profile and the measured data was minimized (see also Section 5.2). This approach yields best fits for a cooling rate of 0.36 °Cyr1 for the measured Mg-profile in plagioclase P_20 (Fig. 8a) and 0.0011 °Cyr1 for plagioclase P_500 (Fig. 8c). For the more convex profile shape of plagioclase P_20, continuing

K. Faak et al. / Geochimica et Cosmochimica Acta 140 (2014) 691–707

cooling in the modeling calculations below Tcrim does not produce any change in calculated profile shapes (Fig. 8b). On the other hand, for the flatter profile of plagioclase P_500, continuing cooling with the assumption of instantaneous equilibration at the surface produces convex profile shapes that do not match the observed profile shape (Fig. 8d). This implies that either the assumption of instantaneous equilibrium with an infinite reservoir at the interface does not hold for this sample (e.g., kinetic barriers, see Section 4.3.4 above) or that there was a rapid increase in cooling rate below Tcrim. (i.e., the assumption of continuous linear cooling is not valid). This aspect will be explored in more detail in a subsequent paper where this sample may be considered along with many others in the context of spatial distribution of cooling history. We note here that it is impossible to obtain such flat profiles by diffusion without violating one of the two assumptions mentioned above. The Tcbulk values obtained above can be evaluated now using alternate approaches for calculating this quantity. Tcmean can be calculated based on the Dodson (1973) formulation using the activation energy E of the diffusion coefficient, the grain size and geometry of the system, and the cooling rate estimated above. For the two examples here, this yields Tcmean = 1034 °C for P_20 and Tcmean = 775 °C for P_500. The similar, but somewhat lower, values compared to those obtained above are expected because (a) the Dodson model cannot account for the effect of XAn-concentration gradient on diffusion

Concentration MgO[ wt%]

(a) 0.14

5.2. Sensitivity of the obtained profile shape to the applied cooling rate The sensitivity of the modeled profile shape to changes in the cooling rate is illustrated in Fig. 9. The strongly curved Mg-profiles of plagioclase P_20, is very sensitive to changes in the cooling rate and changing dT/dt from the best-fit value of 0.36 °Cyr1 (Fig. 9a), to 0.5 °Cyr1 (Fig. 9b) and 0.2 °Cyr1 (Fig. 9c) leads to significant visual misfits of the profile. Fitting of the flat Mg-profiles of plagioclase P_500 is less sensitive to changes in the cooling rate (Fig. 9d-f), but when the cooling rate is changed from 0.0011 °Cyr1 (which gives the best fit; Fig. 9d) to 0.004 °Cyr1 (Fig. 9e) or 0.0001 °Cyr1 (Fig. 9f), the calculated and measured profiles are not in good agreement. 6. PLAGIOCLASE AND CLINOPYROXENE PAIRS FROM OTHER SETTINGS The geospeedometer presented above has a broad range of applications and may be used to constrain the cooling history of plagioclase-clinopyroxene pairs from terrestrial metamorphic and igneous rocks as well as meteorites. However, the method is limited by the minimum Mg-content

(c)

P_20

P_20

0.12 0.10

dT/dt = 3.6 10 °Cyr

0.08

-1

0.06

dT/dt = 5.0 10 °Cyr

-1

-1

Xan

0.7

0.04

0.7

0.6

core 500

-1

Xan

rim

1000

rim

core

1500

500

0.6 0.5

0.4

rim

1000

0.4

rim

10001

500

0.7

0.6

Xan

0.5

0.5

0.4

0.10

500

P_500

0.7

0.6

0.12

core

(f)

P_500

0.7

rim

1500

(e)

P_500

-1

Xan

0.7

0.5

0.4

(d) 0.14

dT/dt = 2.0 10 °Cyr

-1

0.6

0.5

0.02 rim

Concentration MgO[ wt%]

of Mg and (b) the Dodson model averages a Tc-profile that extends infinitesimally close to the interface (this includes a region of steep drop in concentration, and hence, Tc, at the rim).

(b)

P_20

705

0.6

Xan

Xan

0.5

0.4

0.4

0.08

dT/dt = 1.1 10 °Cyr -3

0.06

dT/dt = 4.0 10 °Cyr

-1

-3

dT/dt = 1.0 10 °Cyr

-1

-4

-1

0.04 0.02 rim

core 2004

rim 00

Distance [µm]

600

rim

core 200

rim 400

Distance [µm]

600

rim

core 200

rim 4006

00

Distance [µm]

Fig. 9. Measured (black circles) Mg-concentration profiles in two different plagioclase crystals compared to modeled profiles (black line), which were calculated at different cooling rates. The insets in all panels show the respective XAn profiles, measured along the same traverse as the Mg-profiles. Panels (a)–(c) show the same measured profile in plagioclase P_20, which was fitted with different cooling rates. The best fit is obtained for a cooling rate of dT/dt = 0.36 °Cyr1 (a). If the cooling rate is increased to 0.5 °Cyr1 (b) or decreased to 0.2 °Cyr1 (c), the quality of the fit is visually worse. Panels (d–f) show a comparison of fits at different cooling rates for the measured profile in plagioclase P_500. The best fit is obtained for a cooling rate of dT/dt = 0.0011 °Cyr1 (d). Visual misfits are observed, if the cooling rate is changed to 0.004 °Cyr1 (e) and 0.0001 °Cyr1 (f).

706

K. Faak et al. / Geochimica et Cosmochimica Acta 140 (2014) 691–707

that can be measured in plagioclase with different methods (e.g., electron or ion microprobe) and potentially, by the structural transitions discussed above. The lowest reliable Mg-contents measured with an electron microprobe in this study were around 0.005 wt% MgO, which corresponds to a temperature of 660 °C using the thermometer of Faak et al. (2013) for the conditions of the oceanic crust (i.e., XAn = 0.6, Cpx with 14 wt% MgO and the silica activity buffered by Ol + Opx). The temperatures corresponding to this Mg-content for many different possible scenarios are similar. Therefore, Mg-profiles measured by electron microprobe from rocks with peak temperatures around 700 °C (i.e., most mafic metamorphic rocks containing plagioclase–clinopyroxene, and many meteorites) may be amenable to constrain cooling histories with the Mg-in-plagioclase geospeedometer depending on their cooling rate. For rocks with lower peak temperatures, it will be necessary to use analytical devices that allow measurement of lower Mg-concentrations and consider whether there have been any structural transitions in the feldspar during cooling. 7. CONCLUSIONS A new Mg-in-plagioclase geospeedometer has been presented, based on the diffusive exchange of Mg between plagioclase and clinopyroxene during cooling. The underlying mathematical model and its numerical solution are described and ways to determine initial and boundary conditions for the diffusion model are provided. The new ‘geospeedometer’ and its mathematical model are discussed with regard to the robustness of the obtained cooling rates and the sensitivity of the developed Mg-profile shape in plagioclase on: (i) the distribution of the anorthite-component in plagioclase (since the diffusion of Mg in plagioclase is coupled to XAn, and the partitioning of Mg between clinopyroxene and plagioclase also depends on XAn), (ii) the geometry of the plagioclase crystal, and (iii) the cooling history. For any smooth, continuous cooling history, Mg-profiles produced by diffusive exchange with clinopyroxene tend to be convex upward. However, the detailed shape depends on the shape of the XAn profile. All other things being equal, more homogeneous final Mg-profiles shapes are observed for cooling paths with increasing cooling rates over time. By examining the detailed shape of the zoning profiles within plagioclase this tool has the potential to reveal changes in cooling rates as well as average cooling rates (e.g., Fig. 4 and related discussion). The effect of modeling diffusion of Mg in a 1D vs. a 2D diffusion model was tested for different aspect ratios of plagioclase crystals and it was found that, even though 1D models generally underestimate the cooling rate (i.e., require slower cooling rates than 2D models to obtain the same Mg-concentrations), both models yield the same cooling rate for modeling the short dimension in plagioclase crystals with an aspect ratio P1:2. It has been shown that the measurement and modeling of full zoning profiles is essential for extraction of cooling rates. Additionally, measurement of the zoning of Mg in

the rims of clinopyroxene in contact with the modeled profile in plagioclase is required to test whether clinopyroxene acted as an infinite reservoir (from lack of Mg zoning at rim). Measuring and modeling multiple plagioclase crystals of different sizes within a given sample (or samples that can be shown independently to have had the same thermal history) will allow more confidence that the appropriate initial and boundary conditions are being applied. We conclude that, the Mg-in-plagioclase geospeedometer is a potentially powerful tool for the determining the cooling history of a wide range of rocks with coexisting plagioclase and clinopyroxene. ACKNOWLEDGEMENTS We thank Yan Liang and Philipp Ruprecht for constructive reviews. This study was supported by the Grant CH 166/12-1 from the Deutsche Forschungsgemeinschaft (DFG).

APPENDIX A. SUPPLEMENTARY DATA Supplementary data associated with this article can be found, in the online version, at http://dx.doi.org/10.1016/ j.gca.2014.06.005. REFERENCES Bindeman I. N., Davies A. M. and Drake M. J. (1998) Ion microprobe study of plagioclase-basalt partition experiments at natural concentration levels of trace elements. Geochim. Cosmochim. Acta 62, 1175–1193. Blackburn T., Shimizu N., Bowring S. A., Schoene B. and Mahan K. H. (2012) Zirconium in rutile speedometry: new constraints on lower crustal cooling rates and residence temperatures. Earth Planet. Sci. Lett. 317, 231–240. Blundy J. D. and Wood B. J. (1994) Prediction of crystal-melt partition coefficients from elastic moduli. Nature 372, 452–454. Brown W. L. and Parsons I. (1994) Feldspars in igneous rocks. In Feldspars and Their Reactions (ed. I. Parsons) Nato Adv. Sci. Inst. Ser. C, Math. Phys. Sci. 421, pp. 449–499. Carmichael I. S. E., Nicholls J. and Smith A. L. (1970) Silica activity in igneous rocks. Am. Mineral. 55, 246–263. Carpenter M. A. (1994) Subsolidus phase-relations of the plagioclase feldspar solid-solution. In Feldspars and Their Reactions (ed. I. Parsons) Nato Adv. Sci. Inst. Ser. C, Math. Phys. Sci. 421, pp. 221–269. Chakraborty S. and Ganguly J. (1991) Compositional zoning and cation diffusion in aluminosilicate garnets. In Diffision, Atomic Ordering and Mass Transport, Advances in Physical Geochemistry (ed. J. Ganguly). Springer, Berlin, pp. 120–175 (8). Cherniak D. J. and Watson E. B. (1994) A study of strontium diffusion in plagioclase using Rutherford backscattering spectroscopy. Geochim. Cosmochim. Acta 58, 5179–5190. Coogan L. A., Jenkin G. R. T. and Wilson R. N. (2007) Contrasting cooling rates in the lower oceanic crust at fastand slow-spreading ridges revealed by geospeedometry. J. Petrol. 48, 2211–2231. Costa F. and Chakraborty S. (2004) Decadal time gaps between mafic intrusion and silicic eruption obtained from chemical zoning patterns in olivine. Earth Planet. Sci. Lett. 227, 517–530. Costa F., Chakraborty S. and Dohmen R. (2003) Diffusion coupling between major and trace elements and a model for

K. Faak et al. / Geochimica et Cosmochimica Acta 140 (2014) 691–707 the calculation of magma chamber residence times using plagioclase. Geochim. Cosmochim. Acta 67, 2189–2200. Costa F., Dohmen R. and Chakraborty S. (2008) Time scales of magmatic processes from modeling the zoning patterns of crystals. In Minerals, Inclusions and Volcanic Processes (eds. K. D. Putirka and F. J. Tepley). Rev. Min. Geochem. 69, MSA, pp. 545–594. Dodson M. H. (1973) Closure temperature in cooling geochronological and petrological systems. Contrib. Mineral. Petrol. 40, 259–274. Dodson M. H. (1976) Kinetic processes and thermal history of slowly cooled solids. Nature 259, 551–553. Dodson M. H. (1986) Closure profiles in cooling systems. Mater. Sci. Forum 7, 145–154. Dohmen R. and Chakraborty S. (2003) Mechanism and kinetics of element and isotopic exchange mediated by a fluid phase. Am. Mineral. 88, 1251–1270. Eiler J. M., Baumgartner L. P. and Valley J. W. (1992) Intercrystalline stable isotope diffusion: a fast grain boundary model. Contrib. Mineral. Petrol. 112, 543–557. Faak K., Chakraborty S. and Coogan L. A. (2013) Mg in plagioclase: experimental calibration of a new geothermometer and diffusion coefficients. Geochim. Cosmochim. Acta 123, 195– 217. Ganguly J. (1982) Mg–Fe order-disorder in ferromagnesian silicates II. Thermodynamics, kinetics and geological applications. In Advances in Physical Geochemistry 2 (ed. S. K. Saxena). Springer, Berlin Heidelberg New York. pp. 58–99. Ganguly J. and Tirone M. (1999) Diffusion closure temperature and age of a mineral with arbitrary extent of diffusion: theoretical formulation and application. Earth Planet. Sci. Lett. 170, 131–140. Ghiorso M. S. and Sack R. O. (1995) Chemical mass transfer in magmatic processes IV. A revised and internally consistent thermodynamic model for the interpolations of liquid-solid equilibria in magmatic systems at elevated temperatures and pressures. Contrib. Mineral. Petrol. 119, 197–212. Giletti B. J. and Casserly J. E. D. (1994) Strontium diffusion kinetics in plagioclase feldspars. Geochim. Cosmochim. Acta 58, 3785–3793. Ginibre C., Worner G. and Kronz A. (2007) Crystal zoning as an archive for magma evolution. Elements 3, 261–266. Grove T. L., Baker M. B. and Kinzler R. J. (1984) Coupled CaAl– NaSi diffusion in plagioclase feldspar: experiments and applications to cooling rate speedometry. Geochim. Cosmochim. Acta 48, 2113–2121. Jenkin G. R. T., Farrow C. M., Fallick A. E. and Higgins D. (1994) Oxygen isotope exchange and closure temperature in cooling rocks. J. Metamorphic Geol. 12, 221–235. Karson J. A., Klein E. M., Hurst S. D., Lee C., Rivizzigno P., Curewitz D., Morris A. R. and Party H. D. S. (2002) Structure of uppermost fast-spread oceanic crust exposed at the Hess Deep Rift: implications for subaxial processes at the East Pacific Rise. Geochem. Geophys. Geosyst. 3. http://dx.doi.org/ 10.1029/2001GC000155. Lasaga A. C. (1983) Geospeedometry: an extension of geothermometry. In Kinetics and equilibrium in mineral reactions (ed. S. K. Saxena). Springer–Verlag, New York. pp. 82–114. Lasaga A. C. (1998) Kinetic Theory in the Earth Sciences. Princeton University Press, New Jersey. LaTourrette T. and Wasserburg G. J. (1998) Mg diffusion in anorthite: implications for the formation of early solar system planetesimals. Earth Planet. Sci. Lett. 158, 91–108. Lindstrom R., Viitanen M., Juhanoja J. and Holtta P. (1991) Geospeedometry of metamorphic rocks – examples in the

707

rantasalmi-sulkava and kiuruvesi areas, eastern finland - biotite garnet diffusion couples. J. Metamorph. Geol. 9, 181–190. Liu M. and Yund R. A. (1992) NaSi–CaAl interdiffusion in plagioclase. Am. Mineral. 77, 275–283. Miller S. A., Asimow P. D. and Burnett D. S. (2006) Determination of melt influence on divalent element partitioning between anorthite and CMAS melts. Geochim. Cosmochim. Acta 70, 4258–4274. Morse S. A. (1984) Cation diffusion in plagioclase feldspar. Science 225, 504–505. Mueller T., Dohmen R., Becker H. W., ter Heege J. H. and Chakraborty S. (2013) Fe–Mg interdiffusion rates in clinopyroxene: experimental data and implications for Fe–Mg exchange geothermometers. Contrib. Mineral. Petrol. DOI:10.1007/s00410-013-0941-y. Nekvasil H. (1994) Ternary feldspar/melt equilibria – a review. In Feldspars and Their Reactions (ed. I. Parsons). Nato Adv. Sci. Inst. Ser. C, Math. Phys. Sci. 421, pp. 195–219. Onorato P. I. K., Hopper R. W., Yinnon H., Uhlmann D. R., Taylor L. A., Garrison J. R. and Hunter R. (1981) Solute partitioning under continuous cooling conditions as a cooling rate indicator. J. Geophys. Res. 86, 9511–9518. Ozawa K. (1984) Olivine-spinel geospeedometry: analysis of diffusion-controlled Mg-Fe2+ exchange. Geochim. Cosmochim. Acta 48, 2597–2611. Redfern S. A. T., Henderson C. M. B., Wood B. J., Harrison R. J. and Knight K. S. (1996) Determination of olivine cooling rates from metal cation ordering. Nature 381, 407–409. Ruprecht P. and Cooper K. (2012) Integrating the Uranium-series and elemental diffusion geochronometers in mixed magmas from Volca´n Quizapu, Central Chile. J. Petrol. 53, 841–871. Seifert F. (1977) Reconstruction of rock cooling paths from kinetic data on Fe2+–Mg exchange reaction in anthophyllite. Philos. Trans. R. Soc. Lond. A 286, 303–311. Singer B. S., Pearce T. H., Kolisnik A. M. and Myers J. D. (1993) Plagioclase zoning in modpleistocene lavas from the Seguam Volcanic Center, Central Aleutian Arc, Alaska. Am. Mineral. 78, 143–157. Spear F. S. (1991) On the interpretation of peak metamorphic temperatures in light of garnet diffusion during cooling. J. Metamorph. Geol. 9, 379–388. Vanko D. A. and Laverne C. (1998) Hydrothermal anorthitization of plagioclase within the magmatic/hydrothermal transition at mid-ocean ridges: examples from deep sheeted dikes (Hole 504B, Costa Rica Rift) and a sheeted dike root zone (Oman ophiolite). Earth Planet. Sci. Lett. 162, 27–43. Van Orman J. A., Cherniak D. J. and Kita N. T. (2014) Magnesium diffusion in plagioclase: dependence on composition, and implications for thermal resetting of the 26Al–26Mg early solar system chronometer. Earth Planet. Sci. Lett. 385, 79–88. Wallace G. S. and Bergantz G. W. (2003) Constraints on mingling of crystal populations from off-center zoning profiles: a statistical approach. Am. Mineral. 89, 64–73. Wilding M. C., Webb S. L. and Dingwell D. B. (1995) Evaluation of a relaxation geospeedometer for volcanic glasses. Chem. Geol. 125, 137–148. Zhang X., Ganguly J. and Ito M. (2010) Ca–Mg diffusion in diopside: tracer and chemical inter-diffusion coefficients. Contrib. Mineral. Petrol. 159, 175–186. Associate editor: Wolf Uwe Reimold