Discrimination of thermodynamic and kinetic contributions to the heavy rare earth element patterns in metamorphic garnet

Variations of rare earth element (REE) concentrations in metamorphic garnet are an important source of information of geodynamic and geochemical processes in the deeper Earth. In order to extract this information, the thermodynamic equilibrium and kinetic contributions of the REE uptake in garnet must be distinguished and quantified. Utilizing high‐resolution trace element and μ‐Raman mapping together with combined thermodynamic–geochemical–diffusion models, we demonstrate that the equilibrium and kinetic aspects of the REE uptake in metamorphic garnet can be discriminated by interpreting 2D trace element mapping in a single sample. The heavy (H) REE (Tb to Lu) zoning in the investigated garnet from a high‐pressure blueschist comprises an inner part with an overall decrease from core to inner rim, followed by a concentric zone of HREE enrichment and a drastic HREE decrease towards the outermost rim. The central peak in the garnet core decreases in intensity with decreasing atomic number of the REE. The broad overall shape of this pattern resembles those often observed in metamorphic garnet from different rock types and tectonic settings. Superimposed on this trend is a concentric pattern of minor recurring fluctuations in the HREE concentrations with at least six regularly spaced sets of peaks and troughs along the entire garnet radius. Comparison of the observed inclusion suite, the trace element maps and thermodynamic–geochemical models show that the inner part with decreasing HREE concentrations results from fractional garnet growth in an unchanged mineral assemblage, whereas the REE enrichment zone is caused by the breakdown of titanite. We suggest that the width of the central peak is controlled by the bulk permeability of the interconnected transport matrix and the fraction of matrix minerals that the garnet equilibrates with. The superimposed REE fluctuations result from changing element transport properties of the host rock and mark recurring changes from equilibrium REE uptake to transport‐limited REE uptake in garnet. Such fluctuating element transport properties can be best explained by pulse‐like fluid fluxes that rhythmically change the interconnectivity of the intercrystalline transport matrix. Increasing numbers of published spatially highly resolved REE analyses show that such trace element fluctuations are common in metamorphic garnet indicating that recurring changes in rock permeabilities due to pulsed fluid fluxes are a common phenomenon during metamorphism.

zoning in the investigated garnet from a high-pressure blueschist comprises an inner part with an overall decrease from core to inner rim, followed by a concentric zone of HREE enrichment and a drastic HREE decrease towards the outermost rim. The central peak in the garnet core decreases in intensity with decreasing atomic number of the REE. The broad overall shape of this pattern resembles those often observed in metamorphic garnet from different rock types and tectonic settings. Superimposed on this trend is a concentric pattern of minor recurring fluctuations in the HREE concentrations with at least six regularly spaced sets of peaks and troughs along the entire garnet radius. Comparison of the observed inclusion suite, the trace element maps and thermodynamic-geochemical models show that the inner part with decreasing HREE concentrations results from fractional garnet growth in an unchanged mineral assemblage, whereas the REE enrichment zone is caused by the breakdown of titanite. We suggest that the width of the central peak is controlled by the bulk permeability of the interconnected transport matrix and the fraction of matrix minerals that the garnet equilibrates with. The superimposed REE fluctuations result from changing element transport properties of the host rock and mark recurring changes from equilibrium REE uptake to transport-limited REE uptake in garnet. Such fluctuating element transport properties can be best explained by pulse-like fluid fluxes that rhythmically change the interconnectivity of the intercrystalline transport matrix. Increasing numbers of published spatially highly resolved REE analyses show that such trace element fluctuations are common in metamorphic garnet indicating that recurring changes in rock permeabilities due to pulsed fluid fluxes are a common phenomenon during metamorphism.
K E Y W O R D S garnet, kinetics, REE, thermodynamics 1 | INTRODUCTION Garnet, one of the major phases in many metamorphic rocks, is by far the most versatile mineral with respect to the recording of geodynamic and geochemical processes (e.g., Chernoff & Carlson, 1997;Gaidies et al., 2006Gaidies et al., , 2015Ganguly et al., 1996;George & Gaidies, 2017;Inui, 2006;Konrad-Schmolke et al., 2005;Lasaga & Jiang, 1995;O'Brien, 1997;Spear & Daniel, 2001;Spear & Selverstone, 1983;Yang & Rivers, 2002). It is stable over a wide pressure and temperature range and thus can grow continuously from shallow depths to high pressures, thereby recording processes in the entire lithosphere. Its major and trace element composition is sensitive to subtle chemical changes in the surrounding rock volume (e.g., Jamtveit & Hervig, 1994;Kotkov a & Harley, 2010;McCammon & Kopylova, 2004;Moyen, 2009;Yang & Rivers, 2002), which makes garnet an excellent recorder of mineral reactions, fluid-rock interaction and element transport processes in the host rock (e.g., Konrad-Schmolke, Zack, et al., 2008;Skora et al., 2006). Furthermore, because of the extremely slow diffusive intracrystalline element transport in garnet (e.g., Carlson, 2006Carlson, , 2012Tirone et al., 2005), compositional and isotopic growth zonation is very commonly preserved and can be used to extract information about geological processes from single garnet crystals (e.g., Caddick et al., 2010;Hermann & Rubatto, 2003;Whitehouse & Platt, 2003) and provide information relevant for the analysis of mineral nucleation and growth kinetics (George et al., 2018;George & Gaidies, 2017). This holds especially true for rare earth element (REE) growth zonation, which is commonly preserved even if major element variations are diffusively modified (e.g., Kotkov a & Harley, 2010). As garnet is a major REE-bearing mineral and fractionates strongly between light (L) and heavy (H) REEs, REE concentrations, REE variations and REE distribution patterns extracted from garnet are often interpreted in terms of reflecting geodynamic and geochemical processes in the host rock (e.g., Rubatto et al., 2020). However, it is still heavily debated and unclear to what extent REE uptake in garnet is controlled by thermodynamic and geochemical parameters, such as mineral reactions or changing fluid compositions in the host rock (e.g., Konrad-Schmolke, Zack, et al., 2008), or whether kinetic aspects, such as sluggish element transport in the host rock (e.g., Skora et al., 2006), play the key role during REE uptake. It is evident that both factors, thermodynamic and kinetic, are crucial for the characterization and quantification of geodynamic processes. This is because mineral reactions recorded in the REE patterns reflect the reaction path of the rock, whereas element transport properties can be tied to rock permeabilities and fluid percolation. Furthermore, the distinction between thermodynamic and kinetic contributions to garnet growth is important for garnet geochronology as the distribution of Lu and Sm within the garnet influences the ages derived from Lu/Hf and Sm/Nd geochronology. Unambiguous evidence for the effects of both contributions recorded in a suitable sample is still lacking predominantly due to the lack of sufficiently spatially resolved REE analyses in single garnet crystals.
We employ high-resolution trace element mapping, utilizing laser ablation inductively coupled plasma mass spectrometry (LA-ICP-MS), that illustrates subtle changes in the heavy rare earth element (HREE) composition of a single garnet from a subducted and exhumed high-pressure rock. When combined with μ-Raman mapping that visualizes the distribution of μm-scale mineral inclusions in the garnet and compared with thermodynamic-geochemical models, we are able to distinguish between thermodynamic and kinetic contributions to the REE uptake in garnet. Thus, we obtain detailed information about the reaction history of the host rock as well as its element transport properties during garnet growth. Characteristic REE patterns together with systematic REE variations in our sample reflect the combined effects of mineral reactions as well as recurring fluctuations in rock permeabilities during high-pressure metamorphism.

| Geological framework
The sample investigated in this work is a garnet-bearing blueschist from the high-pressure Punta Balandra unit in the Dominican Republic on the island of Hispaniola. In the north-western part of Hispaniola, on the Samana peninsula, various types of high-pressure basement rocks are exposed (see fig. 1 in Escuder-Viruete & Pérez-Estaún, 2006). Among these are blueschist, eclogites, serpentinites and mélange associations that are formed at different levels of a Late-Jurassic to Eocene subduction complex (Nagle, 1974). The rocks of the Punta Balandra unit, located at the southern shore of the Samana peninsula, record the early stages of southwest-directed subduction of the North American plate beneath the Caribbean island arc (Mann et al., 1991). Within the Punta Balandra unit, mafic high-pressure boudins and blocks are intercalated with mica schists and carbonates.
During subduction, the Punta Balandra unit recorded a prograde pressure-temperature (P-T) evolution along a steep P-T gradient between $0.5 GPa and 350 C and 2.4 GPa and 650 C (Escuder-Viruete & Pérez-Estaún, 2006). Different stages of the subduction process are excellently preserved in various forms of mineral assemblages either as numerous inclusions in large garnet porphyroblasts in different rock types or as metastable matrix assemblages in rocks with different extents of proand retrograde metamorphic overprinting. Several detailed petrological and geochronological works that focussed on the rocks of the Punta Balandra unit enable these rocks to be employed as a natural laboratory for the investigation of subduction processes (Catlos & Sorensen, 2003;Escuder-Viruete et al., 2011, b;Goncalves et al., 2000;Joyce, 1991). Based on thermodynamic modelling and conventional thermobarometry, the metamorphic evolution of the Punta Balandra rocks can be reconstructed in detail and six distinct metamorphic stages have been distinguished based on different co-existing mineral assemblages by Escuder-Viruete and Pérez-Estaún (2006).
Stage 1: In some eclogite facies garnet porphyroblasts, rare relics of the early prograde greenschist to blueschist facies stages are preserved in the form of inclusions of calcic-amphibole, albite, chlorite and phengite. The inferred pressure and temperature conditions are $0.5 GPa and 450 C. Stage 2: Lawsonite inclusions and their pseudomorphic replacement by epidote and white mica together with rutile, phengite, paragonite and clinopyroxene inclusions in garnet are interpreted as relics of the lawsonite blueschist facies stage. Pressure and temperature conditions during that stage were between 0.8 and 1.5 GPa at temperatures between 350 and 450 C.
Stage 3: Prograde blueschist to eclogite facies mineral assemblages of garnet, epidote, phengite, paragonite, clinopyroxene and rutile are also found as inclusions in some eclogite facies garnet porphyroblasts. This stage corresponds to pressure and temperature conditions $1.5-1.8 GPa at 500-550 C. Stage 4: Peak eclogite facies metamorphism of the rocks is reflected in well-preserved garnet + clinopyroxene + phengite + rutile + quartz assemblages. Peak metamorphic conditions were in the order of 2.4 GPa and $650 C. Stage 5: The eclogites are partially hydrated and statically retrogressed to garnet-bearing blueschists with the matrix assemblage sodic amphibole + titanite + epidote + quartz + phengite within which centimetre-sized eclogite facies garnet porphyroblasts can be found. This stage is related to fluid influx and metasomatism between peak metamorphic conditions and 1.3 GPa at 550 C along a cold retrograde exhumation path. Stage 6: Late-stage greenschist facies retrogression is reflected in the breakdown of garnet to form epidote and chlorite as well as in the formation of albite and calcic amphibole from clinopyroxene and sodic amphibole, respectively. This last metamorphic imprint occurred at $0.8 GPa and below 500 C.

| Sample description
The sampled rock is a garnet-bearing blueschist statically equilibrated during early Stage 5 of this subdivision. The blueschist is unfoliated with a medium-grained matrix consisting of up to several millimetre-large, often euhedral, sodic amphibole grains and quartz aggregates within which up to centimetre-sized euhedral garnet porphyroblasts are homogeneously distributed (Figure 1a). Minor phases in the matrix are up to several millimetres in size and include titanites (often euhedral), chlorite that partly replaces garnet and amphibole, white mica, epidote, apatite, zircon and minor iron hydroxide grains that reflect very late weathering of chlorite and sodic amphibole. Sodic amphibole and garnet are inclusion-rich. In sodic amphibole, the most abundant inclusion phase is titanite, most of which has rutile cores or at least small rutile inclusions. Less abundant inclusions are white mica, epidote and omphacite (Figure 1b,c) that is likely relict from the peak metamorphic phase assemblage. Inclusion phases in the cores of the large garnet porphyroblasts are titanite (lacking rutile inclusions), apatite, epidote, minor zircon and quartz. In the garnet rims, rutile is the most abundant inclusion phase, but epidote, apatite, zircon and white mica can also be found in the outermost parts of the grains. It is notable that apatite, epidote and zircon grains are abundant as inclusions throughout the entire porphyroblasts (see Supporting Information) and that almost all of the matrix titanites have rutile cores or inclusions.
Detailed petrological investigations of garnet in different rocks from the Punta Balandra unit demonstrated that garnet records a significant portion of the prograde P-T path between 420 and 575 C and 1.2-2.0 GPa and hence documents the physical and chemical history of high pressure/low temperature (HP/LT) metamorphic processes (Escuder-Viruete & Pérez-Estaún, 2006). The sample in this work was chosen because of its large euhedral garnet grains that, judging from the major element compositional variations and from the inclusion assemblage, record large parts of the prograde history of the host rock, although the matrix mineral assemblage is reflecting a later, blueschist facies stage of the rock's metamorphic evolution.

| Analytical and computational procedures
All major and trace element mineral compositions were measured in a polished thin section that was prepared at the Institute of Earth and Environmental Science at Potsdam University in Germany. In order to choose a garnet crystal that is cut closest to the centre, the largest idiomorphic garnet in the thin section with the highest spessartine component in the core was chosen for further investigation.

| Electron probe microanalyses (EPMA)
Major element mapping and spot analyses of the chosen garnet were carried out with a JEOL JXA-8200 electron probe micro-analyser at the Institute of Earth and Environmental Sciences at Potsdam University. Instrumental conditions during map measurements were set at 15 kV accelerating voltage, 200 nA beam current, 2 μm beam size and a dwell time of 200 ms. For the EPMA spot analyses, the acceleration voltage was 15 kV with a beam current of 15 nA and a beam diameter of 2 μm. Detector counting times were 20 s on peaks and 10 s on background on both sides of the peaks.

| LA-ICP-MS
The two-dimensional (2D) element maps were carried out with the Teledyne Photon Machines Analyte Excite laser ablation system equipped with an aerosol rapid introduction system (ARIS) coupled to an Agilent 7900 Quadruple Inductively-Coupled-Plasma Mass Spectrometer (iCapQs ICP-MS) at Trinity College Dublin, Ireland. The acquisition was performed by defining a selected rectangular area composed of successive line scans (rasters) over the garnet crystal (cf. Ubide et al., 2015). Each raster employed overlapping 12 μm spot ablations with a 100 μm s À1 scan velocity and a 51 Hz repetition rate. During the mapping session, 20 major and trace elements were acquired with the dwell time set at 100 μs, whereas NIST610 standard glass (Jochum et al., 2011) was F I G U R E 1 Petrography of the investigated sample. (a) Thin section photomicrograph of the sample showing several large euhedral garnet grains within a matrix of sodic amphibole (Nam), quartz (Qtz), white mica (WM), chlorite (Chl), titanite (Ttn) and apatite (Ap). The investigated garnet (upper right corner) has a darker core with numerous titanite inclusions (not visible) and a brighter rim with large rutile (Rt) inclusions. The horizontal line within the grain is laser ablation pits. (b) Crossed polarizers photomicrograph of a sodic amphibole with inclusions of relict omphacite (whitish interference colours). (c) Backscattered electron detail image of relict omphacite inclusions in sodic amphibole. The omphacite inclusions indicate the sample eclogite facies metamorphic history, whereas the matrix minerals reflect the blueschist facies retrograde equilibration. measured at the start and end of every mapping sequence. In particular, Ca, Mn and REE + Y elemental maps were measured for garnet investigation and Sr, P and Ti maps to aid detection of mineral inclusions. The total analysis time was 7.5 h for a 3.5 Â 4.5 mm area of the garnet sample.

| Data reduction
Time-resolved laser ablation ICP-MS data were reduced by the Trace_Elements data reduction scheme (DRS) in the Iolite v2.5 data reduction software (Paton et al., 2011), which runs under IGOR Pro 6.3 (WaveMetrics, Inc.). The selected DRS setting employed semi-quantitative normalization relative to the reference material NIST 610 glass. Primary trace element 2D maps were generated using the Iolite CellSpace Image function. Later data processing, such as profile extraction from the trace element mapping as well as detailed selections from the original maps, was performed with MATLAB ® version R2018b. MATLAB ® scripts are available on request from the first author. It is important to mention that in order to avoid misinterpretation of the REE concentrations and REE distribution patterns in garnet by inclusion phases, such as titanite, apatite and epidote, we discarded all analyses from the mapping that contained more than 1000, 10 and 1 ppm of Ti, P and Sr, respectively. Wherever full maps are shown for clarity, this is indicated in the figure captions. Profiles that were extracted from the maps were low-pass frequency filtered to eliminate analytical scatter.

| Confocal μ-Raman spectroscopy
Laser-induced Raman measurements and mapping were carried out at room temperature using a Horiba (Jobin Yvon) LabRam HR Evolution at Gothenburg University. The samples were excited with an air-cooled frequencydoubled 532 nm Nd-YAG laser utilizing an Olympus 50Â objective (numerical aperture = 0.9). The lateral resolution of the unpolarized confocal laser beam was in the order of 1 μm. Spot spectra were generated in the range of 100-4000 cm À1 utilizing a 600 grooves/cm grating and a thermoelectric cooled electron multiplier CCD including a front-illuminated 1600 Â 200 pixel chip. The spectral resolution on polished sample was in the order of 1 cm À1 . The wavenumber calibration was done using the 520.7 cm À1 Raman band on a polished silicon wafer with a wavenumber accuracy usually better than 0.5 cm À1 .
Raman mapping was performed with a spatial resolution of 5 μm with a continuous laser beam and a motorized stage, the speed of which was adjusted to the integration time of the measurements. Data processing of the Raman mapping was performed with the Horiba LabRam ® software version 6. All spectra are background and slope corrected and despiked automatically. The height of the respective peaks was used to produce the maps.

| Secondary ion mass spectrometry
Concentrations of H 2 O and fluid mobile element (FME), namely, Li, B, F and Cl, as well as Mg, Ti, Y and Ce concentrations in garnet, were measured in gold-coated thin sections utilizing a Cameca 7f-GEO secondary ion mass spectrometer with Hyperion RF plasma ion source at the Edinburgh Ion Microprobe Facility (EIMF). A focussed 5.1 nA 16 O À primary beam with an impact energy of 18.5 keV was used to sputter an area of $20 Â 20 μm of the sample surface, generating secondary ions. Positive secondary ions with an energy of 50 ± 20 eV were extracted into a double-focussing mass spectrometer. Energy filtering reduces matrix effects as well as certain molecular interferences. The various isotopes were analysed by peak switching of the magnet and ion signals were counted with an electron multiplier. Prior to analysis, the area to be analysed was precleaned for 180 s by rastering the primary beam over an area of 35 Â 35 μm. This reduces particularly the 1 H background considerably (see below for details). Secondary ion beam centring and mass calibration peak adjustment were undertaken prior to each analysis in an automated routine. During analysis, an effective field aperture of 25 μm was used. The mass resolution used was 600 (M/ΔM) to maximize ion transmission. Count times were 2 s for 26 Mg, 30 Si and 47 Ti; 3 s for 1 H, 7 Li and 89 Y; 5 s for 11 B and 140 Ce; and 10s for 19 F and 35 Cl per cycle (magnet sweep), where each analysis consisted of six cycles.
Count rates were converted into concentrations by comparing Si-normalized count rates of samples with those of calibration standards with known SiO 2 and trace element contents. The following standards were used: a set of five garnet from Reynes et al. (2018) for H 2 O (H 2 O range 0-834 μg/g); ALV519-4-1 MORB glass for F, MgO, TiO 2 , Y and Ce (Jackson et al., 2015;Kumamoto et al., 2017); experimental basaltic glass St8.1.A9 from Lesne et al. (2011) for Cl; and the same glass standard employed for Li and B using concentrations determined independently at EIMF. The measured 1 H signal of the water-free garnet was used as a blank correction for H 2 O calibration. Samples were kept in an airlock at a high vacuum (2.5 Â 10 À8 mbar), whereas measurements took place in a sample chamber at 8.9 Â 10 À10 mbar.
Calibration slopes of H 2 O for garnet and basaltic glass differed by only 11%, which is very close considering much larger fractionations observed for other nominally anhydrous minerals, such as clinopyroxene and orthopyroxene (e.g., Gleeson et al., 2022). 2.1.6 | Thermodynamic-geochemical modelling The thermodynamic-geochemical modelling approach is based on an incremented Gibbs energy minimisation forward model of the phase assemblage along the sample's P-T trajectory and a mass-balanced REE distribution among the stable phases at each increment (e.g., Konrad-Schmolke et al., 2011;Konrad-Schmolke, Zack, et al., 2008). Modal proportions and composition of stable phases were calculated at 200 increments along a highpressure P-T trajectory that starts at 420 C/1.2 GPa and ends at 570 C/1.95 GPa. It was chosen to reflect the pressure and temperature interval for garnet growth (Escuder-Viruete & Pérez-Estaún, 2006).
Phase relations and compositions along the P-T trajectory were calculated assuming water fractionation and fractional garnet crystallization in the chemical system Mn, Na, Ca, K, Fe, Mg; Al, Si, O, H (MnNCKFMASH) utilizing Gibbs energy minimisation performed with the software package Theriak (de Capitani & Brown, 1987), the data extraction procedure published in Lanari et al. (2014) and the  database with recent updates. Solid solution models employed are feldspar (Baldwin et al., 2005), garnet , omphacite (Holland & Powell, 1996), chlorite , white mica (Keller et al., 2005), biotite (White et al., 2007), talc , amphibole (White et al., 2003) and epidote (ideal). A comparison of the results utilizing a different set of thermodynamic standard state data and solid solution models is provided in the Supporting Information. Water as a thermodynamic phase was considered pure H 2 O throughout all modelling. The chemical potential of oxygen was constrained by an ilmenite-magnetite-rutile (IMR) buffer . The chosen bulk rock composition used for the calculations is given in Table 1. For the REE modelling, we consider REEs as behaving thermodynamically passively, that is, without influence on the major element composition and stability of the phases, and the distribution of REEs among the stable phases can be described by bulk partition coefficients. We calculated the REE distribution among the modelled phases at every modelled P-T increment using mass balance constraints and internally consistent partition coefficients (Table 2). Bulk distribution coefficients were calculated at each increment based on the thermodynamically modelled phase assemblage. REEs incorporated in garnet and water were fractionated from the system according to the modelled modal amounts of these phases. The initial trace element concentrations of every element in the bulk rock are given in Table 3, and a more detailed description of the algorithm can be found in the Supporting Information.

| RESULTS
The euhedral garnet grain chosen for detailed analysis has a diameter of $9 mm and clearly exhibits two different growth zones (Figures 1a and 2a). A core region that is slightly darker than the rim comprises 65-70% of the area of the section through the grain. The core contains numerous inclusions of very fine-grained titanite together with lower abundances of apatite (ap), zircon (zrc) and quartz (qtz). The lighter coloured rim overgrowth has fewer inclusions, which consist of apatite, epidote, white mica and relatively large rutile (rt) grains that can be clearly seen as brown inclusion phases ( Figure 2a). The garnet rims show no or very limited signs of post-growth resorption.

| Compositional zoning within the garnet grain
The investigated iron-and calcium-rich garnet is characterized by smooth, continuous and concentric compositional major element transitions from core to rim (Figure 2b,c). The observed major element variations are typical for mafic high-pressure/low temperature (HP/LT) garnet (e.g., , such as increasing Mg and decreasing Ca and Mn contents from core to rim, as well as Fe contents that first increase towards the outer core and then decrease towards the rim. Mn shows a slight increase towards the outer rim, which is associated with a slight change in the slope of the Fe zoning. However, no abrupt or rhythmic changes nor deviations from concentric compositional patterns are visible in the major element zoning ( Figure 2b). In contrast to the major elements, the HREEs and Y show a much more complex compositional pattern in the LA-ICP-MS maps ( Figure 3a). The heaviest REEs (Ho to Lu) and Y have high concentrations in the centre of the crystal and show a continuous decrease across the core through $50% of the garnet radius. In the lighter HREEs (Tb, Dy) this central high is also developed but is less pronounced than in the heavier HREEs. Towards the rim, the HREEs and Y show a broad annular enrichment zone that runs parallel to the euhedral edges of the garnet crystal. Rimward of this annular high, the concentrations of all HREEs greatly decrease towards the outermost rim.
Superimposed on these broad trends are rhythmic fluctuations in the REE concentrations illustrated by several concentric and roughly equally spaced peaks and associated troughs that are visible in all the HREE maps in the outer core and the annular portion of the rim ( Figure 3a). These fluctuations are also visible in profiles ( Figure 3b) that are extracted along the transect arrow shown on the Yb map. Three peaks and four associated troughs are developed in the core of the garnet (peaks denoted by arrows #1, 2 and 3 on the Tm map), two peaks and one trough are within the broad enrichment zone (arrows # 4 and 5) and at least one more peak is visible within the rim of the garnet (arrow #6). These six peaks are also marked with grey bars in Figure 3b, which demonstrates that although relative peak heights differ among the REEs, they occur, with slight variations, at the same position for all elements. In the following section, we focus on these six peaks and their associated troughs.
In addition to the profiles, chondrite-normalized REE distribution patterns ( Figure 4) were extracted from the different growth zones (peaks and troughs) evident in the REE mapping. The slopes and shapes of these REE distribution patterns reflect relative changes in REE concentrations and thus yield information about subtle geochemical changes in the host rock during garnet growth. The extracted REE distribution patterns from the troughs are clearly distinct from those of the peaks. The arrows in each panel mark the relative change of each element between the two sets (peaks and troughs) of the REE distribution patterns. For example, from the centre of the garnet to the first trough the REE distribution patterns are HREE-enriched but straight (grey lines in Figure 4a), whereas within the first trough the REE distribution patterns show that the heaviest REEs (Ho to Lu) are preferentially depleted (black lines in Figure 4a), leading to a slight decrease in the slope and an increase in the curvature of the patterns. This systematic difference becomes T A B L E 2 REE distribution coefficients for the geochemical modelling   Sassi et al., 2000. obvious if all the fluctuations are considered. Notably, peak-to-trough changes in the REE distribution patterns (Figure 4a,c,e,g,i,k) are associated with a fractionation of heavier REEs with respect to lighter REEs, such that the heavier REEs show a more pronounced decrease in the troughs. In contrast, trough-to-peak changes in the REE patterns ( Figure 4b,d,f,hj,l) are characterized by a more or less equal increase in all REEs. These recurrent changes lead to an overall change in the shape of the REE distribution patterns from linear patterns with a steep positive slope in the core towards more curved, convex upward patterns in the broad REE enrichment zone and sinusoidal patterns in the rim. Because of progressive fractionation of heavier REEs into garnet during growth, the crest of the curved REE patterns continuously shifts towards the lighter REEs.

| Raman mapping
A μ-Raman map was undertaken in the same area as the REE maps prior to laser ablation analysis to visualize the distribution of the inclusion phases within the garnet. Figure 5 shows an RGB image of the Raman mapping that illustrates the outline of the titanite and rutile inclusions within the garnet grain. The RGB channels combine the 244-262 cm À1 (green) and 568-607 cm À1 (red) Raman bands, which were selected to constrain the different characteristic bands for titanite and rutile. The μ-Raman mapping shows that within the garnet core titanite is the only Ti-bearing inclusion phase, whereas in the rim rutile is present and titanite absent. Larger titanite grains outside the garnet core can only be found in the matrix and along cracks within the garnet rim and presumably reflect sample retrogression during exhumation. The Raman mapping also shows that there is a small area where rutile and titanite are both stable in a transition zone between the core and rim, which is illustrated on the Raman map by two boundaries that delimit the outermost titanites (white circles defining the white line) and innermost rutile inclusions (blue circles defining the blue line). These mineral stability fields are also plotted on top of the Lu and Ho maps ( Figure 3a) and in the extracted profiles ( Figure 3b) to demonstrate the spatial coincidence between the titanite-rutile transition as indicated by the change in the inclusion assemblage and the REE-enriched zone within the garnet crystal. The first (innermost) appearance of rutile inclusions correlates perfectly with the inner peak (peak #4) of the broad annulus, whereas the outermost occurrence of the titanite inclusions coincides well with the second peak (peak #5) of the enriched zone.
T A B L E 3 Initial REE concentration of the sample as used for the thermodynamic-geochemical model

| SIMS profiles
Concentrations of water and the FMEs Li, B, F and Cl were measured by secondary ion mass spectrometry (SIMS). Mg and Ti were measured to correlate major and trace element data as well as to monitor potential titanite inclusions, respectively. Y was monitored to correlate the FMEs with the REEs and the mapping shown in Figure 3. SIMS measurements were performed at carefully selected spots that were free of inclusions and F I G U R E 3 Rare earth element variations in the investigated garnet grain. (a) HREE + Y maps of the area indicated in Figure 2a. The HREE + Y maps show a complex compositional pattern with an area of high concentrations in the core, a broad enrichment zone at the core-rim transition and six smaller concentric peaks with associated troughs (see Tm map) more or less evenly distributed throughout the entire grain outside the inner core. The transition from titanite (ttn) and rutile (rt) inclusions is shown in the Lu and Ho maps and correlates well with the broad HREE + Y enrichment zone. (b) Core-rim HREE + Y profiles extracted from the maps in (a). All HREE + Y show a broad pattern with decreasing concentrations from the centre of the core to $30% radius followed by an increase peaking at 70% radius and a significant decrease towards the rim. Superimposed on this broad pattern are six evenly distributed minor peaks that appear in all HREE + Y profiles at the same position. Thin black lines are the measured data; bold red lines are low-pass frequency filtered data.
that represent a profile that stretches from the outer part of the central HREE high concentration zone to the annular enrichment zone in the inner garnet rim ( Figure 6). All of the FMEs as well as water show quantifiable concentrations, but only minor compositional variations that do neither correlate with the REE fluctuations (peaks correspond to grey areas in Figure 6), nor to the Y pattern, nor do they correlate with each other. Water, which scatters between 50 and 80 μg/g and F, which varies between 8 and 22 μg/g, show slightly decreasing trends from core to inner rim, whereas Li, apart from a spike about halfway along the profile and B, increase slightly from $1-2 μg/g and 30-70 μg/g, respectively. Cl is relatively constant between 15 and 25 μg/g throughout the profile. Mg, as already indicated in the EPMA measurements (Figure 2), is constantly increasing from core to inner rim, neither correlating with the FMEs nor the REE fluctuations. Ti is constant, apart from two spikes (arrows) that might be due to the influence of titanite inclusions. These Ti spikes also show small F spikes, a slight increase in H 2 O and significant Ce spikes (not shown), consistent with this interpretation.
F I G U R E 4 HREE distribution patterns of peaks and troughs extracted from the maps in Figure 1. From peaks to the associated troughs the convex upward curvature of the HREE distribution patterns is increasing due to preferential fractionation of the heaviest REE, whereas from troughs to the next peak the pattern curvature does not change significantly as all HREE are evenly replenished. See text for further explanation.
F I G U R E 5 Correlation of (RGB) μ-Raman map with the thin section photomicrograph. The μ-Raman map shows that titanite (ttn) inclusions (yellow) are found in the core, whereas rutile (rt) inclusions (red) are restricted to the rim of the garnet grain. Within a narrow zone, both Ti phases coexist. Titanite inclusions in the rim are always associated with cracks and are related to the exhumation of the sample as do the titanites in the matrix. Apatite inclusions (orange) can be found throughout the entire crystal. WM, white mica

| DISCUSSION AND INTERPRETATION
So far, our results demonstrate that the combination of high-resolution analytical methods, such as trace element and μ-Raman mapping, allows for visualization of subtle trace element and mineralogical variations in two dimensions. In the following sections, we demonstrate that in combination with thermodynamic-geochemical modelling, this integrated approach enables distinction of thermodynamic and kinetic contributions to the REE patterns in garnet. With this information, it is possible to trace in great detail geochemical and geodynamic processes by investigating the trace element record stored within garnet crystals. Although shapes and trends in trace element growth zonation in metamorphic garnet are highly variable in detail and are influenced by variations in their trace element bulk rock chemistry or by reactions along the reaction path, REE zoning patterns in particular commonly share some very characteristic features that are also reflected in our sample. It is notable, however, that the size of our sample was $10 Â 10 Â 15 cm, which is not only constraining the bulk rock composition for the thermodynamic and trace element calculations but also defines the maximum extent for which our calculations, for example, regarding rock wide equilibration or diffusion length scales, hold true.

| The inner part of garnet REE pattern
The inner part of the REE zoning profiles in our sample is characterized by a core enriched in HREEs followed by decreasing REE concentrations outwards (Figure 3). This trend is visible in all HREE but with decreasing accentuation towards the lighter REEs as the height of the central peak decreases with decreasing atomic number (Figure 3b). In our sample, this central high is observable for the elements Lu to Ho but almost invisible for the lighter REEs. Such REE patterns with high HREE concentrations in garnet cores together with a decreasing height of the central peak for the lighter REEs are common for garnet in metamorphic rocks (e.g., Skora et al., 2006) and can be simulated with a thermodynamic-geochemical model. Here, we utilized Gibbs energy minimization to model F I G U R E 6 Results of the SIMS measurements of water, fluid mobile elements (FMEs), magnesium and titanium. The positions of the measured spots are along a core-to-inner rim profile and are shown in the Tm map in the upper right panel. The position within the large Tm map is shown in Figure 3. The grey bars in the diagrams indicate the approximate positions of the REE peaks shown in Figure 3. Neither the concentrations of the FMEs nor those of water, MgO and TiO 2 are fluctuating such as the REEs indicated by the lack of correlation with the peak positions. The arrows in the TiO 2 diagram mark spots with a potential contribution of titanite inclusions. mineral modes and major element compositions along a given pressure-temperature trajectory together with a mass-balanced REE distribution among the modelled stable phases to calculate the resulting REE zoning in garnet (e.g., Konrad-Schmolke et al., 2011;Konrad-Schmolke, Zack, et al., 2008). The modelled coexisting mineral assemblage along the chosen P-T trajectory (Figure 7a) during prograde garnet growth is garnet, chlorite, sodic amphibole, epidote, white mica and quartz (Figure 7b), which is in good agreement with the evolution of natural samples from the same outcrop that have similar bulk rock compositions and similar garnet textures (Escuder-Viruete & Pérez-Estaún, 2006). Over the chosen pressure and temperature, interval garnet is growing between 465 and 570 C and pressures between 1.4 and 1.9 GPa by the consumption of chlorite and epidote (Figure 7b). The modelled shape of the major element zoning in garnet ( Figure 7c) resembles that of the natural sample supporting that the chosen P-T trajectory reflects the metamorphic evolution of the sample. Quantitative differences between the modelled and observed compositions are largest in the core, where modelled Fe and Ca concentrations differ by $15% from the observed values and modelled Mn as well as Mg differ by $5% from the natural sample. A parameter study of the thermodynamic modelling that shows the effect of different standard state datasets and solid solution models on the results can be found in the Supporting Information.
The forward-modelled core-to-rim HREE concentration profiles in garnet assuming fractional equilibrium garnet growth and mass-balanced equilibrium REE distribution among the stable phases are shown in Figure 8. The solid lines with and without circles are calculated with and without the titanite-out reaction, respectively. All the modelled profiles are characterized by an overall bell-shaped trend that is pronounced for the heavier REEs but visible in all profiles (note that only one half of the bell is shown). These bell-shaped overall REE patterns reflect the garnet's REE concentration if the crystal is in equilibrium with an unchanged mineral assemblage but with changing modal amounts of the reacting minerals, in this case epidote and chlorite (Figure 7b). The different peak heights reflect the bulk distribution coefficients between garnet and the matrix minerals, whereas the broad core-to-rim decreasing trend in the F I G U R E 7 Summary of the thermodynamic forward modelling results. (a) Calculated P-T pseudosection for the bulk rock composition of our sample (Table 1). The black arrow marks the P-T trajectory along which the rock evolution and garnet growth were modelled. The path follows the results of Escuder-Viruete and Pérez-Estaún (2006) (E + P 06). The grey circles show the proposed garnet growth stages from E + P 06. The grey rectangles mark the peak P-T conditions for the Punta Balandra rocks from E + P 06. The grey arrow shows the proposed retrograde P-T trajectory, and the hatched area marks the P-T region of the matrix assemblage equilibration in our sample. (b) Normalized mineral assemblage coexisting with garnet along the modelled P-T trajectory. Garnet grows from the decomposition of chlorite and epidote. The titanite-out reaction is superimposed with an Avrami-type reaction progress. New garnet and white mica plot on the right axis. (c) The modelled garnet composition reflects well the observed growth zonation in the natural sample.
HREE concentrations reflects element fractionation from the reacting bulk rock volume into the garnet crystal (e.g., Konrad-Schmolke, Zack, et al., 2008;Moore et al., 2013;Otamendi et al., 2002;Skora et al., 2006). Both of these characteristic features indicate that garnet at this stage incorporates REEs in fractional equilibrium (Moynihan & Pattison, 2013) with the matrix minerals (e.g., George & Gaidies, 2017). The continuous fractional equilibration of matrix phases with the growing garnet in our sample is exemplified by the contemporaneous change in the REE concentrations and REE distribution patterns of the titanite inclusions and the garnet crystal. F I G U R E 9 REE compositional correlation of host garnet and titanite inclusions. (a) Eu/Yb ratios in titanite inclusions extracted from the REE mapping. Eu/Yb in the inclusions is increasing from garnet centre to a crest marked by the dash-dotted line. Beyond that crest, the Eu/Yb ratios are decreasing until the titanite-out reaction zone. (b) Eu/Yb ratio in garnet parallels that in the titanite inclusions indicating a continuous equilibration of matrix minerals with the growing garnet. (c) REE distribution patterns in the titanite inclusions. The slope of the REE patterns is decreasing from the garnet core to the garnet mantle followed by an increase towards the reaction zone. Titanites in the reaction zone have also higher REE concentrations caused by the higher availability of REE during titanite breakdown. (d) Map of Eu= ffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi ffi Sm Â Gd p values extracted from the trace element mapping. The Eu= ffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi ffi Sm Â Gd p values are relatively constant $0.9 in the core (yellow) and increase in the garnet rim to values of $1.2 (red). Blue values correspond to pixels with LREE values below detection limit. Notable is the absence of fluctuations in the garnet core such as seen in the HREE mapping. See text for further discussion.
The Eu/Yb ratio in titanite, reflecting the slope of the HREE distribution patterns, is increasing from the inner to the outer core ( Figure 9a) until it reaches a local maximum marked by the dashed-dotted line in Figure 9a,b. The same trend is observed in the garnet crystal ( Figure 9b). The shallowing of the HREE distribution patterns in garnet, indicated by an increasing Eu/Yb ratio, is typical of HREE depletion in the bulk rock because of fractional garnet crystallization (e.g., Rubatto, 2002), and the similarity of the pattern in titanite (change from blue to orange patterns in Figure 9c) indicates that matrix phases, here exemplified by the titanite inclusions, also adjust to that HREE depletion in the reacting rock volume. From the local maximum marked by the dash-dotted lines towards the REE enrichment zone in garnet (stippled lines in Figure 9a,b), the Eu/Yb ratios in titanite and garnet are decreasing and the overall REE concentration in titanite and garnet are increasing (Figure 3 and change from orange to red patterns in Figure 9c) as a result of the enhanced REE liberation during the titanite-out reaction. These overall changes in the REE distribution patterns as well as in the concentrations of REEs in the co-existing titanite suggest a continuous approach to REE reequilibration of matrix phases during garnet growth.
As the mineral assemblage during garnet growth in our example is unchanged until the breakdown of titanite, the calculated REE growth zonation follows almost exactly a Rayleigh fractionation trend in the garnet core (bold stippled line in the Lu profile in Figure 8). This leads to a bell-shaped pattern with a broad central peak and an inflection point at $65% of the garnet radius. However, the central peak in our sample (Figure 3), as well as in many natural garnet (e.g., Corrie & Kohn, 2008;George et al., 2018;Shrestha et al., 2019;Skora et al., 2006), is narrower than expected for a Rayleigh fractionation scenario. This can be explained by the fact that in an unchanged reacting mineral assemblage, that is, chlorite + epidote + sodic amphibole + white mica + titanite, the width of the central peak reflects the rate with which the reacting bulk rock volume is depleted in garnet-compatible elements. The rate of bulk rock depletion is dependent on the rate of newly formed garnet and the element availability in the effective bulk rock volume (e.g., George & Gaidies, 2017). This relation is demonstrated by the exemplary models shown in the Lu plot of Figure 8. The thin solid lines display the modelled Lu zoning in garnet assuming that the effective bulk rock volume is only a fraction of the entire matrix mineral assemblage. The smaller the reacting bulk rock volume relative to the garnet growth rate, the steeper the slope of the HREE profile in the garnet core. Transportlimited REE supply from the reacting rock volume to the garnet surface, as discussed in the next paragraph, will have the same effect and leads also to narrow central peaks with steep concave upward flanks (cf. Skora et al., 2006). It has been shown, however, that garnet growth rates can vary during a rock's metamorphic evolution (e.g., George & Gaidies, 2017) and that deviations from fractional equilibrium major and trace element incorporation during garnet growth occur (e.g., Carlson, 2002), which implies that in order to quantify the extent of thermodynamic equilibrium in a rock volume a representative number of garnet porphyroblasts has to be investigated in detail (e.g., George & Gaidies, 2017;Gaidies et al., 2008).
In general, the widths and slopes of the bell-shaped HREE central peaks in natural garnet are constrained by the endmember scenarios of (1) fractional equilibrium crystallization, if the garnet grows in fractional equilibrium with the entire matrix; (2) partial equilibration with matrix phases if only parts of the matrix phases equilibrate with the growing garnet; and (3) transport-limited REE supply. Natural samples will likely reflect a combination of these scenarios, but shapes and widths of the central peaks can be used to qualitatively constrain the bulk element permeability of the interconnected transport matrix (ITM). Small to moderate amounts of garnet growing in fractional equilibrium with a large reacting rock volume, such as often observed in felsic rocks, will develop a zonation similar to a Rayleigh fractionation profile (e.g., Gieré et al., 2011;Moore et al., 2013;Otamendi et al., 2002). In contrast, large garnet porphyroblasts with limited REE supply from matrix phases, such as in our metamafic sample, will more rapidly deplete the reacting rock volume and develop a narrower and steeper central peak.

| The effect of titanite breakdown on the REE zonation
Towards the rim of the natural garnet crystal, the bellshaped REE trend is reversed and most REEs start to increase and form an annular enrichment zone. Such an enrichment zone following the bell-shaped pattern in the garnet core is also a common feature seen in REE profiles in metamorphic garnet from different rock types and different tectonic settings (e.g., Corrie & Kohn, 2008;George et al., 2018;Gieré et al., 2011;Rubatto et al., 2020;Shrestha et al., 2019;Skora et al., 2006;Usui et al., 2007). As indicated by the correlation between the μ-Raman and REE mapping (Figures 3 and 5), the reversals in the trends of all HREEs that form the broad overall enrichment zone in our sample can be linked to titanite breakdown and the growth of rutile in the sample. In contrast to rutile, which usually contains no or only negligible amounts of REEs (e.g., Meinhold, 2010), titanite is an REE-bearing phase (e.g., Garber et al., 2017) and the REEs liberated by the titanite breakdown are thus incorporated into garnet during this reaction. This clear coincidence is interpreted to indicate that the reaction where titanite breaks down to form rutile is associated with the broad zone of REE enrichment at the core-rim transition in the garnet crystal.
In order to model the effect of titanite breakdown on the REE pattern in garnet, titanite stability was superimposed on the thermodynamic calculations and its breakdown was simulated with an Avrami-type reaction progress between 541 C at 1.72 GPa and 556 C at 1.79 GPa (Figure 6a,b). Superimposed titanite breakdown was preferred over the thermodynamic simulation of Ti phases as the experimental data on the titanite-rutile transition are limited (Angiboust & Harlov, 2017;Manning & Bohlen, 1991). Further, the published data on the titanite-rutile transition suggest a strong dependence on the bulk rock composition (Angiboust & Harlov, 2017). Furthermore, the thermodynamic data on Ti incorporation in major phases, such as amphibole, pyroxene and white mica, are mostly unknown, and the incorporation of Ti into thermodynamic calculations requires knowledge of the oxygen fugacity during the rock's metamorphic evolution because of the interplay of iron oxides and ilmenite. Lastly, we believe that neither titanite nor rutile has a significant influence on the major phase assemblage.
The narrow REE hump in the modelled profiles between $70 and 90% of the garnet radius ( Figure 8) results from titanite breakdown. The effect of the titaniteout reaction on the REE uptake is best demonstrated by a comparison with the modelled curves calculated without titanite breakdown (bold solid lines in Figure 8). Controlled by the REE partition coefficients, the titanite-out reaction leads to an increase in the HREE (Dy to Lu) but not in REEs lighter than Dy. Titanite typically has flat or negatively sloped HREE patterns (Scibiorski et al., 2019;Figure 9) so that the difference in peak heights cannot be related to the REE composition of the titanites, but rather to the combined effect of REE liberation and garnet/bulk rock distribution coefficients.
Several studies have shown that other distinct mineral reactions, such as the breakdown of major phases (amphibole and epidote; Konrad-Schmolke, Zack, et al., 2008), of accessory phases (allanite, monazite, xenotime, ilmenite and titanite; Pyle & Spear, 1999;Gieré et al., 2011;Gaidies et al., 2021;this work) or the resorption of earlier grown garnet (e.g., Kulh anek et al., 2021) can also lead to such an annular increase in the REEs. Characteristic for this stage of equilibrium REE uptake in garnet is that the enrichment peaks for all REEs occur at the same position within the garnet crystal if the REE-liberating phase breakdown occurs over a small garnet growth interval (e.g., Moore et al., 2013). Nevertheless, the enrichment peaks can vary slightly in their positions for different REEs if the REE-liberating reaction is continuous and occurring over a larger garnet growth interval. In that case, the reactants might be able to equilibrate with respect to their REE content during the reaction progress, as demonstrated in Figure 7, and fractionation of HREE over LREE might occur because of the continuous equilibrium REE distribution between reactants and garnet. Hence, in that case, HREEs might be earlier incorporated than LREEs, which are continuously held back in the reactants (e.g., Konrad-Schmolke, Zack, et al., 2008). However, this effect must not be misinterpreted as diffusioncontrolled REE incorporation (cf. Skora et al., 2006), which will be discussed in the following paragraphs.
It should be noted, however, that it might be difficult to unambiguously determine the reactant minerals simply from the reaction-controlled REE enrichment in garnet, if the sample does not preserve such a detailed inclusion assemblage as in our case. This is because the REE concentrations as well as the REE distribution patterns in the reacting minerals control the shape and intensity of the resulting REE increase as well as the resulting REE distribution patterns in the coexisting garnet (cf. Konrad-Schmolke, Zack, et al., 2008). Hence, reactants with similar REE contents, such as apatite, epidote or lawsonite, might leave similar REE signatures in the co-existing garnet.

| The REE fluctuations-Recurring variation in rock permeability
Superimposed on the afore-described broad REE pattern in our sample are distinct recurring smaller peaks at fixed positions relative to the garnet radius that are developed to different extents but are visible in the profiles of all REEs. This is a feature that often occurs in the REE patterns of metamorphic garnet (e.g., Skora et al., 2006;Moore et al., 2013;George et al., 2018; this work) but is often undetected as it requires a high spatial resolution of the analytical spots or even REE mapping to be visualized properly. REE mapping from our sample shows that these recurring smaller peaks (Peaks 1-6) are indeed annular features (Figure 3a) that occur at fixed positions (Figure 3b), which rules out analytical artefacts to be the reason for the slight REE variations. These compositional fluctuations (peaks and troughs) superimposed on the broad REE zoning pattern in our sample offer the possibility to extract further detailed information about variations in the rock element transport properties and permeabilities.
Repetitive compositional variations, such as oscillatory zoning, in metamorphic garnet have been attributed to several geodynamic and kinetic parameters. They have been interpreted, for example, to be the result of rhythmic variations of pressure and temperature during garnet growth (Kohn, 2004;Stowell et al., 2011;Viete et al., 2018) but more likely to reflect garnet growth rate variations (e.g., George & Gaidies, 2017) or changes in the element supply such as by varying fluid fluxes (e.g., Clechenko & Valley, 2003;Dziggel et al., 2009;Jamtveit et al., 1993;Jamtveit & Hervig, 1994;Ranjbar et al., 2016;Zhai et al., 2014). In the following, we will discuss these aspects in detail.

| Pressure and/or temperature variations
Rhythmic variations in pressure and/or temperature are very unlikely reflected in the compositional variations of our garnet crystal. Garnet major element composition is very sensitive to pressure and temperature changes, which makes it a versatile mineral regarding geothermobarometry. Hence, slight variations, especially fluctuating variations in pressure and/or temperature, will result in rhythmic variations in the major element composition of garnet, such as demonstrated in Kohn (2004) and Viete et al. (2018). In our sample, major element patterns (Figure 2) display smooth and continuous changes from core to rim. We take the preservation of significant fluctuations in fast diffusing elements, such as Li (Figure 6), as evidence for the assumption that the smooth major element zonations are actually growth zonations and that extensive diffusional equilibration did not occur during the sample's exhumation. Hence, the absence of major element fluctuation that corresponds with the REE variations ( Figure 6) excludes pressure and/or temperature fluctuations as the cause of the REE fluctuations in our sample.

| Garnet growth rate variations
However, it is more difficult to decide whether garnet growth rate variations or changing fluid fluxes are the cause of the REE fluctuations. If element transport in the ITM is sluggish and REE incorporation at the garnet surface is potentially transport-controlled, growth rate variations will lead to varying REE depletion halos around the growing porphyroblasts (e.g., George et al., 2018). At times of enhanced garnet growth, the higher demand of REEs will, at constantly sluggish element fluxes in the ITM, lead to enhanced depletion of slowly diffusing elements in the ITM surrounding the growing garnet (cf. Carlson, 1989;Chernoff & Carlson, 1997;Meth & Carlson, 2005;Skora et al., 2006). This will in turn lead to decreasing REE concentrations at the garnet surface. At times of slow garnet growth, element transport might keep pace again with the demand at the garnet surface and REE concentration in garnet will increase. On the other hand, the same effect on REE variations in garnet will occur if, at constant garnet growth, the element supply is changing because of fluctuations in the ITM permeability, such as by varying fluid fluxes through the host rock. Hence, the two processes, which might even be interconnected, will have the same effect on the REE pattern in garnet.
Nevertheless, one striking argument against significant growth rate fluctuations during the formation of the garnet core in our sample is the homogeneous distribution of the numerous fine-grained inclusions. Matrix phases are included in porphyroblasts if the advancement of the growing surfaces is preferred over the energetic equilibration of the surface energies between matrix phases and the growing crystal. Therefore, slowly growing phases are more likely to contain fewer inclusions, whereas poikiloblastic grains are interpreted to be the result of rapid crystal growth. Consequently, growth rate variations in our sample should be reflected in a concentrically changing inclusion density pattern, as it is observed if the garnet core and the rim of our sample are compared. Within the core, however, the inclusions are very homogeneously distributed, pointing towards a potentially rapid, but constant advancement of the garnet surfaces throughout the growth of the garnet core.

| Changing fluid chemistry
Especially in settings with high fluid throughput, such as in hydrothermal systems like skarns or, as in the case of our sample, in subduction zones, where chemically contrasting rocks are tectonically juxtaposed, rapidly changing fluid compositions might be expected. Such variations in fluid chemistry might include the concentration of REEs (e.g., Dziggel et al., 2009) and other trace elements (e.g., Jamtveit et al., 1993), the oxygen fugacity (fO 2 ) as well as changes in salinity or pH (e.g., Zhai et al., 2014). Although fluid chemistry changes cannot be ruled out completely to be the reason for the REE fluctuations, some of the observations in our sample rather suggest relatively constant chemical conditions in the fluid. First of all, water concentrations are constantly between 50 and 80 μg/g, which suggests that an aqueous fluid is always present during garnet growth ( Figure 6). The concentrations of the FMEs in this fluid are relatively constant, and none of the minor fluctuations in FME concentration correlates with the observed REE variations (Figure 6). Chlorine is also relatively constant indicating only minor changes in fluid salinity. Hence, even if changes in fluid chemistry occurred and are represented in the minor fluctuations of the FMEs, they did not affect the REE distribution in the garnet. Also, f(O 2 ), a factor that could influence the REE distribution coefficients (e.g., Dalpé & Baker, 2000), was likely constant at least during the growth of the garnet interior. This is indicated by the Eu= ffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi ffi Sm Â Gd p value that is shown in the map in Figure 9d. Eu= ffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi ffi Sm Â Gd p is relatively constant $0.9 (orange colours) in the garnet core until the annular REE enrichment zone. Only the rim of the garnet shows significantly different values that reach up to 1.2 (red colours) indicating lower f(O 2 ) with progressive subduction.

| Changing element transport properties in the host rock
Evidence for changing element transport properties of the host rock during garnet growth results from the investigation of the REE distribution patterns in the garnet crystal. Garnet prefers HREEs over LREEs with a broadly linear relationship on a chondrite-normalized REE distribution diagram, and hence, typical REE distribution patterns in garnet are straight with a steep positive slope when plotted with increasing atomic number (e.g., the grey REE spectra from the garnet core in Figure 4a). Convex upward curved REE distribution patterns, such as in Figure 4l, indicate HREE depletion of the reacting rock as a result of element fractionation into the growing garnet (e.g., Otamendi et al., 2002). This effect is more pronounced with smaller reacting rock volumes; hence, the change in curvature of the REE distribution patterns is indicative of the size of the equilibration volume around the growing garnet.
The REE distribution patterns extracted from the troughs and peaks of the fluctuations differ systematically. The REE distribution patterns in all troughs show a preferred fractionation of heavier versus lighter REEs (Figure 4a,c,e,g,i,k), indicated by stronger concave downward curvature compared with the associated peaks. This difference is clearly indicative of the enhanced removal of HREEs from a smaller reacting bulk rock volume during trough growth compared with the larger equilibration volumes during peak formation. In turn, this can be interpreted as reflecting the presence of a compartmentalized ITM, which might be for example a noninterconnected fluid-filled pore space during the development of the troughs. Slow REE transport in a stagnating fluid column together with ongoing garnet growth and fractional crystallization leads to rapid REE depletion in the reduced ITM volume and the development of concentration gradients around the garnet crystal, which in turn is reflected in the strong fractionation of heavier REEs with respect to lighter REEs during trough development.
In contrast, the curvature of the REE distribution patterns does not change during the development of the fluctuation peaks (Figure 4b,d,f,h,j,l), which indicates a limited fractionation effect and hence an increasing equilibration volume during peak development. It is therefore likely that, at these growth stages, garnet is in equilibrium with a larger rock volume, suggesting that the ITM is interconnected and allows for effective trace element (REE) transport. Hence, during the growth of the 'peak' zones, all REEs are equally replenished, whereas during growth stages with decreasing REE concentrations (the troughs), the supply of heavier REEs could not keep pace with growth resulting in a more convex upward REE distribution pattern. These observations demonstrate that the development of the REE fluctuations superimposed on the broad overall REE pattern in the garnet is likely caused by changing rock permeabilities and therefore changing REE transport properties.

| A conceptual model of garnet growth
In order to support this conceptual model, we use a simplified numerical approach to model the qualitative effect of limited REE supply from a restricted rock volume on the garnet REE zonation. We simulate the development of a concentration gradient around the garnet crystal by superimposing diffusion-limited REE uptake at the garnet surface on the fractional equilibrium model. For this approach, we superimposed a Crank-Nicolson finite difference diffusion model on the incremented thermodynamic-geochemical calculations. For the diffusion simulation, the amount of step increments on the reaction path was reduced to 20 increments for the diffusion-limited model (every 10th increment of the original P-T trajectory). Each of the 20 increments was then again segmented into 100 equally spaced elements on which the Crank-Nicolson finite difference method was applied. As explained in the Supporting Information, we set the time increment to 200 years and the spatial increment to 0.025 mm. Volume diffusion within the garnet crystal was neglected. Intercrystalline REE diffusion was modelled assuming a system size of 10 mm and an Arrhenius-type function for the calculation of the diffusion coefficient D = D 0 exp ÀQ/RT . We used a value of D 0 = 2.8 Â 10 À5 cm 2 /a for the preexponential factor (cf. Skora et al., 2006) and an activation energy of 30 kJ/mol for intercrystalline Lu diffusion. Modelling of the diffusion-controlled REE uptake starts at each calculated increment with the REE concentration defined by the fractional equilibrium REE distribution calculations. Between two increments of the diffusion model, the flux of REE from the surrounding matrix towards the garnet surface is calculated assuming the above-mentioned diffusion parameters. After each calculated diffusion-limited growth interval, the REE concentration at the garnet surface returns to the value calculated for fractional equilibrium crystallization. A detailed description of the algorithm is given in the Supporting Information.
The modelled Lu concentration in garnet assuming fractional equilibrium garnet growth and mass-balanced equilibrium REE distribution among the stable phases is demonstrated by the bold black line in Figure 10a. This is the same profile as in Figure 8 and represents the garnet's Lu concentration for fractional equilibrium with the matrix minerals. A similar curve can also be drawn enveloping the Lu profile of the natural sample (solid black line in Figure 10b). As discussed above, this curve reflects the maximum Lu concentration the garnet can achieve in its respective growth environment. In the modelled case (Figure 10a), this is the fractional equilibrium with the entire rock volume as well as with the reacting titanite during the development of the annular REE increase. In the case of the natural sample (Figure 10b), where the central peak is narrower, hence the equilibration volume might be smaller, this curve can be seen as the maximum Lu concentration the garnet can achieve in a certain background REE flux.
However, in order to explain the REE fluctuations superimposed on this broad pattern, seven segments of the P-T trajectory are shown where periods of diffusionlimited intercrystalline REE transport towards the garnet surface are simulated with our simplified Crank-Nicolson approach (Figure 10a). Diffusion-limited REE supply leads to a decrease of the REE concentration relative to the fractional equilibrium concentration at the garnet surface (thin black line in Figure 10a). After each calculated interval of diffusion-limited REE supply, the REE concentration at the garnet surface returns to the value calculated for fractional equilibrium crystallization leading to the saw-tooth pattern shown by the thin black curves in Figure 10a. The red curve in Figure 10a displays a concentration profile resulting from the fluctuating interplay of fractional equilibrium REE incorporation and diffusion-limited REE uptake, thus reflecting fluctuating rock permeabilities during garnet growth. This modelled curve mimics the REE fluctuations measured in our sample given by the solid red line in Figure 10b.
The different effects of fractional equilibrium and kinetically controlled garnet growth in our sample can be summarized in the conceptual model shown in Figure 11, which involves three fundamental steps: 1. Garnet growth starts in equilibrium with a rock volume that has a background permeability allowing for a certain continuous element flux (Figure 11a, left panel, physical process). REE incorporation into garnet follows fractional equilibrium crystallization controlled by modal changes of matrix minerals in the host rock (Figure 11a, middle panel, chemical process). This background element flux is enabled by an interconnected and potentially fluid-filled pore space in the ITM. The REE concentration in garnet follows the bold black line in Figure 10a and in the right panel of Figure 11a. The equilibration among garnet and (at least fractions of) matrix phases in the natural sample is demonstrated by the simultaneous change in REE concentrations and REE distribution patterns in garnet and titanite ( Figure 6). A rock-wide elemental equilibration at this stage of garnet growth is likely facilitated by the availability of a free fluid phase as an element transport medium at the early stages of metamorphism. 2. The permeability of the rock is decreasing, for example, because of burial-or reaction-induced compaction and/or fluid drainage (e.g., Balashov & Yardley, 1998;Connolly, 2010), which leads to the closure of parts of the pore space and to compartmentalization of the reacting rock volume (Figure 11b, left panel). Decreasing rock permeability shuts off fluid percolation, hinders element equilibration and leads to a strongly transport-controlled REE incorporation into the growing garnet. Garnet growth at this stage causes a rapid REE decrease in the compartmentalized ITM as well as at the garnet surface (Figure 11b, middle panel), which leads to the development of the observed REE troughs (Figure 11b, right panel). The pronounced fractionation of heavier REEs over lighter REEs associated with the trough development ( Figure 4) further supports the arguments for this process. 3. Rock permeability increases until the pore space becomes interconnected again (Figure 11c, left panel), e.g., due to ongoing devolatilization reactions in the host rock, which elevates the pore fluid pressure (e.g., Connolly, 1997;Marti et al., 2021;Oliver & Bons, 2001;Taetz et al., 2018), or due to changing external stresses that might cause various forms of rock dilatancy (e.g., Sibson, 1994;Fusseis et al., 2009;Viete et al., 2018). The increasing permeability allows for increasing fluid percolation and faster element transport in the rock, which facilitates REE transport and replenishes the REE budget available for the growing garnet (Figure 11c, middle panel). REE concentration in garnet returns to that of the initial continuous background element flux (Figure 11c, right panel) without significant fractionation of heavier over lighter REEs as observed during peak development ( Figure 4).
Recurrence of these three process steps together with the titanite-to-rutile transition eventually leads to the development of the observed complex REE zonation patterns and REE distribution pattern in our sample (Figure 11d, right panel). Nevertheless, we point out that the diffusion models are used to demonstrate the feasibility of our conceptual model and that more data on REE diffusion in fluids as well as more sophisticated model approaches are needed to derive quantitative data from this approach. This, however, is beyond the scope of this paper.

| Application to other natural samples
The combination of our natural observations with the above thermodynamic-geochemical model together with the consideration of published examples demonstrate that the broad shape of REE zoning patterns often observed in metamorphic garnet, that is, a bell-shaped REE pattern with decreasing central peak height for decreasing atomic number in the garnet cores, followed by an annular REE enrichment zone, is the result of a change from fractional equilibrium REE uptake during garnet growth in an unchanged matrix mineral assemblage, to a reaction controlled REE uptake during the breakdown of major and/or accessory phases, such as demonstrated in Moore et al. (2013). Hence, such typical REE patterns can be used to infer the reaction path of the host rock (e.g., Moore et al., 2013). In contrast, the kinetic aspects of garnet formation are reflected in the widths and shapes of the central peaks as well as in recurring slight variations in the HREE concentrations that are superimposed on the general trend. These HREE fluctuations reflect the formation of depletion halos around the growing garnet and are therefore the expressions of recurring kinetic barriers that change the length scales of equlibration for the HREE in the ITM. We interpret the fluctuating occurrence of these peaks to be the result of changing interconnectivities of the ITM, hence indicating recurring changes in the rock's element transport properties. As these fluctuations are common in metamorphic garnet, especially in high-pressure samples, it is likely that fluctuations in the element transport properties in metamorphic rocks, especially in subduction zones, are also common. The general reason for that can be that pulselike fluid production and migration is a common process during metamorphism (e.g., Connolly & Podladchikov, 2013;Pollington & Baxter, 2011). 4.6 | Lateral peak shifts: Diffusion-or reaction-controlled?
Another often observed, but yet not fully explained, feature in REE patterns of metamorphic garnet is a significant lateral shift of the REE annuli that are resulting from a reaction-controlled REE uptake towards more outward positions with decreasing atomic number. This feature has first been described by Skora et al. (2006) and is interpreted to be the result of transport-controlled element supply and strongly different diffusion velocities for HREE and LREE. Konrad-Schmolke, Zack, et al. (2008) interpreted this peak shift to be caused by fractional incorporation of REEs into garnet during reactions that occur over a larger pressure and temperature interval, such that heavier REEs are primarily incorporated into garnet at the onset of the reaction, whereas lighter REEs are held back in the reacting phases as dictated by the distribution coefficients of reactants and garnet. This would be the case if reactants and garnet are in thermodynamic equilibrium and the reactant composition is constantly equilibrated, for example, because of recrystallisation, and continuously depleted in HREEs over LREEs. However, in our example, such a peak shift is exemplified by the observation that Peaks 4 and 5 in Figures 3 and 4, which result from the titanite breakdown reaction, are changing in their relative heights, which shifts the peak of the reaction-induced REE enrichment zone towards more outward positions for the lighter REEs. Furthermore, Peaks 5 and 6 in Figures 3 and 4 are also slightly shifted outward with decreasing atomic number of the REE. Regarding the formation of Peak 5, it is very unlikely that the peak shift is caused by different diffusive properties of HREEs and LREEs as the transport distance between reactants (titanite) and product (garnet) is negligible. On the other hand, a reactioninduced fractionation of HREE over LREE is not observed, as shown in Figure 9. An HREE depletion in the reactants and/or a preferential incorporation of LREE in the refractory reactants can therefore not be the case for the shift of Peak 5. Nevertheless, the remaining titanite still seems to be equilibrated with the garnet during the reaction as indicated by the homogeneous replenishment of REEs in titanites from the outer core towards the reaction zone (Figure 9b), a feature that is also observed in the garnet during the development of the REE peaks ( Figure 4). Hence, the observations in our sample cannot unambiguously distinguish between the different scenarios but rather yield arguments for both aspects. It has to be noted, however, that a peak shift in our sample is only observed in the outermost parts of the garnet, when garnet growth likely occurred at conditions with a limited fluid supply and reduced element transport capacities of the host rock, as with ongoing subduction the rock continuously dehydrates and fluid liberation gradually ceases. Furthermore, as indicated by the recurring REE fluctuations, periods of a transportlimited element supply are evident during garnet growth.
Data constraining the intercrystalline REE transport in metamorphic rocks do not exist, but in a few cases, differences in diffusion velocities between LREE and HREE are reported. LREE volume diffusion seems to be significantly slower than HREE volume diffusion in pyroxenes (Van Orman et al., 2001) and, more important for the interpretation of our data, LREE diffuse also slower than HREE during grain boundary diffusion in metal alloys (Higashino et al., 2019;Loewe et al., 2017). Hence, one could argue that the lateral shifts in REE annuli at the outermost parts of metamorphic garnet are developed in a setting with a very limited fluid supply where grain boundary diffusion is the dominant REE transport mechanism. A conceptionally sound model that clarifies the processes that lead to such lateral peak shifts must still be developed as this feature poses an integral question towards the full understanding of REE zonation patterns and their interpretation in terms of geodynamic and geochemical processes during garnet growth in metamorphic rocks.

| CONCLUSIONS
A detailed mineralogical and chemical study of a single garnet grain in a high-pressure metamorphic sample from the Samana peninsula (Dominical Republic) yields insight into the reaction path as well as into the chemical and physical changes during the rock's metamorphic evolution. Based on combined 2D REE and μ-Raman mapping together with thermodynamic-geochemical modelling, a detailed interpretation of the REE variations within the garnet crystal in terms of thermodynamic and kinetic contributions to the REE growth zonation is possible. The overall bell-shaped REE growth zonation observed in the investigated garnet grain reflects fractional garnet crystallization in an unchanged matrix phase assemblage. Differing heights of the central peaks for the different REEs are caused by fractionation of HREE over LREE into the growing garnet dictated by the garnet/bulk rock distribution coefficients. Peak widths and slopes of the peak flanks are controlled by the background REE supply from the matrix assemblage to the garnet surface, which is influenced by the fraction of matrix phases the garnet is in equilibrium with as well as by the overall transport properties of the interconnected transport matrix. A pronounced annular increase in HREE towards the rim of the garnet crystal reversing the overall bell-shaped compositional trend can be attributed to the titanite-out reaction that liberates significant amounts of HREE that cannot be incorporated into newly formed rutile during titanite breakdown. Superimposed on the bell-shaped growth zonation are also recurring REE fluctuations that can be visualized with the trace element mapping. These REE fluctuations might reflect recurring changes in the rock's element transport permeabilities associated with varying fluid fluxes through the host rock during garnet growth. Our results further demonstrate that 2D trace element and μ-Raman mapping yield a plethora of new information that is crucial for future interpretations of complex compositional growth zonation in metamorphic minerals. well as Fred Gaidies and Tetsuo Kawakami for their thorough work and helpful comments and suggestions that improved the manuscript. Katy Evans is thanked for her thorough editorial handling, including important comments that improved the quality of the text. SIMS analyses were supported by a NERC IMF pilot grant to Ralf Halama. David Chew acknowledges past and present support from Science Foundation Ireland (SFI) through research grants 12/IP/1663, 15/IA/3024, 13/RC/2092 and 13/RC/2092_P2 (iCRAG Research Centre). iCRAG is funded under the SFI Research Centres Programme. in calc-pelite and pelite, Gagnon terrane, western Labrador. Geological Materials Research, 4(1), eaaq0234. Zhai, D. G., Liu, J. J., Zhang, H. Y., Wang, J. P., Su, L., Yang, X. A., & Wu, S. H. (2014)

SUPPORTING INFORMATION
Additional supporting information can be found online in the Supporting Information section at the end of this article. Fig. A1: Crossed Nichols photomicrographs of garnet crystals and Raman spectra of epidote inclusions in the investigated sample. a) Few up to 500μm large epidote crystals can be found as inclusions from the inner core throughout the entire garnet crystals. b) Close up image of the epidote grain in the inner core of the garnet crystal. c) Raman spectrum of the epidote grain shown in a) and b) together with a reference epidote spectrum (RRUF database). d) and e) Several epidote inclusions distributed between the outer core to the inclusion-poor rim, demonstrating the occurrence of epidote inclusions throughout the entire inclusion-rich zone in the cores of the garnets in our sample. f) Raman spectra of the epidote grains shown in d) and e). g) Crossed Nichols photomicrograph of two large epidote inclusion in the outer rim of a garnet crystal from the sample. Fig. A2: Intensity mapping (black to green) of the $ 964 cm -1 (PO 4 ) 3stretching mode Raman band overlain over a photomicrograph of the investigated garnet. The wavenumber window was set between 950 and 970 cm -1 . The mapping shows a more or less homogeneous distribution of the apatite inclusions in the garnet crystal. There is no significant change in inclusion density in the area of the REE fluctuations, nor in the titanite-rutile transition zone. Fig. A3: Log f(O2) vs. temperature plot for various buffers that are commonly used to model or constrain the oxygen fugacity in rocks. See text for discussion.  Fig. A4: Modelled garnet composition and modal abundances along the P-T trajectory given in Fig. 7 of the main text utilizing the standard state and solution data from the reference set 2 (Table A1). Fig. A5: Modelled HREE profiles in garnet using the thermodynamic input data from reference set 2. Although the central high and the initial HREE decreasing in the garnet centre are more pronounced as a result of the higher pyroxene abundances in these models, the overall HREE profiles are very similar to those shown in Fig. 8 of the main text. Especially the HREE increase due to titanite breakdown is comparatively affecting the profiles as in the main text calculations. Fig. A6: Results of the thermodynamic modelling utilizing the thermodynamic dataset 1 (Table A1) with the garnet solid solution model of Berman (1990). a) Pseudosection modelled with the bulk rock composition given in table 1 of the main text. b) Modal abundances of the coexisting phases along the modelled P-T trajectory. c) Modelled garnet composition. Fig. A7: Modelled heavy rare earth element zoning patterns utilizing the thermodynamic dataset 1 (Table A1)