Reaction progress of clay minerals and carbonaceous matter in a contact metamorphic aureole (Torres del Paine intrusion, Chile)

This study reports on reaction processes in a transition zone from contact to regional metamorphism by using Raman spectroscopy on carbonaceous matter (RSCM), illite “crystallinity” (Kübler index, KI), chlorite geothermometry, and thermal modeling. The thermal effect due to the emplacement of the Torres del Paine intrusion (TPI, assembly time of ca. 150 kyr) had different consequences for inorganic and organic compounds of the host rock. The thermal alteration of the pre-intrusive regional metamorphosed host rock is documented by elevated RSCM temperatures, high-temperature chlorite generations, and the appearance of epidote and retrograde Fe-rich chlorite. Microprobe analysis on chlorite indicates incomplete re-equilibration as evidenced by various chlorite populations of individual contact metamorphic samples. This study indicates that the maturity of organic matter is the most reliable and unequivocal indicator on timescales of several thousand years to determine the lateral extension of the TPI contact aureole. Raman geothermometry reveals that the lateral extension of the contact-influenced zone expands up to a distance of 1.5 km and, thus, expands to ca. 1.1 km further out than the macroscopically mappable hornfels contact aureole. The best match between measured (Raman geothermometry) and calculated (thermal modeling) ΔT_max values (ΔT=54 ^∘C) is achieved with a total intrusion assembly time of 150 kyr, a magmatic temperature of 800 ^∘C, a two-batch model (batch repose time of 10 kyr) with five pulses per batch, short heating durations (3 kyr), and long pulse repose times (15 kyr).

1 Introduction

Magmatic intrusions emplaced in sedimentary basins produce local heat anomalies which affect both organic and inorganic compounds of the host rock. The heating results in an irreversible structural transformation of organic matter with release of volatiles including CO₂, H₂O, and CH₄, ultimately producing graphite at T>700 ^∘C (e.g., Franklin, 1951; Grew, 1974; Clayton and Bostick, 1986; Buseck and Beyssac, 2014). Maturation of carbonaceous matter in contact with igneous intrusions is well documented in the literature (e.g., Sweeney and Burnham, 1990; Olsson, 1999; Frings et al., 2004; Aarnes et al., 2010). Contact metamorphic processes of inorganic material typically include dehydration, decarbonation, and host-rock melting reactions (e.g., Jamtveit et al., 1992). During these processes, clay minerals become less expandable due to a decrease in the occurrence of smectite interlayers (Pytte and Reynolds, 1989; Esposito and Whitney, 1995; Nadeau and Reynolds, 1981; Abad et al., 2014). Pytte and Reynolds (1989) suggested that the reaction is controlled by kinetic rather than equilibrium factors, meaning that changing physical boundary conditions (i.e., temperature, time, and K⁺ activity) has the potential to disproportionately affect clay mineral compositions. They documented a zone of decreased expandability extending into the country rock up to a distance of approximately twice the intrusion thickness.

Temperature, time, and fluids (i.e., their presence and composition) associated with the thermal anomaly are the most important factors affecting the maturation of organic matter, as well as the crystallinity and chemical composition of clay minerals in a contact metamorphic aureole. The maturity of organic matter has been frequently shown to be foremost a function of maximum temperature (T_max) and is used as a paleo-T_max indicator. In contrast, changes in illite and chlorite “crystallinities” are also dependent on other parameters, such as a definite temperature duration, the activity of cations, and fluid availability (e.g., Merriman, 2005). The nature and effect of these additional parameters and their influence on clay compositions in a contact metamorphic setting has not been fully discerned yet.

Figure 1 (a) Geological map of the Torres del Paine intrusive complex in the Última Esperanza district. Map modified from © Leuthold et al. (2012) and © Süssenberger et al. (2018c). Red stars indicate sampling positions. Black rectangles mark the insets of panels (b) and (c). (b) Sampling positions at the Paso John Gardner profile. (c) Sampling positions at the Monte Almirante profile.

The studied contact aureole of the Torres del Paine intrusion (TPI) is located in the Magallanes–Austral Basin at 51^∘ S (Chile, Fig. 1). The TPI is an isolated Miocene intrusion, located outside the main Cretaceous–Cenozoic Patagonian Batholith, and is therefore insulated from other magmatic sources that could conceal its thermal effects. The TPI was emplaced into the Cretaceous Punta Barrosa and Cerro Toro formations between 12.59 and 12.43 Ma and postdates the main phase of deformation during fold-and-thrust belt development (Halpern, 1973; Michel et al., 2008; Leuthold et al., 2012). The intrusion was rapidly exhumed in an erosion-related exhumation–cooling process. To date, it offers a unique environment to study contact metamorphism due to its excellent outcrop conditions, very good geochronological record of the host rocks, and very detailed studies about the age and incremental growth of the intrusion (Michel et al., 2008; Leuthold et al., 2012; Bodner, 2013). Although the TPI provides excellent outcrop conditions, the transition from the low-temperature outer contact aureole to the pre-intrusive low-temperature regional metamorphism is rarely recognized on an outcrop scale due to the comparatively high degree of background maturation or the regional metamorphic imprint, as well as the relatively dark and fine-grained lithology. The low-temperature metapelitic aureole is typically more extensive than the outcropping hornfels aureole, which is commonly recognized by the first appearance of biotite in the albite–epidote–hornfels facies. A cryptic aureole can extend up to 3 km beyond the outer limit of hornfels observed in outcrops, with an outer margin in the middle anchizone (Abad et al., 2014).

In the present study, we take advantage of a very well dated and studied pluton, located outside of the main batholith and emplaced in a fine-grained sedimentary sequence, with significant organic-matter content. We compare the effect of rapid heating on the evolution of clay mineral reaction progress and the maturity of organic matter with respect to its suitability to reconstruct contact metamorphic overprinting. The analytical data were used, together with the results from thermal modeling, to establish the lateral extension of the Torres del Paine contact metamorphic aureole on the outcrop scale. Samples were collected mostly along three profiles across the margins of the intrusion, from the cordierite-in isograd in the contact aureole to units of the Punta Barrosa and Cerro Toro formations unaffected by contact metamorphism. Although the Punta Barrosa and Cerro Toro formations do not represent the main productive units of the Magallanes–Austral Basin, they are a reservoir and secondary source with considerable organic-matter content. Investigation of the samples involved X-ray diffraction (XRD), scanning electron microscope (SEM) imaging, and electron microprobe (EMP) analyses. Temperatures are derived from Raman spectroscopy on carbonaceous matter (RSCM) and from chlorite geothermometry. The measured T_max values obtained from RSCM are validated against calculated T_max values based on numerical heat-transfer models.

2 Geological setting

The Miocene Torres del Paine intrusion (TPI), situated in the Última Esperanza District (southern Patagonia, 51^∘ S), is genetically linked to the subduction of the Chile rise and the Patagonian Batholith (Michael, 1984, 1991; Baumgartner et al., 2007; Leuthold et al., 2012; Bodner, 2013). The TPI is a shallow laccolithic crustal sill complex with bimodal composition that was assembled between 12.58±0.01 and 12.49±0.01 Ma in at least three batches (Michel et al., 2008; Leuthold et al., 2012) and emplaced into the regional metamorphosed Upper Cretaceous Punta Barrosa and Cerro Toro formations. Here, we use the term batch to describe an intrusive unit with sharp contact with the adjacent units, representing a distinguishable entity in the field. A batch may be constructed of several pulses, i.e., individual magma flows. The general setup of the laccolith comprises a basal mafic sill complex and the overlying Paine granite generations I, II, and III (Michael, 1984; Michel et al., 2008; Leuthold et al., 2012; Michael, 1991). Michel et al. (2008) established a time frame of 90±40 kyr for the emplacement of the Paine granite within multiple pulses. After the first assemblage of granites, three mafic sheets with an overall thickness of ca. 400 m underplated the laccolith between 12.50 and 12.43 Ma (Leuthold et al., 2012). Today, the oldest sheet is located at the top, while the youngest is intruded at the base of the granitic laccolith. The subhorizontal laccolith complex is at its western extremity connected to a vertical feeder system (Baumgartner et al., 2007). Pressure estimates based on contact metamorphic assemblages and granite-solidus thermobarometry suggest that the bimodal felsic–mafic shallow crustal laccolith was emplaced at pressures of ca. 750±250 bars, which corresponds to an emplacement depth of ca. 2–3 km and matches estimates based on hornblende thermobarometry and the stability of prehnite in metacarbonates (Putlitz et al., 2001; Baumgartner et al., 2014; Leuthold et al., 2014). The temperature of the granitic intrusion is constrained to ca. 900–1000 ^∘C by contact metamorphic studies (Bodner, 2013; Baumgartner et al., 2014). Peak metamorphic temperatures between 400–450 ^∘C were derived from the occurrence of prograde prehnite a few meters away from the contact (Putlitz et al., 2001). The ca. 80 km² large intrusion created a direct contact metamorphic aureole, which reaches a thickness of between 150 and 400 m (Putlitz et al., 2001; Bodner, 2013). The thickness of the contact aureole is considerably smaller (ca. 200 m) above the intrusion than below it (Bodner, 2013). Bodner (2013) could show that the intrusion established a relatively small hydrothermal convection cell despite the occurrence of several magmatic pulses.

The maximum ambient temperature of the metasediments prior to the emplacement of the TPI was between 240 and 260 ^∘C and related to regional low-grade metamorphism (Süssenberger et al., 2016, 2018b). The regional anchizonal metamorphic degree is recorded by the maturity of organic matter, chlorite temperatures, and illite “crystallinities” (Süssenberger et al., 2018a, b, c). Generally, a decrease in the metamorphic grade is observed from west to east. The highest metamorphic grade (lower epidote facies) is observed in samples of the Tobífera and Zapata formations and is related to burial and intense deformation during the early Paleogene (e.g. Süssenberger et al., 2018b). More weakly metamorphosed rocks are observed east of the TPI in the Punta Barrosa, Cerro Toro, and Dorotea formations (Süssenberger et al., 2018b). Although the Punta Barrosa and Cerro Toro formations represent the syntectonic basal foreland sequence (Fosdick et al., 2013; Ghiglione et al., 2014), the main thrusting- and shortening-related tectonic uplift took place during the Paleogene (Fosdick et al., 2011; Süssenberger et al., 2018b), preceding the TPI emplacement. Late Miocene to Pliocene denudation is indicated by thermochronometry, which exposed the roof of the pluton and has been assigned to erosion and isostatic rebound (Fosdick et al., 2013).

Table 1 Compilation of results from Raman spectroscopy on carbonaceous matter, chlorite geothermometry, and illite “crystallinity” (Kübler index) for the <2 µm fraction. n.d.: not determined; R: retrograde; RM: regional metamorphic; CM: contact metamorphic.

^a Samples from Tobler (2012). ^b Samples from Bodner (2013). ^c Data from Süssenberger et al. (2018c).

Figure 2 (a) View of the Torres del Paine intrusive complex from the eastern park entrance in the Última Esperanza district. The black rectangle indicates the location of panel (d). (b) Well-exposed Torres del Paine sill which is overlain by the Upper Cretaceous Cerro Toro Formation. The contact between the granite and adjacent metapelite is indicated by the dashed red line. (c) Contact metamorphic hornfels ca. 50 m from the intrusion. Chisel length is ca. 30 cm. (d) Contact between granite and metapelite at the Monte Almirante profile.

3 Material and methods

3.1 Sample material and field observations

Seventeen samples distributed over the contact aureole and the adjacent unaltered host rock were selected for this study (Figs. 1 and 2, Table 1, Supplement 1). The samples were collected along two profiles near Monte Almirante and along one profile at the western limit of the TPI (Figs. 1 and 2). Profiles and individual samples were collected starting from the cordierite-in isograd of the contact aureole moving towards solely regionally affected units of the Punta Barrosa (western Paine mafic complex) and Cerro Toro (Monte Almirante) formations (Fig. 1). Monte Almirante is situated at the eastern end of the laccolith and represents the front end of the magma flow (Baumgartner et al., 2007). Investigated samples from above the intrusion are described in more detail by Bodner (2013). The transition from low-temperature regional to contact metamorphism appears to be obscured by the very fine-grained lithology and is intrinsically difficult to distinguish in the field. The regional metamorphic shale paragenesis consists of albite, chlorite, and illite (dominant clay phase; Süssenberger et al., 2018b). A detailed characterization of the anchizonal to epizonal regional metamorphic host rock, clay mineral compositions, and constraints of fold-and-thrust belt formation are given in Süssenberger et al. (2018a, b, c). The samples were subjected to X-ray diffraction (XRD), scanning electron microscope (SEM) imaging, and electron microprobe (EMP) analysis.

3.2 Illite “crystallinity”

The preparation, sample treatment, and analysis of the XRD patterns were carried out according to the recommendations by Moore and Reynolds (1997). All XRD analyses were performed on a Panalytical Empyrean X-ray diffractometer (installed at the University of Geneva). The analyses were performed in continuous-scan mode using Bragg–Brentano geometry, a step size of 0.013^∘ 2Θ per step, and a counting time or step of 350 s in the range of 4 to 70^∘ (45 kV, 40 mA, Cu K-α). The mineral determination was performed using the software HighScore Plus v. 3.0e. The measurement of the illite “crystallinity” (Kübler index, KI) was performed on glass plates, according to the recommendation of Weber (1972) to prepare “thin” texture compounds using 1.5–2.5 mg cm⁻². The samples were scanned from 6–11^∘ 2Θ, using a step size of 0.0131^∘ 2Θ, and a measuring time of 400 s per step in continuous mode. The illite “crystallinity” was determined using the program Newmod (Reynolds and Reynolds, 1996). KI values may empirically range from 0.060Δ^∘ 2Θ (ideally ordered muscovite) to 1Δ^∘ 2Θ (poorly ordered, I–S mixed layers). The KI values were calibrated against the Crystallinity Index Standard (CIS) introduced by Warr and Rice (1994), and the revised boundary limits of 0.32 and 0.52^∘ 2Θ for the anchizone were applied (Warr and Ferreiro Mählmann, 2015). All measurements were performed in the glycolated state.

3.3 Raman spectroscopy

Carbonaceous matter in metapelites and metasandstones was analyzed by Raman spectroscopy. The degree of graphitization is controlled by the maximum temperature achieved during regional metamorphism (Pasteris and Wopenka, 1991; Beyssac et al., 2002) and contact metamorphism (Aoya et al., 2010). Raman spectra were obtained using a confocal LabRAM spectrometer at the University of Geneva, equipped with a green 532.12 nm Nd:YAG laser coupled to an optical microscope (Olympus BX51, ×100 objective lens). The Raman spectrometer was calibrated with a silica standard. Between 10 and 15 spectra on individual organic-matter particles were obtained for each sample, and the average value was used for further analysis. Temperatures were obtained on thin sections following the measuring conditions and peak-fitting procedure described in the literature for regional low-grade metamorphic samples (Kouketsu et al., 2014) and for contact metamorphic samples (Aoya et al., 2010). The error associated by using the full-width-at-half-maximum (FWHM) D1 is around ±30 ^∘C (Kouketsu et al., 2014), although relative uncertainties may be smaller.

3.4 Chlorite thermometry

Chlorite temperatures are derived from thermodynamic calculations based on the method proposed in Lanari et al. (2014). The calculation is built upon four linearly independent end-members, which are amesite, clinochlore, daphnite, and sudoite, and accounts for Fe–Mg, Tschermaks, and vacancy substitution. The chlorite geothermometer is restricted to chlorites with Si <3 a.p.f.u. and does not require the determination of Fe³⁺ (i.e., it assumes that $\mathrm{\Sigma }\mathrm{Fe}={\mathrm{Fe}}^{\mathrm{2}+}$). The assumption of $\mathrm{\Sigma }\mathrm{Fe}={\mathrm{Fe}}^{\mathrm{2}+}$ is a simplified approach, as the redox state is likely to affect the distribution of other elements in the crystalline structure of chlorite as well as its formation temperature (Lanari et al., 2014). For this study, a semi-empirical equation is employed (equation Chl 2 in Lanari et al., 2014) which can be used over a P–T range of 1–20 kbar and 100–500 ^∘C.

Table 2 Rock properties and initial conditions used in thermal modeling.

^a Data from England et al. (1980). ^b Data from Waples and Waples (2004).

3.5 Thermal modeling

Thermal modeling was used to estimate the lateral extension of the contact metamorphic aureole of the TPI as a function of time and the number of magmatic pulses and batches. Further, it was evaluated whether the TPI provides a reasonable heat source and duration for the carbonaceous matter and clay mineral transformation. The thermal modeling was performed by systematically varying the magmatic temperature, heating duration, number of magmatic pulses, and repose time between pulses in a full-factorial design (Supplement 2). We performed three model runs with one, two, and three batches, respectively. The repose time between batches was chosen to have a total model duration close to the 90 kyr estimated for the felsic magmatism (i.e., 12.50–12.43 Ma) or the 160 kyr estimated for the total magmatic activity, including the mafic underplating (Michael, 1984, 1991; Michel et al., 2008; Leuthold et al., 2012). Thermal modeling was performed using MATLAB v. R2013b with code provided by Gerya (2010), as shown in Eq. (1).

$$\begin{array}{}\text{(1)}& {\displaystyle \frac{\partial T}{\partial t}}={\displaystyle \frac{\mathit{\kappa }}{\mathit{\rho }{c}_p}}\mathrm{\Delta }T,\end{array}$$

where κ is thermal conductivity, T is temperature, t is time, c is the specific heat capacity, and ρ is the density (compare Table 2). A 2D model was preferred over a 1D model to account for the location of the sampled profiles at the eastern extremity of the intrusion. The employed model is simplified by neglecting heat advection by fluid flow and heat production through radioactive decay. This simplification is justified by negligible advective heat transport for fluid expulsion of less than 10 wt % through a low permeable system, such as metapelites (Podladchikov and Wickham, 1994). The effects of latent heat during crystallization and the heat consumption during prograde metamorphic reactions are considered to cancel each other out (Bowers et al., 1990; Bodner, 2013). The area of the domain numerically modeled by a 169×191-point grid is a 5.7 km wide and 5 km deep box containing a uniform metapelitic lithology and a sill of 800 m thickness, which was emplaced in a depth between 2 and 3 km (Supplement 1). In the model setup, one time step reflects 1000 years of thermal changes. The boundary values at the edges of the grid were fixed, allowing heat to leave the system. The top boundary value at the Earth's surface has a fixed temperature of 20 ^∘C. The initial temperature of the host rock at depth is defined by a geothermal gradient of 30 ^∘C km⁻¹. The boundary values at the edge of the penetration were fixed during active intrusion and allowed to vary during cooling. Calculated T_max values were taken along a transect departing with an angle of 45^∘ downwards from the eastern extremity of the intrusion. The model assumes various assembly times, a thermal diffusivity of 10⁻⁶ m² s⁻¹, a granitic temperature between 800 and 1000 ^∘C, and up to 15 individual intrusion pulses. The batch repose time ranges between 10 and 84 kyr for the two-batch model and between 12.5 and 40.5 kyr for the three-batch model. The relatively high intrusion temperature of 1000 ^∘C is based on mineral thermometry (Leuthold et al., 2012). Further parameters for the pelitic host rock are a thermal conductivity of 1.37 $\mathrm{W}{\mathrm{m}}^{-\mathrm{1}}{\mathrm{K}}^{-\mathrm{1}}$ (England et al., 1980) and a specific heat capacity of 910 $\mathrm{J}{\mathrm{K}}^{-\mathrm{1}}{\mathrm{kg}}^{-\mathrm{1}}$ (Waples and Waples, 2004) and for the granitic intrusion are 2.95 $\mathrm{W}{\mathrm{m}}^{-\mathrm{1}}{\mathrm{K}}^{-\mathrm{1}}$ and 850 $\mathrm{J}{\mathrm{K}}^{-\mathrm{1}}{\mathrm{kg}}^{-\mathrm{1}}$ (Waples and Waples, 2004). All parameters are listed in Table 2. The validity of calculated T_max values depends on the quality of the input data, which include, for example, the thermal conductivity of individual lithologies, the heating duration, and the maximum temperatures reached. In order to estimate which of the various possible scenarios of the progression of the intrusion is most probable, the parameter space describing these scenarios must be examined. The root-mean-square error was used to describe the difference between the calculated and Raman T_max values and, thus, provides an indication of the probability of the corresponding scenario. The deviation from measured values, as shown in the optimization plots, was only calculated for the first six measurement points (240, 600, 900, 1100, 1200, 1500 m), as measurement points beyond this distance were not affected by contact metamorphism.

3.6 Data analysis

All statistical data analyses were performed in RStudio v. 1.2.5033 (RStudio Team, 2019), running R v. 3.6.3 (R Core Team, 2020). Temperature data in relation to the distance to the intrusion were fitted using a moving average with a window size of 2. The Raman vs. chlorite temperature was compared using Deming regression (Deming, 1964) in the R package “deming” v. 1.4, taking into account the measurement errors in both values. For the Raman temperature vs. KI comparison, Theil–Sen regression (Sen, 1968; Theil, 1950) from the R package deming v. 1.4 was used instead, since KI has no quantifiable error. Dendrograms were calculated on the Euclidean distances, using an average linkage algorithm.

Figure 3 Carbonaceous matter from the Torres del Paine intrusion in reflected light microscopy. (a–d) Individual samples, labeled in the top right corners, ordered by increasing metamorphic degree. Scale bars equal 20 µm.

Figure 4 Density distribution of temperatures obtained from Raman spectroscopy of carbonaceous matter from the Torres del Paine intrusion, ordered by increasing temperature. Sample name is given in top left corner of each panel; the dashed grey line indicates the mean.

4 Results

4.1 Carbonaceous matter and Raman temperatures

All samples show a high level of maturation, which is indicated by the high degree of reflectance (Fig. 3). Generally, it is difficult to differentiate between solid bitumen (pyrobitumen) and inertinite. The organic matter occurs either within a micrinitic groundmass (i.e., tiny pits of inertinite material) or as larger particles which appear to be coal-like and comprise good-quality telocollinite and vitrinite. Inertinite exhibits classic bogen structures (Fig. 3b). Raman measurements were performed on larger coal-like particles. The Raman temperatures are partly taken from the dataset in Süssenberger et al. (2018c). Temperatures derived from carbonaceous matter range between 242±16 (CPA 16-12) and 480±26 ^∘C (CPA 15-11; Table 1). The spread of individual analyses within one thin section is shown in Fig. 4. The between-sample variation is much larger than the within-sample variation (0.4 %), and, thus, our range of Raman temperatures is very reliable. Temperatures between 242 and 252 ^∘C, as detected for the samples located further from the intrusion – CPA 16-9, CPA 16-10, CPA 16-11, and CPA 16-12 (Fig. 1a, c) – are considered to reflect the regional metamorphic background prior to the emplacement of the TPI (Süssenberger et al., 2018b). Temperatures >260 ^∘C are considered thermally altered by the TPI. An increase in temperature as a function of proximity to the intrusion can be observed for the Monte Almirante profile (i.e., samples CPA 16-15, CPA 16-16, CPA 16-17, CPA 16-18, and 09PR44; Fig. 1) and for the Paso John Gardner profile (i.e., samples CPA 15-10, CPA 15-11, CPA 15-12, CPA 15-13, and CPA 15-15; Fig. 1). The highest temperatures are observed for samples (CPA 15-12 and CPA 15-13) located at the western extremity (Paso John Gardner) of the intrusive complex (Fig. 1, Table 1). At distances exceeding 1500 m from the intrusion, samples show typical regional metamorphic temperatures (Fig. 5).

Figure 5 Temperatures obtained from carbonaceous matter as a function of distance towards the Torres del Paine intrusion. The red line indicates the moving average with a window size of 2.

4.2 Chlorite temperatures and generations

The chemical chlorite analyses used for geothermometry are taken from the dataset in Süssenberger et al. (2018a). T calculations were performed for 13 samples, and results are presented in Table 1. The calculated mean chlorite temperatures range between 220±45 ^∘C (sample CPA 16-10) and 378±43 ^∘C (sample 09PR44). Based on the derived chlorite temperatures, three populations can be distinguished: (i) chlorite which reflects T_max conditions of >260 ^∘C set during contact metamorphism (CM), (ii) chlorite which reflects regional metamorphic temperature conditions of 200–260 ^∘C (RM), and (iii) chlorites which grew during retrograde metamorphism at temperatures <200 ^∘C (R). The temperatures are given as approximate value or ranges, since the expected temperatures for contact and regional metamorphism may vary with the sample location. Regional metamorphic samples, such as CPA 16-10, consist generally of two chlorite populations, RM and R, whereas samples which experienced contact metamorphism, such as CPA 16-18, often feature all three chlorite populations (i.e., CM, RM, and R; Table 1, Figs. 6 and 7). Microprobe analyses on chlorites from sample CPA 16-18 reveal that chlorite generations RM and CM occur in distinct clusters distributed in the matrix (Fig. 7). Retrograde chlorite (R) formed after the contact metamorphism as a replacement for epidote and titanite (Fig. 7).

Figure 6 (a) Chlorite temperature density plot, comparing chlorite generations in a regional metamorphic sample (16-10) and in a contact metamorphic sample (16-18) from the Torres del Paine intrusion. (b–c) Chlorite dendrogram plots indicate potential chlorite generations.

Figure 7 (a) Scanning electron microscope (SEM) image of sample 16-18. Chlorite point analyses are marked by red stars. Label CM indicates chlorites recording contact metamorphic temperatures, and label RM represents chlorites which record regional metamorphic temperatures. (b) Allanite rim surrounding an epidote grain. (c) Optical and backscattered-electron images showing retrograde chlorite generation replacing epidote. (d) Retrograde chlorite generation replacing titanite.

4.3 Kübler index based on XRD analyses

The Kübler index (KI) was determined for 11 metapelites in the <2 µm grain size fraction. Part of the data used here is taken from the datasets in Süssenberger et al. (2018b, c). The values vary between 0.15 and 0.34Δ^∘ 2Θ (Table 1, Fig. 8) and indicate upper anchizonal to epizonal metamorphic conditions. KI values show only a faint trend of decreasing values (i.e., better “crystallinities”) with decreasing distance towards the intrusion. Epizonal KI values are recorded in the western part of the TPI, where samples experienced higher contact metamorphic temperatures. At the eastern side of the TPI, the KI values remain mostly invariant as a function of temperature.

Figure 8 Comparison between T_max obtained from Raman spectroscopy on carbonaceous matter and (a) chlorite T_max and (b) KI values. Note that Raman- and chlorite temperatures have a strong correlation with a small confidence interval and regression slope of approximately 1. In contrast, the uncertainty in the KI is too large for any reliable temperature reconstruction.

Figure 9 Optimization maps for the different model runs showing the deviation between measured and modeled temperatures in the contact aureole of the Torres del Paine intrusion. (a) One-batch model at 800 ^∘C intrusion temperature. (b) One-batch model at 1000 ^∘C intrusion temperature. (c) Two-batch model at 800 ^∘C intrusion temperature. The model providing the best fit is marked with a red margin. (d) Two-batch model at 1000 ^∘C intrusion temperature. (e) Three-batch model at 800 ^∘C intrusion temperature. (f) Three-batch model at 1000 ^∘C intrusion temperature.

Figure 10 Comparison between measured Raman temperatures and modeled temperatures in the contact aureole of the Torres del Paine intrusion, using the best framework parameters (Fig. 9c). We are using the two-batch model at 800 ^∘C intrusion temperature. The best-fitting model is marked in red and represents a model with short heating duration, five pulses, and long repose times (i.e., 15 kyr).

4.4 Thermal modeling

Numerous parameter settings of the thermal model were simulated, followed by comparing corresponding calculated T_max values with measured T_max values (RSCM). Representative examples of thermal modeling results are shown in Figs. 9 and 10, as well as in the Appendix. The growth of the TPI was simulated by a single batch (Fig. 9a, b), a two-batch (Fig. 9c, d), or a three-batch model (Fig. 9e, f), with magmatic temperature varying between 800 and 1000 ^∘C. The best match between measured and calculated T_max values for a one-batch model and an intrusion temperature of 800 ^∘C is achieved by applying 13 pulses, each having a heating duration of 3 kyr, and a repose time of 15 kyr (Fig. 9a; ΔT=54 ^∘C). If an intrusion temperature of 1000 ^∘C is assumed, the best-fitted model corresponds to five pulses, with each having a heating duration of 3 kyr, and a repose time of 15 kyr (Fig. 9b; ΔT=76 ^∘C). The best fit for a two-batch model for an intrusion temperature of both 800 and 1000 ^∘C is achieved by applying a batch repose time of 10 kyr; five pulses per batch, each having a heating duration of 3 kyr; and a total repose time of 15 kyr (Fig. 9c, d, 10; ΔT=54 ^∘C and ΔT=82 ^∘C, respectively). The best fit for a three-batch model and an intrusion temperature of 800 ^∘C can be accomplished by applying a batch repose time of 27.5 kyr; five pulses per batch, each having a heating duration of 5 kyr; and a repose time of 5 kyr (Fig. 9e; ΔT=67 ^∘C). If an intrusion temperature of 1000 ^∘C is assumed, the best-fitted model corresponds to one pulse with a heating duration of 3 kyr and a batch repose time of 37.5 kyr (Fig. 9f; ΔT=81 ^∘C).

Although relative temperature differences are small between models, the overall best fit is achieved with a total assembly time of 150 kyr: a magmatic temperature of 800 ^∘C, a two-batch model (batch repose time of 10 kyr) with five pulses, short heating durations (3 kyr), and long pulse repose times (i.e., 15 kyr; Figs. 9c and 10). The best-fitting one-batch model would have an equal predictive ΔT value but would lead to excessively long magmatic-assembly durations.

5 Discussion

5.1 Constraining the lateral extension of the contact aureole of the Torres del Paine intrusion

The Torres del Paine contact aureole can be divided into three zones. First is the inner contact aureole, which is characterized by the appearance of K-feldspar and the breakdown of muscovite and biotite (Bodner, 2013). Second is the outer contact aureole, which is characterized by the appearance of cordierite and biotite, the breakdown of chlorite, and the modal decrease in muscovite (Bodner, 2013). Bodner (2013) calculated reaction temperatures of 480 and 540 ^∘C for reactions of the cordierite-in and the K-feldspar-in isograds, respectively. Third is the wider contact aureole, which is more difficult to differentiate and will be referred to as the contact-influenced zone (CIZ). The CIZ is characterized by elevated Raman temperatures, indicating anchizonal to epizonal metamorphic conditions. It records only slightly higher temperatures than the regionally metamorphosed host rock. Significant differences in the type and quantity of clay minerals are not expected to occur due to the relatively high degree of regional metamorphism. Besides the regional metamorphic mineral paragenesis (albite, chlorite, and illite), the CIZ is recognized by the appearance of epidote and retrograde Fe-rich chlorite (Fig. 7). Pseudomorphic incomplete replacement for epidote and titanite by chlorite and allanite rims around epidote crystals indicate the circulation of magmatic and hydrothermal fluids (Fig. 7b–d). Although this chlorite cannot be unequivocally ascribed to one of the three observed populations, it is likely that it corresponds to the low-grade retrograde population (R).

The lateral extension of the CIZ is determined by the measured T_max derived from RSCM and by thermal modeling. Based on RSCM, the CIZ can locally be identified at a distance of 1.5 km and, thus, spreads ca. 1.1 km further out than the hornfels aureole (Fig. 5). Temperatures obtained from RSCM are maximum temperatures and are not biased by other low-grade metamorphic organic-matter generations. This observation is in line with our own results (Fig. 4) and with studies performed in contact metamorphic aureoles (Velde and Lanson, 1993; Abad et al., 2014). However, the uncertainty in the exact distance of samples with respect to the subterraneous extension of the intrusion is a potential source of error.

Mathematically, an approximation for the lateral extension of the contact metamorphic aureole can be obtained by employing Eq. (2), which describes the distance x that an isotherm recedes in a certain time t from the contact downwards as the intrusion cools (Walther and Wood, 1984):

$$\begin{array}{}\text{(2)}& x=\sqrt{k\cdot t},\end{array}$$

where k describes the effective thermal diffusivity of around $\mathrm{7}\times {\mathrm{10}}^{-\mathrm{3}}$ cm² s⁻¹ (Hanley et al., 1978). Equation (2) implies that the ca. 450–500 ^∘C isotherm, which was initially located at the contact, would advance ca. 200 m through the host rock within 2 kyr, ca. 470 m within 10 kyr, and ca. 570 m within 15 kyr. This corresponds to an average rate of advance of $\mathrm{1.5}\times {\mathrm{10}}^{-\mathrm{7}}$ to $\mathrm{1.2}\times {\mathrm{10}}^{-\mathrm{7}}$ cm s⁻¹. A heating duration of 2 kyr corresponds roughly with the expectation of observable temperatures for samples CPA 15-11 and 09TP19 (Table 1).

Thermal modeling reveals that one- and two-batch models with magmatic temperatures of 800 ^∘C can best match measured RSCM values with increasing distance from the intrusion (Figs. 9 and 10). While thermal models generally overestimate the first T_max value at a 240 m distance from the intrusion, they do not provide enough heat to match T_max values at distances of >1000 m (Fig. 10, Supplement 1). Multiple-pulse models with short heating durations (i.e., 3 kyr) and rather long repose times (i.e., 15 kyr) are preferred as they are able to fit the T_max values of the CIZ better. These models avoid the heating of sediments above pre-intrusive temperatures by constantly allowing the sediments to cool and laterally equilibrate temperatures. The assumption of there being several discrete magma pulses is corroborated by geochronological, petrological, and field observations (Michel et al., 2008; Leuthold et al., 2012; Bodner, 2013). The overall best fit of ΔT_max values (i.e., 54 ^∘C) is achieved with a total assembly time of 150 kyr, a magmatic temperature of 800 ^∘C, a two-batch model (batch repose time of 10 kyr) with five pulses per batch, short heating durations (3 kyr), and long pulse repose times (15 kyr; Figs. 9c and 10). This observation is in agreement with calculations by Bodner (2013) and Leuthold et al. (2012) for the total assembly duration of the TPI. Although a similarly good fit can be achieved with a one-batch model composed of 13 pulses, each having a heating duration of 3 kyr, and a repose time of 15 kyr, the total assembly duration of 215 kyr in this model setup would be too high and not supported by the geochronological data. A total assembly duration of 90±40 kyr, as established by Michel et al. (2008) for the emplacement of the Paine granite, would be too short to fit our model data. Therefore, the calculated assembly duration of 150 kyr implies that the three mafic sheets which underplated the laccolith between 12.50 and 12.43 Ma (Leuthold et al., 2012) thermally affected the sediments at the sampled profile. On the one hand, this is quite conceivable as the mafic sill complexes were emplaced at higher temperatures (>900 ^∘C; Leuthold et al., 2014) compared to the Paine granite; and on the other hand, this might also imply that their subterraneous extension is probably larger than expected from outcrop mapping (Fig. 1). Principally, larger ΔT_max values at increased distance from the intrusion (between 600 and 1500 m) may be caused by other scenarios. For instance, exothermic hydration reactions in the wider contact aureole could have increased the T_max value of the host rock by ca. 30–40 ^∘C (Bodner, 2013), thus leading to higher observed RSCM temperatures at a greater distance from the intrusion than our models can explain. Additionally, limited fluid circulation, as indicated by stable-isotope measurements, would play an important role in heat transport (Baumgartner and Valley, 2001; Baumgartner et al., 2014; Siron et al., 2017). We note that a general lack of fluid circulation would decrease the heat equilibration within the host rock, while our models underestimate the naturally occurring heat transport. This means that a lack of fluid activity cannot explain a deviation between model and RSCM data in the direction we observe. It could, however, be feasible that channelized fluid transport occurred that would have led to higher heat equilibration within the host rock while at the same time locally indicating a reduced fluid activity if measurements were conducted outside of the fluid channels.

5.2 Reaction progress between clay minerals and organic matter

The rapid temperature increase during the emplacement of the TPI followed by retrogressive alteration caused differences in the reaction progress between carbonaceous matter and clay minerals. Carbonaceous matter has reacted to completion on timescales of several thousand years, whereas clay minerals reacted incompletely upon thermal changes.

The good correlation between T_max values obtained from RSCM and the T_max values from chlorite generations (i.e., CM) indicates that some of the chlorites re-equilibrated or re-crystallized during contact metamorphism (Fig. 8). Other chlorite generations, yielding lower formation temperatures, were excluded for this figure, as they do not record the contact metamorphic event (i.e., these chlorites record regional deformation and/or retrograde metamorphism). Contrary to the chlorites, which record various temperature–time–chemistry conditions, all organic-matter particles seem to have attained a similar state of reaction during the contact metamorphic conditions (Fig. 4). While carbonaceous matter and, to some extent, chlorite re-equilibrated with the contact metamorphic conditions, illite crystals remained apparently unaffected in the CIZ. However, at distances of <500 m from the intrusion, epizonal values of 0.15Δ^∘ 2Θ reflect a thermal influence of the granitic intrusion. The KI data suggest upper anchizonal to epizonal temperature conditions, which broadly overlap with the temperatures obtained from RSCM (Table 1, Fig. 8). Figure 8 demonstrates that Raman and chlorite T_max values have a strong correlation with a small confidence interval and regression slope of approximately 1. In contrast, the uncertainty in the KI is too large for reliable temperature reconstruction and thermal modeling as required for this study. In particular, a faint temperature increase as indicated by Raman temperatures for samples CPA 16-10 to CPA 16-18 is not evident solely based on KI values (Table 1). A similarly discrepancy is observed for sample CPA 15-12 yielding a Raman temperature of ca. 375 ^∘C and an epizonal KI value of 0.29Δ^∘ 2Θ (Table 1). It can be concluded that illite “crystallinity” does not entirely record the contact metamorphic event, and thus, a similarly detailed study to that performed for chlorite would be necessary to reveal different illite and muscovite generations (e.g., using illite polytypes and geochemical data). However, it is difficult to assess if the illite reaction kinetics truly lag behind, as the obtained illite results are potentially limited by the methodology we employed. Contrary to chlorite, no textural or compositional analyses of the illite–muscovite assemblages were undertaken. KI measurements were determined on bulk samples for a given grain size range (<2 µm) and, thus, are potentially a mixture of different illite–mica generations. In contrast, chlorite crystallization temperatures were obtained from point analyses on individual crystals. However, comparable observations of a sluggish reaction rate of illite compared to carbonaceous matter in a contact metamorphic setting were also reported by Velde and Lanson (1993), Olsson (1999), Suchý et al. (2004), and Abad et al. (2014). Possible explanations for the slow reaction kinetics of illite and some chlorite generations (RM) within the CIZ include a short heating time of a few thousand years; lack of induced tectonic stress; and probably a low activity of fluid and of K⁺, Mg⁺, and Fe²⁺, as the transformation from illite to muscovite is significantly inhibited at low water : rock ratios (Whitney, 1990; Velde and Lanson, 1993). Whitney (1990) showed that if the host rock is dry, little recrystallization is expected to occur. Remnants of fluid flow activity, such as veins departing from dikes, are well documented at the eastern part of the laccolith (Siron, 2017). Numerous dikes depart from the intrusion and can be traced at up to several hundred meters into the host rock (Fig. 2d). Fissure-like, millimeter-wide veins at their end document locally limited fluid infiltration into the host rock (Siron, 2017). The absence of a significant oxygen-isotopic infiltration front supports the assumption of there being a limited amount of fluid that has affected the adjacent host rock (Baumgartner and Valley, 2001). Siron (2017) was able to identify minor magmatic fluid infiltration near to the contact (ca. 100 m) by low δD values and concluded that the fluid was successively expelled towards the end of the multiple rapid crystallization events that formed granite I (12.49 Ma; Leuthold et al., 2012). Based on the aforementioned studies, we believe that a sluggish recrystallization of some chlorite and illite generations may have happened due to a lack or low mobility of magmatic fluid and/or low availability of cations. Alternatively, the K⁺ activity is a major component in illite reactions. We resolved equations in Huang et al. (1993) for [K⁺], which allowed us to estimate the likely K concentration in the magmatic fluids, using our estimates for temperature, assembly time, and change in smectite content from 0.07 to 0.02 (Süssenberger et al., 2018b). We found that indeed the [K⁺] of the magmatic fluids must have been as low as 0.37 ppb under these framework conditions.

Thermal maturity reactions for carbonaceous matter show fairly similar activation energies (E_a 9–12 kcal mol⁻¹) for different geological environments such as geothermal systems or contact metamorphic settings (Barker, 1989). However, the rate of reaction can significantly vary between different geological settings, as it is a function of temperature. The difference in reaction rate is also expressed by the time required for carbonaceous matter to stabilize thermal maturation. Barker and Pawlewicz (1986) showed that the duration of heating is not important anymore in hydrothermal systems after 10 kyr and in contact metamorphic systems after 1–10 years. Similar results were obtained by Hill et al. (2004) and Mori et al. (2017). Generally, a system with low free water content and low porosity and permeability would be less likely to respond quickly to a thermal impulse that would require rapid dissolution and precipitation from solution (Velde and Vasseur, 1992).

Our results and similar observations from other studies (e.g., Barker and Pawlewicz, 1986; Sweeney and Burnham, 1990; Hill et al., 2004) indicate that the heating duration is essential for the clay mineral reaction progress and that KI values are mainly controlled by kinetics (Merriman and Frey, 1999). Merriman (2005) concludes that the rate at which clay mineral reactions progress from mature to supermature assemblages can range from 10⁴ to 2×10⁹ years or more. Derived kinetic values (low activation energies and reaction-rate scaling) support the case for the time sensitivity of illite mineral reactions (Velde and Vasseur, 1992). On the contrary, carbonaceous matter maturation is, in the first instance, a function of temperature, and after reaching T_max it ceases to have a significant reaction.

5.3 Textural local equilibrium and chlorite generations

Three chlorite populations can be distinguished based on crystallization temperatures derived from thermodynamic calculations (Figs. 6 and 7). Each chlorite population reflects chlorite re-crystallization and re-equilibration events in different time–temperature–chemistry conditions. The chlorite populations comprise (1) chlorite which records T_max conditions (T>260 ^∘C) set during contact metamorphism (CM), (2) chlorite recording regional metamorphic temperature conditions of 200–260 ^∘C (RM), and (3) chlorite which re-equilibrated during retrograde metamorphism at temperatures <200 ^∘C (R). Petrographic observations indicate that chlorite populations R and CM predominantly formed as replacements for other minerals (i.e., epidote and titanite; Fig. 7) rather than as pervasive pore filling or grain coating. Since epidote is interpreted as a contact metamorphic phase, the replacing chlorite has to be regarded as a retrograde mineral phase. With decreasing distance from the intrusion, the amount of chlorite affected by contact metamorphism gradually increases (Fig. 6). Interestingly, a sample which experienced a relatively high metamorphic temperature (sample CPA 16-18, 296 ^∘C) exhibits all three chlorite populations, indicating that chlorite of regional metamorphism (RM) did survive regardless of contact metamorphism and retrograde alteration (Figs. 6 and 7). Thus, the studied chlorites indicate incomplete re-equilibration throughout the rock history, and in polymetamorphic environments, such as in this study, individual chlorites do not reflect T_max values but rather the local temperature–chemistry conditions controlled by local kinetics on the microscale. Note, that temperatures above the expected Raman T_max value for contact metamorphism, as identified in the chlorite dendrogram plot, may indicate chlorites of detrital origin (Fig. 6).

The textural distribution of different chlorite generations within one sample reveals local clusters of chlorite grains that attained equilibrium in varying time–temperature conditions (Fig. 7). The coexistence of different chlorite generations within a sample and their distribution across clusters supports the argument for the pivotal role of water-saturated pore spaces and the local availability of cations (i.e., Mg⁺ and Fe${}^{\mathrm{2}+/\mathrm{3}+}$) in controlling the reaction progress and re-equilibration and neocrystallization of chlorite. Fluids can be released directly from the intrusion or through continuous and discontinuous metamorphic reactions (e.g., Baxter and Caddick, 2013). Evidence for magmatically induced fluid flow in the inner contact aureole of the TPI, <150 m away from the magmatic body, is given by chlorine and hydrogen stable-isotope data in biotite (Siron et al., 2017; Bodner, 2013). At distances >150 m, chlorine content in biotite is indistinguishable from unaltered biotites (Bodner, 2013; Siron et al., 2017). At even larger distances from the TPI, Bodner (2013) observed hydration reactions, which supports our hypothesis of inhomogeneous fluid distribution and cation activity resulting in the patchy pattern of different chlorite generations. However, pore spaces filled with non-aqueous fluids, such as petroleum or natural gas, impede recrystallization, transport, and dissolution processes (e.g., Whitney, 1990). In the investigated samples, we consider such non-aqueous pore fluids to be unlikely, since anchizonal temperatures of the host rock were already achieved 40 Myr prior to the intrusion and any potentially generated oil or gas would have already migrated away from the source place.

Differing grain sizes, porosities, and textural inhomogeneities in the metapelitic host rock result in kinetically controlled reaction progressions and promote the observed cluster pattern of chlorite (Fig. 7). Local areas with higher porosity contribute to higher permeability and increased fluid infiltration relative to clay-dominated matrix areas. As a result, pore fluids facilitate the local diffusional exchange and buffer mineral reactions in pore spaces, so the chlorite chemistry responds more thoroughly to local changes in composition as a function of temperature. In these areas, chlorites stay connected with the local micrometer-scale environment and may record ongoing reactions at any time if favorable conditions are attained. In clay-rich matrix areas, the permeability is significantly reduced, and the exchange within adjacent areas is increasingly limited. In other words, the overall achieved rock temperature is, in comparison to the availability of interstitial fluids and elements, of minor importance.

6 Conclusions

The thermal effect of the emplacement of the TPI on the enclosing rocks had different consequences for inorganic and organic compounds of the host rock. The impact of the contact thermal alteration of the pre-intrusive regional metamorphosed host rock is documented by elevated Raman temperatures, high-temperature chlorite generations, and the appearance of epidote and retrograde Fe-rich chlorite. However, if and to what extent illite was affected by the contact metamorphism cannot be unequivocally revealed by KI values alone. Although KI values are a well-known tool to decipher regional metamorphic conditions, their applicability to contact metamorphic settings has often been reported to be unreliable in the literature. This observation is confirmed by our study, which implies that KI, in contrast to RSCM, does not allow reliable metamorphic temperature reconstructions as required for thermal modeling. Temperatures obtained by RSCM are the most reliable and unequivocal indicators on timescales of several thousands of years to determine the lateral extension of the TPI contact aureole.

The contact-influenced zone can be identified at a distance of 1.5 km and spreads ca. 1.1 km further out than the outcropping hornfels aureole. The combination of measured T_max values (RSCM) and calculated T_max values (thermal modeling) reveals that the TPI was emplaced by multiple pulses. The best match between measured and calculated ΔT_max values (i.e., 54 ^∘C) is achieved with a total assembly time of 150 kyr, a magmatic temperature of 800 ^∘C, a two-batch model (batch repose time of 10 kyr) with five pulses, short heating durations (3 kyr), and long pulse repose times (15 kyr).