Multiple-reaction thermobarometry

Here we explain how to configure THERMOCALC for multiple-reaction thermobarometry (“calcmode 2”). If you want to know more about THERMOCALC’s other functions and its history, see the overview. If you want to start using the program for the first time, see Get started with THERMOCALC. Documentation and tutorials can be found on the THERMOCALC help page. Note that the SCOlP and COlP geobarometers by Ziberna, Green & Blundy (2017) can be downloaded as a bundle including the necessary documentation. These barometers can be applied to igneous rocks containing clinopyroxene + orthopyroxene + olivine ± spinel, at crustal pressures.

On this page:

Accessing multiple-reaction thermobarometry mode

In the prefs file tc-prefs.txt, set

calcmode 2

to make THERMOCALC run in multiple-reaction thermobarometry mode.

What does it do?

In multiple-reaction thermobarometry mode, the user must have mineral analyses for an assemblage that is believed to represent the equilibrium assemblage, or a subset of the equilibrium assemblage. The analysed compositions are given to THERMOCALC, which converts them to activities using the HPx-eos, and incorporates them into a set of thermodynamic equilibrium constraints. Each of these constraints acts as an independent estimate of pressure or temperature. Using some assumptions about the uncertainties involved and how they are correlated, THERMOCALC infers one optimal estimate of pressure and/or temperature from the set of constraints, along with an uncertainty, a goodness-of-fit value and other diagnostic parameters.

THERMOCALC 3.50 offers two types of multiple-reaction thermobarometry calculation:

  • average pressure (avP): calculate at multiple temperatures and find the best fit to estimate T
  • average temperature (avT): calculate at multiple temperatures and find the best fit to estimate P.

Another possibility is ∆P or ∆T – relative pressures and temperatures. The user has to run these procedures manually.

General operation notes

avP, avT, avPT

AvP, avT and avPT (Powell & Holland, 1994) are used for estimating absolute values of pressure, temperature or both.

  • As usual, you will need a prefs file, scriptfile, and the axfile containing your chosen HPx-eos set.
  • Calcmode 2 is used interactively. The prompts will guide you through the process.
  • You will need to provide your analysed mineral compositions to THERMOCALC, which will then convert them to activities of end-members.
  • THERMOCALC will ask you if you want to exclude any end-members. You should exclude any ordered end-members in your phases – see the headers to the x-eos in axfiles or -it files.
  • You will be given a choice between avP, avT and avPT modes. We tend to be suspicious of avPT. The question is, can the compositional data actually contain enough information to constrain both P and T uniquely, even if the x-eos were perfect? So we suggest calculating avP at multiple temperatures, choosing the temperature at you get the best sigfit/smallest uncertainty, or alternatively calculating avT at multiple pressures. There may not be much pressure or temperature information present, as a result of weaknesses in the x-eos, a shortage of the necessary information in the compositional data, or both.
  • Using avP mode as an example, THERMOCALC:
    • assembles an independent set of reactions among end-members in the various model phases;
    • uses each of those reactions separately to calculate a value for P at the specified T;
    • considers the correlated uncertainties on the positions and slopes of each reaction, which are derived from
      • the correlated uncertainties in the enthalpies of the end-members, obtained from fitting the Holland & Powell dataset,
      • supposedly independent uncertainties in the activities of the end-members, derived from uncertainties in the x-eos and in the compositional analysis;
    • adjusts the positions and slopes of the reactions, respecting the uncertainties and correlations, until all of the reactions coincided at an “average pressure”.
  • THERMOCALC provides a lot of output related to the statistics of this process, and it’s important to analyse it carefully in order to gain geological insight. There is an on-screen explanation of the statistical parameters involved, but for more help, see the references in the Documentation section below, especially Optimal geothermometry and geobarometry.
  • You may want to try to improve your results by, for example, excluding phases that seem to be causing problems – is this telling you that they are not part of the equilibrium assemblage? Or, you might want to remove, say, ferric end-members that have particularly uncertain activities.
  • Never use results that are calculated with an incomplete set of reactions!

∆P, ∆T, ∆PT

∆P, ∆T and ∆PT (Worley & Powell, 2000) can be applied where one wants to find the difference in P (or T) for two (or more) samples with the same assemblage. Much smaller uncertainties are involved in estimating this difference than in estimating absolute values, because many of the uncertainties affect both calculations in the same way, and can therefore be ignored. Taking ∆P on two samples as an example:

  1. Do avP on each sample; in each case, check that the sigfit value is okay, and if so, record the value of avP.
  2. Repeat, but add the script justrelu (“just relative uncertainties”) to the scriptfile; in each case, check that this produces similar P to before, and record the sd(P) values (ignore the horrifying sigfit values).
  3. You can now find ∆P from the difference in pressures in (1), and find the uncertainty on it by adding the standard deviations in (2).

Specifying thermo data

As usual, the user must specify the thermodynamic input they want to use – i.e. HPx-eos family and accompanying version of the Holland & Powell dataset. See the THERMOCALC download guide.

Specialist scripts

In addition to the mineral composition scripts, and scripts related to ∆P and ∆T above, THERMOCALC 3.50 (Dec 2020 onwards) offers numerous scripts for controlling an avP or avT calculation. Key scripts are shown in the table, while additional scripts to facilitate automation of calculations can be found in the sample script file below.

diagramPT <minP maxP minT maxT>Provides the P-T window within which THERMOCALC should look for reactions.
calcT <T(ºC)> or <minT maxT interval>
in an avP calculation, or
calcP <P(kbar)> or <minP maxP interval>
in an avT calculation
The temperature(s) at which to calculate avP, or pressure(s) at which to calculate avT.
ptforaxe <P(kbar) T(ºC)>An initial best guess of P and T, used to make an initial calculation of the activities of end-members.
maxsdp <standard deviation (kbar)>
in an avP calculation, or
maxsdt <standard deviation (ºC)>
in an avT calculation
If no isr block is set (see next entry), THERMOCALC seeks its own independent set of reactions. The maxsdp and maxsdt scripts specify the maximum uncertainties in P or T that a reaction can have, in order to be included in this set.
isr, oisrThe independent set of reactions can be specified via a block of isr scripts: see sample scriptfile below. Alternatively, if THERMOCALC has found its own reactions, oisr instructs it to print the corresponding isr block in the logfile, to be used in future scripts.

Mineral compositions and scripts

Similarly to pseudosection calculation mode, multiple-reaction thermobarometry calculations require the user to set mineral compositions via scripts. However, for pseudosection calculations, these scripted mineral compositions are starting guesses for the non-linear equation solver to use in calculating mineral compositions. By contrast, in calcmode 2, these scripts are the way that the user tells THERMOCALC about the analysed compositions of their minerals . To make this clear, you can use the script xyz <value> instead of xyzguess <value> if you prefer (THERMOCALC will handle both identically).

The user must then convert their analytical mineral compositions into the variables used by THERMOCALC. Start by recalculating your mineral formulae and assigning ions to sites, making appropriate assumptions about ferrous vs ferric iron. Strictly, you could then optimise your analytical site occupancies subject to constraints of charge balance and site totals, as in Powell & Holland (2008); however, Droop’s rules or equivalent are also available. The Teaching Phase Equilibria pages at SERC Carleton College offer some very useful mineral recalculation spreadsheets, contributed by several authors.

THERMOCALC’s compositional variables are written in terms of site fractions of ions on mixing sites. For each x-eos, you can find definitions of the variables, and of the expected oxygen bases and mixing sites, in the -it file or axfile. It may sometimes be tricky to reconcile the mixing sites expected in the x-eos with those in the spreadsheet, in which case feel free to ask for help and we will compile a list of the necessary conversions. You do not need to find values for order variables, which start with Q. These can usually be set to a small value such as 0.001, except in the case of amphibole (hb, gl, act), for which you should use the amphibole Q1Q2 spreadsheet:


All on-screen output is recorded in the file tc-<project>-o.txt.

THERMOCALC generates extensive statistical output, and the user is strongly recommended to engage with it. Some explanation of the output is printed on-screen, but for more help, see the references in the Documentation section, especially Optimal geothermometry and geobarometry.


Sample scriptfile for avP calculations:

Some important references are (if you’re unable to access them, contact a member of the development team):

Documentation from the old Mainz website. This shows a nice example of using the diagnostics to optimise the avPT estimate. (Note that we no longer process analytical data using the AX program, since the user can’t control the thermodynamic data input):

You may also find something useful in the following material, from the 2006 THERMOCALC shortcourse in São Paulo run by RP and Richard White: