To further understand the interactions of CO2-brine-rock at geological time scales, in this study, a 1D reactive transport model of CO2 intrusion into sandstone of the Longtan Formation (P2l) in the Huangqiao area, China, was constructed based on site-specific data. The simulation time is consistent with the retention time of CO2 in the Longtan sandstone Formation and is set to 20 Ma. The reactive transport model is calibrated and revised using the measured data for sandstone samples from Well X3 (i.e., the natural analogue). By comparing the simulation results with measured data for the natural analogue, the long-term geochemical reactions are investigated. The simulation results indicate that the brine-rock interactions induced by CO2 can be roughly divided into two stages. First, susceptible minerals (e.g., chlorite, ankerite, calcite, and feldspar minerals) dissolve rapidly under acidic conditions formed by the dissolution of CO2. The precipitation of siderite is facilitated by the dissolution of ankerite and chlorite. Smectite-Ca and dawsonite precipitate due to the dissolution of anorthite and albite, respectively. Dawsonite begins to convert into smectite-Na when albite is completely dissolved. As the reactions continue, intermediate products (i.e., illite, smectite-Na, and smectite-Ca) generated in the first stage become the reactants and subsequently react with CO2 and brine. These three clay minerals are not stable under acidic conditions and transform into kaolinite and paragenetic quartz in the later stage of reaction. Comparing the simulation results of the Base Case with the measured data for the natural analogue and inspired by previous studies, the scour of kaolinite is supposed to have occurred in this region and is considered in the revised model by introducing a coefficient of the scour of kaolinite (i.e., Case 2). The simulation results of Case 2 fit well with the measured data on mineral assemblage, and the trend of the sandstone porosity growth caused by the CO2-brine-rock reaction is captured by our simulation results. The combination of numerical simulation and natural analogue study indicates that the joint effects of long-term CO2-brine-rock reactions and scour of kaolinite increase the pore space of the host rock and result in an increase in quartz content in the sandstone.

1. Introduction

Capturing carbon dioxide (CO2) from the atmosphere and then injecting it into suitable geologic formations is considered an option to compensate for anthropogenic atmospheric emissions of CO2 [1]. Such formations could include deep saline aquifers, depleted gas and oil reservoirs, and unmineable coal seams. For a long time, deep saline aquifers have appeared to be the most promising storage places for large quantities of CO2 [2]. The wide global distribution of saline aquifers and their occurrence close to the potential CO2 emission sources, such as coal-fired power plants, is particularly appealing [3]. The storage capacity of deep saline aquifers for CO2 is very large and has been estimated to be as low as 320 Gt CO2 [4] and as high as 200,000 Gt CO2 [5].

Pure CO2 is slightly reactive; however, once dissolved in brine, it forms H2CO3, which is a weak acid. This change causes a decrease in the pH of the brine and makes the brine reactive to the minerals it contacts. The impact of acidic fluids on rocks is illustrated by the damage from acid rain to historic buildings [6]. Generally, the solid skeleton of the saline aquifers consists of sandstone, with burial depths between 800 and 2000 m. The salinity values vary between 100,000 (10% NaCl by weight) and 300,000 ppm (30% NaCl by weight) in most saline aquifers [7]. The interactions between saline aquifer pore fluid, CO2, and rock create significant changes in the chemical and mineralogical structures of the aquifer formation [8] and affect the mechanical behavior of reservoir rock [9]. To assess the feasibility and security of current or future CO2 storage projects, the mineral alteration and the evolution of rock physical properties induced by the injected CO2 need to be investigated. Related studies have been conducted by experimental methods [1013], numerical simulations [1420], and natural analogues [2124]. The short duration (several hours to several years) of the laboratory or field experiments, combined with the difficulty of observing mineral reactions and the sluggish nature of these reactions, results in a very limited observational database [25]. Therefore, numerical modeling is often employed for long-term CO2 geological storage assessment, which can provide valuable insights into the chemical and physical consequences of CO2 injection into the subsurface environment. However, such simulations rely on thermodynamic and kinetic databases that have uncertainty [26]. Consequently, the numerical simulation results of long-term CO2-brine-rock interaction should be calibrated or validated by the natural accumulations of CO2 (i.e., natural analogues) which have retained CO2 for millions of years [27, 28].

Numerous natural CO2 reservoirs have been found in the world, and some of them have been used as analogues to provide natural evidence for evaluating the impacts of CO2 geological storage over long time scales [2934]. Ren et al. [35] reported that the Huangqiao CO2 field in eastern China is unique with the coexistence of CO2 and hydrocarbon accumulations and that the CO2-rich fluid induced substantial reformation of the Longtan sandstone reservoir. The core test results show that carbonate and feldspar minerals are completely dissolved in the sandstone and that the main component of the rock is quartz with small amounts of clay minerals (detailed in Section 3). These authors conducted a set of CO2-brine-sandstone interaction experiments from 60°C to 120°C and from 10 MPa to 40 MPa, which correspond to burial depths from 1500 m to 4000 m, to understand the corrosion of CO2-rich fluid on the sandstone reservoir. Experimental results indicated that the corroded mineral is mainly carbonate and that the corrosion reached equilibrium after 72 h. The corrosion of feldspar was not detected by Ren et al. [35] because aluminosilicate mineral alteration is very slow under ambient deep-aquifer conditions and is not amenable to an experimental study lasting a few days [36]. To date, how the sandstone in the Huangqiao CO2 reservoir became quartzose and what were the major reactions among the CO2, brine, and sandstone are still unknown.

To further understand the long-term interaction among the dissolved CO2, brine, and sandstone and the impact of the reactions on the physical and mineralogical characteristics of the strata, in this paper, the natural CO2 reservoir in the Huangqiao area, Southeast China, is selected as the natural analogue for long-term CO2-brine-rock interactions. A one-dimensional geochemical reactive transport model is constructed to reproduce the brine-rock reactions induced by CO2 over long time scales. The simulation results are calibrated by the observations in the Huangqiao natural CO2 reservoir.

2. Geological Background

The research area, the Huangqiao CO2 reservoir, is located in Huangqiao Township, Jiangsu Province, Southeast China, and has been commercially exploited for more than 35 years since 1983. The Huangqiao region is in the southern North Jiangsu Basin. Its tectonic structure is the Huangqiao complex syncline in the northeastern Nanjing Depression (Figure 1(a)).

The Huangqiao region is widely covered by Quaternary deposits. All other strata are basically well preserved, except for Paleogene, Lower Cretaceous, Jurassic, and Upper Triassic strata (detailed stratigraphic information is shown in Figure 2). As presented in Figure 2, eight CO2 reservoirs have been detected in this area with a proven CO2 reserve of [38]. The origin of CO2 has been inferred to be the mantle or magma in previous studies [38, 39], and the deep faults are the main channels for CO2 migration to shallow strata. Cui [40] suggested that the CO2 confined in the deeper strata of P1q migrated upward along the fault (as shown in Figure 1(b)) to the Permian Longtan Formation (P2l), the target strata of our study, in the Neogene period (approximately 20 Ma) because of tectonic movement. The overlying strata (P2d) consisting of ~50 m thick mudstone with siltstone layers inhibit the upward migration of CO2 due to their low permeability and high capillary entry pressure [41]. Consequently, a large amount of CO2 has accumulated in the Longtan Formation to form the CO2 reservoir there. The Longtan Formation is an ~140 m thick unit of interbedded black mudstone and grey-black sandstone, with thin coal interlayers. The sandstone is the main CO2-producing section in Well X3 [38]. According to the local burial history (Figure 3), P2l has undergone slight subsidence since the Neogene Period, and the burial depth of P2l in this period is approximately 1600 m with a temperature of approximately 70°C.

3. Core Mineralogy Test Data Analysis

Several wells were drilled in this region for the exploration or exploitation of oil and CO2 (the distribution of the drillings can be found in Figure 1b in [38]), most of which were drilled through or into the P2l sandstone formation and P2l was sampled. The mineral compositions of the core samples were obtained through X-ray diffraction (XRD) and thin section analysis, which were carried out by the East China branch of SINOPEC. The testing procedures and methods are not the focus of our study and are not discussed here. To meet the format of the simulation code TOUGHREACT (detailed in the next section) for defining mineral abundance, the mass fractions of minerals directly measured by XRD analysis were converted to volume fractions.

A total of 18 core samples of P2l sandstone formation from three wells (i.e., Wells X1, X2, and X3) were tested, and the sampling depth and mineral composition of each sample are listed in Table 1. The lithology of core samples from Well X1 is fine- to medium-grained feldspathic quartzose sandstone, with quartz contents (in volume fraction) of 41.4% to 66%. Illite is the main clay mineral, with the content ranging from 14.7% to 22.2%. Feldspar minerals are dominated by plagioclase with contents of 3.1% to 16%. K-feldspar and carbonate minerals (e.g., calcite, ankerite, and siderite) are found in a portion of the samples from Well X1, with contents ranging from a few percent to more than a dozen percent. Similar to the samples from Well X1, the mineral composition of samples taken from Well X2 is dominated by quartz (44.4%-69.3%), illite (7.3%-34%), and plagioclase (6.7%-19.1%). Kaolinite, K-feldspar, and ankerite are not detected in the samples from Well X2. Calcite and siderite are detected in most of the samples from Well X2, with abundances of several percent. Four samples from Well X3 were tested, and the testing results indicate that the lithology of the formation is medium-grained quartzose sandstone with extremely high contents of quartz (85.7%-89.1%). The contents of clay minerals of samples from Well X3, with volume fractions of 4% to 8%, are relatively low compared to those of samples taken from Wells X1 and X2. In addition, the abundance of kaolinite in samples from Well X3 is significantly higher than those in samples taken from Wells X1 and X2, and chlorite that is commonly found in Wells X1 and X2 is not detected in Well X3. Plagioclase, calcite, and siderite, which are common in samples from Wells X1 and X2, are not found in those from Well X3. K-feldspar and ankerite are not detected in samples from Well X3.

As shown in Figure 1(b), Well X3 is located northwest of Well X2 at a distance of ~1 km. Well X1 is located approximately 4 km southeast of Well X2 (as presented in Figure 1b in [38]). In other words, the three wells are not far apart, so why are the mineral compositions of samples from Well X3 so different from those of samples from Wells X1 and X2? It should be noted that CO2 is yielded only in Well X3 in the Longtan sandstone Formation and no CO2 is found in the same strata in Wells X1 and X2. Ren et al. [35] proposed that the existence of CO2-induced brine-rock interactions in Well X3 is the main reason for the differences between the sandstone reservoirs of Well X3 and Well X2. However, due to the limited duration (120 h) of their experiments, only carbonate minerals were found to be eroded, which could not explain why the sandstone samples from Well X3 do not contain feldspar minerals and why only quartz and a small amount of clay minerals are detected. Therefore, we intend to use a combination of numerical simulation and natural analogy to explore whether the brine-rock interactions induced by CO2 resulted in such large differences in the core sample mineral compositions between Well X3 and Wells X1 and X2. If so, what are the main geochemical reactions, and how long do they take?

4. Numerical Modeling Approach

All present simulations were carried out using the nonisothermal, reactive transport modeling code TOUGHREACT [43], which has been widely used for CO2 geological sequestration [4446]. This code was developed by introducing reactive transport into the existing framework of a nonisothermal multiphase multicomponent fluid and heat flow simulator TOUGH2 [47]. Complicated thermophysical-chemical processes are considered under various thermohydrological and geochemical conditions of pressure, temperature, ionic strength, and water saturation. TOUGHREACT can be applied to one-, two-, or three- dimensional porous and fractured media with physical and chemical heterogeneity. The code can accommodate any number of chemical species present in liquid, gas, and solid phases. Further descriptions about the selection of the kinetic and thermodynamic parameters are provided in Appendix A.

5. Model Setup

5.1. Conceptual Model and Boundary Conditions

As shown in Figure 1(b), due to the sealing effect of the overlying P2d mudstone strata, CO2 transported upward along the fault migrated laterally and accumulated in the Longtan sandstone Formation. Well X3 is located in this CO2 reservoir, and the sandstone samples taken from the well are presumed to have been eroded by the CO2-rich flow. To simulate the effect of long-term CO2-brine-sandstone interaction on the physical properties (e.g., porosity) and mineral composition of the sandstone reservoirs, much specific and detailed information is required, e.g., the thickness and lateral extension of the sandstone strata, the physical properties and mineral composition of the sandstone, and its spatial heterogeneity. Before moving into the all-encompassing investigation, general features and issues related to these processes need to be explored. This exploration can be done by investigating a one-dimensional CO2 source-sandstone system that abstracts the site-specific features and is able to represent characteristics that are common to many such systems.

As presented in Figure 4, a one-dimensional horizontal column is used to represent the CO2 source-sandstone system in the model. Considering the extremely long modeling period (20 Ma), only three grids with different roles are set in the model to reduce the CPU load. The grid on the left is set as the source of the CO2, and an overpressure of 0.5 MPa is assigned to it to ensure the CO2 flow into the sandstone formation that is represented by the grid in the middle. The size of the middle grid is set to to represent the sandstone formation that was invaded by CO2 as revealed by Well X3. The pressure of the left grid is assumed to be constant throughout the simulation (i.e., it represents a Dirichlet boundary), which implies that the CO2 consumed by the reactions will always be replenished. The geochemical reactions stop once the native water is drained by CO2. To ensure that the CO2-brine-rock interactions continue in the simulation, a water- (brine) saturated grid with infinite volume is set in the right to represent the lateral semi-infinite aquifer.

5.2. Initial Conditions

The geochemical simulation should begin as the CO2 invades the Longtan sandstone Formation, and the initial conditions in our model should correspond to the conditions at this moment. As mentioned above, the CO2 invaded the target sandstone strata in the Neogene Period (approximately 20 Ma) along the fracture. Therefore, it is a major challenge to obtain the initial physical parameters and mineralogical composition of the host rock and the initial native water components, because all the measured physical and mineral data are the present-day observations.

As mentioned in Section 3, CO2 is observed only in Well X3 in the target strata, and the sandstone samples taken from the same formation from adjacent Wells X1 and X2 are presumed not to have been eroded by CO2. Moreover, the burial depth of the P2l formation after the Neogene Period changed little, and its temperature and pressure conditions were relatively stable (see Figure 3). Therefore, as previously assumed by Watson et al. [33], the mineralogy and physical properties of the sandstone formation outside the CO2 reservoir (i.e., Wells X1 and X2) are taken as unmodified and can be approximated as the initial conditions for the sandstone formation in Well X3.

5.2.1. Initial Physical Parameters

Most hydrogeological and physical parameters (i.e., porosity, permeability, and density) of the sandstone in our model are set to the average values from test results of the samples from Wells X1 and X2 (listed in Table 2). The burial depth of P2l is approximately 1600 m after the Neogene Period, and the ambient temperature is approximately 70°C (see Figure 3). Therefore, the initial temperature is set to 70°C and is assumed to be constant over time. In addition, the initial pore water pressure of the strata is set to 16 MPa by hydrostatic equilibrium. As stated before, an overpressure of 0.5 MPa is set in the left grid, and the pressures of left and right grids are set to be constant throughout the simulation. The relative permeability and capillary pressure parameters are set empirically to represent the properties of sandstone in the Huangqiao area and are listed in Table 2.

5.2.2. Initial Conditions for Geochemistry

The initial conditions for geochemistry include the initial mineral composition and water chemistry. The average mineral composition of samples taken from Wells X1 and X2 is set as the initial mineral assemblage for sandstone in our model. Plagioclase is treated as a 1 : 1 mixture of endmembers albite and anorthite in the model to simplify the mineralogy [48, 49]. Similarly, smectite is treated as a 1 : 1 mixture of smectite-Na and smectite-Ca (detailed in Table 3). The measured native formation water chemical composition, taken from Ren et al. [35], is used during the initialization of the simulations. The resulting water composition is used as the initial condition for the subsequent modeling (Table 4).

6. Results and Discussion

6.1. Mineral Alternation

We started the first simulation using the parameters given in Tables 24 (called the Base Case simulation). The flooding of CO2 decreases the pH of the solution and makes it acidic which induces the dissolution of some minerals. To understand the long-term mineral alteration pattern in sandstone, the variations in the mineral abundances (expressed in the volume fraction of mineral in the medium) are presented in Figure 5. As shown in the figure, feldspar minerals including anorthite and albite (Figure 5(a)) and K-feldspar (Figure 5(b)) dissolve significantly and disappear after approximately 8000 years, 22,000 years, and 100,000 years, respectively. For clay minerals, chlorite (Figure 5(c)) is the most reactive under acidic conditions formed by the dissolution of CO2 and dissolves completely after 8000 years. Illite (Figure 5(d)) precipitates slightly between 300 and 10,000 years. Then, it begins to dissolve and disappears after approximately 1.2 Ma. K+ to form illite is provided by the dissolution of K-feldspar. The carbonate mineral siderite (Figure 5(e)) precipitates slightly at the expense of chlorite and then begins to dissolve. The dissolution of ankerite (Figure 5(e)) also provides Fe2+ for the precipitation of siderite. Calcite (Figure 5(e)) dissolves rapidly under high CO2 pressure conditions and disappears after approximately 70 years. As a secondary mineral (does not exist initially), dawsonite (Figure 5(e)) is generated after approximately 30 years and precipitates significantly with a maximum abundance of 4.6% after 80,000 years. The formation of dawsonite is facilitated by the dissolution of albite. After 80,000 years, when albite has dissolved completely, the dawsonite begins to dissolve and disappears after approximately 170,000 years. The dissolution of dawsonite provides Na+ for the precipitation of smectite-Na (Figure 5(f)). In addition, Ca2+ needed for the formation of smectite-Ca (Figure 5(f)) is derived from the dissolution of anorthite and calcite. The simulation results indicate that illite, smectite-Na, and smectite-Ca are unstable under acidic conditions and eventually convert to kaolinite and paragenetic quartz. Consequently, kaolinite (Figure 5(g)) and quartz (Figure 5(h)) precipitate significantly throughout the simulation period and become the dominant minerals at the end of the simulation.

The CO2-brine-rock geochemical reaction processes can be roughly divided into two stages. Major reactants in the first stage (probably the first 10,000 to 100,000 years of the simulation) include calcite, ankerite, chlorite, and feldspar minerals. Carbonate minerals siderite and dawsonite and clay minerals smectite-Na, smectite-Ca, and illite are products of reactions in the first stage and become reactants to react with CO2 and water in the second stage when the first reactants (i.e., chlorite and feldspars) have completely dissolved. The precipitation of carbonate minerals could immobilize CO2 (i.e., mineral trapping). Our simulation results indicate that carbonate minerals may dissolve as the reaction continues. However, on the time scale that may be concerned in the CO2 geologic storage project, i.e., 5,000 to 10,000 years, the carbonate minerals (e.g., siderite and dawsonite) precipitate significantly. It should be noted that illite, which is usually supposed to be stable under neutral or alkaline conditions [50, 51], becomes the main reactant in the second stage of the reaction. Our simulation results indicate that illite could convert to kaolinite in the acidic environment formed by CO2. The main products in the second stage are quartz and kaolinite, which are stable under acidic conditions. Finally, the mineralogy of sandstone is dominated by quartz and kaolinite, with volume fractions of 69.1% and 23.9%, respectively. The main CO2-related reactions mentioned above are summarized as follows:

6.2. Porosity Changes

Potential mineral dissolution and precipitation may affect the porosity of the sandstone. Variation in sandstone porosity is simulated in our model and illustrated in Figure 6. The variation in porosity is calculated using changes in the volume fraction of minerals as follows [43]: where is the porosity, superscript is the number of the minerals, is the volume fraction of mineral in the rock, and is the volume fraction of the inactive portion in the rock.

As shown in Figure 6, three peaks appear on the porosity curve over time. The first peak appears after 700 years, which is caused by the significant dissolution of calcite, ankerite, albite, and anorthite. After this time, when calcite is completely dissolved, the porosity of sandstone decreases to a certain extent due to the precipitation of smectite-Ca, dawsonite, and siderite. After approximately 3000 years, the variation in sandstone porosity is dominated by the dissolution of siderite, illite, and feldspar minerals, which results in an increase in porosity. The second wave crest of the porosity curve arises after approximately 50,000 years. Then, the porosity of rock decreases between 50,000 and 100,000 years, dominated by the precipitation of smectite-Na and dawsonite. After 100,000 years, the dissolution of dawsonite and smectite-Na results in an increase in sandstone porosity and the third peak of porosity curve appears after 200,000 years with a value of 13.2%. After this point, the geochemical reactions in the system are dominated by the conversion of smectite-Na and illite to kaolinite and paragenetic quartz. The joint effects of these two reactions reduce the pore space of the host rock, and the porosity of the rock decreases to 3.12% after 1 Ma. The simulation results indicate that the geochemical reaction reaches dynamic equilibrium after 1 Ma and that the mineral contents and rock porosity no longer change.

6.3. Comparison with the Natural CO2 Reservoir

A comparison between modeling results and field observations of the natural CO2 reservoir is instructive, because it helps identify issues and limitations in the models and places where the models need to be improved.

The observed mineral compositions of the natural CO2 sandstone reservoir in the Huangqiao area are listed in Table 1 (i.e., core samples from Well X3) and are characterized by the extremely high abundance of quartz (>85%), small amounts of clay minerals (4-8.1%), and the absence of chlorite, carbonate, and feldspar minerals. The simulated mineral composition in the late stage of reaction in the Base Case shows good agreement with field data, including the absence of chlorite, carbonate, and feldspar minerals; significant dissolution of illite; and overgrowth of quartz. However, the simulated abundance of kaolinite is extremely high (>20%), which is inconsistent with the measured data. In addition, the simulated porosity of sandstone, i.e., 3.12% at the end of the simulation, is much lower than the average value of measured data (15%). The contradictions between modeling results and the measured data in the mineral abundances (mainly in kaolinite content) and rock porosity indicate that there are still some shortcomings in our models. The best knowledge and efforts should be used to eliminate these shortcomings.

Higgs et al. [24] analyzed the sandstone cores from a high CO2 concentration well (Garvoc-1) in the Otway Basin, Australia, and found that kaolinite is the dominant authigenic clay mineral with the bulk fraction distributed between 1.8% and 22.2% (the largest value of kaolinite bulk fraction is very close to our simulation result). They also found that the authigenic kaolinite occurred as very fine, grain-coating clays that may tend to be washed away by fluid flow. The loose sand grains (quartz) that are very common in the Garvoc-1 cuttings suggest that the fine kaolinite has been washed away. This process can explain why the abundance of kaolinite in sandstone cores from the same well is unevenly distributed, with the lowest value of 1.8% and the largest value of 22.2%. Cao et al. [52] found that the authigenic kaolinite in the sand reservoir migrates away instead of precipitating in situ in a strong flow system. The sandstone cores in Well X3 in the Huangqiao area contain loose quartz grains and have developed a large number of granules and intergranular pores [35]. Moreover, kaolinite was identified as the main fine mineral associated with fine migration-related formation damage in previous studies [5356]. Considering the above evidence, we propose that a portion of authigenic kaolinite has been washed away in the Longtan sandstone Formation under the strong flow system caused by the intrusion of CO2, which can explain why the simulated abundance of kaolinite in the Base Case is much higher than the observations from cores. As the scour and formation of kaolinite occur simultaneously, the process of kaolinite washing away must be coupled in reactive transport modeling.

The invoked simulation program TOUGHREACT (ECO2N module) does not consider minerals as fluid components. To describe the erosion process of kaolinite, the new function must be added to the existing program manually. For simplicity, we introduced a coefficient of scour (Sco) for a mineral to describe the erosion process of the mineral. The resulting equation is as follows: where and are the volume fractions of the -th mineral at times and , respectively; is the coefficient of scour of the -th mineral and represents the proportion of mineral that is washed away each year (in this study, it is set to for kaolinite through trial and error); and is the length of the current time step in years.

The scour of minerals affects the variation in porosity. These variations can be evaluated using equation (2).

The updated volume fraction of mineral obtained in equation (3) is used in equation (2) for the calculation of porosity influenced by precipitation/dissolution and possible scour of minerals. Only the scour of kaolinite is assumed to have taken place in the Huangqiao area. To examine its effects on the mineral composition and porosity changes in the Longtan sandstone Formation, an additional simulation (i.e., Case 2) was performed by taking into account the scour of kaolinite in the Base Case.

Variations in the mineral composition and rock porosity of these 2 sets of simulations (i.e., Base Case and Case 2) are shown in Figures 7 and 8, respectively. It should be noted that to maintain consistency with the measured data in Table 1; the mineral abundances presented in Figure 8 are expressed in the volume fraction of rock grain (the ratio between the volume fraction of medium and volume fraction of rock grain is (1-porosity)). Considering that fluid flow in the deep subsurface is generally very slow, the erosion rate of kaolinite would be very small, and it would take thousands of years or more to affect the mineral composition and rock porosity. Here, the coefficient of scour of kaolinite is set equal to heuristically after a number of trial calculations. The scour of kaolinite would decrease the volume of kaolinite, thereby affecting the rock porosity and abundance of minerals.

As shown in Figure 7(c), it takes nearly 100,000 years to show a significant difference between the two cases with and without consideration of the scour of kaolinite when chlorite, feldspar, and carbonate minerals are completely dissolved. The simulation results show that the scour of kaolinite does not affect the abundance of illite (Figure 7(a)) and smectite-Na (Figure 7(b)). Affected by erosion, the abundance of kaolinite (Figure 7(c)) decreases significantly and decreases to zero after approximately 1 Ma. Consequently, the abundance of quartz (Figure 7(d)) increases prominently, and quartz becomes the dominant mineral of the rock.

The pore space becomes larger as the kaolinite is washed away by fluid flow. As shown in Figure 8, the simulated rock porosity in Case 2 increases significantly due to the scour of kaolinite between 50,000 and 600,000 years, with a maximum porosity of 29.3%. After 600,000 years, the joint effects of mineral dissolution/precipitation and scour of kaolinite decrease the pores of the host rock to some extent, and the porosity of the sandstone is approximately 27% after 1 Ma.

As mentioned above, CO2 is proposed to have invaded the Longtan sandstone Formation at approximately 20 Ma. However, the exact duration of the CO2-brine-rock interaction is difficult to determine because it stops once the water has been drained by CO2. By comparing the simulation results with the measured data for samples from Well X3 (i.e., the natural analogue), we find that the mineral assemblage in Case 2 after approximately 500,000 years fits well with the measured data (as shown in Figure 9). The simulated mineral composition at this time point consists of 87.5% quartz, 3.3% kaolinite, 3.7% illite, and 0.02% smectite. The simulated abundance of quartz lies within the observed range (85.7~89.1%) of sandstone samples from Well X3, and that of kaolinite is slightly higher than the measured data (1.7~2.8%). The simulated volume fraction of illite fits well with the measured data while that of smectite is slightly smaller than the measured data. In addition, the simulated porosity in Case 2 after 500,000 years is 29%, while the measured average porosity of sandstone samples from Well X3 is approximately 15%. The simulated porosity in the Base Case is approximately 9% after 500,000 years and decreases as the reaction continues, which is inconsistent with the field observations. Therefore, the simulated changes in porosity give us more confidence that the scour of kaolinite has occurred in this area and affected the mineral composition and porosity variation of the rock, although the simulated rock porosity and abundance of kaolinite in the revised model (i.e., Case 2) are higher than the measured data. The simulated main evolution trends of minerals are consistent with field observations, including the complete dissolution of chlorite, carbonate, and feldspar minerals and the precipitation of kaolinite and paragenetic quartz. Given the complexity of the long-term geochemical reactions and the variety of factors that may impact, we deem that these deviations between simulation results and measured data are tolerable and that the revised model (Case 2) can reproduce the process of long-term CO2-water-rock interaction in the Longtan sandstone Formation in the Huangqiao area.

7. Summary and Conclusions

A one-dimensional reactive transport model of CO2 intrusion into the Longtan sandstone Formation was constructed, and the alteration patterns of minerals and changes in porosity for the sandstone induced by the intrusion of CO2 over a long time scale (20 Ma) were investigated. The simulation results were calibrated using the measured mineralogical and physical data for sandstone cores from Well X3, and the following conclusions can be drawn.

The geochemical processes induced by the intrusion of CO2 can be roughly divided into two stages. First, the flooding of CO2 decreases the pH of the solution, resulting in the dissolution of some susceptible minerals (e.g., chlorite, calcite, ankerite, and feldspars). The dissolution of chlorite and ankerite supplies Fe2+ for the precipitation of siderite. Dawsonite and smectite-Ca form with the dissolution of albite and anorthite, respectively. Dawsonite begins to transform into smectite-Na after albite is completely dissolved. K-feldspar converts into illite with Mg2+ provided by the dissolution of chlorite. Second, as the reactions continue, intermediate products (i.e., illite, smectite-Na, and smectite-Ca) generated in the first stage become the reactants and subsequently react with CO2 and brine. Kaolinite and quartz are the final products of geochemical reactions.

Compared with the field data, the abundance of kaolinite simulated in the Base Case is much higher. By taking into account the scour of kaolinite in the simulation (i.e., Case 2), the simulated volume fractions of quartz and illite after approximately 500,000 years fit well with the measured data for samples taken from Well X3 in the Huangqiao CO2 gas field (i.e., the natural analogue). We also find that the joint effects of CO2-water-rock interactions and scour of kaolinite increase the porosity of the host rock, which is consistent with field observations. Therefore, the scour of kaolinite is likely to have occurred in this area and significantly affected the mineral composition and variation in porosity of the rock.

Further investigation is still required, especially using 2D or 3D models, to fully elucidate the geochemical and physical processes associated with the intrusion of CO2 into the sandstone formation. The range of problems concerning CO2-brine-rock interactions is very broad, and the present simulation results are specific to the parameters and conditions considered. Given the need and the general inadequacy of other reliable methods to predict the long-term interactions of CO2-brine-rock and reservoir responses, the combination of numerical modeling and natural analogue study would greatly facilitate and promote our understanding.


A.Kinetic and Thermodynamic Parameters

Both equilibrium and kinetic mineral dissolution/precipitation can be simulated in TOUGHREACT. Considering the sluggish nature of the CO2-brine-rock reactions, kinetic formulations are used for all the mineral-water reactions in the model. For kinetically controlled mineral dissolution and precipitation, a general form of the rate law [57, 58] is used. where the subscript denotes the kinetic mineral index, is the dissolution/precipitation rate of the mineral (positive values indicate dissolution, and negative values indicate precipitation), is the rate constant that is temperature dependent, is the specific reactive surface area per kg H2O, and is the kinetic mineral saturation ratio. The superscripts and must be determined from experiments; usually, but not always, they are taken as equal to one.

The dissolution and precipitation of minerals are often catalyzed by H+ and OH-. Therefore, for many minerals, the kinetic rate constant includes three mechanisms (i.e., acid, neutral, and base mechanisms) [59, 60] or where superscripts or subscripts , , and indicate neutral, acid, and base mechanisms, respectively; is the rate constant at 25°C; is the activity of the species; is the power term (constant); is the gas constant; and is the absolute temperature. The parameters for the kinetic rate calculation of specific minerals are similar to those used by Xu et al. [19], which were derived from Palandri and Kharaka [59] (see Table 5).

The thermodynamic data for aqueous species, gases, and minerals were mostly taken from the built-in the database of TOUGHREACT, i.e., the EQ3/6 V7.2b database [61], as previously done by Xu et al. [28, 36] and Tian et al. [45].

Data Availability

The underlying data can be provided if necessary by the corresponding author upon request.

Conflicts of Interest

The authors declares that they have no conflicts of interest.


This work was jointly supported by the National Program on Key Research and Development Project (No. 2016YFB0600804), by the National Natural Science Foundation of China (Grant No. 3A417C371425), by the National Major Special Project (No. 2016ZX05016-005), by the 111 Project (No. B16020), and by the Research Grant Scheme for Doctoral Candidates of Jilin University (Grant 10183201841).