Arsenic mobilization in shallow aquifers due to CO2 and brine intrusion from storage reservoirs

We developed an integrated framework of combined batch experiments and reactive transport simulations to quantify water-rock-CO2 interactions and arsenic (As) mobilization responses to CO2 and/or saline water leakage into USDWs. Experimental and simulation results suggest that when CO2 is introduced, pH drops immediately that initiates release of As from clay minerals. Calcite dissolution can increase pH slightly and cause As re-adsorption. Thus, the mineralogy of the USDW is ultimately a determining factor of arsenic fate and transport. Salient results suggest that: (1) As desorption/adsorption from/onto clay minerals is the major reaction controlling its mobilization, and clay minerals could mitigate As mobilization with surface complexation reactions; (2) dissolution of available calcite plays a critical role in buffering pH; (3) high salinity in general hinders As release from minerals; and (4) the magnitude and quantitative uncertainty of As mobilization are predicated on the values of reaction rates and surface area of calcite, adsorption surface areas and equilibrium constants of clay minerals, and cation exchange capacity. Results of this study are intended to improve ability to quantify risks associated with potential leakage of reservoir fluids into shallow aquifers, in particular the possible environmental impacts of As mobilization at carbon sequestration sites.

Geologic CO 2 sequestration (GCS) is considered a promising approach for mitigating CO 2 emissions from centralized sources [1][2][3][4][5][6][7][8][9] . One major concern is the risk of CO 2 and/or brine leakage from deep sequestration reservoirs through highly-permeable zones such as faults and abandoned wells into overlying underground sources of drinking water (USDW) 10 . Carbon dioxide itself is not hazardous to water quality, but increased CO 2 concentrations in shallow groundwater aquifers could reduce pH and enhance geochemical reactions between groundwater and aquifer sediments, resulting in release and mobilization of toxic trace metals 11,12 . An additional risk is the leakage of reservoir brine, which may contain toxic substances, into USDWs 13,14 .
To assess the risk of CO 2 /brine leakage to overlying USDWs and to detect signatures of aquifer quality changes at early stages, various approaches have been conducted with lab-scale experiments 11, 15-17 , short-term field-scale tests [18][19][20] , numerical modeling [21][22][23][24][25][26][27] , and natural analog observations 13,28 . Most of these studies focus on CO 2 leakage and its impacts on groundwater quality, but a limited number of studies have examined the leakage of brine with/without CO 2 into shallow aquifers. Keating et al. 29 observed the upward migration of CO 2 and saline water under natural conditions, which affected the salinity and trace metal concentrations in shallow groundwater. To date, modeling approaches combined with laboratory/field observations are necessary for studies of the geochemical impacts of leaked CO 2 in shallow aquifers, to reduce the uncertainties of modeling itself and to interpret the observation data with appropriate reaction patterns. Yang et al. 27 developed an inverse multicomponent geochemical modeling approach to interpret responses of water chemistry to the introduction of CO 2 into a set of laboratory batch reactors containing carbonate-poor and carbonate-rich potable aquifer sediments. Bacon et al. 21 applied multiphase reactive transport modeling to identify potential trace metal release mechanisms under elevated CO 2 conditions in a carbonate aquifer with both batch and column experiments. However, limited studies have been conducted with such combined approaches for geochemical mechanisms under saline conditions to interpret the responses to brine leakage, which is also significant for quantitative risk analysis of GCS projects.
The primary goal of this paper is to elucidate the mechanisms of trace metal mobilization with an integrated experimental and simulation framework. The case study includes elevated CO 2 conditions at the Chimayo site in northern New Mexico, a natural analog with CO 2 upwelling. Arsenic (As) is relatively rich within the sediments and is a potential source of high As concentration in local water, thus As mobilization is of specific interest. Batch experiments of water-rock-CO 2 interactions are conducted with both fresh groundwater and saline water, to mimic scenarios of CO 2 leakage with and without deeper formation brine. Quantitative interpretation of As mobilization due to leaked CO 2 and brine utilizes an inverse reactive transport modeling approach. A specific objective for this study is to identify the responses of As release under different water salinities and to quantify the key parameters controlling As mobilization processes.

Results
Water chemistry changes. When CO 2 is introduced into the water-rock batch systems, pH decreases immediately because of CO 2 dissolution in both background (BG) and saline (S) samples (Fig. 1a,b). After the sharp drop, pH remains stable at ~5.9, even after CO 2 injection is stopped. The simulated pH also shows a sudden drop within 1 h of CO 2 injection, then a slight increase, and finally reaches a steady-state. This phenomenon was also observed by Shao et al. 17 and Yang et al. 27 . This process is mainly affected by CO 2 diffusion and mineral reactions, and the dissolution of minerals (especially calcite) consumes hydrogen ions and causes the solution pH to increase slightly between 4 and 72 h 17 . For the CO 2 -free (control) experiments, both measured and calculated pH remain stable at ~8.3 for BG and ~7.3 for S (Fig. 1a,b), and good matches with high confidence are achieved between calculated and experimental results.
Synthetic groundwater and saline water were used in our experiments, and the fast dissolving minerals (such as calcite) are not abundant in the system. This makes it hard to track major ion concentration changes in the water samples during the experiments, especially for Mg, Na, Si, Cl and SO 4 in the saline-water reactor when CO 2 is introduced. Only Ca shows a more than 50% concentration increase compared to the CO 2 -free experiments, indicating that calcium minerals (especially calcite) dissolved during the experiments (Fig. 1c,d).
With CO 2 intrusion into the shallow groundwater aquifer, trace metals of environmental concern might be released. Arsenic is of specific interest in our study due to its high concentration in the shallow groundwater of the Chimayo site 13 . Figure 1e-f shows As concentration changes of the batch experiments and simulated results. The simulated results match well with the experimental measurements, which capture the trends of As concentration changes in all four cases. Without CO 2 introduction, it shows an increase of As concentration and reaches to equilibrium after sediment and water mixed for both BG and S systems. When CO 2 is introduced to the reactors, it shows a sharp increase of As concentration at the initial time, suggesting a large amount of As released from the sediments due to CO 2 intrusion. After a few hours, As concentration starts to decrease slowly, and reaches equilibrium after 26 days. For both reactors with and without CO 2 introduction, the BG reactors show higher increases of As concentrations than that in the S reactors (20.1 µg/L vs. 15.7 µg/L for the sharp increase with CO 2 introduction and 17.9 µg/L vs. 3.8 µg/L of the control reactors), which indicates that salinity impacts the behavior of As mobilization. One possible reason is that high aqueous salinity reduces the dissolution of minerals (e.g. calcite and clay minerals), affects the water-rock equilibrium and system pH, which further hinders As release. Although As concentration exceeds the maximum contamination level (MCL) of the U.S. Environmental Protection Agency (EPA) in the beginning when CO 2 is introduced for both BG and S reactors, it drops below the MCL after 8-day exposure, which might not be considered as a long term concern of the USDWs. This As concentration drop was also observed in other studies 17,30,31 . It is notable that the BG reactor without CO 2 injection shows a large increase of aqueous phase As concentration to exceed the MCL. However, the field samples show low concentrations of As (~1 µg/L) without CO 2 and brine exposure 13 . This reveals the limitation of batch experiments in that they overestimate the water-rock and water-rock-CO 2 reactions within well-mixed water-sediment systems and large reaction areas of mineral surface, which is also pointed out by other researchers 18, 27 . Geochemical reaction parameters. The inverse modeling approach is applied for estimating geochemical reaction parameters for each experiment in this study. Estimated mean values of geochemical reaction parameters from the four experiments are shown in Table 1. The coefficient of variation (CV, the ratio of the standard deviation to the mean value) is also shown in the table to suggest the extent of variability of each parameter. Generally, the estimated values are in agreement with those reported by literature 12,25,32,33 , which suggests high confidence of the parameter estimation approach. However, the parameters with large CV values indicate high uncertainties among different specific cases. Especially, mineral reaction rates and mineral surface areas of calcite and clay minerals vary with up to two orders of magnitudes among the four cases, indicating that these parameters should be carefully chosen for simulations, and our estimation can provide a reference for determining the range of these parameters. To verify that the estimated mean values could represent the scenarios of all the cases, the model was simulated with these mean values for all the four experiments. The results indicate that these parameters are reasonable for geochemical reaction simulations especially for the analysis of As behaviors with/without CO 2 introduction under a large range of salinity conditions.
A composite sensitivity analysis of the geochemical parameters was conducted to explore the most sensitive parameters to the experimental results by calculating their composite sensitivity coefficients 34,35 . The most sensitive parameters are listed in Fig. 2. Generally, Sensitivity is increased when CO 2 is introduced, indicating that the related reactions are pH dependent. With CO 2 introduction, the model is highly sensitive to the calcite reaction rates (cal_rkf1 and cal_rkf2), hematite dissolution rate (hem_rkf), surface areas of calcite (cal_amin) and hematite (hem_amin), kaolinite adsorption surface area (kao_soh_ssa), and kaolinite adsorption equilibrium constant (logK_kao). The cation exchange selectivity and CEC showed relatively high sensitivities for all the four cases in our study, indicating that cation exchange reactions are significant for aqueous phase ion concentrations including major ions (Ca, Mg, Na, K) and trace metals (As).
Overall, the inverse modeling approaches provide reasonable geochemical reactions for our study, suggesting that it is possible to predict groundwater chemistry responses with CO 2 intrusion in other cases. However, there are still some limitations with such approaches. For example, the range and initial value of key parameters should be carefully chosen to obtain a reasonable outcome; there might be more than one solution leading to a good representation of the observation data; and uncertainties of the model structure and observation data could also impact the estimations 36 . Therefore, experimental studies of mechanisms of geochemical reactions are meaningful to provide reasonable reaction pathways and the ranges of key parameters under different conditions. Our study using the combined application of batch experiments with geochemical simulation approaches can provide a good example for future research with multiple reaction parameters for CO 2 and/or brine leakage studies.

Discussion
Arsenic mobilization mechanisms. Arsenic mobilization behavior has been widely discussed by researchers because it is one of the major concerns for groundwater quality in the event of CO 2 leakage. Most of the studies consider the behavior of As to be largely related to adsorption/desorption onto/from the surfaces of clay and Fe-rich minerals 30,37,38 , and some researchers believe that As-bond mineral dissolution (i.e. As-carbonate  minerals and Arsenopyrite) is the primary source of As contamination 12,39,40 . The two mechanisms are sometimes combined to explain As mobilization behavior in experimental and modeling approaches 12,25 . It is hard to quantitatively determine the formula and fractions of such As-bond minerals because of their low concentrations, and estimated formula and volume fractions are usually used in simulations 12 . To verify As mobilization mechanisms and to demonstrate the role of As mineral dissolution and surface adsorption/desorption, a small volume fraction (<0.01%) of As-calcite (CaCO 3 ·As 2 O 5 ) was added to the numerical model as an uncertainty parameter to analyze its impact to As behavior. The results show that the sensitivity of As concentrations to such mineral dissolution is far less than that to surface adsorption/desorption process, and As-calcite dissolution shows minimal impact on As concentration in the aqueous phase. The reason is because the mineral dissolution rate is far less than surface complexation processes, and its dissolution is limited during the experiments. According to our simulation results, adsorption/desorption onto/from the surfaces of clays, especially kaolinite, controls As mobilization with water-rock-CO 2 interactions. Arsenic also shows adsorption/desorption to illite, hematite and smectite, but with a few orders of magnitude less than that adsorbed on kaolinite. Figure 3 shows the simulated adsorbed As on kaolinite with CO 2 introduction to both BG and S reactors. With CO 2 introduction, there is an immediate sharp drop of adsorbed As, and it slowly re-adsorbs to local equilibrium in 7-15 days afterwards. The surface complexation reactions of the clay minerals (S represents mineral sites) could be written as 41 : The mass balance for the total numbers of reactive sites (omit other adsorbed metals) is: [SOH] T is related to the adsorption site density and mineral volume fraction, and the surface functional groups are competitive for the available sites. During the experimental time in our study, clay minerals do not show significant volume change, thus we can assume [SOH] T to be a constant value for this discussion. With the sudden drop of pH (increase of H + ) at the initial time, it promotes reaction (2), and SOH 2 + rapidly occupies the surface complexation sites, which results in As release from SAsO 4 2− to the aqueous phase (H 3 AsO 4 ). With the buffering effect of calcite dissolution (CaCO 3 + H + ⇌ Ca 2+ + HCO 3 − ), the pH of the system starts to increase slowly to a certain extent (Fig. 1a,b), and a series of local equilibriums are reached with SAsO 4 2− /H 3 AsO 4 for gradual As adsorption onto clay minerals (Equation (3)). Without CO 2 introduction, because of the low As concentration in aqueous phase, the adsorbed As releases from the clay mineral adsorption sites to reach to a series of local equilibriums, which results in an increase in As concentration in the water (Fig. 1e,f).
Environmental implications. Trace metal contamination, especially As, is always a major concern when considering the potential risks of CO 2 and saline water leakage from GCS sites. Results of this study suggest that, in general, As may be considered an insignificant long-term concern in a CO 2 rich environment because of clay adsorptions. Likewise, in a saline environment, high concentrations of major ions (Ca, Mg, Na et al.) could hinder As release from the clay mineral sites. In such circumstances, if As is present, the reservoir brine might contain low concentrations of As and other trace metals to begin with (i.e. not due to enhanced desorption and mobilization). In many cases, increased salinity of USDW via the leaked saline water may likely be a larger concern than any associated released the trace metals 42 . Observed data from the Chimayo site suggest that the brackish water leaked through faults provides a source of arsenic and other heavy metals 13 . Kirsh et al. 39 and Karamalidis et al. 40 showed As concentration increased in their experiments under both freshwater and saline environments, Scientific RepoRts | 7: 2763 | DOI:10.1038/s41598-017-02849-z indicating that As-rich mineral dissolution becomes dominant in such cases. Arsenic mobilization mechanisms should therefore be treated as a site-specific issue, with risk analyses of target GCS sites guided by preliminary assays of As content of both reservoir and USDW strata.

Materials and Methods
An integrated framework of systematically combined batch experiments and reactive transport simulations has been developed for better understanding As mobilization mechanisms with CO 2 leakage into shallow aquifers. The sediment samples for batch experiments are collected from the Chimayo site, a natural CO 2 analog located at the eastern margin of Espanola basin in northern New Mexico 13 . The shallow drinking water aquifer at this site has been investigated by a series of studies 25,29,43 , and water chemistry was analyzed to evaluate mobilizations of trace metals. Observed data show evidence of upward migration of CO 2 and saline water along the faults, whereas the locations far from the faults are not impacted by CO 2 or brine. This site provides an example for analyzing the impacts of CO 2 /brine leaked from a GCS site on overlying USDW quality.
Batch Experiments. Due to high As concentrations in the Chimayo groundwater exposed to CO 2 and saline water, and the relatively high amount of As present in the sediment samples (up to 147 mg/kg), the behavior of As caused by CO 2 and saline water intrusion is of specific interest. Batch experiments were conducted to mimic the influx of CO 2 and saline water into the aquifer in order to evaluate As mobilization due to the reactions between the sediments and CO 2 (Fig. 4). Details of the experiments were described by Viswanathan et al. 25 . Two sets of experiments were conducted, with "background" (BG, low salinity, TDS ~200 mg/L) and "saline" (S, high salinity, TDS ~4000 mg/L) synthetic groundwater, to represent the conditions with/without saline water intrusion. The major components of each synthetic groundwater are listed in Table 2. With each set of experiments, one reactor was exposed to CO 2 (maintained 1 atm) and one reactor was kept as a steady-state, CO 2 -free control. Before CO 2 injection, the sediment samples were exposed to the synthetic groundwater for ~3 days to reach a steady state. Water samples were then collected over a 26-day experimental period: 14 days for CO 2 injection and 12 days for post injection. Major ion and As concentrations were analyzed subsequently.
Geochemical Modeling Approach. In this study, a reactive transport model considering aqueous species complexation, mineral dissolution/precipitation, adsorption/desorption, and cation exchanges was created, in order to evaluate the mechanisms of As mobilization. All simulations were performed with TOUGHREACT V2 and ECO2N 44,45 .
A total of 54 aqueous complexes that have high impacts on the results were selected for the simulations. The aqueous complexes and their equilibrium constants are listed in Table 3. For all the batch experiments, the reactors were exposed to an oxygen-rich environment (atmosphere). Under this condition, As (V) is considered a dominant form in the aqueous phase, which also corresponds to the conditions of the site (pe > 4). Precipitation/  were included in the model. Precipitation of possible secondary minerals was allowed to constrain major ion chemistry in the system. The reaction rates and surface areas of the primary minerals were treated as uncertainty parameters. Although it is hard to detect trace As-bearing minerals (such as As-carbonate minerals) by XRD or SEM-EDS analysis, it is possible that As-rich minerals exist in the sediments 12 and impact As concentration in the aqueous phase. In this study, the initial amount of trace As minerals were treated as uncertainty parameters. Adsorption/desorption of As from clay/Fe-rich mineral surfaces was considered as an important process for As mobilization with CO 2 and saline water intrusion. Hematite, kaolinite, illite and smectite were chosen as principal adsorbents, because they were relatively abundant in the sediment samples and also widely reported by former studies 12,27,46 . Arsenic aqueous species HAsO 4 2− and H 2 AsO 4 − were chosen as major surface adsorption ions 41 . Adsorption/desorption reactions are controlled by the total amount of reactive sites (product of amount of adsorbent, site density and adsorbent surface area) 20 , which has high uncertainty for different samples 41 . The local equilibrium of adsorbed As is affected by salinity and pH 37,38,46 . Therefore, adsorbent surface area and surface complexation equilibrium constant were selected as uncertainty parameters.
Cation exchange reactions were also considered in the model for major cations, which might affect the response of trace metals and pH to CO 2 and saline water intrusion. The Gaines-Thomas convention was used in this study 47 . The site-specific parameters of cation selectivity coefficient and cation exchange capacity (CEC) were not measured for this site, and they were treated as uncertainty parameters.
To obtain the best estimations for the uncertainty model parameters, the nonlinear parameter estimation program PEST 48 was applied. The 26-day experimental data were used to estimate the uncertainty parameters listed above via inversion by minimizing the objective function J 48 where E i (p) is the sub-objective function from chemical species i, N (=9) is the number of chemical species, w ii is the weighting coefficient for the lth measurement of the ith species, which is computed with the inverse of the standard deviation of the experimental data 49 , and u l i and ũ l i are the simulated and observed concentrations of Ca, Mg, K, Na, Si, Cl, SO 4 , As and pH. A composite sensitivity analysis of the uncertainty parameters was also conducted in order to determine the most sensitive parameters 48 .