Code for the paper: I.N. Bindeman, O.E. Melnik The rises and falls of zirconium isotopes during zircon crystallisation.

Geochem. Persp. Let. (2022) 24, XX–XX | https://doi.org/10.7185/geochemlet.2241

Copy two folders (functions and Results) and one file (main.m) into the same folder. File main.m runs the program by using subroutines (functions) in the functions folder and results of the simulation appear in the Results folder as a text file and a series of graphs (see Fig. S-2 in supplementary materials). Use lines 14–24 in main.m to change input parameters: Zr concentration in the melt, water content (XH2O), size of the cell (length), and ratio of diffusion coefficients of a other trace elements of interest to Zr in the melt. Lines 32 (Tend) and 33(tfin) set up the final temperature of crystallisation (in °C) and time in years. Cooling and crystallisation are set from T of zircon saturation as a linear T-t function along the phase diagram. Phase diagram (Temperature-% crystals) is set using Piwinski (1968) with a concave down function to account for eutectic crystallisation, a much smaller T change for the last 30 % of crystallization. Longer crystallisation durations may require increasing the number of time steps nt = 500 (line 14) to a greater number.

You can choose the trace element by name at line 41 (inside the functions folder the function kdHf.m allows to use partition coefficients of trace elements taken from Melnik and Bindeman, 2018).

File TemperatureHistory.m controls the temperature-time history, which is linear in the present model; users can construct and insert their own, for example the heat propagation equations from Bindeman and Melnik (2016) for natural Temperature-time histories in and around cooling intrusions. File ZrSaturation.m chooses the zircon saturation equation among published models: Watson and Harrison (1983) is used in the present model, but users can uncomment the Boehnke et al. (2013) line to use it instead. DiffusionCoefficient.m parameterises diffusion coefficients of Zr as a function of temperature and water content in the melt from Zhang et al. (2010).