Experimental evidence for the preservation of U-Pb isotope ratios in mantle-recycled crustal zircon grains

Zircon of crustal origin found in mantle-derived rocks is of great interest because of the information it may provide about crust recycling and mantle dynamics. Consideration of this requires understanding of how mantle temperatures, notably higher than zircon crystallization temperatures, affected the recycled zircon grains, particularly their isotopic clocks. Since Pb2+ diffuses faster than U4+ and Th+4, it is generally believed that recycled zircon grains lose all radiogenic Pb after a few million years, thus limiting the time range over which they can be detected. Nonetheless, this might not be the case for zircon included in mantle minerals with low Pb2+ diffusivity and partitioning such as olivine and orthopyroxene because these may act as zircon sealants. Annealing experiments with natural zircon embedded in cristobalite (an effective zircon sealant) show that zircon grains do not lose Pb to their surroundings, although they may lose some Pb to molten inclusions. Diffusion tends to homogenize the Pb concentration in each grain changing the U-Pb and Th-Pb isotope ratios proportionally to the initial 206Pb, 207Pb and 208Pb concentration gradients (no gradient-no change) but in most cases the original age is still recognizable. It seems, therefore, that recycled crustal zircon grains can be detected, and even accurately dated, no matter how long they have dwelled in the mantle.

during transport into the mantle affects the zircon U-Th-Pb isotope ratios. In nature U and Th ions are always tetravalent, but Pb has two oxidation states, 4+ and 2+, the latter forming the most stable compound because of the "inert-pair effect" that affects the heaviest elements of the B subgroups. The diffusivity of Pb 4+ is similar to U 4+ and Th 4+ but Pb 2+ diffuses much faster. Accordingly, it is essential to determine the Pb oxidation state in zircon. Watson et al. 19 proposed that it is Pb 2+ but Kramers et al. 20 , to explain the scarcity of observable strong diffusion effects in pre-magmatic zircons (see discussion in 21 ), suggested that it is tetravalent. The latter, however, found no support from XANES measurements because these revealed that Pb 2+ is far more abundant than Pb 4+ 22 . It seems, therefore, that either radiogenic Pb is mostly generated as Pb 2+ or, if some is generated as Pb 4+ , it is reduced to a divalent state due to the strong reducing environment of the zircon 23 .
Accepting this, Arrhenius equations [24][25][26] indicate that reheating at T > = 1000 °C increases the Pb 2+ diffusion coefficient to values at which it may easily migrate outwards from the zircon crystal (D ≈ 3.37 * 10 −24 m 2 s −1 at 1000 °C; D ≈ 6.2 * 10 −20 m 2 s −1 at 1300 °C). The U 4+ and Th 4+ diffusion coefficients are about 5 to 7 orders of magnitude lower, thus implying that diffusion can substantially disturb the U-Th-Pb isotope ratios of zircon xenocrysts, especially those recycled through the mantle because of the long timescales and high temperatures involved.
The thermally enhanced Pb 2+ diffusion may produce two extreme effects depending on the diffusive behavior of the zircon environment: (i) progressive loss of radiogenic Pb from the zircon to the outside ultimately resetting the U-Th-Pb geochronometers (ii) homogenization of the Pb concentration throughout the crystal with no Pb lost from the zircon grain. Most researchers only consider the first effect. Pilot et al. 1 , for example, proposed that under shallow mantle conditions (~1000 °C) zircon grains would totally lose their radiogenic Pb in less than 150 million years. This would set an upper limit for mantle residence time during which older-than-host zircon could be detected. Moreover, considering that temperature increases with depth and diffusivity increases exponentially with the temperature (Fig. 1) prevalence of the Pb loss effect implies that crustal zircon grains would not be detected if they were carried to greater depths with temperatures ≥1200 °C (90 to 130 km depending on the geotherm). In these conditions Pb would be lost in less than 0.05 Myr from crystals with sizes similar to the grains found in common rocks. Yet, the above estimations might be incorrect because they implicitly assume that the material hosting zircon has much higher Pb 2+ diffusivity and partitioning than the zircon itself, and consequently may behave as a Pb 2+ sink. Certainly, a zircon grain surrounded by a high-diffusivity material, such as silicate melt, will be quickly depleted in Pb (or, more probably, dissolved before diffusion may cause perceptible effects, see Figure 14 in 21 ). However, if zircon xenocrysts are surrounded by, or included in, minerals that hardly accept Pb, i.e., with very low Pb partitioning such as olivine 27,28 , orthopyroxene 29 , etc., the Pb loss from zircon grains might be insignificant and age resetting never occurs.
The objective of this paper is to understand whether crustal zircon xenocrysts extracted from mantle-derived rocks can yield reliable U-Th-Pb ages, which is the essential premise to track their origin and understand how they were recycled in the mantle. To this end, we studied the diffusive behavior of the U-Th-Pb system in natural unaltered zircon grains extracted from two rocks: a Variscan tonalite (SAB51) with chemically zoned but isotopically homogeneous zircon precisely dated at 318 ± 2 Ma; and a Cambrian-Ordovician orthogneiss (SAB50) in which most zircon grains have a ca. 480 Ma-485 Ma rim and a ca. 600 Ma-610 Ma core 17 . The former was selected to consider how Pb diffusion affects isotopically homogeneous but chemically heterogenous zircon grains. The latter was selected to consider how diffusion affects zircon grains with different age components.

Experimental Procedure
To determine the effects of diffusion we heated several sets of zircon grains from the two samples mentioned before at different temperature and over different times. These heated grains were then analysed for U-Th-Pb isotopes with the IBERSIMS SHRIMP ion microprobe. The results were compared to unheated pristine zircon grains from the same rocks and it was assumed that the differences resulted solely from diffusion during the experiments.
Given the long time involved in geological processes and the necessarily short duration of experiments, diffusion in laboratory conditions needs the heating to be at much higher temperatures to produce comparable effects (Fig. 1). Unfortunately, the experimental temperature is limited by the thermal decomposition of zircon, which may hinder, or even prevent, work with natural specimens. At atmospheric pressure pure ZrSiO 4 decomposes to SiO 2 + ZrO 2 at 1673 ± 10 °C 30 , but the transformation occurs at much lower temperature, around 1400 °C, if: zircon contains impurities 31 ; or reacts with a silica-subsaturated phase such as forsteritic olivine to form baddeleyite + orthopyroxene 32 . Váczi et al. 33 determined that natural zircon grains heated in open crucibles at 1400 °C for just 96 hours always partially decompose, and that the crucible material exerts a great influence on the zircon behavior. We have found that, after one week, fresh natural zircon grains embedded in forsterite begin decomposing at around 1300 °C, and at 1400 °C if embedded in orthopyroxene or garnet. They, however, may remain stable for months at temperatures up to 1500 °C if embedded in microcrystalline silica (cristobalite at experiment temperatures), probably because the excess SiO 2 suppresses the reaction zircon = baddeleyite + silica. This, therefore, permits higher experimental temperatures to be reached than when the zircon is embedded in olivine, pyroxene or garnet. The microcrystalline silica additionally provides an environment of very low Pb 2+ partitioning, analogous in this regard -and only in this regard-to mantle olivine or orthopyroxene.
Accordingly, we loaded zircon concentrates of about 100 grains from each rock type embedded in pure microcrystalline silica and heated them in small alumina crucibles in a N 2 atmosphere at: 1300 °C for 30, 90, and 180 days; and at 1500 °C for 30 days.The last conditions being roughly equivalent to 1 yr at 1300 °C ( Fig. 1, see Supplementary Material). The microcrystalline silica was prepared by heating smoked silica in a platinum crucible at 1400 °C for 24 h and then grinding the product, mostly cristobalite, in an agate mortar to a grain size less than 10 µm. This ensured maximum contact between zircon grains and the silica host. After each heating experiment, zircon grains were recovered by dissolving the silica in warm HF, first washing it with a saturated H 3 BO 3 solution and then with ultrapure water to eliminate all traces of the acid. The so released zircons grains were dried, mounted in epoxy, documented with the SEM, and analyzed with the SHRIMP as described in the Samples and Methods section. To check whether some Pb had migrated from zircon to the enclosing cristobalite, one run (1400 °C/ 90 days) was vacuum impregnated in epoxy, cut in several slices, and the silica grains surrounding the zircon analysed with a QP-LACIPMS system. In all cases Pb was below detection limit thus indicating that cristobalite acted as an effective Pb sealant. Figure 2A shows the Wetherill concordia plots of pristine and heated zircon grains from the Variscan tonalite SAB51 ( Table 1 in Supplementary Material). Pristine zircon grains define a neat cluster of concordant ages with a mean of 318.6 ± 1.5 Ma. Zircon grains heated at 1300 °C for 30 days perceptibly spread up and down along concordia, this effect increases with time and temperature to reach a maximum in the 1500 °C/30 days experiment. Figure 3 compares these results to two alternative 3D numerical simulations of Pb diffusion for the same T/t conditions, one considers a diffusive surrounding, similar to a melt, and the other a non-diffusive surrounding about three orders of magnitude less diffusive than zircon. The experimental results do not fit well with the simulation that considers a diffusive surrounding, because this will cause tremendous Pb loss resulting in highly discordant and very low U-Pb ages (Fig. 3A). The experimental results do, however, fit well with the simulation that considers a non-diffusive environment (Fig. 3B) except for a small fraction of the experimental data with younger 206 U/ 238 Pb ages between 250 Ma to 300 Ma peaking at about 285 Ma (Fig. 3C). Remarkably, these slightly discordant ages always correspond to spots analyzed in the region of large inclusions of albitic glass (Fig. 3D) presumably representing molten inclusions of albite (melting point ≈1100 °C at 1 atm 34 ), that are very abundant in the pristine SAB51 zircon grains. The inclusions apparently act as local Pb sinks enhancing diffusion and removing Pb from the zircon so lowering the apparent U-Pb age. Figure 2B shows the Wetherill concordia and age distribution plots for the orthogneiss SAB50 zircon grains ( Table 2 in Supplementary Material). Pristine unheated zircon grains show two well-defined modes: a narrow peak at ca. 480 Ma that corresponds to the rims, and a wider peak at ca. 605 Ma that corresponds to the cores. A few grains contain older cores not shown in the figure. On heating, the 480 Ma mode widens but the peak position remains constant. The 605 Ma mode, in contrast, is replaced by several randomly distributed small peaks so that the 605 Ma age becomes unrecognizable. Comparing the 1500 °C/30 days experiment to 3D numerical simulations shows again that the experimental data do not fit with the simulation that considers a diffusive surrounding because this produces highly discordant rim ages peaking at ca. 380 Ma, far from the true mode (Fig. 4C). On the other hand, the simulation with a non-diffusive surrounding closely matches the experiments, except for the peak labeled "Y" in Fig. 4B, which can also be attributed to some Pb loss by diffusion into small bands of, mostly albitic, glass. These glass bands generally appear at the rim-core interfaces (Fig. 4A) filling microcracks that were most likely caused by the non-uniform thermal expansion of the rim-core systems (see 35 and references therein).

Discussion and Conclusions
Extrapolation of the experimental results to mantle recycling of crustal zircon grains requires understanding first whether or not these dwelled in the mantle surrounded by a diffusive or non-diffusive environment and, second, whether the experimental conditions, especially silica embedding, adequately simulate the environment for zircon grains included in mantle minerals.
Regarding the second question, we must consider the following: The few available data about Pb in natural quartz indicate concentrations at the sub-ppb level (e.g. 36 ). Our own LA-ICPMS measurements, even in quartz crystals coexisting with Pb-rich feldspars are always below detection limit (≈1 ppb). This indicates that the incorporation of Pb 2+ impurities into quartz crystals is very limited, as expected from the inverse relation between impurities uptake and ionic radius show by Chakraborty 37 . Unfortunately, the Pb 2+ diffusivity in silica polymorphs has not been studied. From the varied experimental data summarized in 38 it follows that diffusivity decreases with increasing ionic radius for alkaline elements, for example Ca 2+ (Ø ≈ 1.12 Å) diffusivity is orders of magnitude slower than Na + (Ø ≈ 1.18 Å). It would be logical to expect, therefore, that Pb 2+ (Ø ≈ 1.29 Å) would diffuse much slower than Ca 2+ . The combination of very low Pb uptake and low diffusivity causes silica embedding to be an effective Pb sealant. This is supported by the above mentioned fact that the Pb concentration in silica grains surrounding zircon after a 1400 °C/90 days experiments always had Pb below LA-ICPMS detection limit.
Regarding the first question, it should be noted that mafic and ultramafic melts have elevated Pb 2+ diffusivities 39 but they are extremely subsaturated in zircon 40,41 .They are capable, therefore, of dissolving all zircon grains in contact with them almost instantaneously, in fact much before diffusion might cause perceptible effects (see Fig. 14 in 21 ). We must assume, therefore, that crustal zircon recycled in the mantle could only survive if they are shielded from the melt, i.e., included in mantle minerals. Pb 2+ is highly incompatible in olivine, orthopyroxene and garnet, moderately incompatible in spinel and clinopyroxene, and compatible in plagioclase and apatite [42][43][44] . Compared to zircon, Pb 2+ diffusion is faster in clinopyroxene, plagioclase and apatite [45][46][47] , of the same order in orthopyroxene 47 , no data are available for olivine and aluminosilicate garnets. Nevertheless, the very low Pb 2+ partitioning into olivine, orthopyroxene, and garnet suggests that these minerals will likely seal zircon xenocrysts included in them so preventing the migration of Pb out of the zircon grain. In contrast, plagioclase, apatite and, to some extent, spinel and clinopyroxene will not have the same effect.
Plagioclase is limited to the uppermost mantle layer; spinel and apatite are scarce, so these three minerals are unlikely to host many recycled zircon grains. Therefore, with the probable exception of zircon grains included in clinopyroxene, most zircon xenocrysts would be enclosed by minerals that act as a barrier for the diffusive Pb migration loss, i.e., similar in this regard to our experimental conditions.
If this is the case, our results indicate that most zircon grains do not lose radiogenic Pb and thus never fully reset their isotopic clocks. Lead diffusion may occur within each zircon crystal and tends to homogenize the Pb concentration. Even for highly heterogeneous zircon grains such as from the orthogneiss SAB50 homogenisation occurs very quickly, in about 500 years at 1300 °C, and 3 * 10 5 years at 1100 °C (Fig. 1). This effect modifies the U-Pb and Th-Pb isotope ratios to a level that solely depends on the initial U, Th and radiogenic Pb distributions 21 . This situation, however, can be slightly modified if recycled zircon grains contain large inclusions of glass, albite or other easy-to-melt material. If the temperature is high enough for these to fuse, the so-formed liquid may form streaks or pockets that act as finite Pb sinks, thus depleting the surrounding zircon and causing localized slightly discordant U-Pb ages.
It might be argued that our experimental results do not account for the structural state of recycled zircon grains, which may range from crystalline to metamictic. However, as summarized in 48 , moderately radiation-damaged zircon grains completely reorganize their structure at 900 °C in matter of hours or days, becoming highly crystalline before diffusion may cause any effect. This may causes atom-scale reorganization of Pb 49 and even the formation of metallic Pb nanospheres 50 . These effects may certainly degrade the quality of SIMS measurements, but the effects do not seem to be significant. Heavily damaged or amorphous zircon grains, on the other hand, will segregate into their constituent oxides 48 [and references therein] so that they are not found as inherited zircons.
In summary, recycled crustal zircon grains that were initially homogeneous or nearly homogeneous will yield accurate U-Pb crystallization ages, independent of how long they have remained in the mantle. Initially heterogeneous zircon grains, on the other hand, would show a spread of subconcordant ages along concordia, reverse discordias, and spurious high U-Pb ages as shown in the numerical models of 21 .
Our results also indicate that Pb loss because of diffusion imposes no limit on the depth that recycled zircon grains may reach; depth is only limited by the zircon transformation to the high-density polymorph reidite. According to the experimental data of Ono et al. 51 , this occurs at P (GPa) = 8.5 + 0.0017 × (T − 1200) in the temperature range 1100-1900 K and sets the upper stability limit of zircon in the mantle at ~250-270 km if zircon crystals are shielded from any melt present.

Samples and Methods
Samples. For this work we selected two samples containing abundant unaltered zircon with no, or very little, common lead. The first sample, SAB51, is a Variscan tonalite from the Sanabria zone, in NW Iberia, which contains abundant zircon grains. Most of them are long narrow prisms (150-300 µm) terminated by short pyramids, transparent, colorless to pinkish, and contain abundant inclusions, mostly albite. Previous SHRIMP work (authors' unpublished data) revealed that SAB51 zircon grains are isotopically homogeneous, with a U-Pb age of 318 ± 2 Ma, but chemically heterogeneous, with U ranging from about 230 ppm to 1700 ppm.
The second sample, SAB50, is a Cambro-Ordovician orthogneiss from the Ollo de Sapo formation of the Central Iberian Zone 52 . Rocks from this formation have an abnormally high zircon inheritance so that about 90% of their zircon grains consist of ca. 485 Ma rims around ca. 605 Ma, in some cases older, cores [52][53][54][55] . The sample studied here is medium-grained, with a marked augen structure defined by large K-feldspar crystals, and is slightly metamorphosed under lower greenschists facies. It contains abundant zircon grains. Most of the crystals are short stubby prisms terminated at both ends by long pyramids with total length between 100 and 250 µm. The crystals are transparent, pinkish, with small inclusions and optically noticeable ellipsoidal cores occupying most of the central prismatic body.
In both cases, zircon concentrates were separated from about 15 kg of crushed fresh rock sieved between 300 µm and 50 µm. Separation was done by panning, first in water and then in ethanol. The concentrates were purified first with an Nd magnet, and then by hand picking to remove all traces of other minerals.
Batches of 60-80 grains of each rock were repeatedly analyzed with the SHRIMP before starting the diffusion experiments to ensure they always yielded the same age distribution, i.e., the initial heterogeneity of the zircon loads does not mask the effect of heating. Different batches of the tonalite SAB51 zircon grains always yielded the same age. Different batches of the orthogneiss SAB50 always yielded the same age in the rims but some minor differences in the cores caused by variable, but always small, fraction of older-than-Ediacaran cores. This variation is not relevant for the experiments because only those grains with cores younger than 800 Ma were considered.
Diffusion experiments. Volume diffusion in zircon is practically insensitive to pressure, for this reason the experiments were carried out in open crucibles within a N 2 atmosphere using the following procedure: cylindric alumina crucibles with a diameter of 12 mm and height of 14 mm were filled with a 5 mm thick bed of fine-powered microcrystalline silica. This was mainly cristobalite, obtained by heating smoked silica at 1500 °C for 3 days in a platinum container and grinding the resulting product in an agate mortar. About 100 zircon grains were then transferred on top of the silica layer using a needle. The zircon grains were covered with another 5 mm thick silica layer, and the whole load was then pressed with a steel piston. The crucibles were placed in a flat bottom alumina boat and transferred to the center of a tube furnace capable of sustaining 1600 °C for long periods of time. The alumina tube was closed at either extreme with a water-cooled steel gauge that permitted an inert atmosphere, N 2 in this case, to be maintained inside the tube for the full duration of the experiments. These were carried out at 1300 °C over 30 days, 90 days and 180 days, and at 1500 °C over 30 days, which is roughly equivalent to one year at 1300 °C (Fig. 1). Once the experiments were finished, the crucibles were slowly cooled at room temperature and zircon grains released from the silica envelope by dissolving it in warm HF, decanting the fluid and washing the zircon grains first with a saturated H 3 BO 3 solution and then with ultrapure water to eliminate all traces of acid. The so released zircons grains were dried, mounted in epoxy, documented with the SEM, and analyzed with the SHRIMP as described below. In this way, we recovered between 80% to 40% of the zircon load, depending on how the heated grains behaved: some of them were so broken and cracked that they were useless for SHRIMP analyses.
Once mounted and polished, zircon grains were studied by optical and cathodoluminescent imaging, coated with a 10 nm thick gold layer, and analyzed for U-Th-Pb using a SHRIMP IIe/mc ion microprobe at the IBERSIMS laboratory of the CIC-University of Granada, Spain. The SHRIMP U-Th-Pb analytical method generally followed that described in 41 , and detailed at www.ugr.es/ibersims. Uranium concentration was calibrated using the SL13 reference zircon (U: 238 ppm). U/Pb ratios were calibrated using the TEMORA-II reference zircon (417 Ma 56 ) which was measured every 4 unknowns. When required, common lead was corrected from the Scientific REPORTS | (2018) 8:12904 | DOI:10.1038/s41598-018-30934-4 measured 204 Pb/ 206 Pb, using the model of terrestrial Pb evolution of Cumming and Richards 57 . Data reduction was done with the SHRIMPTOOLS software (downloadable from www.ugr.es/~fbea) using the STATA ™ programming language. The SHRIMP results are included as downloadable files.
Numerical simulations. For these we used COMSOL ™ , a commercial finite element software that links a given geometry with multiple partial differential equations 58 . To simulate the diffusion of U 4+ and Pb 2+ we used Cherniak and Watson 25 Arrhenius equations: D Pb = 0.0776 e (−545000/RT) m 2 s −1 and D U = 1.63 e (−726000/RT) m 2 s −1 and resolved Fick's second law: ∂ c i /∂ t = ∇. (D i ∇ c i ); where c i and D i are the concentration and the diffusion coefficient of the species i, respectively, and t is the time. For the "diffusive" models, we assumed the zircon was surrounded by a material with the same diffusivity as a felsic silicate melt as calculated in 59 : D Pb = exp (−9.08−4.32*w−(28512−7900*w)/T) where w represents the percentage of water. For the "non-diffusive" models we assumed a D Pb three orders of magnitude smaller than in zircon.
Calculations were done for two zircon model geometries that simulate the most common geometries of the tonalite SAB51 and orthogneiss SAB50 zircon grains. For SAB51 we used a 170 × 78 × 56 µm prismatic crystal terminated by a short pyramid at one end and concentrically zoned, with 5 µm wide zones. For SAB50 we used a 160 × 75 × 65 µm prismatic crystal terminated by a short pyramid at either end that contains a large ellipsoidal core with dimensions of 100 × 70 × 56 µm. The core had a concentrically zoned 50 × 35 × 25 µm internal domain also with 5 µm wide zones. In both models the ages were calculated by averaging the concentration of 238 U, 206 Pb and 207 Pb inside a 20 × 17 µm ellipse that simulates the surface analyzed by the SHRIMP ion microprobe with a standard 120 µm Kohler aperture. The ellipse is located in an XY working plane that cuts the crystal in two halves, and randomly moved on that plane within the crystal section always including more than one zone. Details of the procedure can be found in 21 . The COMSOL ™ applications, either as COMSOL, Java, Matlab, or VBA model files, are available from F. Bea upon reasonable request.