The Absolute Age and Origin of the Giant Gypsum Geode of Pulpí (Almería, SE Spain)

Abstract

Subaqueous gypsum (CaSO₄·2H₂O) crystals are relatively common in epithermal systems where sulfide ore deposits are present. The Giant Geode of Pulpí (Almería, SE Spain) hosts some of the largest (up to 2 m in length) subaqueous gypsum crystals discovered to date. Here, we present the first U-series ages of its crystals and reconstruct the oxygen and hydrogen isotopic composition (δ¹⁸O and δ²H) of the Pulpí paleo-aquifer from which the crystals formed by using stable isotopes of gypsum hydration water. We successfully dated the onset of gypsum precipitation in the Geode at 164 ± 15 ka. However, the extremely low U concentration (11 ppb) and relatively high detrital Th content (²³⁰Th/²³²Th 3.2) hinder accurate dating other gypsum samples. The δ¹⁸O and δD values of the paleo-aquifer during the growth of the crystals aligned with the local meteoric water line, suggesting that the sulfate-enriched mother solution consisted of meteoric water that recharged the aquifer during that period. The mean isotopic composition of the Pulpí paleo-aquifer (δ¹⁸O = −6.5 ± 0.1‰ and δ²H = −42.3 ± 0.5‰) during the formation of the crystals was similar to the current groundwater in this area (δ¹⁸O = −6.1 ± 0.8‰, δ²H = −42 ± 6‰). The isotopic differences observed in samples collected from distinct locations and in individual crystals were probably related to changes in the isotopic composition of the aquifer, as a consequence of varying climate that impacted on the isotopic composition of rainwater during thousands of years in this region. Our results indicated that subaqueous selenite crystals may be useful for paleo-hydrological reconstructions. However, improving the current analytical techniques for dating gypsum with low U concentrations will be essential to obtain accurate and reliable records from Quaternary gypsum cave crystals in the future.

1. Introduction

The Giant Geode of Pulpí (Almería, SE Spain) is a cavity where the surfaces are completely covered with prismatic selenite crystals (CaSO₄·2H₂O) up to 2 m long [1,2] (Figure 1). The Geode was discovered in 1999 within the Mina Rica of Pilar de Jaravía and is now an important touristic and scientific attraction due to its spectacular characteristics. Indeed, the Giant Geode of Pulpí is currently the only gypsum geode open to the public worldwide [3].

To date, some investigations have addressed the mechanisms that led to precipitation of the massive selenite crystals [1,2,4]. Fluid inclusion microthermometric analyses revealed that subaqueous gypsum precipitation occurred at temperatures from 10 to 35 °C, with the majority of analyses suggesting 20 to 25 °C [4] from a low-salinity solution (from 0.0 to 0.5 wt % eq NaCl; [2]). These results discard gypsum crystallization from evaporated or unevaporated seawater, in addition to precipitation from a high-temperature hydrothermal solution. The possible presence of thermal water at deeper levels did not significantly influence the crystallization of gypsum, but controlled the previous stages of barite, celestine and Fe-oxide precipitation [4]. Sulfur and oxygen isotope analyses suggested that the source of sulfate for the first stages of gypsum precipitation was the dissolution of microcrystalline anhydrite (CaSO₄) and interlayered Triassic gypsum hosted in the carbonate bedrock [2,4]. The ⁸⁷Sr/⁸⁶Sr values in gypsum and celestine in the Mina Rica mine (0.7104 ± 0.0001 and 0.7106 ± 0.0001, respectively, [5]) are significantly higher than in seawater during the Phanerozoic and at present (⁸⁷Sr/⁸⁶Sr 0.709; [6]). Dissolution of Triassic gypsum (⁸⁷Sr/⁸⁶Sr = 0.707 − 0.708; [7]) and Miocene gypsum (⁸⁷Sr/⁸⁶Sr 0.709; [8]), which are abundant in the surrounding bedrock of the Sierra del Aguilón, cannot explain the strontium isotope values observed in gypsum and celestine. Weathering of silicate minerals in the soils has been proposed as a source of strontium for the relatively high ⁸⁷Sr/⁸⁶Sr values observed in the Giant Gypsum Geode [5], since strontium from the edaphic reservoir commonly has high isotope ratios [9].

U-series analyses of a carbonate coating on a gypsum crystal in the shallower mine passages (~30 m above the Giant Geode level) yielded an age of 60 ka. This indicates that gypsum precipitation, at that level, occurred prior to that time [4]. Despite several attempts being made to date the Geode, an absolute age of its gypsum crystals has not been provided so far. In this study, we performed the first successful U-series analyses of the Giant Gypsum Geode. Furthermore, we studied the oxygen and hydrogen stable isotopes (δ¹⁸O and δ²H) of gypsum hydration water to reconstruct the isotopic composition of the groundwater from which the Geode crystallized, in order to: (1) investigate the characteristics of the paleo-aquifer water (i.e., meteoric water recharge vs. presumed marine intrusion) and (2) to evaluate the potential and drawback of subaqueous selenite crystals as archives to reconstruct the δ¹⁸O and δ²H values of paleo-aquifers.

2. Materials and Methods

2.1. Sampling

Nineteen gypsum samples (~1–5 g each) from the Mina Rica mine were collected in December 2019 at different altitudes, ranging from 93 to 128 m above sea level (asl) In the particular case of the Giant Gypsum Geode (98 m asl), we collected a sample from the base of a crystal (MI-15) that was in contact with the underlying celestine and, therefore, represented the first stage of gypsum precipitation at the Geode level. In addition, a 22 cm long gypsum crystal (crystal GE-1), which was removed from the Geode entrance during habilitation works, was taken for analysis. Sample MI-15 and the base of the crystal GE-1 were analyzed for U-series dating (Table 1). We attempted to date 3 additional gypsum samples but analyses were unsuccessful because of the low U content (2 ppb). In addition, 22 subsamples (~50 mg) for stable isotope analyses were extracted along the main growth axis of the crystal GE-1, every centimeter using a dentist tool (Dremel^®) (Table 2). Overall, we obtained the age of 2 samples by U-series dating and 42 samples were analyzed for oxygen (δ¹⁸O) and hydrogen (δ²H) isotopes of gypsum hydration water.

2.2. Analytical Methods

2.2.1. U–Th Dating

U-series analyses were performed at the Institute of Geology and Mineralogy (University of Cologne, Germany), following the procedure outlined in Obert et al. [10]. Gypsum samples were ground to fine powder that were divided into aliquots (3 for sample GE-1 and 4 for sample MI-15). The samples (~150–200 mg) were dissolved in 10 mL 7 M HNO₃ and placed on a hotplate at 120 °C overnight. A mixed ²²⁹Th-²³³U-²³⁶U tracer was added to the solution, which was again put on the hotplate at 120 °C overnight for sample spike equilibration and subsequently dried down. The ²²⁹Th-²³³U-²³⁶U mixed tracer at Cologne is a combination of the IRMM-3636 U double spike with a nominal ²³³U/²³⁶U ratio of 1.01906 [11] and enriched ²²⁹Th spike, calibrated against pure U-Th metal solutions. Accuracy of the tracer calibrations was further ensured by repeated measurements of silicate rock standards that were close to or near secular equilibrium

After sample-spike equlibration, the samples were taken up again in 10 mL 7 M HNO₃, and 2 mL concentrated HCl and 0.5 mL H₂O₂ were added to remove any organic material. Prior to chemical separation, the samples were dried down again and re-dissolved in 10 mL 7 M HNO₃. Uranium–Th separation from the sample matrix was conducted by ion exchange column chemistry during which t samples were passed through a two column chemical separation protocol. We used BioRad AG1-X8 resin (100–200 mesh) in custom-made columns with a column volume (CV) of 1 mL. The resin was cleaned using 6 M HCl, H₂O and 7 M HNO₃ in several steps and conditioned, before the sample solution was passed through the column in 7 M HNO₃ media in individual 1 mL steps, until the whole sample was loaded onto the column. The matrix was rinsed from the resin with two CV of 7 M HNO₃. Subsequently, the Th fraction was collected in 5 individual CV 6 M HCl, then U was eluted using 3 CV 1 M HBr. Both fractions were dried down, and Th was taken up again in 1 mL 7 M HNO₃ and passed through the same AG1-X8 column to purify the Th cut from any residual Ca. Both the purified Th and the U fraction were dried down, treated with H₂O₂ to remove any organic residue from the resin, and re-dissolved in 0.2 M HCl + 0.01 M HF (Th) and 0.14 M HNO₃ (U) for measurement. The U–Th separation of each aliquot was conducted in independent columns and the concentrations and isotope ratios of each aliquot were measured individually (Table 1), allowing to asses, the external reproducibility of the method.

Concentrations and isotope ratios were measured using a Thermo Scientific Neptune MC-ICPMS with a central SEM equipped with a RPQ. We used an Aridus II desolvator system following a standard sample bracketing procedure. For U, CRM112A reference material was used, and for Th, IRMM35 and IRMM36 standards doped with CRM112A were used. Measured U isotope ratios were corrected for mass bias using the known ratio of the ²³³U/²³⁶U double spike. Thorium samples were doped with the CRM112A U standard, using the known ²³⁸U/²³⁵U ratio (137.88) for mass bias correction. The ‘Faraday Cup-SEM’ yield was corrected using the known ²³⁴U/²³⁸U or ²³⁰Th/²³²Th ratio of the respective bracketing standard. We used a mathematical correction in which the U–Th isotopic composition of the detrital contamination was arbitrary estimated, in a manner resembling the approach by Ludwig and Paces [12]. ²³²Th was used as the index for the correction and assumed a typical crustal Th/U ratio of 3.5 (atomic ratio), with ²³⁴U/²³⁸U and ²³⁰Th/²³⁸U activity ratios near secular equilibrium. In our model, we specifically used the initial activities ratios ²³²Th/²³⁸U = 1.21 ± 50%, ²³⁴U/²³⁸U = 1 ± 10% and ²³⁰Th/²³⁸U = 1 ± 10%. Ages and initial (²³⁴U/²³⁸U) ratios were calculated using an iterative approach using the decay constants of Cheng et al. [13] and Jaffey et al. [14]. The respective age uncertainties of each aliquot were derived from Monte Carlo simulations taking the uncertainties of the activity ratios into account. The ages of all the aliquots of each sample (n = 4 for MI-15 and n = 3 for GE-1) were averaged and the error is given as 1 standard deviation (1σ).

2.2.2. Stable Isotopes of Gypsum Hydration Water

The powdered gypsum samples were dried at 45 °C for 24 h, and then placed under vacuum to a pressure of ~5 × 10⁻³ mbar for at least 3 h at room temperature to remove adsorbed water. This low vacuum pumping is effective at removing adsorbed water with no detectable loss of hydration water [15]. The isotopic analyses of gypsum hydration water were conducted at University of Almería (Almería, Spain) using a Heat Induction Module (IM-CRDS, Picarro©) coupled to a Cavity Ring-Down Spectrometer (CRDS, Picarro© L2140i). The water vapor released after heating the powdered samples (~10 mg) by the IM to a final temperature of 250 °C was transferred on-line to the CRDS analyzer that measured the isotopic values (δ¹⁸O and δ²H) of the vapor with a frequency of 1 s. The Gaussian peak generated from the samples was integrated by the Picarro software that used an approach similar to that described in Bauska et al. [16]. The results were standardized to the Vienna Standard Mean Ocean Water (V-SMOW), Standard Light Antarctic Precipitation (SLAP) and Greenland Ice Sheet Precipitation (GISP) by analysing four gypsum standards before each set of 30–35 samples. The gypsum standards were previously calibrated against liquid water using the cryogenic extraction method by Gázquez et al. [15]. Each sample was analyzed 3–4 consecutive times. The mean analytical precision (1SD) of the 3–4 repetitions was 0.1‰ for δ¹⁸O and 0.9‰ for δ²H. The drift of the CRDS instrument was monitored and corrected (if needed) by measuring one of the gypsum standards every 10–12 samples.

The isotopic composition of the gypsum original solution (called paleo-aquifer in this study) was calculated by using isotope fractionation factor at 25 °C, since gypsum in the Geode crystallized at 20–25 °C [4]; 1.0034 for α¹⁸O_gypsum-water and 0.981 for α²H_gypsum-water [17], where:

$${\alpha }^{18}{O}_{gypsum-water}=\frac{{\delta }^{18}{O}_{gypsum}+1000}{ {\delta }^{18}{O}_{water}+1000}$$

and,

$${\alpha }^2{H}_{gypsum-water}=\frac{{\delta }^2{H}_{gypsum}+1000}{ {\delta }^2{H}_{water}+1000}$$

These fractionation factors, especially α¹⁸O_gypsum-water, are very insensitive to temperature and have been previously used to reconstruct the isotopic composition of ancient water bodies [8,18,19,20].

3. Results and Discussion

3.1. The Age of the Giant Geode

The U concentration in samples MI-15 and GE-1 were ~11 ppb and ~3 ppb, respectively, (Table 1). These concentrations are lower than in subaerial gypsum stalactites analyzed in previous studies (i.e., 94–150 ppb U; [20]), but are in agreement with the U content reported in the subaqueous selenite crystals from the caves of Naica (2 to 10 ppb U, [21]). We attempted to analyze other gypsum samples from the Mina Rica mine, but the U contents were even lower (2 ppb), impeding the obtainment of additional ages. As ²³⁰Th/²³²Th activity ratios in the analyzed samples were 20 (2.0 to 3.2), a correction for ²³⁰Th contamination was needed to obtain accurate ages. The corrected ages were 21–56 ka younger than the uncorrected ages (Table 1). The age of sample MI-15, that represents the early stage of gypsum precipitation at the Geode level, was 164 ± 15 ka (n = 4). This indicated that gypsum precipitation in the Giant Geode started during the Marine Isotope Stage 6 (MIS 6) (ca. 191 to ca. 131 ka; [22]). The age of sample GE-1 was 108 ± 87 ka (n = 3), more imprecise because of the lower U content, but within analytical errors with the age of sample MI-15. These results indicated that gypsum samples with U content 10 ppb can be dated by the method used here, however lower U concentrations result in more imprecise and unreliable ages.

Previous studies reported U-series ages of a calcite coating that covers selenite crystals in the upper mine levels (30 m above the Giant Geode level) [4]. Although these authors did not provide information about the origin of this coating (i.e., subaerial/subaqueous) they interpreted that conditions for calcite precipitation prevailed at 60 ka; this may indicate that gypsum precipitation in the Sierra del Aguilón mountain range, and probably in the Giant Geode, ceased before that time. Thus, we concluded that gypsum precipitation in the Giant Geode extended during a maximum period of 100 ka, from 164 ± 15 ka to sometime before 60 ka.

3.2. Paleoclimatic Significance of Subaqueous Gypsum Crystals

The selenite crystals of the Giant Geode and the Mina Rica mine formed in subaqueous conditions [2]. During crystallization, the two molecules of structurally-bound hydration water of gypsum (CaSO₄·2H₂O) are incorporated from the mother solution (i.e., paleo-aquifer). In consequence, the δ¹⁸O and δ²H values of gypsum hydration water (δ¹⁸O_gypsum and δ²H_gypsum hereafter) reflect those of the mother solution (δ¹⁸O_aquifer and δ²H_aquifer hereafter). The oxygen and hydrogen isotope fractionation factors for the gypsum-solution pair (α¹⁸O_gypsum-water and α²H_gypsum-water) are well-known and very insensitive to temperature [17]. Thus, by using the δ¹⁸O_gypsum and δ²H_gypsum the isotopic composition of the paleo-aquifer at the time of gypsum precipitation (δ¹⁸O_aquifer and δ²H_aquifer) can be calculated (Table 2). Here we use this property of gypsum to reconstruct the isotopic composition of the paleo-aquifer of Pulpí during gypsum crystallization.

The δ¹⁸O_aquifer during the formation of the investigated samples ranged from −5.3 to −7.1‰ (mean ± 1 standard deviation = −6.2 ± 0.4‰) and the δ²H_aquifer ranged from −36.6 to −45.0‰ (−40.5 ± 2.2‰). The d-excess of the paleo-aquifer (d-excess = δ²H_aquifer − 8 × δ¹⁸O_aquifer) ranged from 3.4 to 14.4‰ (9.1 ± 2.3‰). These values are close to those of modern rainwater and groundwater in this region (δ¹⁸O = −6.1 ± 0.8‰, δ²H = −42 ± 6‰, d-excess~6.3 ± 4.8‰; [23,24]) (Figure 2). In particular, the samples collected from the Giant Geode (sample MI-15 and 22 subsamples from the crystal GE-1) crystallized from a solution with a δ¹⁸O_aquifer value of −6.0 ± 0.2‰, δ²H_aquifer value of −39.1 ± 1.0‰ and d-excess of 8.5 ± 2.1‰.

The δ¹⁸O_aquifer and δ²H_aquifer pairs are aligned with modern rainwater values in this region (Figure 2). This indicates that seawater (with δ¹⁸O~0‰ and δ²H~0‰) or fluids from deeper magmatic/metamorphic processes (usually with higher δ¹⁸O values; [25]) did not contribute significantly to the solution from which gypsum precipitated. In contrast, the gypsum crystals recorded the isotopic composition of meteoric water seepage that recharged the paleo-aquifer of Sierra del Aguilón at ca. 164 ka and during MIS 6. Indirectly, these results indicate that during the formation of the Geode the Mediterranean Sea level was at a similar or lower position than at present (Figure 3). Otherwise, seawater would have influenced the solution by increasing salinity and isotopic values, which was not seen in our results of gypsum hydration water, or in previous results of microthermometry of fluid inclusions [2]. A relatively low sea level position during the formation of the Geode is supported by Mediterranean Sea level reconstructions from subaqueous carbonate speleothems that recorded a lowstand of ~−60 m relative to present during MIS 6 [26].

Given the uncertainties in the ages of the analyzed samples and the lack of additional dates as a result of the difficulties we found in the U–Th analyses, robust interpretations about the evolution of the Pulpí paleo-aquifer are not possible. However, the range of isotopic values observed between the various groups of samples from Mina Rica (Giant Geode vs. the rest of the mine levels) can be due to the fact that gypsum formed under different climate conditions. The isotopic variability of the aquifer water was likely controlled by changes in the δ¹⁸O and δ²H values of rainwater, as a result of long-term variations in atmospheric temperature and changes in the source of moisture. Indeed, both parameters govern the modern seasonal variations in δ¹⁸O, δ²H and d-excess values of rainwater in this region [27,28], and also in the past [20,29]. Lower δ¹⁸O and δ²H values are expected for meteoric water seepage during glacial (colder) periods [29,30,31]. In contrast, higher isotopic values of the recharge indicate rainwater infiltration to the aquifer during an interglacial (warmer) period (e.g., modern conditions) [29,30,31].

Studies of stable isotopes in fluid inclusions in subaqueous carbonate speleothems from the central Mediterranean reported a change in δ¹⁸O by 1.7‰ in paleo-groundwater waters from the glacial MIS8 to the interglacial MIS5c [30]. Similarly, a δ¹⁸O difference by 2‰ in paleo-groundwater during the transition from MIS6 to MIS5 was found in a subaqueous flowstone from Devil’s Hole (southwest USA) [31]. Considering that the range of reconstructed δ¹⁸O values of the Pulpí paleo-aquifer ranged by 1.7‰ during crystallization of gypsum, it is possible that the gypsum precipitation extended over a long period of time and that the crystals recorded major climate changes.

In particular, the δ¹⁸O values of the Pulpí paleo-aquifer at 164 ± 15 ka during the onset of gypsum precipitation in the Giant Geode (−6.5‰, sample MI-15) are more positive than during the crystallization of crystals stratigraphically younger inside the geode (maximum δ¹⁸O_aquifer of −5.3‰), but more negative than in the fluids that generated gypsum in other mine levels (minimum δ¹⁸O_aquifer of −7.1‰). A paleoclimate reconstruction from carbonate speleothems (i.e., Gitana Cave [22], ~50 km from Pulpí) suggests that dry/cold conditions during MIS6 prevailed between 175–154 ka and peaked at ~157 ka, while there is evidence for a wetter/warmer interstadial around 151 ka. Sea surface temperatures reconstructed from Alboran Sea sediments (core ODP 077; [32]) show that temperature declined gradually by ~4 °C from ~180 to ~140 ka. Altogether, these paleoclimatic reconstructions demonstrate that despite climatic conditions during MIS6 in southern Iberia being generally cold and dry, there were significant variations in rainfall amount and temperature that probably resulted in changes in the δ¹⁸O of rainwater in this region. It is possible that the gypsum crystals in the Aguilón Mountain Range formed during different climatic phases of MIS6, or even during MIS5 in the case of the younger samples with higher δ¹⁸O and δ²H values of gypsum hydration water.

4. Conclusions

U-series dating of gypsum has revealed that the Giant Geode of Pulpí grew during the Upper Pleistocene, starting at least 164 ± 15 ka and probably until 60 ka. The Pulpí paleo-aquifer was occupied by freshwater at that time, rather than by brackish solution or seawater. The δ¹⁸O and δ²H values of the paleo-aquifer during the formation the gypsum crystals varied by 1.7‰ for δ¹⁸O and by ~9‰ for δ²H. This indicates that gypsum in the Aguilón Mountain range crystallized probably over a long period of time, in which a succession of wetter/drier and colder/warmer climatic conditions took place. Climate variability led to changes in the isotopic composition of rainfall and resulted in secular variations of the isotopic composition of the Pulpí paleo-aquifer that were ultimately recorded by the isotopic composition of the gypsum crystals.

Table 1. U–Th results and calculated ages of samples MI-15 (4 repetitions) and GE-1 (3 repetitions) from the Giant Gypsum Geode. We used a mathematical correction in which the U–Th isotopic composition of the detrital contamination was arbitrarily estimated, in a manner resembling the approach by Ludwig and Paces [12]. ²³²Th was used as index for the correction and assumed a typical crustal Th/U ratio of 3.5 (atomic ratio), with ²³⁴U/²³⁸U and ²³⁰Th/²³⁸U activity ratios near secular equilibrium. In our model, we specifically used the initial activity ratios ²³²Th/²³⁸U = 1.21 ± 50%, ²³⁴U/²³⁸U = 1 ± 10% and ²³⁰Th/²³⁸U = 1 ± 10%.

Sample ID	U (ng/g)	Th (ng/g)	Th/U	(230Th/232Th)	2σ ±	(234U/238U) Corrected	2σ ±	(230Th/238U) Corrected	2σ ±	Uncorrected Age (ka)	2σ ±	Corrected Age (ka)	2σ ±	Initial (234U/238U) Corrected	2σ ±
MI-15-A	8.98	5.76	0.641	3.202	0.124	0.788	0.034	0.598	0.019	212	41	179	36	0.649	0.092
MI-15-B	9.84	6.19	0.630	3.131	0.127	0.810	0.038	0.567	0.021	173	25	143	23	0.716	0.074
MI-15-D	11.73	7.99	0.681	3.174	0.079	0.849	0.010	0.637	0.019	196	14	166	14	0.756	0.021
MI-15-E	11.20	8.01	0.715	3.122	0.059	0.868	0.009	0.661	0.014	200	11	169	11	0.788	0.019
M.15-Mean Age												164	±15
GE-1 A	2.11	2.79	1.321	2.032	0.109	0.956	0.060	0.805	0.076	264	111	208	100	0.922	0.135
GE-1 B	2.04	1.29	0.630	2.225	0.176	1.054	0.054	0.340	0.044	63	8	42	7	1.060	0.060
GE-1 C	2.68	1.93	0.722	2.505	0.204	0.979	0.053	0.486	0.061	101	16	75	15	0.974	0.066
GE-1-Mean Age												108	±87

Table 1. U–Th results and calculated ages of samples MI-15 (4 repetitions) and GE-1 (3 repetitions) from the Giant Gypsum Geode. We used a mathematical correction in which the U–Th isotopic composition of the detrital contamination was arbitrarily estimated, in a manner resembling the approach by Ludwig and Paces [12]. ²³²Th was used as index for the correction and assumed a typical crustal Th/U ratio of 3.5 (atomic ratio), with ²³⁴U/²³⁸U and ²³⁰Th/²³⁸U activity ratios near secular equilibrium. In our model, we specifically used the initial activity ratios ²³²Th/²³⁸U = 1.21 ± 50%, ²³⁴U/²³⁸U = 1 ± 10% and ²³⁰Th/²³⁸U = 1 ± 10%.

Sample ID	U (ng/g)	Th (ng/g)	Th/U	(230Th/232Th)	2σ ±	(234U/238U) Corrected	2σ ±	(230Th/238U) Corrected	2σ ±	Uncorrected Age (ka)	2σ ±	Corrected Age (ka)	2σ ±	Initial (234U/238U) Corrected	2σ ±
MI-15-A	8.98	5.76	0.641	3.202	0.124	0.788	0.034	0.598	0.019	212	41	179	36	0.649	0.092
MI-15-B	9.84	6.19	0.630	3.131	0.127	0.810	0.038	0.567	0.021	173	25	143	23	0.716	0.074
MI-15-D	11.73	7.99	0.681	3.174	0.079	0.849	0.010	0.637	0.019	196	14	166	14	0.756	0.021
MI-15-E	11.20	8.01	0.715	3.122	0.059	0.868	0.009	0.661	0.014	200	11	169	11	0.788	0.019
M.15-Mean Age												164	±15
GE-1 A	2.11	2.79	1.321	2.032	0.109	0.956	0.060	0.805	0.076	264	111	208	100	0.922	0.135
GE-1 B	2.04	1.29	0.630	2.225	0.176	1.054	0.054	0.340	0.044	63	8	42	7	1.060	0.060
GE-1 C	2.68	1.93	0.722	2.505	0.204	0.979	0.053	0.486	0.061	101	16	75	15	0.974	0.066
GE-1-Mean Age												108	±87

Table 2. Oxygen and hydrogen stable isotopes of gypsum hydration water (δ¹⁸O_gyp and δ²H_gyp) in samples from the Mina Rica Mine. The isotopic composition of the paleo-aquifer (δ¹⁸O_aquifer, δ²H_aquifer and d-excess) were calculated from δ¹⁸O_gyp and δ²H_gyp using known isotope fractionation factors at 25 °C [17].

Sample ID	Description	Location (m asl)	δ18Ogyp (‰)	1SD	δ18Oaquifer (‰)	δ2Hgyp (‰)	1SD	δ2Haquifer (‰)	d-Excess (‰)
GE-1	Geode crystal (1 cm from the base)	98	−2.2	0.2	−5.6	−57.1	1.1	−38.9	5.7
GE-2	Geode crystal (2 cm from the base)	98	−2.1	0.1	−5.5	−58.4	0.9	−40.1	3.8
GE-3	Geode crystal (3 cm from the base)	98	−1.9	0.1	−5.3	−57.4	1.8	−39.1	3.5
GE-4	Geode crystal (4 cm from the base)	98	−2.5	0.0	−5.9	−57.4	1.2	−39.1	8.2
GE-5	Geode crystal (5 cm from the base)	98	−2.5	0.2	−5.9	−57.6	1.3	−39.3	8.1
GE-6	Geode crystal (6 cm from the base)	98	−2.4	0.2	−5.8	−56.2	1.5	−37.9	8.2
GE-7	Geode crystal (7 cm from the base)	98	−2.5	0.2	−5.9	−56.3	1.9	−38.1	9.3
GE-8	Geode crystal (8 cm from the base)	98	−2.7	0.2	−6.0	−57.1	0.7	−38.8	9.5
GE-9	Geode crystal (9 cm from the base)	98	−2.6	0.1	−6.0	−58.7	0.3	−40.5	7.3
GE-10	Geode crystal (10 cm from the base)	98	−2.7	0.1	−6.1	−57.8	1.8	−39.6	9.4
GE-11	Geode crystal (11 cm from the base)	98	−2.6	0.1	−6.0	−58.6	0.5	−40.4	7.6
GE-12	Geode crystal (12 cm from the base)	98	−2.6	0.1	−6.0	−58.9	0.4	−40.7	7.5
GE-13	Geode crystal (13 cm from the base)	98	−2.6	0.2	−6.0	−57.3	1.1	−39.0	8.8
GE-14	Geode crystal (14 cm from the base)	98	−2.7	0.1	−6.1	−57.8	0.9	−39.6	9.3
GE-15	Geode crystal (15 cm from the base)	98	−2.8	0.1	−6.2	−57.1	0.2	−38.9	10.5
GE-16	Geode crystal (16 cm from the base)	98	−2.8	0.1	−6.1	−57.3	0.9	−39.0	10.1
GE-17	Geode crystal (17 cm from the base)	98	−2.5	0.2	−5.9	−56.9	0.4	−38.7	8.7
GE-18	Geode crystal (18 cm from the base)	98	−2.7	0.1	−6.1	−55.4	2.6	−37.1	11.6
GE-19	Geode crystal (19 cm from the base)	98	−2.8	0.1	−6.1	−55.3	2.2	−37.1	12.0
GE-20	Geode crystal (20 cm from the base)	98	−2.7	0.2	−6.1	−57.9	0.6	−39.7	9.1
GE-21	Geode crystal (21 cm from the base)	98	−2.8	0.2	−6.2	−58.2	0.7	−40.0	9.4
GE-22	Geode crystal (22 cm from the base)	98	−2.8	0.1	−6.2	−57.7	0.8	−39.5	10.2
MI-1	Main mine gallery	128	−3.2	0.2	−6.6	−62.2	0.4	−44.1	8.8
MI-2	Main mine gallery	128	−2.5	0.2	−5.8	−54.9	0.5	−36.6	10.1
MI-3	Main mine gallery	128	−2.9	0.1	−6.3	−60.4	0.6	−42.2	8.1
MI-5	Main mine gallery	128	−3.2	0.3	−6.6	−60.9	0.7	−42.7	10.3
MI-10	Main mine gallery	128	−3.0	0.1	−6.3	−61.8	1.3	−43.6	7.1
MI-11	Main mine gallery	126	−3.7	0.4	−7.1	−60.5	0.9	−42.3	14.4
MI-13	Main mine gallery	113	−2.5	0.1	−5.9	−60.8	0.2	−42.6	4.7
MI-14	Main mine gallery	110	−2.9	0.2	−6.3	−58.4	0.4	−40.2	10.3
MI-15	Giant Gypsum Geode	97	−3.1	0.0	−6.5	−58.9	0.3	−40.7	11.5
MI-16	Main mine gallery	106	−3.5	0.0	−6.9	−61.6	0.4	−43.4	11.9
MI-18	Main mine gallery	102	−3.7	0.2	−7.1	−62.8	0.5	−44.6	12.4
MI-19	Main mine gallery	93	−2.9	0.2	−6.3	−60.3	0.4	−42.1	8.2
MI-20	Main mine gallery	93	−3.5	0.2	−6.8	−60.4	1.3	−42.2	12.5
MI-21	Main mine gallery	128	−3.3	0.1	−6.7	−63.1	1.3	−45.0	8.6
MI-22	Main mine gallery	128	−3.3	0.1	−6.7	−62.0	0.2	−43.8	9.4
MI-24	Main mine gallery	128	−2.9	0.1	−6.3	−60.1	0.9	−41.9	8.3

Table 2. Oxygen and hydrogen stable isotopes of gypsum hydration water (δ¹⁸O_gyp and δ²H_gyp) in samples from the Mina Rica Mine. The isotopic composition of the paleo-aquifer (δ¹⁸O_aquifer, δ²H_aquifer and d-excess) were calculated from δ¹⁸O_gyp and δ²H_gyp using known isotope fractionation factors at 25 °C [17].

Sample ID	Description	Location (m asl)	δ18Ogyp (‰)	1SD	δ18Oaquifer (‰)	δ2Hgyp (‰)	1SD	δ2Haquifer (‰)	d-Excess (‰)
GE-1	Geode crystal (1 cm from the base)	98	−2.2	0.2	−5.6	−57.1	1.1	−38.9	5.7
GE-2	Geode crystal (2 cm from the base)	98	−2.1	0.1	−5.5	−58.4	0.9	−40.1	3.8
GE-3	Geode crystal (3 cm from the base)	98	−1.9	0.1	−5.3	−57.4	1.8	−39.1	3.5
GE-4	Geode crystal (4 cm from the base)	98	−2.5	0.0	−5.9	−57.4	1.2	−39.1	8.2
GE-5	Geode crystal (5 cm from the base)	98	−2.5	0.2	−5.9	−57.6	1.3	−39.3	8.1
GE-6	Geode crystal (6 cm from the base)	98	−2.4	0.2	−5.8	−56.2	1.5	−37.9	8.2
GE-7	Geode crystal (7 cm from the base)	98	−2.5	0.2	−5.9	−56.3	1.9	−38.1	9.3
GE-8	Geode crystal (8 cm from the base)	98	−2.7	0.2	−6.0	−57.1	0.7	−38.8	9.5
GE-9	Geode crystal (9 cm from the base)	98	−2.6	0.1	−6.0	−58.7	0.3	−40.5	7.3
GE-10	Geode crystal (10 cm from the base)	98	−2.7	0.1	−6.1	−57.8	1.8	−39.6	9.4
GE-11	Geode crystal (11 cm from the base)	98	−2.6	0.1	−6.0	−58.6	0.5	−40.4	7.6
GE-12	Geode crystal (12 cm from the base)	98	−2.6	0.1	−6.0	−58.9	0.4	−40.7	7.5
GE-13	Geode crystal (13 cm from the base)	98	−2.6	0.2	−6.0	−57.3	1.1	−39.0	8.8
GE-14	Geode crystal (14 cm from the base)	98	−2.7	0.1	−6.1	−57.8	0.9	−39.6	9.3
GE-15	Geode crystal (15 cm from the base)	98	−2.8	0.1	−6.2	−57.1	0.2	−38.9	10.5
GE-16	Geode crystal (16 cm from the base)	98	−2.8	0.1	−6.1	−57.3	0.9	−39.0	10.1
GE-17	Geode crystal (17 cm from the base)	98	−2.5	0.2	−5.9	−56.9	0.4	−38.7	8.7
GE-18	Geode crystal (18 cm from the base)	98	−2.7	0.1	−6.1	−55.4	2.6	−37.1	11.6
GE-19	Geode crystal (19 cm from the base)	98	−2.8	0.1	−6.1	−55.3	2.2	−37.1	12.0
GE-20	Geode crystal (20 cm from the base)	98	−2.7	0.2	−6.1	−57.9	0.6	−39.7	9.1
GE-21	Geode crystal (21 cm from the base)	98	−2.8	0.2	−6.2	−58.2	0.7	−40.0	9.4
GE-22	Geode crystal (22 cm from the base)	98	−2.8	0.1	−6.2	−57.7	0.8	−39.5	10.2
MI-1	Main mine gallery	128	−3.2	0.2	−6.6	−62.2	0.4	−44.1	8.8
MI-2	Main mine gallery	128	−2.5	0.2	−5.8	−54.9	0.5	−36.6	10.1
MI-3	Main mine gallery	128	−2.9	0.1	−6.3	−60.4	0.6	−42.2	8.1
MI-5	Main mine gallery	128	−3.2	0.3	−6.6	−60.9	0.7	−42.7	10.3
MI-10	Main mine gallery	128	−3.0	0.1	−6.3	−61.8	1.3	−43.6	7.1
MI-11	Main mine gallery	126	−3.7	0.4	−7.1	−60.5	0.9	−42.3	14.4
MI-13	Main mine gallery	113	−2.5	0.1	−5.9	−60.8	0.2	−42.6	4.7
MI-14	Main mine gallery	110	−2.9	0.2	−6.3	−58.4	0.4	−40.2	10.3
MI-15	Giant Gypsum Geode	97	−3.1	0.0	−6.5	−58.9	0.3	−40.7	11.5
MI-16	Main mine gallery	106	−3.5	0.0	−6.9	−61.6	0.4	−43.4	11.9
MI-18	Main mine gallery	102	−3.7	0.2	−7.1	−62.8	0.5	−44.6	12.4
MI-19	Main mine gallery	93	−2.9	0.2	−6.3	−60.3	0.4	−42.1	8.2
MI-20	Main mine gallery	93	−3.5	0.2	−6.8	−60.4	1.3	−42.2	12.5
MI-21	Main mine gallery	128	−3.3	0.1	−6.7	−63.1	1.3	−45.0	8.6
MI-22	Main mine gallery	128	−3.3	0.1	−6.7	−62.0	0.2	−43.8	9.4
MI-24	Main mine gallery	128	−2.9	0.1	−6.3	−60.1	0.9	−41.9	8.3

Our results suggest that the isotopic composition of subaqueous gypsum crystals can be used for paleo-hydrological reconstructions. However, to achieve this objective, new analytical advances in obtaining accurate and precise ages of hydrothermal gypsum with sub-ppb U levels are required. We anticipate that gypsum crystals formed from unevaporated freshwater will be utilized in the future to reconstruct the isotopic composition of paleo-aquifers or hydrothermal fluids, with important implications for paleoclimate research.