1. Introduction
The Cycladic islands in the central Aegean Sea are remnants of an Alpine orogenic belt that originally stretched from mainland Greece in the west to Turkey in the east. The geological record currently preserved in the Cyclades is discontinuous; the islands often represent different structural levels of the fragmented and extended Alpine orogen. The Eastern Cycladic island of Donoussa (Figure 1) and the adjacent island of Naxos, for example, represent strongly contrasting metamorphic environments and timing. Whereas the andalusite-sillimanite metamorphic sequences of Donoussa reflect Late Cretaceous Buchan-type high-temperature/low-pressure (high T/P) metamorphism, the Naxos metapelite sequence represents Miocene Barrovian kyanite-sillimanite metamorphism that overprinted the Eocene high-pressure metamorphism of the Cycladic Blueschist Unit (CBU).
(a) Geological Sketch Map of the Cycladic massif (modified after Katzir et al. Reference Katzir, Matthews, Garfunkel, Schliestedt and Avigad1996) showing the rock sequences of the Cycladic Blueschist unit (CBU), Miocene granitoids and Upper Tectonic Unit. 1. CBU Eocene high-pressure metamorphic rocks. 2. Eocene high-P metamorphic rocks overprinted by Oligocene/Miocene Barrovian low- to medium-pressure metamorphism. 3. Barrovian low- to medium-pressure Oligocene-Miocene metamorphic rocks. 4. Miocene granitoids. 5. Upper tectonic unit rocks (undifferentiated). Donoussa is separated from the Naxos Cycladic Blueschist Unit by the Naxos Paros Detachment System (NPDS). (b) Map of Greece showing the geotectonic zones and main outcrops of ophiolites. The Internal Hellenides, comprising the Pelagonian to Rhodope Zones, are characterized by a polymetamorphic-magmatic history during the Alpine orogeny, whereas the External Hellenides (which include the island of Crete) underwent continuous sedimentation until the Miocene (Papanikolaou, Reference Papanikolaou2021). The locations of the Asterousia Crystalline Complex (ACC) in the uppermost tectonic unit of southern Crete and Donoussa in the Cyclades (Do) are shown on the map.

The Alpine-Attic-Cycladic orogenic events of the CBU comprise blueschist- to eclogite-facies rocks that formed during the subduction of the Pindos Ocean floor and the subsequent collision of Apulian microplate fragments with the European plate. These rocks were then variably affected by Oligocene/Miocene greenschist–amphibolite facies metamorphism and Miocene igneous intrusion (Altherr et al. Reference Altherr, Schliestedt, Okrusch, Seidel, Kreuzer, Harre, Lenz, Wendt and Wagner1979; Reinecke et al. Reference Reinecke, Altherr, Hartung, Hatzipangiotou, Kreuzer, Harre, Klein, Keller, Geenen and Boeger1982; Matthews and Schliestedt, Reference Matthews and Schliestedt1984; Schliestedt et al. Reference Schliestedt, Altherr, Matthews and Helgeson1987; Okrusch and Bröcker, Reference Okrusch and Bröcker1990; Avigad et al. Reference Avigad, Matthews, Evans and Garfunkel1992; Bröcker et al. Reference Bröcker, Kreuzer, Matthews and Okrusch1993; Bröcker and Franz, Reference Bröcker and Franz1998; Putlitz et al. Reference Putlitz, Cosca and Schumacher2005; Ring et al. Reference Ring, Glodny, Will and Thomson2010; Dragovic et al. Reference Dragovic, Samanta, Baxter and Selverstone2012; Dragovic et al. Reference Dragovic, Baxter and Caddick2015; Roche et al. Reference Roche, Laurent, Giovanni, Cardello, Jolivet and Scaillet2016; Glodny and Ring, Reference Glodny and Ring2022; Lamont et al. Reference Lamont, Roberts, Searle, Gardiner, Gopon, Hsieh, Holdship and White2023 a, b). Exhumation of the Cycladic high P/T metamorphic rocks has taken place during accretion and orogeny followed by Miocene backarc extension and southward-directed slab roll-back that led to the formation of a Cordilleran-type metamorphic core complex of which Naxos forms part (Bonneau, Reference Bonneau1984; Lister et al. Reference Lister, Banga and Feenstra1984; Avigad and Garfunkel, Reference Avigad and Garfunkel1989, Reference Avigad and Garfunkel1991; Avigad et al. Reference Avigad, Garfunkel, Jolivet and Azañón1997; Trotet et al. Reference Trotet, Jolivet and Vidal2001; Brun and Facenna, Reference Brun and Facenna2008; Jolivet and Brun, Reference Jolivet and Brun2010; Ring et al. Reference Ring, Glodny, Will and Thomson2010; Jolivet et al. Reference Jolivet, Faccenna, Huet, Labrousse, Le Pourhiet, Lacombe, Lecomte, Burov, Denèle, Brun and Philippon2013, Reference Jolivet, Menant, Sternai, Rabillard, Arbaret, Augier, Laurent, Beaudoin, Grasemann, Huet and Labrousse2015, Reference Jolivet, Menant, Clerc, Sternai, Bellahsen, Leroy, Pik, Stab, Faccenna and Gorini2018; Lamont et al. Reference Lamont, Searle, Gopon, Roberts, Wade, Palin and Waters2020; Searle and Lamont, Reference Searle and Lamont2020).
In a review of worldwide blueschist belts, Maruyama et al. (Reference Maruyama, Liou and Terabayashi1996) noted that whereas most Cordilleran-type blueschist belts are associated with well-developed arcs and high T/P metamorphic rocks, collision-type blueschists frequently feature poorly developed calc-alkaline magmatic arcs. Such is the case of the Attic-Cycladic metamorphism wherein, in contrast to the wealth and depth of knowledge on the CBU, much less is understood about the metamorphic evolution in the upper levels of the subduction zone. The dismembered Upper Tectonic Unit of the Cyclades is interpreted to represent structurally higher levels located in the upper crust during CBU subduction (e.g., Bonneau, Reference Bonneau1984; Altherr et al. Reference Altherr, Kreuzer, Lenz, Wendt, Harre and Dürr1994). The tectonic contact between the Upper Tectonic Unit and the underlying CBU unit is everywhere characterized by low-angle extensional faults that omitted a thick part of the crust between the two units (Avigad and Garfunkel, Reference Avigad and Garfunkel1989, Reference Avigad and Garfunkel1991; Jolivet et al. Reference Jolivet and Brun2010, Reference Jolivet, Menant, Sternai, Rabillard, Arbaret, Augier, Laurent, Beaudoin, Grasemann, Huet and Labrousse2015). The Upper Tectonic Unit on the eastern Cycladic islands of Donoussa, Nikouria and Anafi (Figure 1) features Late Cretaceous metamorphic rocks that coexist with calc-alkaline igneous intrusives (Reinecke et al. Reference Reinecke, Altherr, Hartung, Hatzipangiotou, Kreuzer, Harre, Klein, Keller, Geenen and Boeger1982; Altherr et al. Reference Altherr, Kreuzer, Lenz, Wendt, Harre and Dürr1994; Beeri-Shlevin et al. Reference Be’eri-Shlevin, Avigad and Matthews2009; Martha et al. Reference Martha, Dörr, Gerdes, Petschick, Schastok, Xypolias and Zulauf2016; Koutsovitis et al. Reference Koutsovitis, Soukis, Voudouris, Lozios, Ntaflos, Stouraiti and Koukouzas2022; Martha et al. Reference Martha, Xypolias, Cheng, Dörr, Gerdes, Hezel, Kutzschbach, Millonig, Schmeling, Marschall, Müller and Zulauf2025).
The Late Cretaceous metamorphism is widely identified in Upper Tectonic Unit sequences throughout the Aegean Sea area. It is found on the central and northern Cyclades islands of Syros, Tinos and Ikaria; the eastern island of Samos (e.g., Altherr et al. Reference Altherr, Kreuzer, Lenz, Wendt, Harre and Dürr1994; Katzir et al. Reference Katzir, Matthews, Garfunkel, Schliestedt and Avigad1996; Patzak et al. Reference Patzak, Okrusch and Kreuzer1994; Zeffren et al. Reference Zeffren, Avigad, Heimann and Gvirtzman2005; Pe-Piper and Photiades, Reference Pe-Piper and Photiades2006; Ring et al. Reference Ring, Okrusch and Will2007; Souklis and Stockli, Reference Soukis and Stockli2013; Mavrogonatos et al, Reference Mavrogonatos, Magganas, Kati, Bröcker and Voudouris2021; Broecker and Berndt, Reference Bröcker and Berndt2024) and has also been recognized on Skyros Island in the northern Aegean (Boundi et al. Reference Boundi, Papanikolaou, Bosio and Montemagni2024). South of the Cyclades, in the uppermost tectonic unit of Crete, the Asteroussia mountains and other outcrops form a similarly aged belt of Late Cretaceous high T/P rocks and calc-alkaline intrusive granitoids generically referred to as the Asteroussia Crystalline Complex (Bonneau, Reference Bonneau1972; Seidel et al. Reference Seidel, Seidel, Okrusch, Kreuzer, Raschka and Harre1976; Langosch et al. Reference Langosch, Seidel, Stosch and Okrusch2000; Martha et al. Reference Martha, Dörr, Gerdes, Krahl, Linckens and Zulauf2017, Reference Martha, Zulauf, Doerr, Binck, Nowara and Xypolias2019; Zulauf et al. Reference Zulauf, Dorr, Albert, Martha and Xypolias2024). The different occurrences may not be directly comparable, as there are indications of a different geological and P–T evolution (e.g., lack of andalusite and sillimanite in metapelites, no Cretaceous intrusive rocks or dikes and an ongoing discussion about the presence or absence of metamorphic soles). The focus of this work is a detailed petrological analysis of the P–T path of the Late Cretaceous Cycladic Upper Tectonic Unit metamorphism on Donoussa Island. The study develops upon previous petrological, microstructural and dating studies (Altherr et al. Reference Altherr, Kreuzer, Lenz, Wendt, Harre and Dürr1994; Martha et al. Reference Martha, Xypolias, Cheng, Dörr, Gerdes, Hezel, Kutzschbach, Millonig, Schmeling, Marschall, Müller and Zulauf2025).
Garnet is a key mineral in quantifying metamorphic PT conditions in metapelites, both through its divalent cation chemistry and frequent preservation of growth and retrograde zoning. This study combines classical cation exchange thermobarometry of metapelites and calcic amphibolites with phase diagram P–T pseudosections (equilibrium assemblage diagrams for a fixed bulk composition) calculated in the MnO-Na2O-CaO-K2O-FeO-MgO-Al2O3-SiO2-TiO2-O2 (MnNCKFMASHTO) and NCKFMASHTOCr systems, respectively, using up-to-date THERMOCALC software (Holland and Powell, Reference Holland and Powell2011; White et al. Reference White, Powell, Holland, Johnson and Green2014a; Holland and Powell, Reference Holland, Green and Powell2018; Green et al. Reference Green, Holland, Powell, Weller and Riel2025). This combined approach allows us to determine the P–T stability field of garnet + Al-silicate + biotite metapelite assemblages up to their initial melting, particularly distinguishing the critical role that increased bulk-rock MnO plays in stabilizing garnet at lower P–T than would otherwise occur. The recent petrological analysis of metapelite rocks by Pattison and Forshaw (Reference Pattison and Forshaw2025) has pointed out the importance of understanding the thermodynamics of the spessartine component of garnet. The metapelite studies are complemented by measurements of the T conditions and phase diagram of the cogenetic amphibolite assemblages. Additionally, we broaden the scope of the thermometry to include refractory accessory mineral (RAM) oxygen isotope thermometry (Valley et al. Reference Valley2001) of boudinaged quartz–andalusite–sillimanite lenses and veins within the metapelite schists. RAM thermometry has successfully been applied to Barrovian syn-metamorphic quartz–Al-silicate veins on Naxos (Putlitz et al. Reference Putlitz, Valley, Matthews and Katzir2002).
2. Geological and mineralogical setting
Donoussa is unique in the Cyclades in being one of only two islands that are entirely assembled from Late Cretaceous high T/P metamorphic metalpelite and amphibolite rocks and marbles, with contemporary intrusive peraluminous granitoids (Altherr et al. Reference Altherr, Kreuzer, Lenz, Wendt, Harre and Dürr1994; Figure 2a). Nikouria, a small Upper Tectonic Unit island separated from Amorgos by a high-angle normal fault, features a similar makeup and timing (Durr et al. Reference Dürr, Seidel, Kreuzer and Harre1978 a, b; Dürr, Reference Dürr1985; Rosenbaum et al. Reference Rosenbaum, Ring and Kühn2007). The eastern Cycladic Island of Anafi, by contrast (Figure 1), features a far more complex Upper Tectonic Unit sequence of meta-ophiolites, amphibolite and calc-silicates, garnet-biotite paragneiss, quartzites and marbles, ophiolites and igneous granitoids and diorites (Reinecke et al. Reference Reinecke, Altherr, Hartung, Hatzipangiotou, Kreuzer, Harre, Klein, Keller, Geenen and Boeger1982; Be’eri-Shlevin et al. Reference Be’eri-Shlevin, Avigad and Matthews2009; Martha et al. Reference Martha, Dörr, Gerdes, Petschick, Schastok, Xypolias and Zulauf2016; Koutsovitis et al. Reference Koutsovitis, Soukis, Voudouris, Lozios, Ntaflos, Stouraiti and Koukouzas2022). Ophiolites are a common feature of the Hellenides (Figure1b), where they provide evidence of past oceanic basins caught up in orogenic subduction processes. The Cycladic Upper Tectonic Unit sequences feature numerous ophiolite tectonic slices (Katzir et al. Reference Katzir, Garfunkel, Avigad, Matthews, Tymaz, Yilmaz and Dilek2007, Reference Katzir, Matthews, Garfunkel, Schliestedt and Avigad1996).
Photos of field exposures and petrographic thin sections on Donoussa. (a) General oblique view of Donoussa island from the south showing the approximate boundaries of the metamorphic complex (Metam. Complex) consisting of interlayered pelite gneisses, amphibolites, calc-silicates and intrusive granitoids enclosed from above and below between thick marble (Mbl.) horizons. Q = Pleistocene deposits as defined by Dürr, Reference Dürr1985. The Donoussa island topographic image is generated using Google 3D software. (b) Amphibolite exposure featuring thin alternating laminae of C(I) calcic-hornblende amphibolite rock and C(II) Mg-Fe orthoamphibolite rock. (c) Quartz-andalusite nodule enclosed within metapelite host rock. The quartz + andalusite ± sillimanite lenses, nodules and boudins are rotated into alignment with and enclosed by the S2 foliation of the enclosing metapelite host rock. Individual crystals of rose-coloured andalusite 0.5∼1cm in size can be seen in the bottom right of the nodule. (d) Mesoscale recumbent D2 fold in metapelite. The fold axis (yellow double-ended arrow) is NNE directed (20/10°) and is part of a major antiformal structure involving both marble and schists. (e) Cross-polarized petrographic image of metapelite rock showing two large andalusite crystal grains enclosed by the quartz, muscovite and biotite matrix which define the S2 schistosity. The rotated andalusite grain on the right of the photo appears to have grown prior to the D2 deformation and therefore is pre- or syn-kinematic. Biotite inclusions within this andalusite are orientated at an angle to the S2 schistosity and exhibit a preferred orientation. The preferred orientation of the inclusions (S1) is evidence for an earlier deformation event (D1) that occurred pre- or syn-andalusite growth. (f) Plane-polarized image showing orientated biotite crystals that define the S1 schistosity. The andalusite encloses biotite and quartz crystals but is wrapped around by the S2 schistosity, indicating that the development of the main metamorphic D2 mineral assemblage postdated the andalusite growth. (g) Aggregate of fibrous sillimanite grains defining the S2 lineation in a sillimanite-biotite-quartz rock. Aggregates of small quartz and biotite grains also define the S2 orientation. (h) Scanning electron microscope image of a calc-silicate thin section. The actinolite grain (act) clearly overgrows the main calc-silicate assemblage comprising clinopyroxene (cpx), garnet (g), calcite (cc) and epidote (ep).

Several studies have proposed a correlation of the high T/P metamorphism with a Late Cretaceous volcanic arc of the Pelagonian unit on mainland Greece or with igneous arcs in northern Turkey (Beer-Shlevin et al. Reference Be’eri-Shlevin, Avigad and Matthews2009; Martha et al. Reference Martha, Dörr, Gerdes, Petschick, Schastok, Xypolias and Zulauf2016: Koutsovitis et al. Reference Koutsovitis, Soukis, Voudouris, Lozios, Ntaflos, Stouraiti and Koukouzas2022; Boundi et al. Reference Boundi, Papanikolaou, Bosio and Montemagni2024). On Anafi, Koutsovatis et al. (Reference Koutsovitis, Soukis, Voudouris, Lozios, Ntaflos, Stouraiti and Koukouzas2022) proposed that Late Cretaceous granitoids intruded a high T/P metamorphic sequence in a subducting slab beneath the southern external margin of the Pelagonian microcontinent. The granitoid magma rise to shallow crustal levels was facilitated by upward-directed ‘corner flow’ intrusion mechanisms. However, the observation of Late Permian and Late Cretaceous arc-type granitoids on Crete led Zulauf et al. (Reference Zulauf, Dorr, Albert, Martha and Xypolias2024) to infer that the Asteroussia Crystalline Complex was derived from the Late Permian/Late Cretaceous igneous arc located further north in the Strandja-Rhodope region of northern Greece and Bulgaria. This volcanic arc source was adopted for Donoussa by Martha et al. (Reference Martha, Xypolias, Cheng, Dörr, Gerdes, Hezel, Kutzschbach, Millonig, Schmeling, Marschall, Müller and Zulauf2025), who propose rapid (<2.6 Ma) quartzite forearc deposition, high T/P metamorphism and post-peak-temperature igneous dyke intrusion during rapidly retreating southward-directed subduction slab rollback.
The 1:50,000 scale geological mapping of Donoussa was done by Dürr (Reference Dürr1985). The metamorphic succession is dominated by thick marble sequences intercalated with alternating series of amphibolite-facies metapelite schists and gneisses, amphibolites, marbles and calc-silicates (Figures 2a and 3). The presence of volcanic and ultramafic detritus and detrital zircon U–Pb ages led Martha et al. (Reference Martha, Xypolias, Cheng, Dörr, Gerdes, Hezel, Kutzschbach, Millonig, Schmeling, Marschall, Müller and Zulauf2025) to propose that the protoliths of the metasedimentary rocks were deposited in a Campanian forearc basin.
Simplified geological map of Donoussa showing main geological units (modified after Dürr, Reference Dürr1985). The locations of the main sampling sites for the EPMA and geochemical studies are indicated in the figure, along with the mineralogical rock type sampled at each site. Location coordinates are given in Tables S1 and S2. The intrusive granitoid dykes and stocks are mainly restricted to the northwest and southeast of the island. They comprise mildly foliated muscovite aplites; two mica granites; and rare dioritic dykes (Altherr et al. Reference Altherr, Kreuzer, Lenz, Wendt, Harre and Dürr1994). The mineral abbreviations in the text and figures are gt (g in phase diagrams) = garnet and = andalusite, sill = sillimanite, pl = plagioclase, ksp = k-feldspar, bi = biotite, mu = muscovite, q = quartz, tour = tourmaline, ilm = ilmenite, st = staurolite, chl = chlorite, hbl = Mg-hornblende, tsch = tschermakite, ep = epidote, cpx = clinopyroxene, act = actinolite, cc = calcite, anth = anthophyllite, ged = gedrite, sph = sphene.

Rocks and Mineral Assemblages. The principal metamorphic rocks and mineral assemblages defined in this study are:
-
A. Marbles. The layered marbles (both massive and interlayered with metapelites, amphibolites and calc-silicates) consist of medium- to coarse-grained calcite crystals with thin microscopic-scale quartz laminae.
-
B. Metapelites. Critical metapelite mineral assemblages on Donoussa recognized here and by Altherr et al. (Reference Altherr, Kreuzer, Lenz, Wendt, Harre and Dürr1994) feature two main parageneses (Figure 3):
-
B(I): Garnet + sillimanite + biotite + plagioclase + quartz + muscovite + ilmenite + tourmaline.
-
B(II): Sillimanite ± andalusite + muscovite + biotite + quartz + tourmaline ± ilmenite.
-
The less common garnet-bearing B(I) assemblages mostly occur in thin cm-scale bands interlayered on an outcrop scale within the mineralogically dominant B(II) assemblages. Separate andalusite and sillimanite metamorphic zones could not be identified in the field, though staurolite was found as relict inclusions within andalusite.
-
C. Amphibolites. Amphibolites outcrop as both massive and schistose Ca-hornblende types, often interbedded with paler calc-silicate layers. An additional type of banded amphibolite rock consists of thin layers of Ca-amphibolite interlayered with Mg-Fe amphibolite (anthophyllite-gedrite, Figure 2b). Amphibolite mineralogy is thus characterized as two main types (Figure 3).
-
C(I): Hornblende-tschermakite ± plagioclase (bytownite/labardorite) + biotite + epidote + quartz + clinopyroxene.
-
C(II): Anthophyllite-gedrite + plagioclase (oligoclase) + biotite + garnet + quartz.
-
The C(I) assemblages commonly feature sphene (titanite) as an accessory mineral (Altherr et al. Reference Altherr, Kreuzer, Lenz, Wendt, Harre and Dürr1994).
-
D. Calc-silicates. Calc-silicate assemblages consist of calc-silicate layers at the contacts between marbles and amphibolite or interbedded within amphibolite. The mineral assemblages include:
-
D(I): Garnet + clinopyroxene + hornblende ± actinolite ± plagioclase ± epidote + quartz + calcite.
-
The metamorphic series are locally intruded by small intrusions of peraluminous granites and diorite dikes (Altherr et al. Reference Altherr, Kreuzer, Lenz, Wendt, Harre and Dürr1994). Compositionally, the granites have a peraluminous calc-alkaline chemistry strongly matching that of Anafi Island granitoids (data of Altherr et al. Reference Altherr, Kreuzer, Lenz, Wendt, Harre and Dürr1994 replotted by Koutsovitis et al. Reference Koutsovitis, Soukis, Voudouris, Lozios, Ntaflos, Stouraiti and Koukouzas2022), with mineralogy comprising coarse-grained equigranular quartz, albite and microcline, together with biotite ± muscovite ± garnet ± tourmaline.
Abbreviations used in the text to describe mineral assemblages (based on Holland and Powell, Reference Holland and Powell1998) are listed in the caption to Figure 3.
Structures and petrography. Kinematic indicators on the deformation history on Donoussa are limited in number, as is generally true of the Late Cretaceous Upper Tectonic Unit on Anafi (Reinecke et al. Reference Reinecke, Altherr, Hartung, Hatzipangiotou, Kreuzer, Harre, Klein, Keller, Geenen and Boeger1982). The detailed studies of Martha et al. (Reference Martha, Xypolias, Cheng, Dörr, Gerdes, Hezel, Kutzschbach, Millonig, Schmeling, Marschall, Müller and Zulauf2025) identify three sequential deformation-metamorphic transitions during the high T/P event: D1 (staurolite grade: ∼550°C, 4 kbar); D2 (sillimanite grade: 660°C, 3 kbar); D3 retrograde greenschist-facies metamorphism (∼540°C, 1.5 kbar). The main penetrative regional S2 schistosity is associated with the D2 phase and is well-defined in high-grade B(II) sill + bi + mu rocks (Figure 2e, f, g). Fold axes in schistose rocks trend NNE (Figure 2d). Lineation directions, mainly evident in the fibrolite sillimanite in metapelites (Figure 2g), are also N to NNE directed. The D2 structural elements are compatible with NE-directed shearing (Martha et al. Reference Martha, Xypolias, Cheng, Dörr, Gerdes, Hezel, Kutzschbach, Millonig, Schmeling, Marschall, Müller and Zulauf2025). The third, non-penetrative retrograde D3 greenschist-facies overprint of the S2 foliation recognized by these authors is possibly represented here by the boudinage texture of sillimanite in Figure 2g.
Andalusite porphyroblasts frequently occur in metapelites but are overgrown by the peak metamorphic B(II) sillimanite assemblage (Figure 2e, f). The andalusite is chronologically the earlier formed Al-silicate, but an andalusite to sillimanite isograd could not be identified in the field. Staurolite inclusions were identified within the andalusite porphyroblasts by both Altherr et al. (Reference Altherr, Kreuzer, Lenz, Wendt, Harre and Dürr1994) and Martha et al. (Reference Martha, Xypolias, Cheng, Dörr, Gerdes, Hezel, Kutzschbach, Millonig, Schmeling, Marschall, Müller and Zulauf2025). The andalusite porphyroblasts contain numerous inclusions of earlier quartz and biotite that exhibit a linear foliation texture partially rotated into the main S2 schistosity (Figure 2e), indicating that andalusite porphyroblasts were rotated together with their foliated biotite inclusions during the D2 deformation phase. This observation matches the aligned biotite-quartz inclusions found in andalusite and garnet [Martha et al. (Reference Martha, Xypolias, Cheng, Dörr, Gerdes, Hezel, Kutzschbach, Millonig, Schmeling, Marschall, Müller and Zulauf2025)] and supports their proposal that the rocks were affected by the earlier tectono-metamorphic event D1 featuring generating an S1 schistosity characterized by staurolite-biotite-quartz assemblages. Andalusite porphyroblast growth continued into the D2 phase, where it preceded the formation of sillimanite. Staurolite was not found in the B(I) and B(II) mineral assemblages associated with the D2 tectonometamorphic phase.
The metapelite exposures also feature boudins, lenses and veins composed of coarse-grained quartz, muscovite andalusite with fibrous sillimanite overgrowths (Figure 2c). Rotation of these boudins and lenses into parallelism with the S2 schistosity and overgrowth and rotation of andalusite porphyroblasts in metapelites indicates that they grew prior to or during the major NE-directed shear movement along the S2 plane. Three quartz-andalusite samples are used in our oxygen isotope geothermometric study.
Geochronology. Mineral K–Ar and Rb–Sr cooling ages determined on hornblende, biotite and muscovite from amphibolites and pelitic gneisses on Donoussa fall in the Late Cretaceous at 65–59 Ma (Altherr et al. Reference Altherr, Kreuzer, Lenz, Wendt, Harre and Dürr1994) and are similar to closure ages on Anafi (∼70 Ma) and southern and eastern Crete (76–66) Ma (Seidel et al. Reference Seidel, Seidel, Okrusch, Kreuzer, Raschka and Harre1976; Reinecke et al. Reference Reinecke, Altherr, Hartung, Hatzipangiotou, Kreuzer, Harre, Klein, Keller, Geenen and Boeger1982; Martha et al. Reference Martha, Dörr, Gerdes, Petschick, Schastok, Xypolias and Zulauf2016, Reference Martha, Dörr, Gerdes, Krahl, Linckens and Zulauf2017, Reference Martha, Zulauf, Doerr, Binck, Nowara and Xypolias2019). The close correspondence between hornblende and biotite ages on Crete led Seidel et al. (Reference Seidel, Seidel, Okrusch, Kreuzer, Raschka and Harre1976) and Martha et al. (Reference Martha, Dörr, Gerdes, Krahl, Linckens and Zulauf2017) to propose rapid cooling from peak temperatures.
Most intrusive granitoids on Donoussa gave the same range of cooling mineral ages (64–59 Ma) as the metamorphic rocks, suggesting a coeval decompression history. A whole rock Rb–Sr isochron age of 81.9 ± 2.2 Ma for garnet two-mica granitoids was taken to reflect a slowly evolving ∼20 Ma P–T cycle (Altherr et al. Reference Altherr, Kreuzer, Lenz, Wendt, Harre and Dürr1994). In direct contrast to this extended timescale, U–Pb age determinations on zircons in quartzite were interpreted to reflect that sedimentary deposition, high T/P metamorphism and igneous intrusion occurred within a short 2.6 Ma timespan beginning with forearc sedimentary deposition at 75.8 Ma and terminating with 73.2 Ma monzogranite dike intrusion (Martha et al. Reference Martha, Xypolias, Cheng, Dörr, Gerdes, Hezel, Kutzschbach, Millonig, Schmeling, Marschall, Müller and Zulauf2025). This was followed by cooling below 600°C with sphene (titanite) formation at 71.3 Ma, D3 folding and greenschist-facies metamorphism and a younger generation of granitoid dyke intrusion at 64-60 Ma. Following on the earlier work of Zulauf et al. (Reference Zulauf, Linckens, Beranoaguirre, Gerdes, Krahl, Marschall, Millonig, Neuwirth, Petschick and Xypolias2023, Reference Zulauf, Dorr, Albert, Martha and Xypolias2024), the rapid sequence of events was interpreted as reflecting fast trenchward forearc retreat of a late Campanian subduction accretion system.
3. Methods
Mineral and whole rock chemistry. The main locations at which samples were taken for the EPMA and whole rock geochemical study are plotted in Figure 3. The metamorphic mineral assemblages were determined by optical microscopy in thin section and by back-scattered electron imaging (BSE) on a scanning electron microscope. Mineral chemistries for P–T determination using cation exchange thermobarometry were determined using a Jeol JXA 8600 electron microprobe. Formula reduction to cations and the calculation of thermodynamic activities were made using the program AX62 (link from https://www.esc.cam.ac.uk/directory/tim-holland). Details of EPMA analytical procedures are given in Be’eri-Shlevin et al. (Reference Be’eri-Shlevin, Avigad and Matthews2009). A total of 5 garnet-bearing metapelite rocks, 10 amphibolites (both Ca & Fe-Mg hornblende types), one calc-silicate and one granitoid were analyzed by EPMA. The full list of individual analyses is given in Supplementary Table S1. Cation mole fractions (e.g., XMngt) are expressed relative to the total number of cations in the formula site.
The trace and major element whole rock chemistry in ppm and wt.%, respectively, of 10 samples of metapelites and amphibolites (Supplementary Table S2) was determined at the Geological Survey of Israel using ICP-AES and ICPMS. Details of the methods and errors of the whole rock analyses may be found in Benkovitz et al. (Reference Benkovitz, Matthews, Teutsch, Poulton, Bar-Matthews and Almogi-Labin2020). Duplicate measurements on metapelite rock KD284 (analyses KD284A and KD284B in Table S2) indicate the high reproducibility of the whole-rock values. Sample coordinates are given in Tables S1 and S2.
Oxygen isotopes. Quartz-andalusite vein samples from Donoussa were sliced and crushed, and mineral separates were hand-picked. Three different samples were analyzed (Table 2): two from boudinaged lenses (KKD134,149) and one from a vein (VD298). The oxygen isotope analyses were made with the laser fluorination extraction systems at the Hebrew University of Jerusalem, which used a Merchantek EO 30W CO2 laser ablation station, BrF5 reagent and a Micromass 602ES mass spectrometer. Additional analyses were made using the laser fluorination system at the University of Wisconsin-Madison using a CO2 laser, BrF5 reagent and a Finnigan-MAT 251. Analytical procedures follow those described by Valley et al. (Reference Valley, Kitchen, Kohn, Niendorf and Spicuzza1995). All results are reported in δ values as ‰ relative to the SMOW standard and are standardized to UWG-2 garnet = 5.80 ‰ VSMOW, which is tied to NIST 28 quartz = 9.59‰ (Valley et al. Reference Valley, Kitchen, Kohn, Niendorf and Spicuzza1995).
Mineral thermobarometry data for (B(I)) garnet-sillimanite-biotite-muscovite-plagioclase rocks

Refractory accessory mineral (RAM) thermometry data

α = (1000 + δ18Oquartz/(1000 + δ18Oandalusite).
Phase diagram calculations. The P–T pseudosection calculations use a combination of MAGEMin software (version 1.4.7) and thermocalc software (Holland and Powell, Reference Holland and Powell2011: version 3.50, database tcd-633). thermocalc databases are based on continuously updated experimental data. The advantage of the pseudosection calculations that is utilized here is their ability to explore mineral composition variations in complex bulk systems. (A review of their benefits over standard phase diagrams or petrogenetic grids can be found in https://serc.carleton.edu). The metapelite and amphibolite diagrams were calculated with thermocalc but using MAGEMin initially to ensure finding all stable assemblages. MAGEMin is a newly developed open-source parallel function software (Riel et al. Reference Riel, Kaus, Green and Berlie2022, https://github.com/ComputationalThermodynamics/MAGEMin) that minimizes the Gibbs free energy of multiphase and multicomponent systems to provide a stable, consistent and rapid phase equilibrium prediction routine. It has proved very useful in guiding calculations using the tcd-633 dataset. Theriak-Domino software (de Capitani and Brown, Reference de Capitani and Brown1987; de Capitani and Petrakakis, Reference de Capitani and Petrakakis2010), which also calculates multiphase equilibria using Gibbs free energy minimization, does not work with recent thermocalc databases. All the phase diagram calculations are made using molar quantities (e.g., bulk-rock MnO mol %) normalized to the whole-rock chemical analyses.
4. Results
4.a. Electron microprobe mineral chemistry
4.a.1. Metapelites
Feldspars are mainly oligoclase (XCafsp= 0.26–0.32); biotite is essentially a phlogopite-annite-eastonite solid solution (XMgbi = 0.44–0.51), and neither mineral shows core-to-rim compositional zoning (Table S1). In contrast, with one exception (sample KD139), the garnets show marked core-to-rim zoning of Fe, Mg and Mn and high Mn contents (XMngt = 0.12–0.32), with rims generally showing higher XMngt values than cores (Table 1; Table S1).
The detailed shape of the core–rim garnet zoning is revealed by an EPMA profile of a 200 μm radius garnet from the B(I) gt–sill–bi–mu sample KD156 (Figure 4). Core profiles maintain near-constant values of XFegt, XMggt and XMngt, but as the rim exterior is approached (∼40 μm from the edge), marked shifts occur characterized by a rise in XMngt and decreases in XFegt, XMggt and Mg/(Mg+Fe) mol ratios. The XCagt EPMA values are relatively constant throughout the grain (0.28 ± 0.2) but are consistent with the equilibrium 0.32 to 0.26 range at 650°C to 550°C calculated later (section 5.b.1). These basic features are evident in the diffusion profiles for Fe, Mg, Mn and Ca described by Carlson (Reference Carlson2006).
Core–rim EPMA profile of garnet sample KD156. The core shows no significant change, but marked compositional shifts occur in the rim (ca <40μm from the grain edge), with XMggt and XFe2+gt decreasing but XMngt increasing. The Mg/(Mg+Fe2+) ratio also decreases from core to rim – consistent with retrograde cooling exchange. Refer to Sections 4.b.1 and 5.a.2 for the description of the XCagt core-rim profile.

The core-to-rim Fe-Mg-Mn zoning in garnets found in this study contrasts with the lack of zoning observed in garnets by Martha et al. (Reference Martha, Xypolias, Cheng, Dörr, Gerdes, Hezel, Kutzschbach, Millonig, Schmeling, Marschall, Müller and Zulauf2025). The study of Fe-Mg interdiffusion shows that the diffusion coefficient DFe-Mg is independent of garnet composition but enhanced by water presence (Zhang et al. Reference Zhang, Bowen, Zhao and Yang2019). Using their data at a nominal temperature of 723°C (see thermometry below), a simple one-dimensional diffusion calculation for dry conditions (DFe-Mg = 10−23 m2/sec) would require 1.25 Ma to produce a 40 mm wide diffusional zone. The presence of >100 ppm water would enhance DFe-Mg by two orders of magnitude, giving a much shorter timescale of 1250 a.
The formation of garnet in the metapelite rocks is consistent with studies showing that the presence of bulk-rock Mn above typical pelitic values for the Scottish Dalradian Barrovian (0.01 to 0.03 wt.%; Atherton and Brotherton, Reference Atherton and Brotherton1982) enlarges the garnet P–T stability field (Symmes and Ferry, Reference Symmes and Ferry1992; Droop and Harte, Reference Droop and Harte1995; Maher et al. Reference Mahar, Baker, Powell, Holland and Howell1997; White et al. Reference White, Powell and Johnson2014b). The Donoussa whole rock metapelite samples have Mn = 0.04–0.07 wt.% (Table S2), but the presence of these garnets in thin cm-scale lamellae within the metapelite matrix suggests that local Mn rock contents may be even higher within these bands. As a feature of the protolith clastic pelite sediments, such Mn-enrichments most likely reflect a redox front whereby pore solution Mn2+ became oxidized to Mn-oxides (Thomson et al. Reference Thomson, Higgs, Croudace, Colley and Hydes1993). In this respect, the metapelite garnets preserve a Mn redox geochemical feature of their marine sediment origin.
4.a.2. Amphibolite-calcsilicate rocks
The EPMA data of Ca hornblendes from the amphibolite and one calc-silicate rock (KD44) is plotted on a diagram of Mg/(Mg+Fe2+) vs. Si cation pfu (Figure 5a). The classification diagram is after Leake et al. (Reference Leake, Woolley, Arps, Birch, Gilbert, Grice, Hawthorne, Kato, Kisch, Krivovichev, Linthout, Mandarino, Maresch, Nickel, Rock, Schumacher, Smith, Stephenson, Ungaretti, Whittaker and Youzhi1997). Compositional zoning was not observed for the Ca-amphiboles, and they are classified as either magnesiohornblende-actinolite or tschermakite-ferrotschermakite solid solutions: hereafter ‘Mg-hornblende’ and ‘Tschermakite’ types (Figure 5a; Table S1). The two types of Ca amphibole occur together only in one sample (KD 158), but even here they occur as separate minerals not in contact, rather than showing evidence of unmixing or overgrowth of one by the other. The formation of two separate amphibole types in the same rock samples is a common occurrence in amphibolite facies rocks (Veblen and Ribbe Reference Veblen and Ribbe1981) and could suggest the existence of a miscibility gap. Some studies, however, suggest that there is a compositional gap related to disequilibrium processes rather than a true solvus phenomenon (e.g., Grapes and Graham, Reference Grapes and Graham1978). In one tschermakite-bearing rock (KD 127), the Ca amphibole is surrounded by a cummingtonite rim. This textural type has been interpreted to reflect lower pressure exchange during the retrograde reaction: tschermakite + SiO2 = cummingtonite + anorthite + fluid (Yardley, Reference Yardley1989, pp. 101–2).
Calcic amphibole chemistry in amphibolites and calc-silicate rocks. (a) Molar Mg/(Mg+Fe) vs. Si cations classification plot (modified after Leake et al. Reference Leake, Woolley, Arps, Birch, Gilbert, Grice, Hawthorne, Kato, Kisch, Krivovichev, Linthout, Mandarino, Maresch, Nickel, Rock, Schumacher, Smith, Stephenson, Ungaretti, Whittaker and Youzhi1997). Sample numbers are indicated in the legend. The analyses fall into two separate groups: 1) the magnesiohornblende (actinolite) group and 2) the tschermakite (ferrotschermakite – ferrohornblende) group. (b) Plot of hornblende-plagioclase ed-tr and ed-ri temperatures (Holland and Blundy, Reference Holland and Blundy1994) versus Si cations. The temperatures given by the Tschermakite group identified in (a) are higher than those given by the Magnesiohornblende group. Actinolites give the lowest temperatures (∼550°C), suggesting that the latter group reflect retrograde cooling. The orange star in Figure 5a refers to the calculated hornblende composition in section 5.4.

EPMA data for orthorhombic and monoclinic Mg-Fe amphibole analyses are given in Fig. S3 and discussed in section 5.c. below.
Petrographically, calc-silicate rocks comprise subhedral intergrowths of garnet, diopside, Ca-amphibole, epidote and calcite. These minerals are overgrown by euhedral grains of later actinolite (Figure 2h). The hornblende of these rocks belongs to the ‘Tschermakite’ type rather than the ‘Mg-hornblende’ type (Figure 5a).
4.b. Refractory accessory mineral (RAM) oxygen isotope data
The three rocks analysed in this study are almost entirely composed of quartz with conspicuous isolated large (cm-size) porphyroblasts of pink andalusite (Figure 2c). Fibrolite sillimanite overgrowths indicate that the rocks entered the sillimanite P–T field during their formation. As noted in section 2 above, the host rocks are metapelites. Several hand-picked separates were analysed in each sample to check for isotopic homogeneity. The relatively narrow standard deviation of mineral analyses given by each rock sample is close to the analytical uncertainty of the laser fluorination procedures and suggests that there is no isotopic zoning in the samples. Such a lack of zoning in the RAM mineral was also found in an ion microprobe study of δ18O in corundum on Naxos (Turnier et al. Reference Turnier, Katzir, Kitajima, Orland, Spicuzza and Valley2020).
The oxygen isotope analyses were determined on both quartz and andalusite separates from lenses and veins that were rotated into parallelism with the main S2 foliation in the schists. Mean quartz-andalusite fractionations expressed as 1000 ln α(q-and) range between 2.44 ± 0.13‰ and 2.78 ± 0.12‰ (1SE) (Table 2).
5. Discussion
5.a. P–T–X conditions of garnet-bearing metapelite assemblages
5.a.1. Fe/Mg exchange thermometry between garnet and biotite
This exchange reaction has been a mainstay of the thermometry of garnet-biotite-bearing metapelites since the publication of the experimental calibration by Ferry and Spear (Reference Ferry and Spear1978). The calibration (labelled here FS78) can be applied using ideal mixing activity models, which were considered valid to within (Ca+Mn)/total cations <0.2. Most garnets in this study plot close to or above this limit. Several alternative activity-composition models have also been suggested for this reaction (e.g., Bhattacheria et al. Reference Bhattacharya, Mohanty, Maji, Sen and Raith1992), but we retain use of the original ideal solution formulation as the model for this thermometer. The second thermometer calibration used in this work (labelled here HP11/AX62) uses the activity-composition models AX62 (https://www.esc.cam.ac.uk/directory/tim-holland), which are compatible with the thermodynamic data of Holland and Powell (Reference Holland and Powell2011). The calculations used the standard equation ΔG°(T, P) = −RTlnK + PΔV, with the enthalpy, entropy and volume data needed to calculate ΔG° taken from Holland and Powell (Reference Holland and Powell2011). The calculations were made using constant ΔVs (also for the GASP barometry below). The thermometric data are given in Table 1 for individual core and rim analyses. Similar overall ranges are observed with each calibration: cores 530–670°C, rims 480–590°C (FS78); cores 520–650°C, rims 460–560°C (HP11/AX62). Consistent with the core-rim profiles in Figure 4, cores give higher temperatures than rims using both calibrations, i.e., the core-to-rim temperatures record a retrograde cooling history from peak or near peak metamorphic temperatures.
5.a.2. Garnet-aluminium silicate-plagioclase (GASP) barometry
The barometry is based on the partitioning of Ca between garnet and plagioclase in the presence of Al-silicate and quartz (Ghent, Reference Ghent1976). As previously noted (Section 4.b.1) above, the relatively flat core-rim profile for Ca is within the range anticipated for diffusional exchange at equilibrium. An alternative model for Mg/(Mg+Fe) and Mn enrichment in amphibolite-grade metapelites envisages a retrograde net-transfer reaction due to the dissolution of garnet rims and Fe-Mg exchange with biotite (Kohn and Spear, Reference Kohn and Spear2000). In this case, cation exchange P–T methods resulted in errors of hundreds of degrees and 3–6 kbar compared with thermobarometry of nearby rocks and petrogenetic grids. Net reaction was proposed to bring about local change in the chemistry of a biotite matrix and produce errors in the thermobarometry. Such large P–T errors are clearly not the case in our study, and we found no variation in the core-to-rim compositions of biotite and plagioclase, whose modal abundance largely buffers Fe-Mg and Ca exchange with garnet, respectively. Accordingly, we have adopted the retrograde diffusive exchange mechanism to model the garnet zoning.
Here, we use the equations of Hodges and Spear (Reference Hodges and Spear1982) (labelled HS82) with FS78 temperatures to calculate pressures using a non-ideal mixing model for Ca-Mg exchange and a fixed activity coefficient γAn = 2.0 for the anorthite component of plagioclase. In the second calculation, we estimated pressures using HP11/AX62 temperatures. The results of these calculations are summarized in Table 1. Pressure estimates show variations: cores 2.0–3.2 kbar, rims 1.1–2.8 kbar (FS78/HS82); cores 3.0–3.9 kbar, rims 1.7–3.5 kb (HP11/AX62). The P–T data are plotted in Figure 6a on a diagram of the Al-silicate polymorphs. Errors for single points are stated as ±40°C for temperature and ±0.5 kbar for pressure. The higher garnet core pressures given using the HP11/AX62 calibration relative to FS78/HS82 are evident, although it can be noted that a relatively small adjustment of the HS82 γAn value to 1.8 would give similar P values from both calculations. Rim temperatures and pressures are lower than cores using both calibrations, and all rim values are within the andalusite stability field (Figure 6a). As is common, the peak metamorphic sillimanite remains in the metamorphic assemblage even though the retrograde P–T conditions move into the andalusite field. The lower core temperatures of sample KD139 (500–520°C) will be discussed in section 5.c.2. An additional trend evident from the EPMA data is the inverse relationship of XMngt to temperature (Figure 6b).
Metapelite P–T graphs deduced from the gt + sill + bi + plag EPMA data in Table 1 and Table S1. (a) P–T diagram of temperatures and pressures deduced by cation exchange thermobarometry (Sections 5.a.1 and 5.a.2). The red-coloured P–T datapoints are calculated using the Holland and Powell (Reference Holland and Powell2011) thermocalc data set (HP11/AX62). The blue P–T datapoints are calculated using the Ferry and Spear (Reference Ferry and Spear1978) and Hodges and Spear (Reference Hodges and Spear1982) equations (FS78/HS82). ΔVs was assumed constant in each calculation. Similar core and rim temperature ranges are obtained from both calibration sets, but lower pressure estimates are given by the FS78/HS82 combination. Each calibration sets show lower P–T values in rims compared to cores, consistent with retrograde cooling and decompression. Errors for individual datapoints are set as ± 40°C and ± 0.5 kbar, respectively. The aluminium silicate phase diagram is calculated using the Holland and Powell (Reference Holland and Powell2011) dataset. (b) A plot of HP11/AX62 temperatures vs. mol fraction of Mn in garnet (XMngt data in Table 1) showing the general increase in XMngt with decreasing temperature and lower rim values compared to garnet cores. The black line is a linear fit to the data.

The most closely related groups of rocks in terms of rock types, petrology and age are found 250 km further south in the Asteroussia Crystalline Complex of south-central Crete (Figure1), where sillimanite ± andalusite + garnet + biotite + cordierite rocks are identified (Seidel et al. Reference Seidel, Seidel, Okrusch, Kreuzer, Raschka and Harre1976, Reference Seidel, Okrusch, Kreuzer, Raschka and Harre1981; Koepke and Seidel, Reference Koepke and Seidel1984; Langosch et al. Reference Langosch, Seidel, Stosch and Okrusch2000; Martha et al. Reference Martha, Dörr, Gerdes, Krahl, Linckens and Zulauf2017). Peak P–T conditions of 650–700°C and 4–6 kbar were constrained using garnet-cordierite Mg/Fe thermobarometry and petrogenetic grid modelling (Seidel et al. Reference Seidel, Seidel, Okrusch, Kreuzer, Raschka and Harre1976; Koepke and Seidel Reference Koepke and Seidel1984).
5.a.3. P–T–X equilibrium phase diagrams for metapelite rock KD284
The main question addressed here concerns the P–T–X conditions that stabilize the B(I) garnet + biotite + sillimanite + muscovite assemblage instead of the more commonly observed B(II) sillimanite ± andalusite ± biotite ± muscovite assemblage. Bulk-rock chemical analyses of the Donoussa metapelites (Table S2) indicate a shale protolith (Kolodner, Reference Kolodner1999). White et al. (Reference White, Powell and Johnson2014b) utilized new activity composition relations for manganese-bearing minerals to calculate the P–T pseudosections in the metapelite system MnNCKFMASHTO system using thermocalc and the HP2011 dataset. These calculations show that an increase of Mn in a metapelite bulk-rock composition at amphibolite facies conditions leads to marked enlargement of the garnet P–T stability field. We follow the same calculation method for bulk-rock chemical analysis of metapelite rock KD284, which features B(I) laminae within a B(II) rock matrix (chemical analysis 284A in Table S2; in the following account, KD284 will refer to the rock, and KD284A refers to the chemical analysis used in phase diagram calculations). The pseudosections were constructed for three bulk-rock MnO compositions of 0.03, 0.07 and 0.10 mol%. The first value of 0.03 mol% represents the average Dalradian metapelite (Atherton and Brotherton, Reference Atherton and Brotherton1982). The second value is the bulk-rock analysis KD284A. The final value of 0.10 mol % is the MnO value used by White et al. (Reference White, Powell and Johnson2014b) and is better expressive of the garnet-bearing Mn-enriched lamellae within bulk rock KD284 and other garnet-bearing rocks.
The intensive variable parameters applied to the equilibrium calculations are defined in Table 3: Fe3+/FeT) = 0.15, XH2O = 1 (water at saturation) and f O2 = QFM+2. The bulk composition is adjusted so that there is just enough H2O to saturate all hydrous phases at the beginning of melting. Unlike the metapelite bulk-rock compositions of White et al. (Reference White, Powell and Johnson2014b) and Pattison and Goldsmith (Reference Pattison and Goldsmith2022) for the Buchan assemblages of the Scottish Dalradian, magnetite does not appear in phase diagrams for KD284, which accords with our petrographic findings.
thermocalc parameters for P–T pseudosection calculations

The P–T pseudosections show that increasing bulk-rock Mn enlarges the garnet + biotite + Al-silicate + muscovite P–T fields (outlined in red colour in Figure 7). At bulk-rock Mn = 0.03 mol%, these garnet assemblages are not observed in the diagrams, and the major mineral parageneses with Al-silicates are the B(II) sill + bi + mu & and + bi + mu (+ q +pl + ilm) assemblages (Figure 7a). Staurolite-bearing assemblages are restricted to the lower temperature–higher pressure segments of the phase diagram. However, since staurolite is only recognized in the D1 inclusion assemblages but not in syn-D2 assemblages (Section 2 above), it follows that temperatures and pressures of the B(II) assemblages fall below the staurolite-in reactions in Figure 7a, i.e., ∼4.5–6 kbar at T > 600°C.
thermocalc P–T pseudo sections in the MnNCKFMASHTO system for Donoussa metapelite rock KD 284 (Table 2; analysis KD284A) illustrating the progressive replacement of the Al-silicate + biotite+ muscovite assemblages (and + bi + mu & sill + bi + mu) by the garnet-bearing assemblages (g + and + bi + mu & g +sill + bi + mu) with increasing bulk-rock MnO. The sub-solidus garnet-bearing fields are outlined in red font, and the Al-silicate + biotite+ muscovite subsolidus fields are outlined in white font. (a) MnO = 0.03 mol%. (b) MnO = 0.07 mol%. (c) MnO = 0.10 mol%. Staurolite is stable on the higher-pressure side of its univariant boundaries with the Al-silicate +biotite assemblages and at higher temperatures than chlorite-bearing assemblages. Magnetite (mt), K-feldspar (ksp) and cordierite (cd) are restricted to pressures lower than 2.5 kbar and temperatures greater than 575°C.

Garnet becomes an additional phase to the and/sill + bi + mu P–T fields with an increase in bulk rock MnO. At MnO = 0.07 mol%, this takes the form of the B(I) gt + sill + bi + mu assemblage at T > 600°C (Figure 7b), whilst the pseudosection at bulk-rock MnO = 0.10 mol% shows a significant enlargement of garnet stability (Figure 7c). The upper temperature limit of all solid assemblages is the solidus at T > ∼650–660°C. Higher bulk-rock MnO mol% values clearly enlarge the P–T field at which garnet becomes a stable additional phase.
An alternative aspect of the role that bulk-rock Mn exerts on equilibrium assemblages is illustrated by the P–bulk-rock MnO mol% diagram at 600°C (Figure 8a). At bulk-rock MnO ∼ 0.10 mol% (vertical red dotted line) and P < 4.5 kbar, the assemblage gt + sill + bi + mu is stable until P ∼ 3.7 kbar, where the and + bi + mu field is entered. As pressures decrease below 3.5 kbar, the stability of the gt + and + bi + mu assemblage requires increasing bulk-rock MnO until cordierite and K-feldspar appear in the diagram at P ∼ 2 kbar. Correspondingly, the calculated XMngt contours increase from 0.30 to 0.44 (Figure 8a). These values are higher than those given by the garnet EPMA analyses and therefore are unrealistic for the Donoussa rock. At the lower bulk-rock MnO values (vertical red dotted lines for 0.07 and 0.03 mol%), bi + mu + and/sill are the stable assemblage.
thermocalc P–T pseudosections in the MnNCKFMASHTO system for the Donoussa metapelite sample KD284. (a) P–bulk-rock MnO mol% relations at T = 600°C. The contours on the right side of the figure show the calculated values of XMngt in mineral assemblages in which garnet is stable at 600°C. The three vertical red lines indicate the values of bulk-rock MnO mol% used in Figure 7 (from left to right: 0.03, 0.07, 0.10 mol%). (b) Calculated contours of log fO2 as a function of temperature for bulk-rock MnO = 0.10 mol%. Relative to the QFM oxygen buffer at 3 kbar, the contours correspond to logfO2 = QFM + 2.3.

The P–T and P–bulk-rock MnO pseudosections indicate that low bulk-rock MnO and higher pressures stabilize garnet-free assemblages in the high T/P metamorphism, whereas higher bulk-rock MnO and lower pressures expand the stability field of garnet (Figures 7 and 8a).
5.a.4. Oxygen fugacity of the metamorphic rocks
The oxygen fugacity value taken in the P–T pseudosection calculations is log fO2 = QFM+2, where the bulk composition O2 content was adjusted to the fO2 determined for Barrovian metapelites from Glen Clova, Dalradian, Scotland (Ague et al. Reference Ague, Baxter and Eckert2001) and is also compatible with wet chemical analyses at that locality, which gave Fe 3+/Fe T = 0.14 (Chinner, Reference Chinner1960). Indeed, it has been argued that Fe2O3 is an excellent measure of redox (e.g., Diener and Powell, Reference Diener and Powell2010; Forshaw and Pattison, Reference Forshaw and Pattison2023). An independent calculation of oxygen fugacity in Donoussa metapelites is made here using the equilibrium reaction:
$\begin{align}8/&3\;quartz + 4/3\;sillimanite + 2\;haematite \\&= 4/3\;almandine + {O_2}\end{align}$
with free energies taken from the HP2011 dataset and calculated activities of haematite (in ilmenite) and almandine (in garnet). At lower P and T, andalusite and annite are used in place of sillimanite and almandine. The log fO2 values vary from −16 to −19.5 (Figure 8b) and relative to the QFM line calculated at 3 kbar (https://fo2.rses.anu.edu.au) plot at log10 fO2 = QFM + 2.3. This is close to the value of QFM +2 that was used in the present phase diagram metapelite calculations (Ague et al. Reference Ague, Baxter and Eckert2001) and widens the potential use of this oxygen buffer value to Buchan high T/P metamorphism.
5.b. Petrogenetic implications of the metapelite data
5.b.1. Thermobarometric record of metapelites
The cation thermobarometry P–T estimates from the gt + sill + bi + plagioclase + quartz core and rim samples reflect a retrograde P–Tpath that best encompasses the g-bi-and/sill fields in the MnNCKFMASHTO system calculated at bulk-rock MnO = 0.10 mol% (Figures 7c and 9a). The thermocalc calculated XMngt contours in the gt + and/sill + bi + + mu fields increase from 14 to 34 mol% with decreasing P–T (Figure 9b). Thus, the EPMA core-to-rim XMngt increase with decreasing T evident in Figure 6b clearly cannot represent a prograde path through the chl + gt + bi + mu field, which would feature a decrease in XMngt. Similarly, a prograde P–T path through the st + gt + bi + mu field would not lead to the observed range of Mn zoning in the garnets. The Mn zoning, however, adequately represents garnet retrograde P–T decrease in the gt + and/sill + bi + mu field. Thus, both the calculated and measured XMngt show a consistent increase that reflects the cooling-decompression P–T path. The observed increase in XMngt with lower P and T provides independent confirmation of the retrograde P–T path.
Diagrams illustrating projections of the garnet metapelite thermobarometry data and EPMA XMngt values onto the thermocalc P–T pseudosections of sample KD284 at bulk-rock MnO = 0.10 mol%. Data sources for the P–T calculations are given in Table 1. P–T fields of mineral + liquid assemblages above the solidus are plotted on the pseudosections. (a) Mean thermobarometry estimates using HP2011/AX62. Red points = garnet core data, white = garnet rim data. (b) P–T trend of EPMA-measured XMngt values plotted on the thermocalc pseudosection together with calculated XMngt contours (expressed as XMngt × 100) in the garnet P–T stability fields. The EPMA-measured XMngt data (Table 1) are plotted on the mid-point of the calculated XMngt × 100 contours in the gt + and + bi + mu & gt + sill +bi +mu fields. Red data points are the garnet core analyses; white data points are rim analyses. Both diagrams reveal well-defined down P–T trends.

Petrographically, sample KD139 is a B(II) sill + bi + mu rock containing sporadic garnet grains. The HP2011/AX62 core temperatures given by this garnet are unusually low (520–530°C, Table 1), i.e., within the chl + bi + mu field (Figure 9a), whereas the XMngt is the highest measured in our samples (31–33 mol%). The garnet also shows a small andradite component (Table S1). Pressure values show a considerable variation (3.2, 1.7 kbar, mean = 2.4 ± 0.7 kbar). These P–T values possibly reflect the D3 retrograde overprint recognized by Martha et al. (Reference Martha, Xypolias, Cheng, Dörr, Gerdes, Hezel, Kutzschbach, Millonig, Schmeling, Marschall, Müller and Zulauf2025).
5.b.2. Anatectic melting in metapelites?
Studies of Miocene igneous intrusions in the Cyclades reveal that garnet-tourmaline two-mica leucogranite sills and dykes have a crustal source, i.e., can be classified as S-type granites (Pe-Piper et al. Reference Pe-Piper, Kotopouli and Piper1997; Matthews et al. Reference Matthews, Putlitz, Hamiel and Hervig2003; Pe-Piper and Piper, Reference Pe-Piper and Piper2005; Lamont et al. Reference Lamont, Roberts, Searle, Gardiner, Gopon, Hsieh, Holdship and White2023a). Similar igneous mineralogy is found in the Late Cretaceous peraluminous garnet-tourmaline two-mica granitoids on Donoussa and were assigned to an I-type to transitional I/S type calc-alkaline affinity (Altherr et al. Reference Altherr, Kreuzer, Lenz, Wendt, Harre and Dürr1994). The comparison of Anafi granitoid mineralogy and geochemistry with data on Donoussa and Eastern Crete led Koutsovitis et al. (Reference Koutsovitis, Soukis, Voudouris, Lozios, Ntaflos, Stouraiti and Koukouzas2022) to interpret all three localities to represent regionally widespread I-type/minor S-type granitoids with calc-alkaline affinities. The presence of minor S-type granitoids with a metasedimentary protolith clearly is a feature of both the Miocene and Cretaceous magmatic rocks.
The Donoussa granitic rock sample analyzed by EPMA (KD266: Figure 3, Table S1) is a small pegmatitic granitoid dike intruding amphibolite rock, composed of pl + ksp + q ± bi ± mu ± gt ± tour. The garnet is spessartine-rich; the plagioclase is albite (Table S1), as previously noted by Altherr et al. (Reference Altherr, Kreuzer, Lenz, Wendt, Harre and Dürr1994) for garnet-tourmaline two-mica granitoids. The calculated P–T solidus for KD284 indicates that melting begins at temperatures above 650°C–660°C in the pure system (Figures 7 and 9). The proximity of peak metamorphic conditions to the solidus suggests that melting of the metapelites at greater depth could also have generated granitoid magmatic liquids. The assemblages immediately above the melting curves are, respectively, bi + sill + mu + liquid and gt + bi + sill + mu + liquid (Figure 9a). At pressures of 3–4 kbar, the k-feldspar-bearing assemblage bi + ksp + sill + liquid becomes stable at temperatures above the solidus.
The calculated liquid composition using thermocalc in such melting indicates a hydrous alkali alumino-silicate melt (Table 4). Compared to early and late fractionated garnet-tourmaline two-mica granitoids and the aplites (Table 4), this liquid composition is closest to the late fractionated granitoid composition. Thus, it is feasible that part of the liquid involved in cogenetic granitoid formation could be derived from anatectic melting of metamorphic pelitic gneisses and schists. It should be emphasized, though, that the phase diagram calculation assumes that the fluid is exhausted immediately on liquid coming in, and the amount of melt generated would be small. Local fluxing by additional water released at the metamorphic peak or magmatic fluid would be essential to generate larger melt amounts.
Calculated liquid composition and representative granitoid compositions from Donoussa

5.b.3. Context with Buchan-type metamorphism
This section explores the perspectives gained from the classic high T/P metamorphism recorded by the Buchan Zone metamorphosed rocks, particularly aspects related to staurolite breakdown in metapelites. Although they are not specifically included in the database, the metapelite rocks of Donoussa belong to the andalusite-sillimanite Mineral Assemblage Sequence MAS2, which essentially reflects an isobaric rise in temperature (Pattison and Forshaw, Reference Pattison and Forshaw2025). The MAS2 assemblages in this database include important Buchan sequences of the Scottish Dalradian. Staurolite formation associated with the earlier S1 schistosity (Section 2 above) led Martha et al. (Reference Martha, Xypolias, Cheng, Dörr, Gerdes, Hezel, Kutzschbach, Millonig, Schmeling, Marschall, Müller and Zulauf2025) to interpret the D1 event as the initial prograde phase of an overall clockwise P–T–t path.
The occurrence of staurolite relics within andalusite poikiloblasts in Donoussa metapelites led Altherr et al. (Reference Altherr, Kreuzer, Lenz, Wendt, Harre and Dürr1994) to infer that the prograde P–T path crossed the univariant reaction in KFMASH:
Altherr et al. (Reference Altherr, Kreuzer, Lenz, Wendt, Harre and Dürr1994) based their prograde path on the experimental calibration of the reaction for an intermediate composition of magnesian staurolite (Hoschek, Reference Hoschek1969). At PH2O = 3 kbar, the reaction is crossed at T = 600°C. The equivalent reaction in the MnNCKFMASHTO system is the conversion of the st + bi + mu assemblage to form bi + mu + and and/or gt + and + bi + mu (Figure 7). These reactions resemble the andalusite overprint of an earlier staurolite assemblage in the Buchan zone high T/P metamorphism, where staurolite is replaced by andalusite in the Boyndie Bay coastal section and the Huntly to Aberchirder section (Pattison and Goldsmith, Reference Pattison and Goldsmith2022; Pattison and Forshaw, Reference Pattison and Forshaw2025). Here too, staurolite is found as relics within andalusite porphyroblasts, and manganiferous garnet is also found within staurolite-andalusite rocks (Hudson, Reference Hudson1975; Hudson and Johnson, Reference Hudson and Johnson2015). Pattison and Goldsmith (Reference Pattison and Goldsmith2022) suggested the reaction of staurolite (higher P) to cordierite + andalusite (lower P) and related this reaction to two metamorphic heating cycles separated by an interval of cooling and decompression of ∼1 kbar. Similar cooling and decompression between two metamorphic cycles (staurolite D1 to sillimanite-grade D2) could also be relevant for Donoussa.
The Buchan zone metamorphic rocks notably developed in a backarc volcanic environment with high T/P heating being driven by contemporary calc-alkaline igneous intrusion (Johnson et al. Reference Johnson, Kirkland, Viete, Fischer, Reddy, Evans and McDonald2017). Polymetamorphism is a feature of the Alpine orogenesis in the Hellenides (Papanikolaou, Reference Papanikolaou2021). Based on the analogy with the Buchan zone polymetamorphism given above, an additional tectonic hypothesis can be proposed for Donoussa. Tectonic reversal accompanying switching from active subduction and volcanism (D1 staurolite stage) to slab rollback (Lister and Forster, Reference Lister and Forster2009, Johnson et al. Reference Johnson, Kirkland, Reddy and Fischer2015, Reference Johnson, Kirkland, Viete, Fischer, Reddy, Evans and McDonald2017) would bring the backarc volcanism and associated high T/P metamorphism directly above the subduction zone (D2 sillimanite stage). This Late Cretaceous magmatism and metamorphism would then serve as a prelude to major southward transport toward the Cyclades during the large-scale slab retreat that eventually brought the Upper Tectonic Unit on Donoussa in tectonic juxtaposition with the Cyclades Blueschist Unit on Naxos.
5.c. Petrology of Donoussa amphibolites
5.c.1. Geochemical provenance of amphibolites
Geochemical whole rock studies on the granitoid rocks of Anafi, Donoussa and the Asteroussia units of Crete indicate their calc-alkaline magmatic arc origin, characterized mainly by I-type granites (Reinecke et al. Reference Reinecke, Altherr, Hartung, Hatzipangiotou, Kreuzer, Harre, Klein, Keller, Geenen and Boeger1982; Altherr et al. Reference Altherr, Kreuzer, Lenz, Wendt, Harre and Dürr1994; Langosch et al. Reference Langosch, Seidel, Stosch and Okrusch2000; Martha et al. Reference Martha, Dörr, Gerdes, Petschick, Schastok, Xypolias and Zulauf2016, Reference Martha, Dörr, Gerdes, Krahl, Linckens and Zulauf2017; Koutsovitis et al. Reference Koutsovitis, Soukis, Voudouris, Lozios, Ntaflos, Stouraiti and Koukouzas2022; Kneuker et al. Reference Kneuker, Dörr, Petschick and Zulauf2015). Total alkali silica (TAS) major element plots on five Donoussa amphibolite samples (Table S2) show that both the Ca- and Fe-Mg types are basaltic to basaltic andesite in composition (Supplementary Fig. S1). Trace element chemical analyses also show calc-alkaline affinity on REE/chrondrite, rock/MORB spider and Ti-Zr and Nb-Zr-Y trace element discrimination diagrams (Supplementary Fig. S2), comparable with the calc-alkaline igneous affinity found for granitoids (Altherr et al. Reference Altherr, Kreuzer, Lenz, Wendt, Harre and Dürr1994). The Ca- and Fe-Mg amphibolite types do not show marked differences on the trace element diagrams. The amphibolite geochemistry is thus compatible with a continental igneous arc origin evident for the granitoid rocks of Anafi and Donoussa (Altherr et al. Reference Altherr, Kreuzer, Lenz, Wendt, Harre and Dürr1994; Be’eri-Shlevin et al. Reference Be’eri-Shlevin, Avigad and Matthews2009; Koutsovitis et al. Reference Koutsovitis, Soukis, Voudouris, Lozios, Ntaflos, Stouraiti and Koukouzas2022).
5.c.2. Amphibolite/calc-silicate thermometry
Amphibolite and calc-silicate temperatures are calculated using the online version (https://www.esc.cam.ac.uk/directory/tim-holland) of the hornblende-plagioclase (hb-pl) thermometers of Holland and Blundy (Reference Holland and Blundy1994), which includes the amphibole Fe3+ calculation used in AX. Since quartz is commonly found in amphibolite C(I) assemblages, temperatures were calculated using both the edenite-tremolite (ed-tr) and the edenite-riebeckite (ed-ri) pairs. Temperatures average 693°C ± 42°C for ‘Tschermakite’ type amphiboles, whereas ‘Mg-hornblende’ temperatures average at 619°C ± 41°C (Figure 5b; Table S1). This significant temperature difference could reflect errors in the thermometer model calculation or retrograde cation equilibration to lower temperatures in the ‘Mg-hornblendes’. The temperatures of amphibolites on Anafi Island (80 km S of Donoussa, Figure 1) calculated using the data of Be’eri (Reference Be’eri2003) and Be’eri-Shlevin et al. (Reference Be’eri-Shlevin, Avigad and Matthews2009) yield average values = 687°C ± 18°C and 655°C ± 23°C, respectively, for the two calibrations. Their better correspondence with the Donoussa ‘Tschermakite’ temperatures therefore suggests that the lower ‘Mg-hornblende’ temperatures on Donoussa reflect retrograde cooling to lower temperatures. This notion is supported by the fact that the lowest temperatures (<600°C) are given by actinolite hornblendes (Figure 5b). The peak metamorphic temperatures given by the ‘Tschermakite’ hornblendes are ∼30°C higher (but within errors) than the maximum values of ∼ 660°C indicated for sub-solidus assemblages and melting in metapelites (Figures 6a and 7).
‘Tschermakite’ analyses in calc-silicate rock KD44 give mean ed-tr temperatures of 702°C ± 18°C. Such temperatures fall within the hb-pl range for Ca-amphibolites (Figure 5b). Late-formed actinolite in KD44 gives a lower temperature of 550°C, consistent with its formation during cooling.
5.c.3. Phase diagram amphibolite P–T constraints
Amphibolite pseudosection calculations in the system NCKFMASHTO were made for C(I) ‘Tschermakite’ type sample KD149, a calcic amphibolite from the banded amphibolites (Table S1). In the absence of a direct measurement, Fe3+/ΣFe is taken as 10 mol% as is commonly done for MORB-type compositions. Such an assumption is necessary because natural rocks impose f O2 via PT and assemblage for a particular bulk composition. In a first calculation, MAGEMin software (Riel et al. Reference Riel, Kaus, Green and Berlie2022) was used to rapidly determine the P–T stable assemblages using the 6.33 dataset. The final pseudosections were then calculated using thermocalc using the same dataset (Figure 10). The two calculation methods give identical results.
Thermocalc P–T calculation of the equilibrium reactions in the system NCKFMASHTOCr for the ‘Tschermakite’ type amphibolite sample KD149 (Table S2). The orange arrows indicate the approximate cooling-decompression path deduced by the metapelite thermobarometry. The red highlighted line indicates the retrograde sphene (titanite)-in reaction. Cummingtonite becomes part of the assemblage at T < 700°C and P < 3 kbar (lower right corner in plot), but this boundary can shift to lower T for a slightly more reducing rock composition (see text in section 5.c.2).

The highest-grade assemblages for KD149 plot in P–T fields consisting of hb + bi + fsp + cpx + ilm + q. During cooling and decompression in this field (indicated for metapelites by the orange arrow), sphene comes in as a later phase at 540–560°C at the expense of clinopyroxene and ilmenite (Figure 10). The retrograde incoming of sphene is consistent with the U-Pb zircon dating of sphene (a.k.a. titanite) in calc-silicate rock to cooling below 600°C at 71.3 Ma (Martha et al., Reference Martha, Xypolias, Cheng, Dörr, Gerdes, Hezel, Kutzschbach, Millonig, Schmeling, Marschall, Müller and Zulauf2025). The calculated mineralogy compares well with the C(I) assemblage of the hbl + pl + cpx +q + bi +sphene rocks. All predicted assemblages are sub-solidus. Garnet is not predicted at the metamorphic pressures of this study but is present in some rocks (Table S1). The amphibolite garnet is relatively enriched in Mn (XMngt > 0.25) and possibly reflects the garnet stability-enhancing effect of Mn recognized in metapelites. The thermocalc-calculated composition of the hornblende mineral in KD149 at 620°C and 3.5 kbar is Si cations = 6.45 and Mg/(Mg+Fe) = 0.259. This value (orange star in Figure 5a) falls within the range of measured ferrotschermakite EPMA analyses for sample KD149.
Retrograde cummingtonite is not predicted by the diagram except at high temperature outside the metapelite cooling P–T path (Figure 10). An independent calculation shows that the cummingtonite-in boundary would shift down by 100°C (from 673°C down to 573°C at 3 kbar) if the bulk composition is more reducing, higher in Na and lower in Ti. The retrograde P–T path would then intercept this boundary.
5.c.3. Fe-Mg amphiboles
The stability of the garnet + anthophyllite + gedrite assemblage in the Mg-Fe amphiboles is most probably due to the rocks being more iron-rich than the more common anthophyllite-cordierite rocks (Spear, Reference Spear1980). Spear (Reference Spear1980) proposed a solvus cresting at 620°C ± 25°C, but it is difficult to accurately know the formation temperature of these samples relative to this solvus. On one hand, they seem to form a continuous series on the Mg/(Mg+Fe) vs. Si diagram (Fig. S3), suggesting temperatures above the solvus. On the other hand, anthophyllite and gedrite also form as separate minerals in different rocks and in the same rock (e.g., KD57), suggesting proximity to the solvus. High temperature conditions of T > 600°C are compatible with those deduced for the Ca-amphibolites.
5.c.4. Summary of amphibolite petrology
Donoussa amphibolites are interpreted as volcanoclastic sediments of calc-alkaline origin whose provenance is situated within the source continental igneous arc. The amphibolite and calc-silicate temperature data and P–T pseudosections calculated for sample KD149 confirm the high-temperature metamorphism evident in the metapelite geothermometric and phase diagram data (Section 5a above). Furthermore, the amphibolite data provides additional support for retrograde cooling from peak metamorphic temperatures of 660–690°C with sphene coming in as a retrograde mineral at T∼550°C.
5.d. Refractory accessory mineral oxygen isotope thermometry
Quartz-rich lenses and veins in metapelites potentially provide the basis for refractory accessory mineral (RAM) oxygen isotope thermometry of the high-grade metamorphic rocks. The modal dominance of quartz means that it will be the dominant oxygen reservoir for closed-system (at hand-sample scale) oxygen isotope exchange, and its δ18O value will not significantly change during retrograde cooling (Ghent and Valley, Reference Ghent and Valley1998). Thus, the oxygen isotope fractionation between quartz and refractory (i.e., slow diffusion rate of oxygen) accessory minerals such as Al-silicates potentially provides an accurate record of the metamorphic event that formed the veins that is unaffected by retrograde exchange (Valley et al. Reference Valley2001).
Temperatures are calculated using the petrologically calibrated quartz-Al-silicate equation: 1000 ln α = A × 106/T2 = 2.25 × 106/T2 (Sharp, Reference Sharp1995). Such RAM thermometry on quartz-kyanite veins gave temperatures of 635–690°C in high-grade rocks of the amphibolite-facies Upper Series (Koronos) unit on Naxos (Putlitz et al. Reference Putlitz, Valley, Matthews and Katzir2002). These temperatures are consistent with the Miocene Barrovian metamorphism estimates for the Naxos kyanite-sillimanite zone rocks (e.g., Buick and Holland, Reference Buick, Holland, Daly, Cliff and Yardley1989; Baker and Matthews, Reference Baker and Matthews1995; Lamont et al. Reference Lamont, Searle, Gopon, Roberts, Wade, Palin and Waters2020, Reference Lamont, Roberts, Searle, Gardiner, Gopon, Hsieh, Holdship and White2023a, Reference Lamont, Smye, Roberts, Searle, Waters and Whiteb). The mean temperatures deduced for the three Donoussa rock samples are 626°C ± 42°C, 669°C ± 35°C, and 687°C ± 56°C (average 660°C ± 45°C; 2SD) (Table 2). Errors in temperature due to an uncertainty on each mineral analysis of 0.1‰ (2SD) are ± 25°C and only partially account for the range of temperatures given by the three samples. The temperatures are compatible with the other high-temperature estimates for D2 metamorphism on Donoussa.
The main Al-silicate mineral in the lenses is andalusite (Figure 2c). Although such temperatures could represent the retrograde history evident in metapelites, several factors argue against this interpretation. At the 3–4 kbar pressures defined by garnet cores, these temperatures are within the P–T field of sillimanite rather than andalusite (Figure 6a). Also, the petrographic history given in section 2 indicates the boudins and lenses are rotated into parallelism with the S2 schistosity. Together with fibrous sillimanite overgrowth and the rotation of andalusite porphyroblasts in metapelites, the petrography indicates that andalusite grew prior to or during shearing movement along the S2 plane. A retrograde temperature history of andalusite growth does not appear viable.
Several field-based studies suggest that there are no equilibrium fractionation factor differences among Al-silicate polymorphs for oxygen isotopes (Cavoisie et al. Reference Cavosie, Sharp and Selverstone2002; Putlitz et al. Reference Putlitz, Valley, Matthews and Katzir2002; Allaz et al. Reference Allaz, Maeder, Vannay and Steck2005). Thus, although the Sharp (Reference Sharp1995) equation is primarily based on q-ky-sill fractionation data in metamorphic assemblages, its calibration should apply to andalusite-bearing rocks. Quartz-andalusite calibrations that substantially utilize semi-empirical incremental theory calculations (Zheng, Reference Zheng1993; Hofbauer et al. Reference Hofbauer, Hoernes and Fiorentini1994; Vho et al. Reference Vho, Lanari and Rubatto2019) propose measurable fractionations among Al-polymorphs. However, the calcite-corundum study on Naxos (Turnier et al. Reference Turnier, Katzir, Kitajima, Orland, Spicuzza and Valley2020) and the q-ky-sill studies (Cavoise et al. Reference Cavosie, Sharp and Selverstone2002; Putlitz et al. Reference Putlitz, Valley, Matthews and Katzir2002; Larson and Sharp, Reference Larson and Sharp2005) give petrologically sound temperatures.
A possible explanation of the high RAM temperatures is that prograde P–T conditions moved into the sillimanite field with oxygen isotope exchange but without significant mineralogical change of andalusite, apart from the observed fibrous sillimanite overgrowth in one sample (KKD140). Such kinetic overshoot has been recognized in petrological studies on the prograde andalusite-sillimanite transition (e.g., Waters and Lovegrove, Reference Waters and Lovegrove2002; Pattison and Spear, Reference Pattison and Spear2018; Pattison and Forshaw, Reference Pattison and Forshaw2025). The presence of small amounts of Fe3+ and Mn3+ in the aluminosilicate lattice may also lead to an increase of up to 25°C–50°C in the stability of iron-saturated andalusite (Grambling and Williams, Reference Grambling and Williams1985). Such cation substitutions in andalusite are supported by the relatively high fO2 and MnO2 values estimated here.
5.e. Metamorphic P–T paths
The thermobarometric P–T and phase diagram calculations show that the sillimanite-grade metamorphism occurred at temperatures above 600°C with a peak temperature of 660°C–690°C. Retrograde cooling and decompression in the metapelites is indicated by the lower garnet rim Mg/Fe values and by the increase in the measured and calculated XMngt contour values with lower P–T (Figures 6b and 9b). Sphene appears in Ca amphibolite rocks during retrograde cooling below 540°C–560°C.
A schematic clockwise P–T path for the garnet-bearing D2 metapelites is given by the red arrows in Figure 11a. It is drawn to be compatible with the upper boundary of the garnet-bearing B(I) metapelite assemblage at bulk-rock MnO = 0.10 mol% (i.e., below the staurolite-in reactions) and the garnet core and rim P–T values for the high-temperature cooling paths. The P–T path commences within the andalusite field but below staurolite equilibria. A second diagram (Figure 11b) projects the same P–T path for the B(I) garnet assemblages onto the P–T pseudosection at bulk-rock MnO = 0.03 mol%. This latter diagram represents P–T conditions where the garnet-absent metapelite B(II) assemblages develop. Notably, both the B(I) and B(II) metapelite assemblages can form at the same P–T conditions in the same rock, depending on bulk-rock MnO. Their formation does not necessarily require increasing metamorphic temperature.
P–T path of Donoussa metamorphism is compatible with phase diagram calculations for sample analysis KD284A and cation exchange thermobarometry of metapelites, with the principal constraint that prograde conditions did not enter the staurolite field. (a) P–T path for the garnet-bearing rocks at bulk-rock MnO = 0.10 mol%. (b) P–T path shown in Figure 11a plotted for lower Mn abundance rocks (bulk-rock MnO = 0.03 mol%) in which and + bi + mu & sill + bi + mu assemblages are stable without garnet. The metapelite P–T path (black arrowed line labelled Martha25) proposed by Martha et al. (Reference Martha, Xypolias, Cheng, Dörr, Gerdes, Hezel, Kutzschbach, Millonig, Schmeling, Marschall, Müller and Zulauf2025) and the placement of their three tectonometamorphic events (D1, D2, D3) are shown for comparison in each figure.

The P–T path inferred by Martha et al. (Reference Martha, Xypolias, Cheng, Dörr, Gerdes, Hezel, Kutzschbach, Millonig, Schmeling, Marschall, Müller and Zulauf2025) is shown by the blue clockwise loops (labelled Martha25) in Figure 11a,b. The path was based on an equilibrium assemblage diagram in the MnNCKFMASHT system (Pattison and Spear, Reference Pattison and Spear2018) for the Nelson contact metamorphic aureole of British Columbia (Pattison and Tinkham, Reference Pattison and Tinkham2009). The Martha25path envisages prograde decompression during the D1 to D2 transition, which. successively passes with increasing T through the equilibrium fields (q + pl + mu + ilm): chl + bi; gt + chl + bi; gt + st + bi; gt + sill + bi; sill + bi. The order of sequential index mineral formation with increasing grade of metamorphism is thus: chlorite>garnet>staurolite>sillimanite. Constraints on this prograde path come from the P–T pseudosections in this study. The XMn zoning patterns of the garnets from sample KD284 show no evidence for prograde metamorphism in the gt + chl + bi + mu and gt + st + bi + mu fields but are consistent with decreasing temperature and pressure in the gt+ sill + bi + mu field (Figures 8a and 11a). It follows that either garnet growth in KD284 was restricted to the peak D2 phase or that Mn zoning in garnet from earlier staurolite-grade metamorphism was reset during D2 metamorphism. An additional point is the difference in peak metamorphic temperatures deduced in the two studies. Martha et al. Reference Martha, Xypolias, Cheng, Dörr, Gerdes, Hezel, Kutzschbach, Millonig, Schmeling, Marschall, Müller and Zulauf2025 identify a peak metamorphism at P < 3 kbar based on the assemblage kfs + bi + sill; our assemblages only identify kfs + sill + bi + liquid at P > 3 kbar and T > 650°C (Figure 8a).
The clockwise P–T path inferred here is not a unique solution. High T/P metamorphism is frequently characterized by counterclockwise P–T cycles with peak temperatures and retrograde cooling being approached from the low-pressure side. A counterclockwise path (e.g., reversing the direction of the red arrows in Figure11) would see garnet forming at peak metamorphic P–T and continuing to exchange cations during cooling and decompression. The prograde path would develop in the bi + and/sill + mu fields. Such a loop, though, would not be compatible with the early D1 staurolite assemblage unless the breakdown reaction to andalusite involved clear decompression, as inferred for the Buchan metamorphism (Pattison and Goldsmith, Reference Pattison and Goldsmith2022; section 5.b.3 above).
6. Conclusions
The combination of classical thermodynamic calculations with Thermocalc P–T pseudosections for metapelites and amphibolites and oxygen isotope thermometry allows a valid reconstruction of the P–T path of sillimanite/amphibolite high T/P metamorphism on Donoussa Island, Cyclades. Although metamorphic reactions are identified in mineral assemblages, there is no clear evidence for isograd reactions in the field, and the thermombarometry record is dominated by peak and post-peak metamorphism. The high T/P garnet-sillimanite and sillimanite/andalusite metapelite assemblages can be interpreted in terms of Mn heterogeneity at an outcrop scale, which in turn relates to redox processes during the original sedimentation.
The excellent exposures on Donoussa have made it a key island for studying upper plate high T/P metamorphism and granitoid magmatism accompanying the complex subduction-accretion events of the Alpine Hellenide orogen in the Eastern Mediterranean area. Chief among the questions is the location of the original magmatic arc that generated the high T/P metamorphism. Among the several different tectonometamorphic models that are discussed in this work (including tectonic switching during subduction), the recent studies of Zulauf et al. (Reference Zulauf, Linckens, Beranoaguirre, Gerdes, Krahl, Marschall, Millonig, Neuwirth, Petschick and Xypolias2023, Reference Zulauf, Dorr, Albert, Martha and Xypolias2024) and Martha et al. (Reference Martha, Xypolias, Cheng, Dörr, Gerdes, Hezel, Kutzschbach, Millonig, Schmeling, Marschall, Müller and Zulauf2025) stand out for the extensive rapid slab rollback extending from northernmost Greece to the Asterousia Crystalline complex of Crete. Nevertheless, the slab roll-back driven 2.6 Ma timescale on Donoussa proposed by Martha et al. (Reference Martha, Xypolias, Cheng, Dörr, Gerdes, Hezel, Kutzschbach, Millonig, Schmeling, Marschall, Müller and Zulauf2025) for forearc deposition, high T/P metamorphism and post-peak metamorphism monzogranite intrusion presents a challenge for tectonic models. At a subduction angle of 45 degrees and burial to 4 kbar, this timescale would require a convergence rate of >20mm/a, even assuming zero time for sedimentation and heating. A slightly longer (∼76 to 71.5 Ma) high T/P metamorphic interval that includes sphene formation in amphibolites and calc-silicates during rapid retrograde cooling and decompression accompanying slab roll-back could extend the overall P–T timescale to greater than 4.5 Ma. An additional possibility is that the intrusion of the monzogranite dyke at 73.2 Ma could have occurred in older sedimentary host rocks, followed by erosion, recycling and redeposition, which would include detritus from eroded dykes that were older than the one dated by Martha et al. (Reference Martha, Xypolias, Cheng, Dörr, Gerdes, Hezel, Kutzschbach, Millonig, Schmeling, Marschall, Müller and Zulauf2025). Further zircon dating on different metamorphic rocks will undoubtedly broaden our understanding of the Late Cretaceous high T/P metamorphism on Donoussa and its connection to the Cycladic Blueschist Unit above which it is juxtaposed.
Supplementary material
The supplementary material for this article can be found at https://doi.org/10.1017/S0016756826100661
Acknowledgements
The research for this paper was supported in part by grant # 94–128 from the United States – Israel Binational Science Foundation (BSF), Jerusalem. The thorough critical reviews of an unknown reviewer and DRM Pattison lead to substantial modification of an earlier version of the manuscript.
Competing interests
The authors declare none.



