Geofluids
Volume 2018 (2018), Article ID 4834601, 19 pages
https://doi.org/10.1155/2018/4834601

Research Article

CO₂ Leakage-Induced Contamination in Shallow Potable Aquifer and Associated Health Risk Assessment

Chan Yeong Kim,¹ Weon Shik Han,¹ Eungyu Park,² Jina Jeong,² and Tianfu Xu³

¹Department of Earth System Sciences, Yonsei University, Seoul, Republic of Korea
²College of Earth System Science, Kyungpook National University, Daegu, Republic of Korea
³Key Laboratory of Groundwater Resources and Environment, Ministry of Education, Jilin University, Jilin, China

Correspondence should be addressed to Weon Shik Han; hanw@yonsei.ac.kr

Received 29 November 2017; Accepted 31 January 2018; Published 5 April 2018

Academic Editor: Liangping Li

Copyright © 2018 Chan Yeong Kim et al. This is an open access article distributed under the Creative Commons Attribution License, which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.

Abstract.

Leakage of stored CO₂ from a designated deep reservoir could contaminate overlying shallow potable aquifers by dissolution of arsenic-bearing minerals. To elucidate CO₂ leakage-induced arsenic contamination, 2D multispecies reactive transport models were developed and CO₂ leakage processes were simulated in the shallow groundwater aquifer. Throughout a series of numerical simulations, it was revealed that the movement of leaked CO₂ was primarily governed by local flow fields within the shallow potable aquifer. The induced low-pH plume caused dissolution of aquifer minerals and sequentially increased permeabilities of the aquifer; in particular, the most drastic increase in permeability appeared at the rear margin of CO₂ plume where two different types of groundwater mixed. The distribution of total arsenic (As) plume was similar to the one for the arsenopyrite dissolution. The breakthrough curve of As monitored at the municipal well was utilized to quantify the human health risk. In addition, sensitivity studies were conducted with different sorption rates of arsenic species, CO₂ leakage rates, and horizontal permeability in the aquifer. In conclusion, the human health risk was influenced by the shape of As plume, which was, in turn, affected by the characteristics of CO₂ plume behavior such as horizontal permeability and CO₂ leakage rate.

1. Introduction

Carbon capture and storage are considered to be one of the mitigating strategies for reducing CO₂ emissions to the atmosphere [1–3]. Among various carbon capture and storage technologies, CO₂ can be injected into geologically stable formations, which typically have large storage capacities and are capped by low-permeability sealing formations. However, during CO₂ injection activity, injection-induced pressure builds up within the storage formation [4–6]. If any unwanted pathways exist within the sealing formation, CO₂ is able to migrate to the shallow aquifer through these pathways while experiencing phase change from supercritical to gaseous CO₂ [7–10]. Especially in this study, the leakage of gaseous CO₂ into a shallow confined aquifer was considered. Leaked gaseous CO₂ dissolves in the potable groundwater and develops a low-pH plume [11, 12], which induces secondary contamination within the aquifer by enhancing the mobility of toxic heavy metals [13–16].

Released toxic heavy metals are able to migrate with the ambient groundwater induced by natural hydraulic gradient or influenced by enforced gradient due to any active municipal wells nearby. If the contaminated groundwater produced from the municipal well is directly distributed without proper treatment (disinfection or chlorination processes) to residents who use it daily for the purpose of drinking, bathing, cleaning, or other household uses, these residents can be exposed to adverse carcinogenic health risks. The World Health Organization’s (WHO) International Agency for Research on Cancer classified various dissolved heavy metals (e.g., arsenic and lead) and other radioactive metals (e.g., uranium and cesium) as toxic substances hazardous to human health [17]. Among these carcinogenic heavy metals, this study focused on arsenic [18]. According to the US Agency for Toxic Substances and Disease Registry, who ranked hazardous substances based on their occurrence, toxicity, and potential for human exposure, arsenic was ranked the first in their Substance Priority List in 2013 and 2015 [19]. Moreover, arsenic contamination and its detrimental impacts have been reported by various countries such as China [20], Bangladesh [21], Vietnam [22], and India [23] in recent years.

Mobilization of arsenic in a shallow groundwater aquifer due to CO₂ leakage has been investigated at one of the natural analog sites, Chimayó, New Mexico, where CO₂-saturated brackish water was leaked into the shallow aquifer along the fault [7, 13, 24, 25]. At this site, decreased pH and the resulting mobilization of trace metals, including arsenic, were observed. Even if their adverse effects had been alleviated due to the high buffering capacity of the local groundwater aquifer, Keating et al. [7] reported significantly elevated trace metal concentrations at a number of local wells due to the influx of brackish waters. Later, both Keating et al. [13] and Viswanathan et al. [24] integrated the field dataset into a multiphase reactive transport model to understand the behavior of arsenic, since some wells in Chimayó exceeded the maximum contamination level (MCL). In addition to studies targeting natural CO₂ release sites, several experiments have been conducted at field-scale CO₂ injection sites to determine secondary contamination caused by the injected CO₂ [11, 26–28]. Decreased pH, increased concentrations, and subsequent changes in groundwater chemistry such as increased Fe²⁺, Mn²⁺, Mg²⁺, and Ca²⁺ concentrations were observed at both Frio-I brine pilot injection [11, 27] and ZERT field tests [12, 28].

In addition to field-oriented research, several studies have focused on the implementation of numerical simulations to evaluate geochemical behaviors associated with arsenic contamination. For example, Zhang et al. [29] and Xiao et al. [30] utilized the reactive transport simulation and investigated the complex chemical change induced by CO₂ leaked into a shallow aquifer. In contrast to these studies which presented detailed geochemical interactions, Siirila et al. [31] simplified the geochemical processes by solving the advection-dispersion equation with a linear sorption. Without full assessment of geochemical reactions, they were able to account for the movement of toxic elements in complex 3D heterogeneous systems as well as for groundwater contamination-induced carcinogenic health risks with a probabilistic approach. Later, Navarre-Sitchler et al. [32] utilized PFLOTRAN and simulated the mobilization of lead in a complex heterogeneous system by assuming that released gaseous CO₂ instantaneously dissolved in the groundwater aquifer. The proposed work in this study was built upon the framework of the previous studies mentioned above. Certain approaches accounted for the detailed geochemical behavior of toxic heavy metals (e.g., complexation, sorption, mineral dissolution, and precipitation) in the groundwater aquifer, whereas other approaches relied on either simplified chemical reactions (e.g., advection-dispersion equation) or multiphase fluid migration (e.g., gaseous CO₂ leakage into the groundwater aquifer) while accounting for carcinogenic health risks or complex 3D heterogeneous systems. Therefore, the goal of this study was to integrate these two approaches and delineate the multiphase behavior of leaked gaseous CO₂ to a shallow potable aquifer. Moreover, leaked CO₂-induced geochemical changes, such as evolution of water quality and mobility of toxic trace metals (e.g., arsenic), were characterized by adapting the multispecies reactive transport model. Finally, the simulated concentration of arsenic species observed from the assigned municipal well was used to quantify the carcinogenic health risk for chronically exposed humans.

2. Behavior of Arsenic in Subsurface Environment

Arsenic in natural water typically originates from arsenic-bearing minerals that frequently possess sulfur, oxygen, and iron [33]. Generally, these naturally occurring arsenic-bearing minerals include arsenopyrite (FeAsS), realgar (AsS), enargite (CuAsS₄), scorodite (FeAsO₄·2H₂O), and tennantite (Cu₆[Cu₄(Fe,Zn)₂]As₄S₁₃) [29, 34, 35]. When these minerals dissolve, various forms of arsenic species such as , , , , , and can be released into groundwater. As seen in Figure 1, arsenic species can be present in several valence states (−3, 0, +3, and +5), but in natural groundwater, it is mostly found in oxyanions of trivalent arsenite (As(III)) or pentavalent arsenate (As(V)), depending on reducing or oxidizing conditions, respectively. Typically, the mobility and toxicity of As(III) are considered to be much higher than those of As(V) [36–38]. Under the reducing condition, uncharged As(III) species, such as , is dominant below pH 9.2 (Figure 1). However, under the oxidizing condition, is dominant at pH < 6.9, while is dominant at higher pH (pH > 6.9). Moreover, and could exist under extremely acidic and alkaline conditions, respectively. In this study, when injected CO₂ leaked into the shallow potable groundwater aquifer, the dissolution of leaked CO₂ induced the acidification of the ambient groundwater to 4 < pH < 6 [11], where and were dominant under reducing and oxidizing conditions, respectively [33, 39]. Here, the shallow aquifer was considered to be under the reducing condition, and thus, was chosen as the most dominant species.

Figure 1: Eh-pH diagram for arsenic species under H-O-As system at 1 bar and 25°C.

3. Work Flow: Numerical to Probabilistic Quantification

The following study comprised two major parts, namely, numerical prediction of CO₂ transport from a leakage point to a municipal well (numerical simulation) and probabilistic quantification of health risks to humans who have been chronically exposed to a certain toxic heavy metal (health risk assessment), in sequence (Figure 2). Processes for both numerical simulation and probabilistic risk quantification were designed in three steps: “Data,” “Process,” and “Result.” During the “Data” step, input parameters were chosen selectively, and the sampling scheme was determined either deterministically or randomly. Subsequently, a series of calculations using either numerical or probabilistic approaches were conducted in the “Process” step. Finally, in the “Result” step, simulation outputs, such as CO₂ plume distribution, maximum concentrations of selected toxic heavy metals, and human health risk, were analyzed. The link between “numerical simulation” and “health risk assessment” was the profile of arsenic concentration monitored at the municipal well; throughout complex spatial and temporal movements of the CO₂ plume predicted from numerical simulations, a time series of arsenic concentration was observed at the municipal well, which was then used to assess carcinogenic human health risk.

Figure 2: Flow chart explaining the sequence of both numerical simulation and probabilistic health risk assessment. Numerical simulation composed of multiphase and multispecies reactive transport modeling required both physical and chemical parameters within the shallow aquifer. Subsequently, the output of numerical simulation was utilized for calculating probabilistic health risk for carcinogenic effect.

4. Numerical Approach

The multiphase and multicomponent reactive transport simulator, TOUGHREACT, was used to simulate secondary contamination processes in a shallow confined aquifer induced by CO₂ leakage [44] in conjunction with the ECO2N module, which was used to predict the fluid properties of H₂O, NaCl, and gaseous CO₂ [45]. For this study, the original thermodynamic database incorporated in TOUGHREACT was not adequate for delineating the complex chemical reactions associated with arsenic species. Therefore, the revised EQ3/6V7.2b database was specifically adopted to account for chemical reactions of arsenic and related chemical species [43]. In detail, the modified thermodynamic database adopted arsenite () as a primary species and incorporated associated aqueous complexes (e.g., , , H₃AsO₄(aq), HAsO₂(aq), , and HAsS₂(aq)), while accounting for their activity coefficients from the extended Debye-Hückel equation [46]. Moreover, solubility products of arsenic-bearing minerals such as arsenopyrite were included in this thermodynamic database.

Precipitation and dissolution of minerals were simulated kinetically by following the rate law, which was coupled with the equation representing the kinetic rate constant (see (1)). In this equation, the kinetic rate constant is dependent not only on temperature but also on the pH as shown below [47]:

where  (J/mol) is the activation energy and  (mol/m²/s) is the rate constant at 25°C with superscripts nu, , and representing neutral, acid, and base mechanisms, respectively.  (J/mol/K) is the gas constant,  (K) is temperature, is the activity of dissolved species under acid or base conditions, and is the power term. is the specific reactive surface area (cm²/g), and is the kinetic mineral saturation ratio of the th mineral. Finally, both and are assumed to be unity.

5. Conceptual Model

5.1. Model Description

The model designed for this study delineates the CO₂ leakage process throughout undetected and unexpected pathways in the sealing formation. In order to understand such processes, the 2D cross-sectional potable aquifer model was designed as shown in Figure 3. The potable aquifer was assumed to be relatively deep (300 m with a reducing condition), where the municipal or high-capacity wells produced a large amount of groundwater, which was then distributed for residential purposes [48, 49]. In addition, the width and thickness of the potable aquifer were assigned as 200 m and 40 m, respectively, with the size of grid blocks of 2 m; the total number of grid blocks was 2,000 (100 × 20). Assuming that the aquifer is located at a depth of 300 m, the initial formation pressure and temperature were assigned as 3 MPa and 25°C, respectively. Top and bottom boundaries of the model were assigned as no-flow conditions assuming that the upper and lower confining formations act as nearly impermeable sealing units. Lateral boundaries (purple-colored grid blocks) were assigned as the Dirichlet condition, where the left and right pressure were 3.2 MPa and 2.8 MPa, respectively. The difference of pressure in lateral boundaries in addition to aquifer properties ( m²) results in 1.94 cm/day of the ambient groundwater flow in the aquifer. The leakage point for gaseous CO₂ was located at a 50 m distance from the left boundary where CO₂ was leaked at a rate of 0.05 kg/s (Figure 3). The CO₂ leakage was maintained for 1 year. Furthermore, the municipal well with a continuous pumping rate of 0.5 kg/s was located 100 m away from the CO₂ leakage point and penetrated a depth of 30 m from the upper confining seal. Finally, the simulation was conducted for 100 years.

Figure 3: Conceptual model representing the shallow potable aquifer with a length and thickness of 200 m and 40 m, respectively. Blue arrows represent ambient groundwater flow with the rate of 2 cm/day. The grid blocks with the purple color stand for constant pressure boundaries; left and right boundary are 3.2 MPa and 2.8 MPa, respectively. Both upper and lower boundaries are assumed to be the no-flow boundary. CO₂ leakage point (red) is located at 50 m away from the left boundary. In addition, the municipal well (green) located at 50 m away from the right boundary penetrates 30 m deep from the upper seal while the pumping activity was only conducted through the screened section.

5.2. Physical and Chemical Parameters

The aquifer was assumed to be homogeneous with its horizontal permeability () and porosity () to be 10–13 m² and 0.2, respectively (Table 1); the anisotropy ratio (/) was assumed to be 0.1. In addition, delineation of gaseous CO₂ transport through the groundwater-saturated aquifer required constitutive relations such as relative permeability and capillary pressure. In this study, Van Genuchten’s functions were adapted, and the relevant parameters are shown in Table 1 [40].

Table 1: Input parameters for shallow groundwater aquifer materials. Relative permeability and capillary pressure follow equations developed by van Genuchten [40].

Aquifer characteristic

Permeability

    Horizontal
10⁻¹³

    Vertical
10⁻¹⁴

Porosity
0.2

Density (kg/m³)
2600

Tortuosity
0.3

Compressibility
10⁻⁸

Diffusivity (m²/s)
10⁻⁹

Relative permeability

: irreducible water saturation
0.20

: irreducible gas saturation
0.05

: exponent
0.457

Capillary pressure

: irreducible water saturation
0.20

: exponent
0.457

: strength coefficient
19.59 kPa

The mineralogical composition of the designated aquifer materials was assumed to be sandstone, which is the typical hosting formation for geological CO₂ sequestration [50, 51] and natural analog CO₂ leakage sites [52, 53]. As an example, Kampman et al. [54] analyzed the mineralogical assemblage of Navajo Sandstone, which was considered to be the primary sourcing aquifer for CO₂-charged brine [55]. Chemical analyses of Navajo Sandstone fluids collected either from adjacent geysers or from springs revealed elevated concentrations of both arsenic and other toxic heavy metals [54, 56–58]. Due to this reason, mineralogical composition in this study was adopted from the composition of Navajo Sandstone while assuming that arsenopyrite is the primary source of arsenic species in the ambient groundwater. A total of 7 primary minerals, that is, quartz, kaolinite, illite, K-feldspar, arsenopyrite, calcite, and magnesite, were chosen, and another 5 secondary minerals (chlorite, dolomite, goethite, oligoclase, and smectite-Ca) were expected to be precipitated (Table 2).

Table 2: Mineral volume fractions for shallow groundwater aquifer materials. A total of 7 primary and 5 secondary minerals were selected.

Primary mineral
Vol. fraction

Quartz
SiO₂

Kaolinite
Al₂Si₂O₅(OH)₄

Illite
K0.6Mg0.25Al1.8(Al0.5Si3.5O₁₀)(OH)₂

K-feldspar
KAlSi₃O₈

Arsenopyrite
FeAsS

Calcite
CaCO₃

Magnesite
MgCO₃

Secondary mineral
Vol. fraction

Chlorite
Mg2.5Fe2.5Al₂Si₃O₁₀(OH)₈
0.00

Dolomite
(CaMg(CO₃)₂)
0.00

Goethite
FeO(OH)
0.00

Oligoclase
CaNa₄Al₆Si₁₄O₄₀
0.00

Smectite-Ca
Ca0.145Mg0.26Al1.77Si3.97O₁₀(OH)
0.00

Mineral volume fraction of each primary mineral was chosen from a similar quantity of Navajo Sandstone with the addition of 1% arsenopyrite; quartz was predominant (81%), and kaolinite accounted for the second-largest amount (16%). Due to the small proportion of carbonate minerals, fluid chemistry possessed less buffering capacity against CO₂ dissolution. With the chosen mineral composition, the batch reaction was conducted to determine a list of primary species and their initial concentrations. As a result, , Ca²⁺, Cl⁻, Fe²⁺, H⁺, H₂O, H₃AsO₃(aq), , , HS⁻, K⁺, Mg²⁺, Na⁺, O₂(aq), SiO₂(aq), and were selected as primary species, and associated aqueous complexes were chosen as secondary species (Table 3).

Table 3: The 16 primary species and 39 secondary species are selected after the batch-reaction simulation.

Primary species (16)

H⁺

Na⁺

Ca²⁺
H₂O
HS⁻
O₂(aq)

Cl⁻
H₃AsO₃(aq)
K⁺
SiO₂(aq)

Fe²⁺

Mg²⁺

Secondary species (39)

CO₂(aq)
H₂S(aq)

Al(OH)₃(aq)

H₃AsO₄(aq)
MgSO₄(aq)

Al³⁺
Fe³⁺

NaCl(aq)

AlOH²⁺
FeCl⁺
HAlO₂(aq)

CaCl₂(aq)

HAsO₄(aq)
NaHCO₃(aq)

CaCO₃(aq)

NaOH(aq)

HAsS₂(aq)

CaOH⁺
H₂(aq)
HS⁻
OH⁻

CaSO₄(aq)

SO₂(aq)

CH₄(aq)

MgCl⁺

As described by Keating et al. [7] and Zheng et al. [43], the concentration of arsenic species in the aquifer media was significantly controlled by both adsorption and desorption processes. Due to this reason, many researchers experimentally measured sorption values () for arsenic species under various conditions (e.g., different minerals and pH) and revealed that of As(III) and As(V) varied from 1.47 L/kg to 275 L/kg within geologic media including subsurface soil and aquifers [39, 59, 60]. In this study, a shallow sandstone aquifer was targeted, composed of over 90% quartz and kaolinite (Table 2). Previous experimental measurements revealed that for quartz and kaolinite was measured at 2 and 19 L/kg, respectively [39]. Therefore, the value was assigned as 10 L/kg.

The multiple parameters required to address the kinetic rates for mineral reactions following (1) are listed in Table 4. In addition, the calculation of reactive surface area for minerals was followed by both Xu et al. [27] and Sonnenthal et al. [61], who assumed that a mineral is a cubic array of truncated spheres, in which the radius of the sphere is assumed to be 0.001 m. In this study, the surface roughness-based area predicted from the spherical radius was reduced by two orders to reasonably represent the reactive surface area. Typically, chemical reactions only occur at selected sites on the mineral surface, and furthermore, only a small fraction of mineral surface is involved in this reaction due to grain coating and armoring. Therefore, the reactive surface areas of most silicate and carbonate minerals were chosen to be approximately 10 cm²/g, similar to those chosen by Knauss et al. [62] and Zerai et al. [63]. Finally, reactive surface areas of clay minerals such as kaolinite, illite, chlorite, and smectite were selected for larger values because of their smaller grain sizes [43].

Table 4: Kinetic data for selected primary and secondary minerals. Only calcite was treated similar to the equilibrium reaction, and the rest of the minerals were kinetically treated with consideration of neutral, acid, and base mechanisms.

Mineral

(cm²/g)
Neutral mechanism
Acid mechanism
Base mechanism
References

(mol/m²·s)

(kJ/mol)

(mol/m²·s)

(kJ/mol)

(H⁺)

(mol/m²·s)

(kJ/mol)

(H⁺)

Quartz
9.8

87.7

[41]

Kaolinite

22.2

65.9
0.78

17.9
−0.47
[42]

Illite

35.0

23.6
0.34

85.9
−0.40
Set to smectites

K-feldspar
9.8

38.0

51.7
0.50

94.1
−0.82
[42]

Arsenopyrite
12.9

62.8

[43]

Calcite
Equilibrium reaction

Magnesite
9.8

23.5

14.4
1.00

[42]

Chlorite
9.8

88.0

88.0
0.50

[42]

Dolomite
12.9

52.2

43.5
0.50

[27]

Goethite
9.8

62.8

[43]

Oligoclase
9.8

69.8

65.0
0.46

[42]

Smectite-Ca

35.0

12.6
0.34

58.9
−0.40
[42]

Porosity was associated with changes in volume due to mineral dissolution or precipitation [64]. The porosity in this study was simulated by the following equation:

where is the number of minerals and and are the volume fractions of th mineral in the rock and nonreactive rock, respectively. Finally, permeability change was calculated using the porosity changes with cubic law [65]:

where and are permeability and porosity, respectively, with subscript representing the initial value.

Even if the proposed approach was capable of simulating water-rock interaction, it would still have some limitations. Firstly, in this study, the sorption effect of arsenic onto the surface of clay minerals was simulated with a linear approach instead of surface complexation. Although the linear approach had limitations in terms of delineating chemical heterogeneities on both temporal and spatial scales, this approach was effective for application to a large-scale simulation by reducing computation [66–68]. Secondly, arsenopyrite is a solid solution of pyrite (FeS₂), in which the ratio of arsenic and sulfur varies depending on their mole fractions. However, in this study, for simplification in predicting the thermodynamic properties, it was assumed that 1 mole fraction of arsenic replaced the sulfur (FeAsS). Finally, it was assumed that arsenopyrite oxidatively dissolves in the presence of common geologic oxidants such as dissolved O₂, , and Fe³⁺. Under acidic conditions, Fe³⁺ quickly oxidizes, and sequential dissolution of arsenopyrite releases arsenic species ( and ) by following (4). Therefore, the stoichiometric reaction for arsenopyrite dissolution, which was used in this numerical simulation, involved the reduction of Fe³⁺ to Fe²⁺ [43].

5.3. Probabilistic Health Risk Approach

Carcinogenic effects on humans who have been chronically exposed to arsenic species through multiple pathways were probabilistically quantified based on individual exposure rate and toxicity, which were suggested by the “Guidelines for Carcinogenic Risk Assessment” of the US Environmental Protection Agency (EPA) [69], Siirila et al. [31], and the EPA Superfund Risk Assessment [70]. In general, three different uptake pathways, namely, ingestion, dermal sorption, and inhalation, were considered. Among these, ingestion and dermal sorption were considered to be major pathways, because humans are often exposed to risk by drinking dissolved toxic species in tap water or showering. Uptake through inhalation was not considered in this study, as it is unlikely that the concentration of vaporized trace metals was high enough to cause carcinogenic effects indoors under normal conditions. Therefore, following Siirila et al. [31], only two exposure pathways (ingestion and dermal sorption) were considered for quantifying exposure rate and toxicity of arsenic.

Arsenic toxicity was predicted from the product of cancer potency factor (CPF) (kg·day/mg) and average daily dose (ADD) (mg/kg·day). Typically, the CPF is different at the individual pathway even for the same toxic metal; in this study, the CPF for ingestion and dermal sorption of arsenic were assumed to be 1.5 and 1.6 kg·day/mg by following IRIS [71] and EPA [72], respectively. Subsequently, ADD was used to assess the individual exposure rate with the following equation:

where is the individual intake per body weight (L/kg/day), AT is the average lifetime, which was assumed to be 25,550 days (70 years), and EF is the standard exposure frequency during 1 year, which was assumed to be 350 (days/year) [73]. The most significant term in (6) is , which is the maximum average of arsenic concentration (mg/L) monitored at the municipal well during the exposure duration (ED) (years). In this study, the ED was assumed to be 30 years following the EPA guideline [74]. Accordingly, , calculated from concentration profiles () of arsenic species at the municipal well in the numerical simulation can be represented as follows:

Specific to ADD, exposure pathways through both ingestion and dermal sorption are defined as follows:

where is the skin surface area per body weight (m²/kg), is the dermal permeability coefficient for arsenic in water, designated as 1.0 × 10⁻⁵ m/hour, is the fraction of skin in contact with water (−), is exposure time of shower per day (hours/day), and CF is the unit conversion factor (0.001 L/m³). Finally, the total risk is the summation of individual risks representing each exposure pathway as shown below:

In this study, the calculation of the total risk was repeated 100,000 times to account for variability, which stands for probabilistic quantification in human health risk. Except for (the concentration of arsenic species obtained from the result of numerical simulation), individual parameters considered in the quantification of total risk were randomly sampled within each intrinsic distribution (Table 5) [31, 75]. Finally, the calculated total health risks were plotted as a cumulative density function, which enables the estimation of the probability exceeding the risk level of concern (10⁻⁴) [76, 77].

Table 5: Parameters for health risk assessment. Each parameter is randomly sampled within the designated value range and distributions to represent individual characteristics of potential victims.

Parameter
Symbol
Unit
Distribution
Values

Exposure duration
ED
[yr]
Constant
30

Exposure frequency
EF
[d/yr]
Constant
350

Averaging time
AT
[d]
Constant
25,550

Ingestion rate per unit body weight
IR/BW
[L·kg/d]
Lognormal

Skin surface area per unit body weight
SA/BW
[m²/kg]
Lognormal

Fraction skin in contact with water

[-]
Uniform

Shower exposure duration

[h/d]
Lognormal

Unit conversion factor

[L/m³]
Constant

Cancer potency factor

[kg·day/mg]
Constant
1.5

1.6

6. Model Scenarios

Case 1 (base case) was designed to delineate CO₂ leakage processes associated water-rock interactions and secondary contamination caused by arsenic species (Table 6). Subsequently, sensitivity studies were conducted with different , CO₂ leakage rates, and horizontal permeability of the aquifer (). First, the degree of sorption intensity was evaluated by varying between 25, 50, and 100 L/kg (Cases 2–4). Second, in Cases 5–7, the effect of CO₂ leakage rate was evaluated by varying its rate between 0.020, 0.025, and 0.030 kg/s. Different CO₂ leakage rates could induce the development of CO₂ plumes with different sizes. For example, as the size of a CO₂ plume increases, larger areas expect to experience water-rock interactions and more dissolution of arsenic species. Consequently, the municipal well captures dissolved arsenic species more when the size of the CO₂ plume is greater, which eventually increases the carcinogenic health risk on humans. Finally, in Cases 8–10, varied ( = 0.2 × 10⁻¹³, 1.0 × 10⁻¹³, and 5.0 × 10⁻¹³ m²) while maintaining the vertical permeability (1.0 × 10⁻¹⁴ m²); increased accelerates the horizontal velocity of ambient groundwater while reducing buoyancy forces on the CO₂ plume.

Table 6: Different parameters for sensitivity studies. Three groups except Case 1 (base case) are designed for sorption rate () (Cases 2–4), gaseous CO₂ leakage rate () (Cases 5–7), and horizontal permeability () (Cases 8–10), respectively.

Case name
Sorption rate 
()
CO₂ leakage rate 
()
Horizontal permeability 
()

Case 1 (base case)
10
0.050 kg/s
 m²

Case 2
25
0.050 kg/s
 m²

Case 3
50
0.050 kg/s
 m²

Case 4
100
0.050 kg/s
 m²

Case 5
50
0.020 kg/s
 m²

Case 6
50
0.025 kg/s
 m²

Case 7
50
0.030 kg/s
 m²

Case 8
50
0.050 kg/s
 m²

Case 9
50
0.050 kg/s
 m²

Case 10
50
0.050 kg/s
 m²

7. Results and Discussion

7.1. Base Case

7.1.1. Migration of Leaked CO₂ Plume within the Shallow Potable Aquifer

Figures 4(a)–4(e) represent the evolution of the leaked CO₂ plume at designated times of 120, 240, 360, 480, and 600 days. The mass centers of the CO₂ plume, shown as red, black, and yellow circles, were calculated and plotted every 120 days; black circles represent the present time of the mass center, and red and yellow circles represent past and future times, respectively. From the leakage point, CO₂ continuously leaked at a rate of 0.05 kg/s only until 365 days. During this period, three flow systems induced by ambient groundwater, CO₂ leakage, and pumping activity interacted with each other (Figure 4(c)), which developed two mixing zones for geochemically different types of groundwater at both the front and the rear margins of the CO₂ plume. The ambient groundwater flow (orange arrows) was developed from the left to the right boundaries at an approximate rate of 1.94 cm/day. Additionally, the CO₂ plume gradually expanded from the leakage point where CO₂ saturation remained at 0.3; the rate of CO₂ flux was approximately 2.1 × 10⁻⁵ kg/(s·m²) (black arrows) adjacent to the leakage point. At the rear margin of the CO₂ plume, two chemically different types of groundwater (ambient groundwater and CO₂-dissolved groundwater) flowed in opposite directions, inducing the development of a vigorous geochemical mixing zone. Subsequently, the CO₂ plume migrated together with the ambient groundwater until it was captured by the municipal well (Figure 4(c)).

(a) 120 days

(b) 240 days

(c) 360 days

(d) 480 days

(e) 600 days

Figure 4: Saturation of leaked CO₂ at different times (120, 240, 360, 480, and 600 days). Red, black, and yellow circles with labels indicate the locations of mass center for CO₂ plume. Red circles in (d) and (e) stand for the evolution of the mass center in multiple times; yellow circles in (e) express the mass center from 600 to 1,080 days.

After CO₂ leakage had stopped at 365 days, the size of the CO₂ plume gradually decreased due to both dissolution to ambient groundwater and extraction from the municipal well (Figures 4(d) and 4(e)); CO₂ solubility predicted by Duan and Sun [78] was 0.83 mol/kg water in this aquifer (3 MPa and 25°C). After movable CO₂ was captured by the municipal well, residually trapped CO₂ governed by irreducible CO₂ saturation (shown in Table 1) remained until complete dissolution to the ambient groundwater; residually trapped CO₂ eventually vanished approximately 5 years after CO₂ leakage had stopped.

The calculated mass center was located close to the plume center while the CO₂ plume migrated toward the municipal well. The migration rate of the CO₂ plume estimated from the mass center locations was approximately 12.5 cm/day until 480 days (Figures 4(a)–4(d)), indicating that the migration rate of the CO₂ plume was faster than the ambient groundwater flow (1.94 cm/day). This is because the pumping activity developed an additional head gradient, which was greater than ambient groundwater flow. After 480 days, the mass center slightly moved back until 600 days (Figure 4(e)), implying that all movable CO₂, which was weighted CO₂ mass at the plume front, was pumped out. After 600 days, continuous CO₂ dissolution occurred at the plume rear, and the pumping activity at the front induced the movement of the mass center slowly.

7.1.2. Induced Geochemical Reactions

The dissolution of gaseous CO₂ into ambient groundwater increased concentration from 0.2 to 1 mol/L while decreasing the pH from 8.0 to 5.5 within the CO₂ plume (Figures 5(a) and 5(d)). Concurrently, both carbonate and silicate minerals were either dissolved or precipitated. In particular, dramatic changes in minerals were localized at both rear and front margins of the CO₂ plume where the vigorous advective mixing of chemically different types of groundwater occurred.

Figure 5: Distribution of pH, carbonate minerals (calcite and magnesite), and associated dissolved species (, Ca²⁺, and Mg²⁺) after 365 days.

Dissolutions in calcite (CaCO₃) and magnesite (MgCO₃) were distinct (Figures 5(b) and 5(c)). Such dissolution released into the ambient groundwater and subsequently induced positive feedback to lower pH. Calculation of the saturation index (SI) using the initial concentration of ambient groundwater indicated that calcite () was more saturated than magnesite (); the initial concentrations of Ca²⁺ and Mg²⁺ were 3.3 × 10⁻³ mol/L and 4.2 × 10–12 mol/L, respectively, in the ambient groundwater. Due to this reason, when CO₂ was leaked, dissolution of magnesite (−6 mol/m³) was greater than that of calcite (−4 mol/m³) within the CO₂ plume. However, even if overall magnesite dissolution was greater than that of calcite, the greatest change in mineral dissolution appeared in calcite (−8 mol/m³), focusing at the rear margin of the CO₂ plume (Figure 5(b)). This implies that the localized dissolution of calcite was primarily induced by the mixing of two chemically different types of groundwater, such as the ambient and CO₂-dissolved groundwater. Initially, the Ca²⁺ concentration in the ambient groundwater was 3.3 × 10⁻³ mol/L (Figure 5(e)). Inside the CO₂ plume, calcite was dissolved and increased Ca²⁺ concentration 10-fold to 2.5 × 10⁻² mol/L. At the rear margin of the CO₂ plume, Ca²⁺ concentration was increased even more (to 3.1 × 10⁻² mol/L). Overall, the distribution of Ca²⁺ was similar to that of calcite (Figures 5(b) and 5(e)), and distributions of both Mg²⁺ and were similar to that of magnesite (Figures 5(c), 5(d), and 5(f)).

Patterns of dissolution and precipitation in silicate minerals were more complex than those of carbonate minerals (Figures 6(a)–6(c)). CO₂ leakage primarily induced the dissolution of K-feldspar (KAlSi₃O₈) (Figure 6(a)); the greatest dissolution (3.5 × 10⁻² mol/m³) occurred at the rear of the CO₂ plume, and the degree of dissolution gradually decreased as the plume approached the municipal well. Dissolution of K-feldspar increased concentrations of K⁺, SiO₂(aq), and in the groundwater (Figures 6(d)–6(f)). Nevertheless, distributions of such species did not imitate the dissolution pattern of K-feldspar. Rather, SiO₂(aq) and K⁺ showed the highest concentrations at the rear and front margins of the CO₂ plume with values of 1.8 × 10⁻⁴ mol/L and 3.5 × 10⁻⁴ mol/L, respectively (Figures 6(d) and 6(e)), but the decrease in concentration occurred uniformly throughout the CO₂ plume (Figure 6(f)). The discrepancy in patterns between K-feldspar and other dissolved species was presumably caused by a combination of both dissolution and precipitation among various silicate minerals such as illite, kaolinite, and chlorite as described below.

Figure 6: Distribution of silicate minerals (K-feldspar, illite, chlorite, and kaolinite) and associated dissolved species (SiO₂(aq), K⁺, and AlO₂) after 365 days.

For the distribution of illite (K0.6Mg0.25Al1.8(Al0.5Si3.5O₁₀)(OH)₂), a small amount (3.3 × 10⁻² mol/m³) was precipitated throughout the CO₂ plume (Figure 6(b)). However, at the rear margin and immediately adjacent to the municipal well, a relatively large degree of illite dissolution was predicted (−9.1 × 10⁻³ and −7.4 × 10⁻² mol/m³, resp.). In contrast to illite, a small amount (−2.9 × 10⁻² mol/m³) of kaolinite was dissolved within the plume, and reversely, a small amount of precipitation was predicted at the rear and front margins (1.9 × 10⁻² and 4.8 × 10⁻² mol/m³, resp.) (Figure 6(c)). Presumably, the dissolution and precipitation of illite and kaolinite would influence the distribution of SiO₂(aq) and K⁺ in addition to K-feldspar dissolution (Figures 6(d) and 6(e)). Finally, the behavior of individual silicate mineral influenced the distribution of , the concentration of which within the plume was lower than that outside (Figure 6(f)). Overall, CO₂ leakage induced dissolution or precipitation of both carbonate and silicate minerals and ultimately changed both the porosity and the permeability of the shallow aquifer. The dissolution of carbonate minerals primarily caused an increase in permeability; permeability increased to 0.43% () within the CO₂ plume, and the most drastic increase (0.77%) occurred at the rear margin of the CO₂ plume (Figure 7).

Figure 7: Distribution of altered permeability (%) after 365 days. The changes in permeability were predicted from the following equation () because its variation was too small. Here, and indicate the initial permeability and altered permeability at designated time, respectively.

Dissolution of arsenopyrite, which was the primary reaction for predicting carcinogenic health risk, occurred only within the CO₂ plume with a dissolved amount of 4.84 × 10⁻⁵ mol/m³ (Figure 8(a)). Similar to that of carbonate minerals, the greatest amount of arsenopyrite dissolution (8.57 × 10⁻⁵ mol/m³) occurred at the rear margin of the plume due to the vigorous mixing of two chemically different types of groundwater. Following (4), oxidative dissolution of arsenopyrite consumed 0.75 moles of O₂(aq) and 1 mole of H⁺ while increasing the concentrations of total arsenic (As), Fe²⁺, and HS⁻ (Figures 8(b)–8(e)). In this study, As represents the summation of primary species such as arsenite ((aq)), which is the by-product of arsenopyrite dissolution, as well as other arsenic species such as , H₃AsO₄(aq), , , , HAsO₂(aq), and HAsS₂(aq). Distribution of As concentration mimics that of arsenopyrite (Figures 8(a) and 8(b)); generally, As concentration within the CO₂ plume was greater than that outside, while the rear margin revealed the highest concentration. However, the other associated species such as Fe²⁺, HS⁻, and O₂(aq) revealed relatively uniform distribution (Figures 8(c)–8(e)). The difference between arsenopyrite-produced species such as As and other associated species (e.g., Fe²⁺, HS⁻, and O₂(aq)) presumably occurred due to the sorption effect, which was accounted with the linear approach. was designated for arsenic species only, and thus, as shown in Figure 8(b), enrichment of As concentration occurred at the rear of the CO₂ plume. In summary, the greatest concentration of As was 4.9 × 10⁻⁷ mol/L at the rear margin of the plume and the average concentration of As within the plume was 2.9 mol/L (Figure 8(b)). The average concentrations of other species such as Fe²⁺ and HS⁻ were 3.9 × 10−7 mol/L and 3.6 × 10⁻⁷ mol/L, respectively (Figures 8(c) and 8(d)).

(a) Arsenopyrite (FeAsS) [mol/m³]

(b)As [mol/L]

(c) Fe²⁺ [mol/L]

(d) HS⁻ [mol/L]

(e) O₂(aq) [mol/L]

Figure 8: Distribution of arsenopyrite dissolution and associated dissolved species (As, Fe²⁺, HS⁻, and O₂(aq)) after 365 days.

7.1.3. Health Risk Assessment of Carcinogenic Effect

To account for carcinogenic health risk, selected species such as pH, As, and arsenite ((aq)) concentrations were monitored at the municipal well for 100 years (Figure 9(a)). Dramatic changes in gaseous CO₂ saturation (), pH, and mass fraction of CO₂ dissolved in groundwater () predicted during 10 years were magnified at the small window (Figure 9(a)). Depending on the profiles of both pH and dissolved arsenic species, two stages (Stages I and II) were characterized. During Stage I (0–6.3 years), the leaked CO₂ plume, which existed as either gaseous CO₂ (, black dotted line) or dissolved CO₂ (, purple dotted line), arrived at the municipal well approximately after 360 days. Once the CO₂ plume arrived at the well, immediate reduction of pH from 8.8 to 5.2 was observed, while As and (aq) concentrations sharply increased to 2.93 × 10⁻⁷ and 1.32 × 10⁻⁷ mol/L, respectively, exceeding the maximum contaminant level (MCL = 1.33 × 10⁻⁷ mol/L, red dotted line) [79]. After 2.1 years, decreased to 0 at the municipal well, implying that all movable gaseous CO₂ was pumped out. Even after all gaseous CO₂ had been diminished due to the pumping activity, residually trapped CO₂ remained within the pores while dissolving into the groundwater. Due to the dissolution of residually trapped CO₂, the mass fraction of dissolved CO₂ () was invariant at 0.035 until 5 years; in this shallow aquifer, flow caused by both ambient fresh groundwater and pumping activity accelerated CO₂ dissolution. Therefore, complete dissolution of residually trapped CO₂ appeared at 6.3 years when became 0 and pH returned to 8.2. In addition, concentrations of As and (aq) reached 3.77 × 10⁻⁷ and 1.70 × 10−7 mol/L, respectively. Immediately beginning with Stage II (6.3–100 years), elevated pH (8.2) of groundwater inhibited the dissolution of arsenopyrite, the reaction of which requires the consumption of H⁺ (see (4)). Due to decreased arsenopyrite dissolution, concentrations of both As and (aq) were stabilized and gradually decreased due to sorption on aquifer media. The pH continuously decreased until 90 years. However, concentrations of As and (aq) reached the background level (6.60 × 10⁻⁸ and 3.01 × 10⁻⁸ mol/L, resp.) at 62 years.

(a)

(b)

Figure 9: (a) Changes in pH, total arsenic (As), and arsenite () concentrations monitored at the municipal well. The maximum contaminant level (MCL) of arsenic species is 1.33 × 10⁻⁷ mol/L. Gaseous CO₂ saturation () and pH for 10 years were also plotted at the magnified window. Two stages were characterized based on the behavior of CO₂ plume. (b) Frequency of predicted health risk accounting for both As and and cumulative density functions representing probabilistic health risk. The bin size of both histograms was designated to be 10⁻⁴, and thus the number of bins was 10 and 20, respectively. Finally, the threshold of risk or the risk level of concern was 10⁻⁴ (red dotted line).

Based on the simulated profiles of both As and (aq), , which is the maximum average concentration calculated from (7), was predicted, and the carcinogenic health risk was quantified following the method described in Section 5.3. Figure 9(b) shows histograms representing calculated frequencies of carcinogenic risk for both As (blue bar) and (aq) (green bar). The carcinogenic risk predicted from (aq) profile, which revealed relatively low concentrations, showed a mean, median, and standard deviation of 4.00 × 10⁻⁴, 4.30 × 10⁻⁴, and 1.69 × 10⁻⁴, respectively. For risk predicted from As concentration profile, the mean, median, and standard deviation were 8.94 × 10⁻⁴, 9.60 × 10⁻⁴, and 3.77 × 10⁻⁴, respectively. Cumulative density functions (blue and green lines) were also plotted together with the risk level of concern (10⁻⁴) [77, 80]. From the cumulative density functions, the risk, which exceeds the risk level of concern, can be considered to have carcinogenic potential after chronic exposure to arsenic-contaminated groundwater. As shown, the risk level predicted from both As and (aq) exceeded the risk level (red dotted line) of concern.

7.2. Sensitivity Studies

7.2.1. Effect of Sorption Intensity () (Cases 2–4)

In this sensitivity study, the intensity of the sorption effect ( of 25, 50, and 100 L/kg) on As was evaluated while gaseous CO₂ was leaked into the shallow aquifer (Table 6). As shown in Figure 4, once gaseous CO₂ was leaked from unidentifiable pathways, it migrated with the ambient groundwater. Here, characteristics of gaseous CO₂ plume such as its size, shape, migrating velocity, and gas saturation () were influenced by multiphase parameters such as capillary pressure and relative permeability (Table 1); the role of in gaseous CO₂ transport was minimal because the sorption typically accounted for the movement of dissolved species within geologic media. Due to this reason, regardless of the variation in , the distribution and behavior of gaseous CO₂ plume were not affected, and therefore the CO₂ plume remained essentially the same in all cases as shown in Figure 4. However, within the gaseous CO₂ plume, various geochemical processes occurred, including reductions in pH, dissolution of arsenopyrite, and interactions between dissolved species. In particular, variation in was anticipated to affect the behavior of dissolved arsenic species after all gaseous CO₂ was pumped out or dissolved into the groundwater.

In Stage I, the gaseous CO₂ plume arrived at the municipal well after 350 days (0.95 years) for all cases (Figure 10(a)). The elevated concentration (2.9 × 10⁻⁷ mol/L) of As occurred simultaneously in all cases because the source of arsenic species was the dissolution of arsenopyrite, which resulted from CO₂ dissolution. Despite the variation in , the arrival times of As at the municipal well were the same because the sorption did not affect the migration of the gaseous CO₂ plume. Therefore, the As profiles evolved similarly until approximately 30 years (the middle of Stage II) when its concentration reached a maximum (3.9 × 10⁻⁷ mol/L). After the arrival of the As peak, the differences between As profiles were amplified until the As concentration reached the background level (6.7 × 10⁻⁸ mol/L). In detail, differences in the slopes of As profiles were small immediately after the arrival of the As peak (30–40 years), but the discrepancy was amplified from 40 years while the slopes for As profiles sharply dropped. Differences in As profiles were attributed to the degree of , which determined the amount of arsenic adsorbed to aquifer media, especially at the rear margin of the CO₂ plume; at this location, the highest As concentration occurred due to the mixing of two chemically different types of groundwater as shown in Figure 8(b).

(a)

(b)

Figure 10: (a) As profiles with different of 25, 50, and 100 L/kg⁻¹ (Cases 2–4) and (b) probabilistic health risk plotted as cumulative density functions. The red dotted line represents the risk level of concern (10⁻⁴).

In these simulations, with increasing , more arsenic was adsorbed to the aquifer media, and the migration of As was therefore retarded. In other words, stronger retardation caused As concentration to be maintained higher and longer in the aquifer, and therefore the arrival of As concentration at the background level was delayed. For example, in Case 2 ( = 25 L/kg), As concentration revealed the earliest recovery (69.3 years) at the municipal well (Figure 10(a)). As increased to 50 and 100 L/kg, the recovery time was delayed to 72.1 and 73.9 years, respectively.

Figure 10(b) represents the predicted probabilistic health risk for Cases 2–4. While calculating the health risk, , the peak calculated from the moving average of As concentration by adopting the designated interval of ED (30 years), influenced the health risk most significantly based on (8), (9), and (10). Since As profiles revealed similar patterns with the same peak values while the only difference being the recovery time, the calculated values for Cases 2–4 were almost the same (3.9 × 10⁻⁷ mol/L). Consequently, predicted carcinogenic health risk for humans was almost identical to variation in . This result implies that variation in was a less influential parameter for assessing health risk for arsenic species. This was because variation in did not affect the characteristics of the CO₂ plume such as size, shape, and migration velocity, which determined the dissolved amount or reaction rate for arsenic sources such as arsenopyrite. Due to this reason, additional simulations were conducted and described in the following section after varying parameters (CO₂ leakage rate and aquifer permeability), which directly affected the size of the CO₂ plume and its migration rate.

7.2.2. Effect of CO₂ Leakage Rate () (Cases 5–7)

In Cases 5–7, the effect of CO₂ leakage rate () on the quantification of carcinogenic health risk was evaluated. Figures 11(a)–11(c) show the distribution of the gaseous CO₂ plume after CO₂ leakage was stopped at 1 year, and Figures 11(d)–11(f) present the distribution of As concentration after 20 years. As increased from 0.020 to 0.025 to 0.030 kg/s, the CO₂ plume approached the municipal well more closely; the calculated mass center of the CO₂ plume, plotted as black circles at a 120-day interval, moved rapidly as increased. In addition, adjacent to the CO₂ leakage point, CO₂ saturation was elevated from 0.27 to 0.30, implying occurrence of active dissolution of gaseous CO₂. As described before, although the movable gaseous CO₂ plume was pumped out after approximately 2 years, residually trapped CO₂ still remained in the pores, concurrently dissolving into the ambient groundwater until 8 years. Increases in dissolved CO₂ concentration due to dissolution of residually trapped CO₂ and resulting low-pH groundwater influenced the distribution of As concentration even after all gaseous CO₂ was pumped out (Figures 11(d)–11(f)). For example, the intensity of governed the size of the CO₂ plume where the active dissolution of arsenopyrite, the source of As, primarily occurred. Therefore, as increased from 0.020 to 0.025 to 0.030 kg/s, the size of As plume expanded at 20 years. Additionally, the effects of dispersion, diffusion, and sorption were amplified while the As plume migrated in the ambient groundwater.

Figure 11: Saturation of leaked CO₂ at different leakage rates after 1 year (365 days). (a) Case 5, (b) Case 6, and (c) Case 7. Black circles with labels indicate the locations of mass center for CO₂ plume. Distribution of As concentration induced by leaked CO₂ at different leakage intensities () after 20 years. (d) Case 5, (e) Case 6, and (f) Case 7.

The As profiles monitored at the municipal well revealed the drastic increase in As concentration immediately after the arrival of the As plume at 1.31, 1.07, and 0.98 years for Cases 5, 6, and 7, respectively (Figure 12(a)). The arrival time of the As plume coincided with that of the gaseous CO₂ plume, implying that multiphase migration of CO₂ governed the movement of dissolved As. The size of the As peaks, all of which were over the MCL, increased with (1.56 × 10⁻⁷, 2.75 × 10⁻⁷, and 4.07 × 10⁻⁷ mol/L for Cases 5, 6, and 7, resp.). In addition, its peak was maintained for a longer time with increased (13.5, 16.6, and 23.0 years for Cases 5, 6, and 7, resp.). Figure 12(b) shows the calculated probabilistic health risk. While calculating the health risk using (9) and (10), an important parameter was , which was the maximum average of arsenic concentration monitored at the municipal well (see (7)). The size of varied with as shown in Figure 12(a), which directly affected the risk prediction. Therefore, as the As concentration increased, the health risk for humans increased accordingly. In particular, the medians of Cases 5, 6, and 7, located at the half percentile in variability, were 3.39 × 10⁻⁴, 6.53 × 10⁻⁴, and 9.56 × 10⁻⁴, respectively, all of which exceeded the risk level of concern (10⁻⁴).

(a)

(b)

Figure 12: (a) As profiles with different leakage intensities () of 0.01, 0.02, and 0.03 kg/s (Cases 5–7) and (b) probabilistic health risk for Cases 5–7 plotted as cumulative density functions (CDF). The red dotted line represents the risk level of concern (10⁻⁴).

7.2.3. Effect of Horizontal Permeability () (Cases 8–10)

From two previous sensitivity studies, it was revealed that the driving force of the As plume was essentially the movement of the gaseous CO₂ plume; depending on the size of the CO₂ plume, the amount of dissolved As was determined. In this study, the magnitude of horizontal permeability (), which altered the velocity of ambient groundwater flow, varied between 0.2 × 10⁻¹³, 1.0 × 10⁻¹³, and 5.0 × 10–13 m² (Cases 8–10, resp.). As increased, the ambient groundwater rate increased from 0.39 to 1.94 to 9.71 cm/day, which induced a change in the shape of the CO₂ plume from oval to flat (Figures 13(a)–13(c)).

Figure 13: Saturation of leaked CO₂ at different horizontal permeabilities () after 1 year (365 days). (a) Case 8, (b) Case 9, and Case 10. Black circles with labels indicate the locations of mass center for CO₂ plume. Distribution of As Concentration induced by leaked CO₂ at different horizontal permeabilities () after 20 years: Cases 8, 9, and 10.

Typically, the shape and movement of the gaseous CO₂ plume within the groundwater aquifer are governed by the balance between buoyancy and viscous forces due to the density contrast of these two fluids. Even the density contrast is amplified as CO₂ leaks into the shallower aquifer, because CO₂ density radically decreases while approaching the surface. The magnitude of the buoyancy number () reflects the change in CO₂ plume shape, where is gravitational acceleration, and are the length and height of the model domain, respectively, is the ambient groundwater flow rate, and is dynamic viscosity of CO₂ (1.47 × 10−5 Pa·s). Finally, is the density difference between groundwater (1,000 kg/m³) and CO₂ (1.842 kg/m³) [81]. The calculated of Cases 8, 9, and 10 was 740.30, 148.06, and 29.61, respectively, implying that the buoyancy force acting on CO₂ plume increased as decreased.

Even if the velocity of the CO₂ plume was the lowest in Case 8 (or the largest ), the CO₂ plume swept the largest area, covering the entire thickness of the aquifer; the calculated mass center moved to the middle height of the aquifer (Figure 13(a)). Due to large coverage of the CO₂ plume, the size of the corresponding As plume was the largest in Case 8 (Figure 13(d)). In contrast, Case 10 with high accelerated the advective groundwater flow and dispersion, resulting in the flattened shape of the CO₂ plume sinking to the bottom of the aquifer (Figure 13(c)). Due to the poor sweeping efficiency of the CO₂ plume, the As plume only migrated beneath the municipal well.

The As concentrations at the municipal well for Cases 8–10 are plotted in Figure 14(a). The arrival time of As concentration was the latest (4.95 years) in Case 8 due to a low . However, due to large coverage of the As plume, the As concentration at the municipal well continuously increased to 3.77 × 10⁻⁷ mol/L until 190 years. In contrast, for Case 10, the arrival time of the As plume was the shortest (0.95 years), reaching a peak of 1.73 × 10⁻⁷ mol/L. However, due to dominance of high advective flow, the As plume was flattened below the municipal well. Therefore, the amount of As captured from the municipal well decreased soon, reaching the background level (0.67 × 10⁻⁷ mol/L) only after 25 years. These results imply that the location of the wellbore (e.g., fully or partially penetrating well, the location of screen interval) and the size of the capture zone (e.g., the pumping capacity) are important characteristics for governing As concentration at the well.

(a)

(b)

Figure 14: (a) As profiles with different horizontal permeabilities () of 0.2 10⁻¹³, 1.0 × 10⁻¹³, and 5.0 × 10–13 m² (Cases 8–10) and (b) probabilistic health risk for Cases 8–10 plotted as cumulative density functions (CDF). The red dotted line represents the risk level of concern (10⁻⁴).

The calculated health risk for each case is plotted in Figure 14(b). As expected, Case 10, which showed the smallest breakthrough of the As concentration (e.g., the smallest ), revealed the lowest risk. Both Cases 8 and 9 showed almost equivalent high-risk prediction even if the profile of the As concentration appeared differently at the municipal well (Case 8: ~80 years; Case 9: ~180 years). While accounting for the risk assessment, the exposure duration (ED) was chosen to be 30 years in this work [74]. Following (7), calculated from the As concentration profile is typically dependent on the duration of ED [31, 82, 83]. For example, can decrease when the breakthrough of As concentration is shorter than the ED. However, when the breakthrough of the As concentration is sufficiently longer than the chosen ED, does not change. Similarly, in both Cases 8 and 9, the profiles of As concentration at the municipal well were sufficiently longer than the chosen ED (Figure 14(a)). Due to this reason, for these two cases were similar to each other, and subsequently the predicted risk levels only showed a slight difference.

8. Conclusion

Even if there is no direct evidence showing that the stored CO₂ has leaked to the shallow aquifers from any major CO₂ injection demonstration sites [50, 84, 85], there exist a few natural analog sites indicating that naturally stored CO₂ has leaked through preexisting fault systems [25, 86–88]. These natural sites where CO₂ leakage is primarily driven by geothermal or tectonic activities are typically less populated with lack of concern in potable groundwater resources. However, as the number of CO₂ injection demonstration activities is growing, the secondary contamination of leaked CO₂ to the shallow potable aquifer becomes an important issue. In this study, with the presence of arsenic-bearing minerals in the aquifer, it is suggested that water-rock interactions induced by CO₂ leakage could mobilize arsenic species to the shallow potable aquifer. Therefore, without proper treatments, any residents who continuously utilize these groundwater resources could have great probability of exposure to cancer-related diseases.

Throughout this study, we developed a 2D confined aquifer model where gaseous CO₂ was leaked, and a nearby municipal well concurrently pumped out both leaked CO₂ and groundwater. Immediately after a leaked CO₂ plume arrived at the municipal well, concentrations of As species increased together, indicating occurrence of active arsenopyrite dissolution in the aquifer media. Subsequently, by analyzing As concentration from the municipal well, the carcinogenic health risk was quantified. The processes characterizing the movement of gaseous CO₂ and associated CO₂-water-rock reactions were simulated with the multiphase, multispecies reactive transport model, and subsequent carcinogenic health risks were predicted with probabilistic approach. The simulation results revealed that the movement of leaked CO₂ plume was governed by local flow fields within the shallow potable aquifer; three driving forces, that is, ambient groundwater flow, CO₂ leakage-driven flow, and pumping-driven flow, were characterized. This complex flow field governed chemical reactions, resulting in the most drastic increase (0.77%) in permeability occurring at the rear margin of the CO₂ plume where the vigorous mixing between the ambient groundwater and CO₂-dissolved fluid accelerated dissolution of the minerals. Additionally, sensitivity studies were conducted while varying the sorption intensity, leakage rate of CO₂, and horizontal permeability.

In summary, key factors that exacerbate the secondary contamination of arsenic species at the municipal well were physical characteristics of CO₂ plume such as shape, size, and migration velocity of CO₂ plume; these physical characteristics govern the area where arsenopyrite dissolves, subsequently affecting As concentration. Furthermore, the size of capture zone (e.g., screen interval, pump capacity) also governed the As concentration in the municipal well. Therefore, when the secondary contamination occurs within the shallow potable aquifer, it is suggested that the aquifer characteristics as well as the amount of leaked CO₂ and its plume size should be evaluated to develop a proper remediation protocol. Additionally, the prompt shutdown of any nearby municipal wells could minimize any potential hazards.

Disclosure

A primitive version of this research was presented at the 2016 Fall General Assembly of the American Geophysical Union.

Conflicts of Interest

The authors declare that there are no conflicts of interest regarding the publication of this paper.

Acknowledgments

The research was partially supported by multiple programs including R&D Project on Environmental Management of Geologic CO2 Storage offered by the Korea Environmental Industry & Technology Institute (Project no. 201700180004), the R&D Convergence Program of National Research Council of Science & Technology in Korea (Grant no. CAP-15-07-KICT), and, finally, the Yonsei University Future-Leading Research Initiative of 2017-22-0044.

References

D. Reichle, J. Houghton, B. Kane, and J. Ekmann, Carbon Sequestration Research and Development, Oak Ridge National Lab., Oak Ridge, Tenn, USA; National Energy Technology Lab., Pittsburgh, Pa, USA; National Energy Technology Lab., Morgantown, WV, USA, 1999.

B. Metz, O. Davidson, H. de Coninck, M. Loos, and L. Meyer, IPCC Special Report on Carbon Dioxide Capture and Storage, Intergovernmental Panel on Climate Change, 2005.

NETL, United States 2012 Carbon Utilization and Storage Atlas, National Energy Technology Laboratory, 2013.

J. T. Birkholzer, Q. Zhou, and C.-F. Tsang, “Large-scale impact of CO₂ storage in deep saline aquifers: a sensitivity study on pressure response in stratified systems,” International Journal of Greenhouse Gas Control, vol. 3, no. 2, pp. 181–194, 2009.

C. W. Strandli and S. M. Benson, “Identifying diagnostics for reservoir structure and CO₂ plume migration from multilevel pressure measurements,” Water Resources Research, vol. 49, no. 6, pp. 3462–3475, 2013.

E. Guyant, W. S. Han, K.-Y. Kim, E. Park, and S.-T. Yun, “Leakage and pressurization risk assessment of CO₂ reservoirs: a metamodeling modeling approach,” International Journal of Greenhouse Gas Control, vol. 54, part 1, pp. 345–361, 2016.

E. H. Keating, J. Fessenden, N. Kanjorski, D. J. Koning, and R. Pawar, “The impact of CO₂ on shallow groundwater chemistry: observations at a natural analog site and implications for carbon sequestration,” Environmental Earth Sciences, vol. 60, no. 3, pp. 521–536, 2010.

K. Pruess, “Leakage of CO₂ from geologic storage: role of secondary accumulation at shallow depth,” International Journal of Greenhouse Gas Control, vol. 2, no. 1, pp. 37–46, 2008.

C. M. Oldenburg, C. Doughty, C. A. Peters, and P. F. Dobson, “Simulations of long-column flow experiments related to geologic carbon sequestration: Effects of outer wall boundary condition on upward flow and formation of liquid CO₂,” Greenhouse Gases: Science and Technology, vol. 2, no. 4, pp. 279–303, 2012.

S. Bachu, “Legal and regulatory challenges in the implementation of CO₂ geological storage: an Alberta and Canadian perspective,” International Journal of Greenhouse Gas Control, vol. 2, no. 2, pp. 259–273, 2008.

Y. K. Kharaka, D. R. Cole, S. D. Hovorka, W. D. Gunter, K. G. Knauss, and B. M. Freifeld, “Gas-water-rock interactions in Frio Formation following CO₂ injection: implications for the storage of greenhouse gases in sedimentary basins,” Geology, vol. 34, no. 7, pp. 577–580, 2006.

J. A. Apps, L. Zheng, N. Spycher et al., “Transient changes in shallow groundwater chemistry during the MSU ZERT CO₂ injection experiment,” Energy Procedia, vol. 4, pp. 3231–3238, 2011.

E. H. Keating, J. A. Hakala, H. Viswanathan et al., “CO₂ leakage impacts on shallow groundwater: field-scale reactive-transport simulations informed by observations at a natural analog site,” Applied Geochemistry, vol. 30, pp. 136–147, 2013.

S. A. Carroll, E. Keating, K. Mansoor et al., “Key factors for determining groundwater impacts due to leakage from geologic carbon sequestration reservoirs,” International Journal of Greenhouse Gas Control, vol. 29, pp. 153–168, 2014.

S. Wang and P. R. Jaffe, “Dissolution of a mineral phase in potable aquifers due to CO₂ releases from deep formations; effect of dissolution kinetics,” Energy Conversion and Management, vol. 45, no. 18-19, pp. 2833–2848, 2004.

L. Zheng, N. Spycher, C. Varadharajan et al., “On the mobilization of metals by CO₂ leakage into shallow aquifers: exploring release mechanisms by modeling field and laboratory experiments,” Greenhouse Gases: Science and Technology, vol. 5, no. 4, pp. 403–418, 2015.

WHO, Chemical Agents and Related Occupations—A Review of Human Carcinogens, International Agency for Research on Cancer, Lyon, France, 2012.

J. Matschullat, “Arsenic in the geosphere—a review,” Science of the Total Environment, vol. 249, no. 1–3, pp. 297–312, 2000.

ATSDR, Support Document to the 2015 Priority List of Hazardous Substances That Will Be Candidates for Toxicological Profiles, Agency for Toxic Substances and Disease Registry, Atlanta, Ga, USA, 2015.

L. Rodríguez-Lado, G. Sun, M. Berg et al., “Groundwater arsenic contamination throughout China,” Science, vol. 341, no. 6148, pp. 866–868, 2013.

U. K. Chowdhury, B. K. Biswas, T. R. Chowdhury et al., “Groundwater arsenic contamination in Bangladesh and West Bengal, India,” Environmental Health Perspectives, vol. 108, no. 5, pp. 393–397, 2000.

A. A. Meharg and M. M. Rahman, “Arsenic contamination of Bangladesh paddy field soils: implications for rice contribution to arsenic consumption,” Environmental Science & Technology, vol. 37, no. 2, pp. 229–234, 2003.

J. M. Hebbar Annapoorna, “Arsenic contamination in groundwater of the areas surrounding Ingaldhal, Chitradurga district, Karnataka state,” International Journal of Geology, Earth & Environmental Science, vol. 6, pp. 1–7, 2016.

H. Viswanathan, Z. Dai, C. Lopano et al., “Developing a robust geochemical and reactive transport model to evaluate possible sources of arsenic at the CO₂ sequestration natural analog site in Chimayo, New Mexico,” International Journal of Greenhouse Gas Control, vol. 10, pp. 199–214, 2012.

E. H. Keating, D. L. Newell, H. Viswanathan, J. W. Carey, G. Zyvoloski, and R. Pawar, “CO₂/brine transport into shallow aquifers along fault zones,” Environmental Science & Technology, vol. 47, no. 1, pp. 290–297, 2013.

R. C. Trautz, J. D. Pugh, C. Varadharajan et al., “Effect of dissolved CO2 on a shallow groundwater system: a controlled release field experiment,” Environmental Science & Technology, vol. 47, no. 1, pp. 298–305, 2013.

T. Xu, Y. K. Kharaka, C. Doughty, B. M. Freifeld, and T. M. Daley, “Reactive transport modeling to study changes in water chemistry induced by CO₂ injection at the Frio-I Brine Pilot,” Chemical Geology, vol. 271, pp. 153–164, 2010.

Y. K. Kharaka, J. J. Thordsen, E. Kakouros et al., “Changes in the chemistry of shallow groundwater related to the 2008 injection of CO₂ at the ZERT field site, Bozeman, Montana,” Environmental Earth Sciences, vol. 60, no. 2, pp. 273–284, 2010.

L. Zhang, H. Parthasarathy, and A. Karamalidis, “Investigation on arsenopyrite dissolution and As (III) migration under geologic carbon storage conditions: A numerical simulation approach,” Greenhouse Gases: Science and Technology, vol. 7, no. 3, pp. 460–473, 2017.

T. Xiao, Z. Dai, H. Viswanathan et al., “Arsenic mobilization in shallow aquifers due to CO₂ and brine intrusion from storage reservoirs,” Scientific Reports, vol. 7, no. 1, article 2763, 2017.

E. R. Siirila, A. K. Navarre-Sitchler, R. M. Maxwell, and J. E. McCray, “A quantitative methodology to assess the risks to human health from CO₂ leakage into groundwater,” Advances in Water Resources, vol. 36, pp. 146–164, 2012.

A. K. Navarre-Sitchler, R. M. Maxwell, E. R. Siirila, G. E. Hammond, and P. C. Lichtner, “Elucidating geochemical response of shallow heterogeneous aquifers to CO₂ leakage using high-performance computing: Implications for monitoring of CO₂ sequestration,” Advances in Water Resources, vol. 53, pp. 45–55, 2013.

P. Lu and C. Zhu, “Arsenic Eh-pH diagrams at 25°C and 1 bar,” Environmental Earth Sciences, vol. 62, no. 8, pp. 1673–1683, 2011.

P. Bhattacharya, D. Chatterjee, and G. Jacks, “Occurrence of arsenic-contaminated groundwater in alluvial aquifers from delta plains, eastern India: options for safe drinking water supply,” International Journal of Water Resources Development, vol. 13, no. 1, pp. 79–92, 1997.

P. Bose and A. Sharma, “Role of iron in controlling speciation and mobilization of arsenic in subsurface environment,” Water Research, vol. 36, no. 19, pp. 4916–4926, 2002.

C. K. Jain and I. Ali, “Arsenic: occurrence, toxicity and speciation techniques,” Water Research, vol. 34, no. 17, pp. 4304–4312, 2000.

G. A. Cutter, “Kinetic controls on metalloid speciation in seawater,” Marine Chemistry, vol. 40, no. 1-2, pp. 65–80, 1992.

E. Berman, Toxic Metals and Their Analysis, John Wiley & Sons, 1980.

P. L. Smedley and D. G. Kinniburgh, “A review of the source, behaviour and distribution of arsenic in natural waters,” Applied Geochemistry, vol. 17, no. 5, pp. 517–568, 2002.

M. T. van Genuchten, “A closed-form equation for predicting the hydraulic conductivity of unsaturated soils,” Soil Science Society of America Journal, vol. 44, no. 5, pp. 892–898, 1980.

J. W. Tester, W. G. Worley, B. A. Robinson, C. O. Grigsby, and J. L. Feerer, “Correlating quartz dissolution kinetics in pure water from 25 to 625°C,” Geochimica et Cosmochimica Acta, vol. 58, no. 11, pp. 2407–2420, 1994.

J. L. Palandri and Y. K. Kharaka, “A compilation of rate parameters of water-mineral interaction kinetics for application to geochemical modeling,” DTIC Document, 2004.

L. Zheng, J. A. Apps, Y. Zhang, T. Xu, and J. T. Birkholzer, “On mobilization of lead and arsenic in groundwater in response to CO₂ leakage from deep geological storage,” Chemical Geology, vol. 268, no. 3-4, pp. 281–297, 2009.

T. Xu, E. Sonnenthal, N. Spycher, and K. Pruess, “TOUGHREACT—a simulation program for non-isothermal multiphase reactive geochemical transport in variably saturated geologic media: applications to geothermal injectivity and CO₂ geological sequestration,” Computers & Geosciences, vol. 32, pp. 145–165, 2006.

K. Pruess and N. Spycher, “ECO2N—a fluid property module for the TOUGH2 code for studies of CO₂ storage in saline aquifers,” Energy Conversion and Management, vol. 48, no. 6, pp. 1761–1767, 2007.

H. C. Helgeson, D. H. Kirkham, and G. C. Flowers, “Theoretical prediction of the thermodynamic behavior of aqueous electrolytes at high pressures and temperatures: IV. Calculation of activity coefficients, osmotic coefficients, and apparent molal and standard and relative partial molal properties to 600oC and 5 kb,” American Journal of Science, vol. 291, no. 10, pp. 1249–1516, 1981.

L. André, M. Azaroual, C. Bernstone, and A. Wittek, “Modeling the geochemical impact of an injection of CO₂ and associated reactive impurities (SO₂ and O₂) into a saline reservoir,” Transport in Porous Media, vol. 108, no. 1, pp. 185–205, 2015.

R. A. Renken, K. J. Cunningham, M. R. Zygnerski et al., “Assessing the vulnerability of a municipal well field to contamination in a karst aquifer,” Environmental and Engineering Geoscience, vol. 11, no. 4, pp. 319–331, 2005.

D. Larsen, R. W. Gentry, and D. K. Solomon, “The geochemistry and mixing of leakage in a semi-confined aquifer at a municipal well field, Memphis, Tennessee, USA,” Applied Geochemistry, vol. 18, no. 7, pp. 1043–1063, 2003.

R. J. Finley, S. M. Frailey, H. E. Leetaru, O. Senel, M. L. Couëslan, and M. Scott, “Early operational experience at a one-million tonne CCS demonstration project, Decatur, Illinois, USA,” Energy Procedia, vol. 37, pp. 6149–6155, 2013.

C. Doughty, B. M. Freifeld, and R. C. Trautz, “Site characterization for CO₂ geologic storage and vice versa: the Frio brine pilot, Texas, USA as a case study,” Environmental Geology, vol. 54, no. 8, pp. 1635–1656, 2008.

B. Dockrill and Z. K. Shipton, “Structural controls on leakage from a natural CO₂ geologic storage site: Central Utah, U.S.A,” Journal of Structural Geology, vol. 32, no. 11, pp. 1768–1782, 2010.

Z. K. Shipton, J. P. Evans, D. Kirschner, P. T. Kolesar, A. P. Williams, and J. Heath, “Analysis of CO₂ leakage through ‘low-permeability’ faults from natural reservoirs in the Colorado Plateau, east-central Utah,” Geological Society, London, Special Publications, vol. 233, pp. 43–58, 2004.

N. Kampman, M. J. Bickle, A. Maskell et al., “Drilling and sampling a natural CO₂ reservoir: implications for fluid flow and CO₂-fluid-rock reactions during CO₂ migration through the overburden,” Chemical Geology, vol. 369, pp. 51–82, 2014.

W. T. Parry, C. B. Forster, J. P. Evans, B. B. Bowen, and M. A. Chan, “Geochemistry of CO₂ sequestration in the Jurassic Navajo Sandstone, Colorado Plateau, Utah,” Environmental Geosciences, vol. 14, no. 2, pp. 91–109, 2007.

N. Kampman, M. Bickle, J. Becker, N. Assayag, and H. Chapman, “Feldspar dissolution kinetics and Gibbs free energy dependence in a CO₂-enriched groundwater system, Green River, Utah,” Earth and Planetary Science Letters, vol. 284, no. 3-4, pp. 473–488, 2009.

W. S. Han, Z. T. Watson, N. Kampman, T. Grundl, J. P. Graham, and E. H. Keating, “Periodic changes in effluent chemistry at cold-water geyser: Crystal geyser in Utah,” Journal of Hydrology, vol. 550, pp. 54–64, 2017.

C. I. A. Carruthers, Metal Mobility in Sandstones and the Potential Environmental Impacts of Offshore Geological CO₂ Storage, The University of Edinburgh, 2016.

M. Sakata, “Relationship between adsorption of arsenic(III) and boron by soil and soil properties,” Environmental Science & Technology, vol. 21, no. 11, pp. 1126–1130, 1987.

Q. H. Hu, G. X. Sun, X. B. Gao, and Y. G. Zhu, “Conversion, sorption, and transport of arsenic species in geological media,” Applied Geochemistry, vol. 27, no. 11, pp. 2197–2203, 2012.

E. Sonnenthal, A. Ito, N. Spycher et al., “Approaches to modeling coupled thermal, hydrological, and chemical processes in the drift scale heater test at Yucca Mountain,” International Journal of Rock Mechanics and Mining Sciences, vol. 42, no. 5-6, pp. 698–719, 2005.

K. G. Knauss, J. W. Johnson, and C. I. Steefel, “Evaluation of the impact of CO₂, co-contaminant gas, aqueous fluid and reservoir rock interactions on the geologic sequestration of CO₂,” Chemical Geology, vol. 217, no. 3-4, pp. 339–350, 2005.

B. Zerai, B. Z. Saylor, and G. Matisoff, “Computer simulation of CO₂ trapped through mineral precipitation in the Rose Run Sandstone, Ohio,” Applied Geochemistry, vol. 21, no. 2, pp. 223–240, 2006.

T. Xu, TOUGHREACT User's Guide: A Simulation Program for Non-isothermal Multiphase Reactive Geochemical Transport in Variably Saturated Geologic Media, V1.2.1, 2008.

C. I. Steefel and A. C. Lasaga, “A coupled model for transport of multiple chemical species and kinetic precipitation/dissolution reactions with application to reactive flow in single phase hydrothermal systems,” American Journal of Science, vol. 294, no. 5, pp. 529–592, 1994.

C. M. Bethke and P. V. Brady, “How the Kd approach undermines ground water cleanup,” Groundwater, vol. 38, no. 3, pp. 435–443, 2000.

A. J. Valocchi, “Describing the transport of ion‐exchanging contaminants using an effective Kd approach,” Water Resources Research, vol. 20, no. 4, pp. 499–503, 1984.

J. A. Davis, J. A. Coston, D. B. Kent, and C. C. Fuller, “Application of the surface complexation concept to complex mineral assemblages,” Environmental Science & Technology, vol. 32, no. 19, pp. 2820–2828, 1998.

U. S. EPA, Guidelines for Carcinogen Risk Assessment, Risk Assessment Forum, Washington, DC, USA, 2005.

B. Means, Risk-Assessment Guidance for Superfund. Vol 1. Human Health Evaluation Manual. (Part A Interim Report (Final)), Environmental Protection Agency, Office of Solid Waste and Emergency Response, Washington, DC, USA, 1989.

I. R. I. S. IRIS, “Chemical assessment summary, arsenic, inorganic,” CASRN 7440-38-2, U.S. Environmental Protection Agency, 1991.

U. S. EPA, Risk Assessment Guidance for Superfund, U.S. Environmental Protection Agency, Washington, DC, USA, 2004.

S. Fahrner, D. Schäfer, F. Dethlefsen, and A. Dahmke, “Reactive modelling of CO₂ intrusion into freshwater aquifers: current requirements, approaches and limitations to account for temperature and pressure effects,” Environmental Earth Sciences, vol. 67, no. 8, pp. 2269–2283, 2012.

U. S. EPA, Risk Assessment Guidance for Superfund, U.S. Environmental Protection Agency, Washington, DC, USA, 1989.

R. M. Maxwell, S. D. Pelmulder, A. F. B. Tompson, and W. E. Kastenberg, “On the development of a new methodology for groundwaterdriven health risk assessment,” Water Resources Research, vol. 34, no. 4, pp. 833–847, 1998.

U. S. EPA, Risk Assessment Guidance for Superfund, U.S Environmental Protection Agency, Washington, DC, USA, 2001.

K. E. McSlarrow, National Oil and Hazardous Substances Pollution Contingency Plan, United States Environmental Protection Agency, JSTOR, 1991.

Z. Duan and R. Sun, “An improved model calculating CO₂ solubility in pure water and aqueous NaCl solutions from 273 to 533 K and from 0 to 2000 bar,” Chemical Geology, vol. 193, no. 3-4, pp. 257–271, 2003.

U. S. EPA, National Primary Drinking Water Regulations, United States Environmental Protection Agency, 2008.

U. S. EPA, Risk Assessment Guidance for Superfund, Vol. 3 (Part A Process for Conduction of Probabilistic Risk Assessment), EPA, Washington, DC, USA, 2001.

R. Span and W. Wagner, “A new equation of state for carbon dioxide covering the fluid region from the triple-point temperature to 1100 K at pressures up to 800 MPa,” Journal of Physical and Chemical Reference Data, vol. 25, no. 6, pp. 1509–1596, 1996.

R. M. Maxwell, S. F. Carle, and A. F. B. Tompson, “Contamination, risk, and heterogeneity: On the effectiveness of aquifer remediation,” Environmental Geology, vol. 54, no. 8, pp. 1771–1786, 2008.

E. R. Siirila and R. M. Maxwell, “A new perspective on human health risk assessment: Development of a time dependent methodology and the effect of varying exposure durations,” Science of the Total Environment, vol. 431, pp. 221–232, 2012.

T. A. Torp and J. Gale, “Demonstrating storage of CO₂ in geological reservoirs: the Sleipner and SACS projects,” Energy, vol. 29, no. 9-10, pp. 1361–1369, 2004.

A. Mathieson, J. Midgely, I. Wright, N. Saoula, and P. Ringrose, “In Salah CO₂ storage JIP: CO₂ sequestration monitoring and verification technologies applied at Krechba, Algeria,” Energy Procedia, vol. 4, pp. 3596–3603, 2011.

W. S. Han, M. Lu, and B. J. Mcpherson, “Characteristics of CO₂-driven cold-water geyser, Crystal Geyser in Utah: Experimental observation and mechanism analyses,” Geofluids, vol. 13, no. 3, pp. 283–297, 2013.

C. D. Farrar, M. L. Sorey, W. C. Evans et al., “Forest-killing diffuse CO₂ emission at mammoth mountain as a sign of magmatic unrest,” Nature, vol. 376, no. 6542, pp. 675–678, 1995.

G. W. Kling, M. A. Clark, H. R. Compton et al., “The 1986 Lake Nyos gas disaster in Cameroon, West Africa,” Science, vol. 236, no. 4798, pp. 169–175, 1987.