Journal of Metamorphic Geology

The thermodynamic properties of 154 mineral end‐members, 13 silicate liquid end‐members and 22 aqueous fluid species are presented in a revised and updated data set. The use of a temperature‐dependent thermal expansion and bulk modulus, and the use of high‐pressure equations of state for solids and fluids, allows calculation of mineral–fluid equilibria to 100 kbar pressure or higher. A pressure‐dependent Landau model for order–disorder permits extension of disordering transitions to high pressures, and, in particular, allows the alpha–beta quartz transition to be handled more satisfactorily. Several melt end‐members have been included to enable calculation of simple phase equilibria and as a first stage in developing melt mixing models in NCKFMASH. The simple aqueous species density model has been extended to enable speciation calculations and mineral solubility determination involving minerals and aqueous species at high temperatures and pressures. The data set has also been improved by incorporation of many new phase equilibrium constraints, calorimetric studies and new measurements of molar volume, thermal expansion and compressibility. This has led to a significant improvement in the level of agreement with the available experimental phase equilibria, and to greater flexibility in calculation of complex mineral equilibria. It is also shown that there is very good agreement between the data set and the most recent available calorimetric data.

Protolith zircon in high‐grade metagranitoids from Queensland, Australia, partially recrystallized during granulite‐grade metamorphism. We describe the zircon in detail using integrated cathodoluminescence, U–Pb isotope, trace element and electron backscatter diffraction pattern (EBSP) analyses. Primary igneous oscillatory zoning is partially modified or obliterated in areas within single crystals, but is well preserved in other areas. A variety of secondary internal structures are observed, with large areas of transgressive recrystallized zircon usually dominant. Associated with these areas are recrystallization margins, interpreted to be recrystallization fronts, that have conformable boundaries with transgressive recrystallized areas, but contrasting cathodoluminescence and trace element chemistry. Trace element analyses of primary and secondary structures provide compelling evidence for closed‐system solid‐state recrystallization. By this process, trace elements in the protolith zircon are purged during recrystallization and partitioned between the enriched recrystallization front and depleted recrystallized areas. However, recrystallization is not always efficient, often leaving a ‘memory’ of the protolith trace element and isotopic composition. This results in the measurement of ‘mixed’ U–Pb isotope ages. Nonetheless, the age of metamorphism has been determined. A correlation between apparent age and Th/U ratio is indicative of incomplete re‐setting by partial recrystallization. Recrystallization is shown to probably not significantly affect Lu–Hf ages. Recrystallization has been determined by textural and trace element analysis and EBSP data not to have proceeded by sub‐grain rotation or local dissolution/re‐precipitation, but probably by grain‐boundary migration and defect diffusion. The formation of metamorphic zircon by solid‐state recrystallization is probably common to high‐grade terranes worldwide. The recognition of this process of formation is essential for correct interpretation of zircon‐derived U–Pb ages and subsequent tectonic models.

Abstract This paper provides methods and a description of a Pascal computer program, thermocalc, for various thermodynamic calculations using the thermodynamic dataset presented in earlier papers in this series (Holland & Powell, 1985; Powell & Holland, 1985). The dataset involves uncertainties on the thermodynamic parameters and therefore allows uncertainties to be calculated on results, for example in geothermometry and geobarometry. Recommendations are made for the uncertainties on activities to be used in calculations on rocks, particular emphasis being placed on preventing underestimates of these uncertainties at small mole fractions. Apposite examples of phase diagram and rock calculations are presented with ouput from thermocalc, demonstrating the utility of the program. Of the rock calculations, the most valuable are considered to be those involving simultaneous combination ‘least squares’of calculated conditions for a set of reactions applicable to a rock. This set of reactions involves the independent reactions which can be written between the end‐members in the minerals in a rock and in the thermodynamic dataset. In contrast to an approach based on specific geothermometers and geobarometers, this approach maximizes the benefit of having an internally consistent thermodynamic dataset. thermocalc is available in IBM PC and Mac versions, from Roger Powell for A$25 or Tim Holland for £10 per version.

We present, as a progress report, a revised and much enlarged version of the thermodynamic dataset given earlier (Holland & Powell, 1985). This new set includes data for 123 mineral and fluid end‐members made consistent with over 200 P–T–X_CO2–f_O2 phase equilibrium experiments. Several improvements and advances have been made, in addition to the increased coverage of mineral phases: the data are now presented in three groups ranked according to reliability; a large number of iron‐bearing phases has been included through experimental and, in some cases, natural Fe:Mg partitioning data; H₂O and CO₂ contents of cordierites are accounted for with the solution model of Kurepin (1985); simple Landau theory is used to model lambda anomalies in heat capacity and the Al/Si order–disorder behaviour in some silicates, and Tschermak‐substituted end‐members have been derived for iron and magnesium end‐members of chlorite, talc, muscovite, biotite, pyroxene and amphibole.

For the subset of data which overlap those of Berman (1988), it is encouraging to find both (1) very substantial agreement between the two sets of thermodynamic data and (2) that the two sets reproduce the phase equilibrium experimental brackets to a very similar degree of accuracy. The main differences in the two datasets involve size (123 as compared to 67 end‐members), the methods used in data reduction (least squares as compared to linear programming), and the provision for estimation of uncertainties with this dataset. For calculations on mineral assemblages in rocks, we aim to maximize the information available from the dataset, by combining the equilibria from all the reactions which can be written between the end‐members in the minerals. For phase diagram calculations, we calculate the compositions of complex solid solutions (together with P and T) involved in invariant, univariant and divariant assemblages. Moreover we strongly believe in attempting to assess the probable uncertainties in calculated equilibria and hence provide a framework for performing simple error propagation in all calculations in thermocalc, the computer program we offer for an effective use of the dataset and the calculation methods we advocate.

Metasedimentary rocks generally contain carbonaceous material (CM) deriving from the evolution of organic matter originally present in the host sedimentary rock. During metamorphic processes, this organic matter is progressively transformed into graphite s.s. and the degree of organisation of CM is known as a reliable indicator of metamorphic grade. In this study, the degree of organisation of CM was systematically characterised by Raman microspectroscopy across several Mesozoic and Cenozoic reference metamorphic belts. This degree of organisation, including within‐sample heterogeneity, was quantified by the relative area of the defect band (R2 ratio). The results from the Schistes Lustrés (Western Alps) and Sanbagawa (Japan) cross‐sections show that (1) even through simple visual inspection, changes in the CM Raman spectrum appear sensitive to variations of metamorphic grade, (2) there is an excellent agreement between the R2 values calculated for the two sections when considering samples with an equivalent metamorphic grade, and (3) the evolution of the R2 ratio with metamorphic grade is controlled by temperature (T). Along the Tinos cross‐section (Greece), which is characterised by a strong gradient of greenschist facies overprint on eclogite facies rocks, the R2 ratio is nearly constant. Consequently, the degree of organisation of CM is not affected by the retrogression and records peak metamorphic conditions. More generally, analysis of 54 samples representative of high‐temperature, low‐pressure to high‐pressure, low‐temperature metamorphic gradients shows that there is a linear correlation between the R2 ratio and the peak temperature [T(°C) = −445 R2 + 641], whatever the metamorphic gradient and, probably, the organic precursor. The Raman spectrum of CM can therefore be used as a geothermometer of the maximum temperature conditions reached during regional metamorphism. Temperature can be estimated to ± 50 °C in the range 330–650 °C. A few technical indications are given for optimal application.

Improved activity–composition relationships for biotite, garnet and silicate liquid are used to construct updated P–T grids and pseudosections for high‐grade metapelites. The biotite model involves Ti charge‐balanced by hydrogen deprotonation on the hydroxyl site, following the substitution , where HD represents the hydroxyl site. Relative to equivalent biotite‐breakdown melting reactions in P–T grids in K₂O–FeO–MgO–Al₂O₃–SiO₂–H₂O (KFMASH), those in K₂O–FeO–MgO–Al₂O₃–SiO₂–H₂O–TiO₂–O₂ (KFMASHTO) occur at temperatures close to 50 °C higher. A further consequence of the updated activity models is that spinel‐bearing equilibria occur to higher temperature and higher pressure. In contrast, the addition of Na₂O and CaO to KFMASH to make the Na₂O–CaO–K₂O–FeO–MgO–Al₂O₃–SiO₂–H₂O (NCKFMASH) system lowers key biotite‐breakdown melting reactions in P–T space relative to KFMASH. Combination of the KFMASHTO and NCKFMASH systems to make Na₂O–CaO–K₂O–FeO–MgO–Al₂O₃–SiO₂–H₂O–TiO₂–O₂ (NCKFMASHTO) results in key biotite‐breakdown melting reactions occurring at temperatures intermediate between those in KFMASHTO and those in NCKFMASH. Given such differences, the choice of model system will be critical to inferred P–T conditions in the application of mineral equilibria modelling to rocks. Further, pseudosections constructed in KFMASH, NCKFMASH and NCKFMASHTO for several representative rock compositions show substantial differences not only in the P–T conditions of key metamorphic assemblages but also overall topology, with the calculations in NCKFMASHTO more reliably reflecting equilibria in rocks. Application of mineral equilibria modelling to rocks should be undertaken in the most comprehensive system possible, if reliable quantitative P–T information is to be derived.

Mineral equilibria calculations in the system K₂O–FeO–MgO–Al₂O₃–SiO₂–H₂O–TiO₂–Fe₂O₃ (KFMASHTO) using thermocalc and its internally consistent thermodynamic dataset constrain the effect of TiO₂ and Fe₂O₃ on greenschist and amphibolite facies mineral equilibria in metapelites. The end‐member data and activity–composition relationships for biotite and chloritoid, calibrated with natural rock data, and activity–composition data for garnet, calibrated using experimental data, provide new constraints on the effects of TiO₂ and Fe₂O₃ on the stability of these minerals. Thermodynamic models for ilmenite–hematite and magnetite–ulvospinel solid solutions accounting for order–disorder in these phases allow the distribution of TiO₂ and Fe₂O₃ between oxide minerals and silicate minerals to be calculated. The calculations indicate that small to moderate amounts of TiO₂ and Fe₂O₃ in typical metapelitic bulk compositions have little effect on silicate mineral equilibria in metapelites at greenschist to amphibolite facies, compared with those calculated in KFMASH. The addition of large amounts of TiO₂ to typical pelitic bulk compositions has little effect on the stability of silicate assemblages; in contrast, rocks rich in Fe₂O₃ develop a markedly different metamorphic succession from that of common Barrovian sequences. In particular, Fe₂O₃‐rich metapelites show a marked reduction in the stability fields of staurolite and garnet to higher pressures, in comparison to those predicted by KFMASH grids.

New activity–composition (a–x) relations for minerals commonly occurring in metapelites are presented for use with the internally consistent thermodynamic dataset of Holland & Powell (, Journal of Metamorphic Geology, 29, 333–383). The a–x relations include a broader consideration of Fe₂O₃ in minerals, changes to the formalism of several phases and order–disorder in all ferromagnesian minerals where Fe–Mg mixing occurs on multiple sites. The a–x relations for chlorite, biotite, garnet, chloritoid, staurolite, cordierite, orthopyroxene, muscovite, paragonite and margarite have been substantially reparameterized using the approach outlined in the companion paper in this issue. For the first time, the entire set of a–x relations for the common ferromagnesian minerals in metapelitic rocks is parameterized simultaneously, with attention paid to ensuring that they can be used together to calculate phase diagrams of geologically appropriate topology. The a–x relations developed are for use in the Na₂O–CaO–K₂O–FeO–MgO–Al₂O₃–SiO₂–H₂O–TiO₂–O₂ (NCKFMASHTO) system for both subsolidus and suprasolidus conditions. Petrogenetic grids in KFMASH and KFMASHTO are similar in topology to those produced with earlier end‐member datasets and a–x relations, but with some notable differences. In particular, in subsolidus equilibria, the FeO/(FeO + MgO) of garnet is now greater than in coexisting staurolite, bringing a number of key staurolite‐bearing equilibria into better agreement with inferences from field and petrographic observations. Furthermore, the addition of Fe³⁺ and Ti to a number of silicate phases allows more plausible equilibria to be calculated in relevant systems. Pseudosections calculated with the new a–x relations are also topologically similar to equivalent diagrams using earlier a–x relations, although with many low variance fields shifting in P–T space to somewhat lower pressure conditions.

Phase diagrams involving solid solutions are calculated by solving sets of non‐linear equations. In calculating P–T  projections and compatibility diagrams, the equations used for each equilibrium are the equilibrium relationships for an independent set of reactions between the end‐members of the phases in the equilibrium. Invariant points and univariant lines in P–T  projections can be calculated directly, as can coordinates in compatibility diagrams. In calculating P–T  and T–x/P–x pseudosections – diagrams drawn for particular bulk compositions – the equilibrium relationship equations are augmented by mass balance equations. Lines in pseudosections, where the mode of one phase in the lower variance equilibrium is zero, and points, where the modes of two phases are zero, can then be calculated directly. The software, THERMOCALC, allows the calculation of these and a range of other types of phase diagram. Examples of phase diagrams and phase diagram movies, with instructions for their production, along with the THERMOCALC input and output files, and the Mathematica^TM functions for assembling them, are presented in this paper, partly in hard copy and partly on the JMG web sites (http://www.gly.bris.ac.uk/www/jmg/jmg.html, or equivalent Australian or USA sites).